J. Phys. Chem. B 1997, 101, 3209-3220
3209
Ion Hydration: Thermodynamic and Structural Analysis with an Integral Equation Theory of Liquids Song-Ho Chong Department of Chemistry, Faculty of Science, Kyoto UniVersity, Kyoto 606, Japan
Fumio Hirata* Institute for Molecular Science, Myodaiji, Okazaki, Aichi 444, Japan ReceiVed: January 13, 1995X
We present a theoretical study for ion hydration based on an integral equation method referred to as the extended reference interaction site method (ex-RISM). We analyze the thermodynamic functions of solvation, especially the partial molar volumes of individual ions at infinite dilution. Special attention is paid to information contained in the partial molar volumes and to the question of whether the partial molar volumes of individual ions reflect the true nature of ion-water interactions. Our results suggest, contrary to the previous work given by Kusalik and Patey (J. Chem. Phys. 1988, 89, 5843), that the partial molar volumes do reflect the nature of ion-water interactions. Concerning the microscopic description of the ion hydration, we revisit the earlier model proposed by Samoilov by defining the activation energy ∆Ei in his model in terms of the ion-water potential of mean force. The theoretical results are in good accord with the earlier model in terms of the classification of ions into the “positive” and “negative” hydrations. We also discuss the structural changes of water due to the presence of an ion utilizing the density derivatives of the solvent distribution functions. Qualitative differences of the density derivatives between the “positively” and “negatively” hydrated ions were observed and found to be consistent with the analysis of the potential of mean force.
1. Introduction The ion hydration is one of the most fundamental physicochemical processes, which bears close connection with many research fields in chemistry including chemical reactions and stability of biomolecules.1,2 It is not surprising that so much effort has been devoted to characterize the properties of the ion hydration experimentally as well as theoretically.3-6 Two naive models were proposed in the earlier stages of the study of ion hydration. These are in sharp contrast with respect to the physical description depending on how they treat ionwater interactions: the dielectric continuum model7 and the “solventberg” model.3 The dielectric continuum model sees a solvent as a structureless continuum and characterizes it with just a single parameter, the macroscopic dielectric constant. It is rather surprising that such a naive model is still prevalent in use in the literature concerning solvation free energy, such as in the Marcus theory for the electron-transfer reaction8 and protein solvation.2 The solventberg model sees a hydrated ion as an ion-water complex due to the chemical bonds and characterizes the property of the complex primarily in terms of the hydration number. A weakness that is common to these earlier models is the disregard of the distinctive feature of the water structure that can be casually summarized as the tetrahedrally coordinated structure resulting from hydrogen bonding.9 Many experiments suggest that this structure in water plays an essential role for the ionic processes in aqueous solution. In 1957 two closely related models for the ion hydration were proposed in the Discussion Faraday Society by Frank and Wen10 and by Samoilov,11 which have provided conceptual guides to analyze microscopically the experimental observation related to ionic processes in aqueous solution. According to the model proposed by Frank and Wen,10 a small ion is surrounded by * Author to whom correspondence should be addressed. X Abstract published in AdVance ACS Abstracts, March 15, 1997.
S1089-5647(96)00878-4 CCC: $14.00
three concentric regions: the innermost “region A” is where water molecules are immobilized due to the electric field produced by the ion, the second “region B” is that in which water is less “icelike”, i.e., more random in organization than “normal” water, and the third “region C” contains normal water. They presumed that the decrease in “icelikeness” or so-called “structure-breaking” effect is caused by the approximate balance in the region B between two competing orientational effects that act on any given water molecule: one of those is the normal structure-orienting influence of neighboring water molecules, while the other is the orienting influence due to the ionic field. This presumption was coincident with that proposed by Gurney.4 Samoilov11 proposed a dynamical characterization of the ion hydration on the basis of an examination of the action of ions on the translational motion of water molecules in the immediate vicinity of the ions. The translational motion of water molecules around an ion can be regarded as the exchange of water molecules that are closest to the ion. If the exchange occurs infrequently relative to that of water molecules in the bulk, the hydration of the ion is considered strong. On the contrary, if the exchange of water molecules around the ion occurs more frequently than that in the bulk, the hydration of the ion can be considered weak. A quantity that characterizes the ionic hydrations in the Samoilov model is the change, due to the influence of the ion, of the activation energy for the exchange process of the closest water molecules, ∆Ei. Samoilov has shown experimentally, from the self-diffusion in water and the temperature coefficients of the ionic mobilities in solutions, that ∆Ei > 0 for Li+ and Na+ and that ∆Ei < 0 for K+, Cs+, Cl-, etc. In the latter case, the water molecules near the ions become more mobile than in pure water. This phenomenon was called “negative hydration”, while the former case, ∆Ei > 0, was called “positive hydration”. The microscopic processes of the ion hydration described in terms of the “negative” and “positive” hydrations by Samoilov and in terms of “region A” and “region © 1997 American Chemical Society
3210 J. Phys. Chem. B, Vol. 101, No. 16, 1997
Chong and Hirata
B” by Frank and Wen are considered the same processes expressed in different languages. Appearance of the negative hydration is intimately related to the existence of tetrahedrally coordinated structure in water. The models proposed by Frank-Wen and Samoilov have been highly regarded by solution chemists, and they still work as key concepts in interpreting experiments. It is our expectation that these models will play even more important roles when the limitation of the earlier continuum and solventberg models is recognized by people in those fields related to solvation. From such a viewpoint, we revisit this rather old problem in physical chemistry with the aid of the extended reference interaction site method (ex-RISM),12-14 the statistical mechanical theory for molecular liquids and solutions. Our primary goal here is to provide microscopic characterization of the structural modification of water in the vicinity of an ion in terms of the site-site pair correlation function of water. The study inevitably touches the problem of the three-body correlation involving the ionwater-water triplet, of which solution has not been solved by means of the integral equation methods except for some limited cases. Thereby, we employ an alternative method developed by Yu and Karplus,15 which characterizes the modification of water structure near ion in terms of the ion-density derivatives of the pair correlation functions for water. In the following section, we examine a conceptual relation between the threebody correlation detected directly by the molecular simulation16 and that obtained from the integral equation method. In section 2 of this paper, we first briefly review the main features of ex-RISM and its relation with the solvation thermodynamics. An equation for the partial molar volumes of individual ions at infinite dilution is derived from the Kirkwood-Buff and ex-RISM theories, the detail of which is given in the Appendix. We then make some conceptual clarification concerning the theoretical description of the structural modification of water in the vicinity of an ion. In section 3, we apply the theory to a variety of ion-water systems and discuss the microscopic characteristics of ion hydration and its relation to solvation thermodynamics. There, we also present a microscopic interpretation of the Samoilov model of the ion hydration in terms of the pair correlation function between ion and water. 2. Theory In what follows, we consider a solution composed of a single solute and a single solvent species for notational simplicity. Extensions of the formalism to more complex solutions are straightforward. We also assume that the solute is monatomic. We denote species with superscripts and interaction sites of molecules with subscripts unless otherwise specified. 2.1. Ex-RISM and Solvation Thermodynamics. In this part of the section, we briefly outline the main features of exRISM,12-14 the derivatives of solvent correlation functions with respect to a solute density,15 and thermodynamic functions of solvation. Most of the theory given here has been described in the literature, except for the derivation of the formula for the partial molar volume (eq 10), which is presented in the Appendix. We begin with the RISM integral equation for a solution that consists of a molecular solvent and a monatomic solute at the infinite dilution limit vv
v
vv
v
v
vv
vv
h ) w *c *w + w *c *Fvh huv ) cuv*wv + cuv*Fvhvv
(1) (2)
where h, c, and w are site-labeled matrices, Fv is the solvent
density, and “*” denotes a matrix convolution product. hvv and huv are the matrices of the solvent(v)-solvent(v) and solute(u)-solvent(v) total correlation functions, respectively, and cvv and cuv are the matrices of the v-v and u-v direct correlation functions, respectively. wv is a solvent intramolecular correlation matrix. The closure used in this work is an analog of the hypernetted chain (HNC) closure
cij ) exp[-βuij + hij - cij] - hij + cij - 1
(3)
where hij and cij are either the v-v or u-v pair correlation functions. β-1 ) kBT, where kB is the Boltzmann constant and T is the absolute temperature. uij is the spherically symmetric pair potential between an ij pair of interacting sites. In the case where each site carries a charge, uij consists of a short-ranged potential u* ij and a long-ranged Coulomb potential, and it is most convenient to cast eqs 1-3 in renormalized form defining the modified direct correlation functions
c*ij ≡ cij + β(uij - u* ij)
(4)
The details of this procedure have been given elsewhere.12-14 Next, we briefly follow the work of Yu and Karplus15 to obtain the derivatives of the v-v correlation functions with respect to a solute density. These functions are required to calculate the solvation energy at infinite dilution (see below). In the next part of this section, it is shown that these functions also give information about the structural changes induced by the presence of a solute. To extract a first-order correction to the solvent correlation functions, Yu and Karplus expanded the RISM equation around the limit Fu ) 0 retaining terms to first order in Fu, where Fu denotes a solute density, and obtained the following equation
δFuhvv ) [wv + Fvhvv]*δFuc*vv*[wv + Fvhvv] + hvu*huv (5) where δFuhvv and δFuc*vv denote the site-labeled matrices of the density derivatives of the v-v total and modified direct correlation functions, respectively. The closure to solve eq 5 can be obtained by expanding eq 3 to first order in Fu as vv vv VV δFuhvv ij ) (δFuhij - δFucij )(hij + 1)
(6)
Thermodynamic functions of solvation can be described in terms of the correlation functions given above. The excess chemical potential, also referred to as the solvation free energy, of ion R at infinite dilution has the form17,18
β∆µ0R ) 4πFv
{
∑s ∫ r2 dr
1 uv 1 uv 2 uv hRs(r) - cuv Rs(r) - hRs(r) cRs(r) 2 2
}
(7)
where the superscript 0 indicates the infinite dilution limit and the Greek and Roman subscripts refer to the interaction sites of the solute and solvent molecules, respectively. The solvation energy at infinite dilution can be written as15
∆0R ) 4πFv
∑s ∫ r2 dr uuvRs(r) guvRs(r) + F2v
vv 4πr2 dr uvv ∫ ∑ ss′(r) δFuhss′(r) 2 s,s′
(8)
The solvation entropy at infinite dilution is given by
T∆s0R ) ∆0R - ∆µ0R
(9)
Theoretical Study of Ion Hydration
J. Phys. Chem. B, Vol. 101, No. 16, 1997 3211
on the basis of the knowledge of ∆µ0R and ∆0R. Derivation of another important thermodynamic function of solvation, the partial molar volume of an ion at infinitely dilute solution, is presented in the Appendix. Here, we only cite the resultant ex-RISM-based expression for the partial molar volume (see eqs A18 and A19)
V h 0* R
)
kBTχ0T
(1 - Fv
∑s
c˜ *uv Rs (0))
We take the excess energy, also called the solvation energy, as a guidepost to explore the relation between U- and Fsolutions. The excess energy per solute for U- and F-solutions are defined as follows
∆U )
(10)
E(V, Nv, Nu) - E(V, Nv, 0) Nu
∆F ) E(V, Nv, 1) - E(V, Nv, 0)
where χ0T denotes the isothermal compressibility of the pure solvent and the tilde denotes the function in Fourier k-space. This equation gives the values for the partial molar volumes of individual ions that correspond to those determined from experiments based on the ultrasonic vibration potentials.19,20 Detailed discussion about the partial molar volumes of individual ions is given in the Appendix. In section 3.3, we are interested in the solvent-mediated potential of mean force (PMF) between an ij pair of interacting sites. This quantity, denoted as Wij, can be obtained from the following equation21
∆0 ) lim
E(V, Nv, Nu) - E(V, Nv, 0) Nu
(15)
) lim E(V′, N′v, 1) - E(V′, N′v, 0)
(16)
) lim E Nuf0
(
) (
(14)
)
(11)
where gij(r) ≡ hij(r) + 1 and denotes a site-site radial distribution function between the ij pair of interaction sites. 2.2. Structural Information in the Density Derivatives. We are interested in the changes in the water structure in the vicinity of an ion, which are essentially dependent on the distance from the ion, and correlation functions that describe these structural changes become the three-particle (u-v-v) correlation functions. Common simulation studies for solutions treat systems that consist of one solute and many solvent molecules and are therefore suited to investigate these threeparticle correlations.16 The integral equation method, however, assumes the translational invariance of the system; i.e., the densities (one-particle distributions) of solute and solvent are assumed to be uniform throughout the system. Under this assumption, the infinite dilution limit is taken, leading to the derivation of eqs 1-3. The solute does not affect the solvent structure in the infinitely dilution limit, and we cannot obtain the information about the solvent structural changes induced by the presence of the solute by solving only eqs 1-3. It is expected that the density derivatives of the solvent correlation functions, δFuhvv’s, may give the information about the structural changes of the solvent due to the presence of solute, but eqs 5 and 6 still assume the translational invariance of the system. The resulting δFuhvv is independent from the location of the solute; i.e., δFuhvv is a two-particle (v-v) correlation function. In this part of the section, we examine the relation between the three-particle (u-v-v) correlation function mentioned above and δFuhvv, which can be obtained from the integral equation method. In what follows, we assume that both solute and solvent molecules are spherically symmetric for notational simplicity, and the subscripts will be omitted. We define U-solution and F-solution as follows. U-solution consists of solvent and solute molecules that are assumed, on the average, to be distributed uniformly. It is therefore translationally and rotationally invariant. F-solution consists of solvent molecules and one solute molecule fixed at the origin of the coordinate system. It therefore loses the translational invariance but keeps the rotational invariance with respect to the origin since we have assumed a spherically symmetric solute. The integral equation methods treats U-solution, while common simulations concerning the infinite dilution work well for F-solution.
(13)
where V denotes the volume of the system and Nu and Nv refer to the numbers of solute and solvent molecules, respectively. Next, we take the infinite dilution limit of U- and F-solutions, where infinite dilution limit is taken as Nu f 0 in U-solution and as Nv f ∞ in F-solution (see below). If U- and F-solutions have a common solvent density, eqs 12 and 13 for the excess energy at infinite dilution become equivalent in the thermodynamic limit. This can be shown as follows
Nuf0
Wij(r) ) -kBT log gij(r)
(12)
V Nv V Nv , ,1 -E , ,0 Nu Nu Nu Nu
V ′f∞ Nv′f∞
where the superscript 0 indicates the infinite dilution result. We have defined V ′ ≡ V/Nu and Nv′ ≡ Nv/Nu, and the final limit is taken such that Fv ) Nv/V ) Nv′/V ′ ) constant. Equation 14 gives the solvation energy of U-solution at infinite dilution (∆0U), while eq 16 gives that of F-solution (∆0F). Corresponding equations for the excess chemical potential and the excess entropy can be derived in a similar way. Thus, thermodynamic functions in U- and F- solutions at infinite dilution are equivalent if they have a common solvent density. Utilizing this, we next study the relation of structural quantities between U- and F-solutions. The expression for the excess energy of U-solution in terms of the correlation functions is as follows (see eq 8):
∆0U ) Fv
∫ 4πr2 dr uuv(r) guv(r) + F2v 2
∫ 4πr2 dr uvv(r) δFuhvv(r)
(17)
The corresponding expression for F-solution is22
∫ dr uuv (|r|) Fv(1) F (r) + 1 (r, r′) - Fvv(2) (r, r′)] ∫ dr ∫ dr′ uvv (|r - r′|)[Fvv(2) F 0 2
∆0F )
(18)
vv(2) In the above equation, Fv(1) (r, r′) are one- and F (r) and FF two-particle solvent distribution functions in F-solution, respectively, and Fvv(2) (r, r′) is a two-particle solvent distribution 0 vv(2) function in a pure solvent. Note that Fv(1) (r, r′) F (r) and FF are actually two-particle (u-v) and three-particle (u-v-v) distribution functions, respectively, since the solute molecule is fixed at the origin in F-solution. As we have assumed that the solute is spherically symmetric, Fv(1) F (r) depends only on r ) |r|, and Fvv(2) (r, r′) depends on a, θ, and R in Figure 1 where F R ≡ |r - r′|. Fvv(2) (r, r) of course depends only on |r - r′| 0
3212 J. Phys. Chem. B, Vol. 101, No. 16, 1997
Chong and Hirata TABLE 1: Parameters for the Solute-Solvent LJ Potentialsa
u*(r) ) 4[(σ/r)1/2 - (σ/r)6] ion (X)
σOX
σHX
OX ) HX
Li+
2.28 2.72 3.16 3.55 2.95 3.55 3.70 3.70
0.87 1.31 1.75 2.14 1.54 2.14 2.29 2.29
1.0 × 10-14 9.3 × 10-15 9.1 × 10-15 9.1 × 10-15 2.0 × 10-14 2.5 × 10-14 1.25 × 10-14 1.25 × 10-14
Na+ K+ K′+ FClG+ G2+ a
σ in angstroms; in ergs.
TABLE 2: Thermodynamic Functions of Solvationa Figure 1. Three-particle (u-v-v) correlation.
and can be written in terms of the correlation functions of the pure solvent as
Fvv(2) (r, r′) ) F2vgvv(|r - r′|) 0
(19) a
We define the correlation functions in F-solution as follows uv Fv(1) F (r) ) FvgF (r)
(20)
Fvv(2) (r, r′) ) F2vgvv F F (a, θ, R)
(21)
Then, eq 18 can be rewritten in the following form:
∫ r2 dr uuv(r) guvF (r) + 4π2 F2v ∫ r2 dr uvv(r) ∫ a2 da ∫ sin θ dθ ×
∆0F ) 4πFv
vv [gvv F (a, θ, r) - g (r)] (22)
where the notation R in eq 21 has changed to r in eq 22. From the equivalence of ∆0U and ∆0F, valid when we take the same solvent density Fv, the relation of the correlation functions between U- and F-solutions can be obtained by comparing eqs 17 and 22. We first obtain the trivial relation
guv(r) ) guv F (r)
(23)
by comparing the first terms in eqs 17 and 22 and next obtain the desired relation
δFuhvv(r) ) 2π
∫ a2 da ∫ sin θ dθ [gvvF (a, θ, r) - gvv(r)]
(24)
by comparing the second terms in eqs 17 and 22. In the righthand side of eq 24, contribution from the solvent molecules in the region far from the solute will be small, since gvv F (a, θ, r) gvv(r) will be small in that region: the solvent structure will not be changed significantly from that of the pure solvent where solvent molecules are far from the solute. Thus, δFuhvv(r) gives structural changes of the solvent in the vicinity of the solute averaged over the u-v distance and orientations. 3. Results and Discussion 3.1. Models. We first describe the potential functions which will be used in the specific applications of the theory of the previous section. The site-site interactions consist of the
ion (X)
∆µ0
∆0
T∆s0
Li+ Na+ K+ FClG+ G2+
-104.3 -80.04 -62.01 -148.1 -86.59 -42.61 -218.7
-114.9 -89.05 -71.49 -168.8 -101.5 -55.55 -240.1
-10.58 -9.019 -9.483 -20.69 -14.95 -12.94 -21.42
In kcal/mol.
Lennard-Jones (LJ) interaction and a Coulomb interaction. For the solvent model, we employ the TIP3P water23 modified by Pettitt and Rossky,24 at temperature of T ) 300 K and the number density of Fv ) 0.033 34 molecules Å-3. Parameters for the solute-solvent LJ potentials are listed in Table 1. Parameters for alkali cations and halide anions have been taken from the literature.25 Other ions we added in Table 1, denoted as G+ and G2+, are artificial ions which were used by Geiger in his molecular dynamics simulation study of the negative hydration effect in aqueous electrolyte solutions.26 It is known from the results of Geiger that G+ ion is a negatively hydrated ion while G2+ ion is a positively hydrated ion. These ions are used to check the validity of our idea for the structural information that the density derivatives of the solvent distribution functions have. 3.2. Thermodynamic Functions of Solvation. Table 2 lists the thermodynamic functions calculated for a variety of infinitely dilute aqueous solutions. It is seen that all ions are stabilized in water due primarily to large energetic stabilization, though the ions disfavor aqueous environment entropically. The size and charge dependence of the solvation free energy and energy is apparent from Table 2: the smaller the size and the larger the solute charge, the more stabilized are ions. For ions of nearly same size with opposite charges (Na+/F-), it is seen that an anion is more stabilized than a cation. This reveals that, in aqueous solutions, the degree of solvation of an anion is greater than that of a cation. This feature is intimately related to the asymmetry of the molecular structure of water,27,28 which was first known almost 60 years ago.29 The entropy of solvation has a peculiar character (Table 2). The interpretation of the trend in the solvation entropy in terms of the solute size or the magnitude of its charge is not straightforward. To investigate the origin of this peculiar behavior of the solvation entropy, we decomposed it into two parts, the entropy change that stems from the solvation process of a LJ solute, i.e., a solute without its charge, and the one which stems from the process of charging up that LJ solute. Table 3 lists the results of this decomposition. From the entropy change corresponding to each process, we can provide clear explanation for each value. The entropy change at the solvation process of
Theoretical Study of Ion Hydration
J. Phys. Chem. B, Vol. 101, No. 16, 1997 3213
TABLE 3: Decomposition of the Solvation Entropya
a
TABLE 4: Partial Molar Volumes of Individual Ionsa
ion (X)
T∆sLJ
T∆scharge
T∆stotal
ion
V h 0R b
V h 0LJ c
∆V hd
V h 0R e (experiment)
Li+ Na+ K+ FClG+ G2+
-1.461 -2.857 -5.006 -4.712 -9.743 -9.630 -9.630
-9.119 -6.162 -4.477 -15.98 -5.207 -3.310 -11.79
-10.58 -9.019 -9.483 -20.69 -14.95 -12.94 -21.42
Li+ Na+ K+ FCl-
-9.1 3.4 17.9 3.5 23.3
4.8 11.2 21.6 14.4 30.3
-13.9 -7.8 -3.7 -10.9 -7.0
-11.2 -7.4 3.4 3.3 23.7
In kcal/mol.
the LJ solute, denoted as T∆sLJ, can be interpreted reasonably well in terms of the size of the solute: the larger the solute size, the more negative the entropy change, and vice versa. The origin of this behavior is largely because of the volume exclusion effect of the solute, although the structural changes of the solvent, especially in case of water, may have some effects. The entropy change at the process of charging up that LJ solute, denoted as T∆scharge, can be well interpreted in terms of the degree of the polarization of the medium induced by the electrostatic field from the solute. The smaller the size and the larger the solute charge, the more the dipole moments of the solvent molecules around the solute are aligned along the radial direction from the solute, and this causes large reduction in the entropy. The behavior of T∆scharge’s in Table 3 can be understood in this manner. The peculiar character of the solvation entropy given in Table 2 is a result of interplay of these two competing effects. Strictly speaking, other thermodynamic functions of solvation should be considered in a similar way as for the solvation entropy. However, since the electrostatic part in the ion-water interaction energy is significant, its contribution dominates the values in the solvation free energy as well as the energy, and contribution from the LJ part can be considered to be small. In what follows, we will see that the partial molar volumes of individual ions should be interpreted in a similar way as we have done for the solvation entropy. As mentioned above, the entropy change at the solvation process of the LJ solute is dominated by its steric effect, i.e., translational contribution, which is apparent from the size dependence of T∆sLJ’s in Table 3. The entropy change at the process of charging up the LJ solute is dominated by an orientational contribution of the solvent, which is also apparent from the size and charge dependence, i.e., surface charge dependence, of T∆scharge’s in Table 3. But these explanations are not the whole story. There will be an orientational contribution to T∆sLJ, in conjunction with the hydrophobic hydrations, and there will be a translational contribution to T∆scharge. There may be a contribution to T∆scharge by breaking solvent structure around the ion due to the negative hydrations. Thus, the separation of translational and orientational contributions is desired to understand the behavior of the solvation entropy in more detail. Such a study, exploiting the GreenWallace expansion method,30 is currently underway in our group, and results will be reported in a subsequent paper. We next discuss another important thermodynamic function of solvation, the partial molar volumes, which are sometimes used to extract information about the ion-solvent interactions.31,32 It is a nontrivial problem to separate an experimentally observed value for a salt, which includes contributions from cations and anions, into contributions from individual ions. The method based on the ultrasonic vibration potential (UVP) has been used to separate the partial molar volumes.19,20 In the Appendix, we give the derivation of an ex-RISM-based formula for the partial molar volumes of individual ions that can be identified as those determined from UVP.
a In cm3/mol. b Theoretical results calculated from eq 10. c Partial molar volumes of corresponding LJ solutes. d The changes in the partial molar volume due to the electrostriction. e Experimental values determined from the ultrasonic vibration potentials (ref 20).
Two factors that determine the value of the partial molar volume of an ion are, as in the solvation entropy case, the magnitudes of its size and charge. The magnitude of the ionic size is manifested in a exclusion volume effect that results in a positive contribution to the partial molar volume, while that of the ionic charge is apparent in an electrostriction effect that gives a negative contribution.33 It has been believed that the latter effect is intimately related to the degree of solvation of an ion. Kusalik and Patey34 (KP), however, suspected that this idea may be false and that the individual partial molar volumes may not in general be used to deduce accurate information about the nature of the ion-solvent interaction. We show below, contrary to KP, that our calculation based on ex-RISM supports the idea that the individual partial molar volumes do reflect the nature of the ion-solvent interaction. Table 4 gives the individual ion partial molar volumes calculated via eq 10. Isothermal compressibility of the pure water that is required in the calculation of eq 10 is calculated from ex-RISM as 5.42 × 10-5 bar for TIP3P water, which is close to the experimental value, 4.57 × 10-5 bar.9 From Table 4, it is seen that our calculated results are in fair agreement with the experimental results based on UVP,20 although significant difference between the observed and theoretical values is seen especially in the result for K+. Two causes are conceivable for the discrepancy. The energy parameters for ion-water interaction employed in the theoretical calculation might not have been well-tuned for comparing with the experiment, or the experimental value for the partial molar volume of K+ is not correct. It is not known at this moment which causes the problem. To extract the exclusion volume and electrostriction effects of ions, the LJ solutes of the corresponding ions are used as references, and we assume that the partial molar volumes of the LJ solutes are exclusively due to the exclusion volume effects. Strictly speaking, this is not true, but the dominant contribution to V h 0LJ comes from the exclusion volume effect. We then regard the difference of the partial molar volume of an ion from that of the corresponding LJ solute as the contribution due to the electrostriction effect. The values of contributions from these effects are also included in Table 4. A number of observations can be made concerning the theoretical results given in Table 4. We take the Na+/F- pair as an representative example. It is of interest to relate the results of Table 4 with those of Table 2 in an attempt to examine how the individual ion partial molar volumes are related to the ionwater interactions. From the thermodynamic consideration given in Table 2, F- ion is more firmly hydrated than Na+ ion although they have almost equal size. (See Table 1. F- ion is slightly larger in size than Na+ ion in terms of σ in LJ parameters.) This is intimately related to the asymmetry of the water molecule, and in general, the solvation of an anion is stronger than that of a cation of the same size in water. This thermodynamic consideration reveals that, for ions of almost equal size, the degree of the electrostriction of the solvent in
3214 J. Phys. Chem. B, Vol. 101, No. 16, 1997
Chong and Hirata
Figure 3. ∆Ei’s for various ions.
Figure 2. Potential of mean force (PMF) curves. Solid line, PMF for O(water)-H(water); dashed line, PMF for O(water)-Li+; dot-dashed line, PMF for O(water)-K+.
TABLE 5: Values of ∆Ei for Various Ionsa ion
∆Eib (ex-RISM)
∆Eic (experiment)
Li+ Na+ K+ FCl-
1.35 0.10 -0.60 2.17 -0.20
0.73 0.25 -0.25 -0.27
a
In kcal/mol. b The results of the present work. See text for details. c The values reported by Samoilov (ref 11).
the vicinity of an anion is greater than that of a cation. From apparent observation that the partial molar volume of Na+ ion is less than that of F- ion, however, one may conclude that Na+ ion is more firmly hydrated than F- ion. This idea is, of course, false because it neglects the fact that F- ion is slightly larger in size than Na+ ion. What should be related to the degree of solvation is the contribution from the electrostriction effect included in the partial molar volume, not the partial molar volume itself. From Table 4, one can see that the degree of the electrostriction effect due to F- ion is greater than that due to Na+ ion, which is consistent with the thermodynamic consideration given above. The calculated results for other ions given in Table 4 can be interpreted in a similar way. As mentioned above, factors that determine the partial molar volume of an ion are the magnitude of its size and charge. The former is related to the steric effect while the latter to electrostriction effect. The bigger the size of an ion, the larger the steric effect. This feature can be seen from V h 0LJ’s in Table 4. The electrostriction effect reflected in ∆V h ’s of Table 4 is consistent with the thermodynamic consideration given in Table 2, where it is clear that the stronger the solvation, the more negative in ∆V h , and vice versa. In particular, the behavior of ∆V h ’s in Table 4 is in accord with that of T∆scharge’s in Table 3.35 The results imply that the contribution from ∆V h to the individual ion partial molar volume can be a good measure of characterizing the ion-water interaction. The reference hypernetted-chain (RHNC) results for the partial molar volumes by KP have exhibited some inconsistency with other thermodynamic functions of solvation concerning the effect of ions with opposite charges.34 Their results imply that their model for the ion-water interaction does not properly account for the asymmetry of the charge distribution of water molecules as far as the partial molar volume is concerned. 3.3. Revisiting the Samoilov Model for the Ion Hydration. Here, we revisit the Samoilov model for the ion hydration by means of the ex-RISM theory. The Samoilov model focuses on the exchange of water molecules closest to an ion with those
Figure 4. Radial distribution functions and their temperature derivatives for pure water. (a) Solid line, gOO; dashed line, δThOO. (b) Solid line, gOH; dashed line, δThOH.
in the bulk.11 The exchange process is viewed as an activated jump, and the ion hydration is characterized by the difference of the activation energies ∆Ei ≡ Ei - E0 associated with those processes. The significance of the activation energies Ei and E0 will now be discussed in detail. Let us consider the translational motion of a water molecule in a dilute aqueous electrolyte solution. The translational motion of a water molecule in the immediate vicinity of an ion and in the bulk water is different because ion-water interactionsare not equivalent to water-water interactions. This difference in the translational motion of a water molecule can be expressed in terms of the mean residence time defined as follows. Let a water molecule spend an average time τ0 in the immediate vicinity of another water molecule. The quantity τ0 is the mean time during which two solvent molecules stay in the immediate vicinity of each other. Similarly, let τi be a mean time for a water molecule to stay in the immediate vicinity of an ion.
Theoretical Study of Ion Hydration
J. Phys. Chem. B, Vol. 101, No. 16, 1997 3215
Figure 5. Solvent radial distribution functions and their derivatives with respect to the density of G+. (a) Solid line, gOO of pure water; dashed lne, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
Samoilov classified the hydration of ions as “positive hydration” in the case τ0 < τi, and in the opposite case, τ0 > τi, as “negative hydration”. Thus, the quantity τi/τ0 bears primary importance in characterizing the ionic hydration in the Samoilov model. This quantity can be directly calculated using molecular dynamic simulations in a straightforward way.36 From the theoretical view point, however, it is desired to reduce the problem of determining τi/τ0 to a more theoretically tractable problem with the aid of the simple reaction rate theory. Let E0 be the activation energy of the process for a water molecule to transfer from the first coordination shell of another water molecule to the next coordination shell. The reaction rate theory in its simplest form tells us the relation between τ0 and E0 as
1 kBT exp(-βE0) ) τ0 h
(25)
The quantity τi is, in a similar way, related to Ei, the activation energy of a water molecule to transfer from the immediate vicinity of an ion to the positions next to the first coordination shell. If we define ∆Ei ≡ Ei - E0, it follows that
τi ) exp(β∆Ei) τ0
(26)
Now, the description of the ion hydration is reduced to the analysis of ∆Ei. The sign of ∆Ei is of special interest: ∆Ei > 0 corresponds to the positive hydration and ∆Ei < 0 to the negative hydration.
Figure 6. Solvent radial distribution functions and their derivatives with respect to the density of G2+. (a) Solid line, gOO of pure water; dashed line, δFuhOO. (b) Solid line, gOH or pure water; dashed line, δFuhOH.
Samoilov extracted ∆Ei from the experimental data for the transport coefficients using the Nernst equation for a model of an activated jump diffusion. It is not necessarily obvious what kind of quantities best corresponds theoretically to ∆Ei obtained by Samoilov, because the earlier method includes an ambiguity largely depending on the diffusion model. It may be most natural to assume that the activated jump takes place along the reaction path concerning the center-to-center distance of an ion and a water molecule. The free energy profile of the reaction is nothing but the potential of mean force curve between the ion and a water molecule, which is averaged over the orientation of the water molecule of interest as well as over all configurations of other water molecules. Then, the activation energy Ei for the transfer of a water molecule from the immediate vicinity of an ion to the next coordination shell can be defined as the difference in the potential of mean force between the first minimum and the first maximum. The activated jump process is always associated with breaking of a “hydrogen bond” between the ion and a water molecule as a dominating factor, along other changes including the reorganization of water molecules. Therefore, it is also sensible to define Ei in terms of the first two extrema in the site-site potential of mean force between an ion and an atom in a water molecule that participate in the “hydrogen bonding”: when cations are concerned, the ion-oxygen (water) pair makes the “hydrogen bond” while the ion-hydrogen (water) pair is responsible for the anion-water “hydrogen bond”. The activation energy E0 for transferring a water molecule from the immediate vicinity of another water molecule to the next coordination shell can be defined in a
3216 J. Phys. Chem. B, Vol. 101, No. 16, 1997
Chong and Hirata
Figure 7. Solvent radial distribution functions and their derivatives with respect to the density of Li+. (a) Solid line, gOO of pure water; dashed line, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
Figure 8. Solvent radial distribution functions and their derivatives with respect to the density of Na+. (a) Solid line, gOO of pure water; dashed line, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
similar way by considering the O-H site pair (or the hydrogen bond pair) between a pair of water molecules. The potential of mean force curves for Li+-O and K+-O are plotted in Figure 2 along with the curve for O-H. From this figure and similar figures for other ions, we can obtain ∆Ei’s for various ions. The results are listed in Table 5 and plotted in Figure 3. From Table 5 and Figure 3, it is seen that Li+, Na+, and F- ions are the positively hydrated ions while K+ and Cl- ions are the negatively hydrated ions on the basis of our definition for the activation energies. Our classification of ions shown in Table 5 are in accord with the original one of Samoilov,11 in which ∆Ei’s are extracted from the experimental quantities. This agreement suggests that our theoretical definition of ∆Ei based on the site-site potential of mean force is a good measure of characterizing the ion hydration in terms of the Samoilov model. 3.4. Analysis of the Density Derivatives. Upon insertion of an ionic solute into water, the structure of water changes in several respects due to the perturbation produced by the ion. In this part of the section, we investigate how the correlations of water in the vicinity of an ion are affected due to the presence of the ion in conjunction with the earlier models of hydrations. We first give a brief outline for the well-known structural features of the pure water manifested in the pair correlation functions and then describe in terms of the density derivatives of the pair correlation functions how the structural features are modified by ions. Although it is known that the theoretical results for water correlation functions have quantitative defects in several respects compared to those obtained from the
simulations or the neutron-scattering experiments,24 they exhibit well-regarded qualitative features of water structure, and the essential points of the present analysis would not be altered. The site-site radial distribution functions, gvv’s, for pure water obtained from the ex-RISM calculations are shown by solid lines in Figure 4. The water structure can be summarized as the icelike tetrahedral structure due to the hydrogen bonding that appears and disappears in a few picoseconds. The primary cause of the tetrahedral coordination is the hydrogen bond between an O-H pair, which shows up in the O-H pair correlation function as the sharp peak located around r ) 1.5 Å in Figure 4b. The main feature is also manifested in an O-O radial distribution function in two distinct manners, the narrowing of the first peak compared to the simple liquids such as neon, which represents the fewer coordination numbers reflecting the tetrahedral coordination, and the appearance of the second peak around r ) 4.3 Å in Figure 4a, which reflects a triplet correlation of water molecules characteristic of the tetrahedral coordination.24 These structural features are affected by several conditions, such as increasing temperature, insertion of an ion into water, etc. We first examine the temperature derivatives, δThvv’s, of the solvent distribution functions28 denoted as dashed lines in Figure 4. These temperature derivatives represent the changes in the site-site distribution functions induced by an increase in temperature. In the figure, δThvv’s are multiplied by T ) 300 K to facilitate comparisons with gvv’s. Both the peak representing the tetrahedral coordination, i.e., the second peak, in gOO and the hydrogen-bonding peak, i.e., the first peak, in gOH are
Theoretical Study of Ion Hydration
Figure 9. Solvent radial distribution functions and their derivatives with respect to the density of K+. (a) Solid line, gOO of pure water; dashed line, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
diminished. Thus, at room temperature, the main features of water structure are decreased by increasing temperature. The innermost peak in δThOO is a consequence of the increased sampling probability of the repulsive configurations at elevated temperatures. Now we investigate the structural changes induced by the presence of an ion. Among those structural changes, the changes that reflect the qualitative differences of the positive and negative hydrations of ions are of particular interest. As pointed out in section 2.2, the derivatives of the solvent distribution function with respect to the solute density, δFuhvv(r)’s, give the structural changes of the solvent in the vicinity of the solute averaged over the distance and orientations. The water structural changes are analyzed using these density derivatives. G+ and G2+ ions are taken as guideposts since, from the simulation study by Geiger,26 G+ is known to be a negatively hydrated ion while G2+ is known to be a positively hydrated ion. The structural changes of water induced by positive density variations of G+ and G2+ ions are shown by the dashed lines in Figures 5 and 6, respectively. In the figures, δFuhvv’s are multiplied by FH2O/2 for ease of comparison with the unperturbed values of gvv’s. The tendency of decrease in the height of the tetrahedral peak, i.e., the second peak, in gOO and of the hydrogen bonding peak, i.e., the first peak, in gOH can be seen in both Figures 5 and 6. Thus, both G+ and G2+ ions tend to diminish the intrinsic structure of the pure water that is formed due to the hydrogen bonding. The feature which distinguishes the hydrations of G+ and G2+ ions is seen in the first peak
J. Phys. Chem. B, Vol. 101, No. 16, 1997 3217
Figure 10. Solvent radial distribution functions and their derivatives with respect to the density of F-. (a) Solid line, gOO of pure water; dashed line, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
position in δFuhOO: in Figure 5(a), the first peak in δFuhOO shifts the first peak in gOO inward, while in Figure 6(a), the first peak in δFuhOO shifts the first peak in gOO outward. Why does the position of the first peak in gOO tend to shift inward due to the perturbation from G+ ion? A similar feature has been seen in the temperature derivative of the O-O distribution function, δThOO in Figure 4a. As has been discussed above, the reason why the first peak in gOO tends to shift inward by an increase in temperature is that a pair of water molecules tends to approach closer to each other due to the increased thermal motion overcoming the core repulsion. From an analogy to this, it can be said that “local temperature” around G+ ions is higher and that the water molecules in the vicinity of G+ ion are more mobile than bulk water molecules. Thus, the reduced position of the first peak in gOO is supposed to be a manifestation of the “higher local temperature” of water molecules in the vicinity of G+ ion. Then, why is it not the case in G2+ ion? It is likely that the water molecules in the vicinity of G2+ ion are highly aligned along the radial direction of the ion due to the strong electrostatic field produced from it. The outward shift of the first peak is gOO can be explained in terms of the increased electrostatic repulsion between a pair of water molecules next to the ion as follows. In the bulk water, the O-O distance is shortened due to the hydrogen bondng compared with the case in which there was no hydrogen atom between the two oxygen atoms. If the water structure is disrupted by aligning their dipole moments along the large electrostatic field of the ion, the O-O distance between the water molecules around the solute will be enlarged due to the electrostatic repulsion between them, which is increased by the absence of the hydrogen bonding. Thus,
3218 J. Phys. Chem. B, Vol. 101, No. 16, 1997
Chong and Hirata
Figure 11. Solvent radial distribution functions and their derivatives with respect to the density of Cl-. (a) Solid line, gOO or pure water; dashed line, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
Figure 12. Solvent radial distribution functions and their derivatives with respect to the density of K′+. (a) Solid line, gOO of pure water; dashed line, δFuhOO. (b) Solid line, gOH of pure water; dashed line, δFuhOH.
the enlarged O-O distance in Figure 6a is a manifestation of the increased alignment of the water molecules along the large electrostatic field of the ion, and it can be said that the solvent molecules are less mobile around G2+ ion. G2+ ion can therefore be regarded as positively hydrated. The explanations given above are consistent with Geiger’s simulation result that has shown that G+ ion is negatively hydrated and G2+ ion is positively hydrated. The water molecules in the immediate vicinity of the negatively hydrated ion are more mobile than bulk water molecules, while the water molecules in the immediate vicinity of the positively hydrated ion have more aligned structure, thus are less mobile. Figures 7-11 give the density derivatives of solvent distribution functions for various ions. By comparing these figures with Figures 5 and 6, it can be said that Li+, Na+, and F- are positively hydrated ions and that Cl- is negatively hydrated. Figure 9 for K+ ion has intermediate features. It is known from our experience that the density derivatives of the solvent correlation functions are sensitive to the u-v potential parameters, especially to σ in the LJ parameters, as they should be. The results presented in Figure 9 suggest, considering the commonly accepted view that K+ is a negatively hydrated ion, that σ in the LJ parameters for K+ ion given in Table 1 may be rather small. Therefore, we added Figure 12 for a modified potassium ion, referred to as K′+ in Table 1, which has a larger σ of the LJ parameters than original K+ to enhance the negative hydration effect by weakening the electrostatic effect. (Note that σ in the LJ parameters for K′+ is the same as that for Clion.) From Figure 12, it is seen that K′+ ion becomes a negatively hydrated ion and behaves as we have anticipated.
4. Concluding Remarks In this paper, we studied the hydration of ions in several aspects on the basis of the ex-RISM theory. We analyzed the thermodynamic functions of solvation, especially the partial molar volumes of individual ions at infinite dilution. Special attention is paid to the information that can be extracted from the partial molar volumes of individual ions at infinite dilution. It can be concluded that the contribution from the electrostriction effect, ∆V h , to the individual ion partial molar volume can be a good measure of characterizing the ion-water interactions. We also presented a microscopic description for ∆Ei, a characterization of the ion hydration by Samoilov, in terms of the ionwater potential of mean force, and found that our classification of ions into the “positive” and “negative” hydrations is in accord with the original one of Samoilov. We then discussed the structural changes of water due to the presence of an ion based on the density derivatives of the solvent distribution functions. Qualitative differences of the density derivatives between the “positively” and “negatively” hydrated ions were observed and found to be consistent with the analysis of the potential mean force. Acknowledgment. S.-H.C. gratefully acknowledges a research fellowship from the Japan Society for the Promotion of Science for Young Scientists (JSPS). We express our thanks to Dr. D. H. Peapus for reading the manuscript carefully. This work is supported by grants-in-aid for the Scientific Research from the Ministry of Education, Science and Culture of Japan. We are also grateful to Research and Development Corporation of Japan (JRDC) for the financial support.
Theoretical Study of Ion Hydration
J. Phys. Chem. B, Vol. 101, No. 16, 1997 3219
Appendix. Partial Molar Volumes of Individual Ions at Infinite Dilution In this Appendix, we consider the partial molar volumes of individual ions at infinite dilution. The partial molar volumes are of particular interest since these quantities are sometimes used to extract information about the ion-solvent interactions.31,32 Our formulation given here closely follows the works of Kusalik and Patey34,37 but is based on the ex-RISM theory. Experimentally, the partial molar volumes at infinite dilution are, in general, first determined for a salt, denoted as V h 02, and 0 then V h 2 is split into individual ion contributions on the basis of a variety of physical assumptions. Theoretically, the partial molar volume of a salt is a well-defined quantity that extrapolates smoothly to its limiting value V h 02 at infinite dilution. However, the partial molar volume of an individual ion is not so well defined theoretically. Since the experimental determination of the partial molar volume of an individual ion depends upon the underlying physical assumptions, it is expected that all the results reported from the experiments cannot be described by a single microscopic formula. There exists an experimental technique that does give the absolute values for the partial molar volumes of individual ions. The method is based on the ultrasonic vibration potentials.19,20 It is therefore quite reasonable to seek a microscopic formula for the partial molar volume of an individual ion that can be identified as that determined from the ultrasonic vibration potentials. Kusalik and Patey (KP)34,37 studied this problem which will be summarized below. The Kirkwood-Buff theory38 gives the exact expression for the partial molar volumes of an n-component system, and KP extended that formalism to mixtures of electrolytes in molecular solvents.34,37 KP also argued what is meant by the partial molar volume of an individual ion and what theoretical quantities should be related to the experimentally determined individual ion contribution. We first briefly outline the formalism of KP for the individual ion partial molar volume at infinite dilution. It is useful to recall here the key definitions in the KirkwoodBuff theory 000 GRβ ≡ h˜ 00;Rβ (0) ) 4π 000 (0) ) 4π CRβ ≡ c˜ 00;Rβ 000 (k) h˜ 00;Rβ
000 (r) dr ∫ r2h 00;Rβ
(A1)
000 (r) dr ∫ r2 00;Rβ
(A2)
000 c˜ 00;Rβ (k)
and are, respectively, the Fourier where transforms of the radial pair and direct correlation functions associated with species R and β in the system. The subscript and superscript zeros on the h’s and c’s denote that these functions are averaged over angular coordinates. In terms of GRβ’s defined above, the Kirkwood-Buff theory yields the partial molar volume of an ion R as
V h 0R )
(A3)
where the superscript 0 indicates the infinite dilution result and subscripts R and V denote the ionic and solvent species, respectively. Note that, from eq A1 to eq A7, subscripts refer to species. The final equality in eq A3 comes from the wellknown relationship
1 + FvG 0vv ) (1 - FvC 0vv)-1 )FvkBTχ 0T χ 0T
where solvent.
G 0Rv ) FvkBTχ 0T (C 0Rv + T 0Rv)
(A4)
(A5)
Then, eq A3 can be written as
V h 0R ) kBTχ 0T (1 - FvC 0Rv - FvT 0Rv)
(A6)
The term T 0Rv is a product of the solute charge and the pure solvent properies, such as the dielectric constant. For ions whose charges are equal in magnitude but opposite in sign, the relation T 0+v ) - T 0-v holds and the T 0Rv’s do not contribute to the observed partial molar volume of the salt, V h 02. Most 0 experiments determine V h 2 and separate it into individual ion contributions on the basis of various assumptions that do not appear to take the terms T 0Rv’s into account. One experimental technique, which is based on the ultrasonic vibration potentials, allows the determination of the individual ion contributions on the basis of the knowledge of V h 02. KP obtained the expression for the individual ion partial molar volume that can be identified as that determined by the ultrasonic vibration potentials as follows: 0 0 V h 0* R ) kBTχ T(1 - FvC Rv)
(A7)
This formula can be regarded as a theoretical microscopic formula for the partial molar volume of an individual ion at infinite dilution. This expression can be formally obtained from eqs A3 and A5 by omitting the term T 0Rv. Now we seek an expression for the partial molar volumes of individual ions at infinite dilution based on the ex-RISM theory following the KP procedure: we first obtain the expression that corresponds to eq A5 and then omit the term that corresponds to T 0Rv to obtain the expression for the individual ion partial molar volume in terms of the site-site correlation functions. At this point, however, it should be noted that all the equations from eq A1 to eq A7 contain subscripts which specify species, not the particular sites on a molecule. In light of eq A7, the resultant expression for the partial molar volume based on the site-site formalism should not be dependent on the particular choice of the site on the solvent species. Indeed, it has been shown that this feature holds for the partial molar volumes of neutral molecular solutes39 and salts.40 As we will see below, the ex-RISM-based formula for the partial molar volumes does satisfy this condition if we closely follow the KP procedure. Let us investigate the small k behavior of site-site correlation functions. We rewrite eq 2, the RISM equation for the u-v correlation functions at infinite dilution in k-space as
˜ v(k) + Fvh˜ vv(k)] h˜ uv(k) ) c˜ uv(k) [w
1 + G 0vv - G 0Rv Fu
) kBTχ 0T - G 0Rv
Expanding h˜ 000 00;Rv(k) in powers of k analytically, KP determined G 0Rv as follows:
(A8)
Note that the small k expansion of the correlation functions of the solvent can be written as41
˜ v(0) + k2w ˜ v(2) + ... w ˜ v(k) ) w
(A9)
h˜ vv(k) ) h˜ vv(0) + k2h˜ vv(2) + ...
(A10)
For charged sites, the longest part of the site-site interaction is the Coulomb potential resulting in a small k divergence; hence, as in eq 4, we express the u-v direct correlation functions as
c˜ uv(k) ) c˜ *uv(k) - 4πβzRZv/k2
(A11)
denotes the isothermal compressibility of the pure where zR is the charge of the ion R and Zv denotes the row
3220 J. Phys. Chem. B, Vol. 101, No. 16, 1997
Chong and Hirata
vector that consists of the partial charges on the solvent sites. For example, in the case where the solvent is a three-site water model, Zv becomes
can also be identified as that determined from the ultrasonic vibration potentials.
Zv ) (zO zH zH)
(1) Marcus, R. A. Annu. ReV. Phys. Chem. 1964, 15, 155. (2) Gilson, M. K.; Honig, B. Nature 1987, 330, 84. (3) Robinson, R. A.; Stokes, R. H. Electrolyte Solutions; Butterworths: London, 1965. (4) Gurney, R. W. Ionic Processes in Solution; McGraw Hill: New York, 1953. (5) Samoilov, O. Y. Structure of Aqueous Electrolyte Solutions and the Hydration of Ions; Consultants Bureau, Enterpr. Inc.: New York, 1965. (6) Franks, F., Ed. Water, a ComprehensiVe Treatise; Plenum: New York, 1973; Vol. III. (7) Born, M. Z. Phys. 1920, 1, 45. (8) Najbar, J.; Tachiya, M. J. Phys. Chem. 1994, 98, 199. (9) Eisenberg, D.; Kauzmann, W. The Structure and Properties of Water; Oxford: London, 1969. (10) Frank, H. S.; Wen, W.-Y. Discuss. Faraday Soc. 1957, 44, 133. (11) Samoilov, O. Y. Discuss. Faraday Soc. 1957, 44, 141. (12) Hirata, F.; Rossky, P. J. Chem. Phys. Lett. 1981, 83, 329. (13) Hirata, F.; Pettitt, B. M.; Rossky, P. J. Chem. Phys. 1982, 77, 509. (14) Hirata, F.; Rossky, P. J.; Pettitt, B. M. J. Chem. Phys. 1983, 78, 4133. (15) Yu, H.-A.; Karplus, M. J. Chem. Phys. 1988, 89, 2366. (16) Matubayasi, N. J. Am. Chem. Soc. 1994, 116, 1450. (17) Singer, S. J.; Chandler, D. Mol. Phys. 1985, 55, 622. (18) Zichi, D. A.; Rossky, P. J. J. Chem. Phys. 1986, 84, 1712. (19) Debye, P. J. Chem. Phys. 1933, 1, 13. (20) Zana, R. F.; Yeager, E. J. Phys. Chem. 1967, 71, 521. (21) McQuarrie, D. A. Statistical Mechanics; Harper and Row: New York, 1976. (22) Matubayasi, N.; Reed, L. H.; Levy, R. M. J. Phys. Chem. 1994, 98, 10640. (23) Jorgensen, W. L. J. Am. Chem. Soc. 1981, 103, 335. (24) Pettitt, B. M.; Rossky, P. J. J. Chem. Phys. 1982, 77, 1451. (25) Pettitt, B. M.; Rossky, P. J. J. Chem. Phys. 1986, 84, 5836. (26) Geiger, A. Ber. BunsenGes. Phys. Chem. 1981, 85, 52. (27) Hirata, F.; Redfern, P.; Levy, R. Int. J. Quantum Chem. 1988, 15, 179. (28) Yu, H.-A.; Roux, B.; Karplus, M. J. Chem. Phys. 1990, 92, 5020. (29) Latimer, W. M.; Pitzer, K. S.; Slansky, C. M. J. Chem. Phys. 1939, 7, 108. (30) Lazaridis, T.; Paulaitis, M. E. J. Phys. Chem. 1992, 96, 3847. (31) Millero, F. J. Chem. ReV. 1971, 71, 147. (32) Marcus, Y. J. Chem. Soc. Faraday Trans. 1993, 89, 713. (33) Hirata, F.; Arakawa, K. Bull. Chem. Soc. Jpn. 1973, 46, 3367. (34) Kusalik, P. G.; Patey, G. N. J. Chem. Phys. 1988, 89, 5843. (35) The referee made the following interesting comment in this connection. In the Born theory of solvation one can eliminate the ion radius by taking ratios of thermodynamic quantities. In particular, the Born theory predicts
(A12)
where zO and zH denote the partial charges of the oxygen and hydrogen sites of the water molecule, respectively. Now that the divergent term is removed, the small k expansion for the modified direct correlation functions defined in eq A11 can be41
c˜ *uv(k) ) c˜ *uv(0) + k2c˜ *uv(2) + ...
(A13)
Using eq A11, we rewrite eq A8 as
˜ v(k) + Fvh˜ vv(k)] h˜ uv(k) ) c˜ *uv(k) [w 4πβzRZv[w ˜ v(k) + Fvh˜ vv(k)]/k2 (A14) Utilizing the well-known compressibility relation,42 it follows
˜ v(0) w ˜ v(0) + Fvh˜ vv(0) ) [1 + Fvh˜ vv(0)] w ) FvkBTχ 0Tw ˜ v(0)
(A15)
where h˜ (0)vv does not depend on the particular choice of the solvent site. Note also that
˜ v(0) ) 0 Zvw
(A16)
due to the charge neutrality of the solvent molecule. Then, from eq A14, we have
˜ v(0) h˜ uv(0) ) FvkBTχ 0Tc˜ *uv(0) w 4πβzRZv [w ˜ v(2) + Fvh˜ vv(2)] (A17) Note that the matrix element of the first term in eq A17 is independent of any particular solvent site label; rather it is a sum over all the solvent sites, and we define
C 0* Rv ≡
∑s c˜ *uv Rs (0)
(A18)
References and Notes
T∆scharge ∆V h
The second term in eq A17 has a form expressed as a product of the solute ionic charge and the pure solvent quantities, Zv and w ˜ v(2) + Fvh˜ vv(2). The latter quantity is related to the dielectric constant of the solvent.41,43,44 The species-labeled GRv in eq A1 should be related to the 40 By comparing eq A5 with site-labeled h˜ uv Rs(0) in some way. eqs A17 and A18, one can see that the first term in eq A17 0 (C 0* Rv) corresponds to C Rv in eq A5 and that the second term in eq A17 is intimately related to T 0Rs in eq A5. Then, in light of the discussion to derive eq A7 from eq A6, we finally have the expression for the partial molar volume of an individual ion 0 0* V h 0* R ) kBTχ T (1 - FvC Rv)
(A19)
In view of eq A18, this equation does not include any particular site on the solvent species, rather a sum over all the interaction sites on the solvent, and can be regarded as the interaction-site formula for the partial molar volume of an individual ion, which
() ∂p
)T
∂T
κ
where κ is the dielectric constant of the solvent. The experimental value of T(∂p/∂T)κ for water at room temperature and atmospheric pressure is 0.714 [kcal/cm3]. From Tables 3 and 4 in the present paper, one finds values of T∆scharge/∆V h as 0.656, 0.79, 1.2, 1.47, and 0.74 for Li+, Na+, K+, F-, and Cl-, respectively. Agreement with the Born theory prediction is surprisingly good for Li+, Na+, and Cl-. These results suggest that solvent electrostriction and molecular reorientation actually affect T∆scharge and ∆V h in much the same way, since that is what is implied by the above relation. Although the answer to the question of why K+ and F- do not show such agreement is not known at this moment, the comment may be suggesting a way of studying the limitation of the continuum model for ion solvation. (36) Impey, R. W.; Madden, P. A.; McDonald, I. R. J. Phys. Chem. 1983, 87, 5071. (37) Kusalik, P. G.; Patey, G. N. J. Chem. Phys. 1987, 86, 5110. (38) Kirkwood, J. G.; Buff, F. P. J. Chem. Phys. 1951, 19, 774. (39) Ohba, M.; Kawaizumi, F.; Nomura, H. J. Phys. Chem. 1992, 96, 5129. (40) Perkyns, J.; Pettitt, B. M. J. Chem. Phys. 1992, 97, 7656. (41) Sullivan, D. E.; Gray, C. G. Mol. Phys. 1981, 42, 443. (42) Lowden, L. J.; Chandler, D. J. Chem. Phys. 1973, 59, 6587. (43) Hoye, S.; Stell, G. J. Chem. Phys. 1976, 65, 18. (44) Chandler, D. J. Chem. Phys. 1977, 67, 1113.