Scientific paper A Microscopic Theory of Solvation of Monoions Tomaž Urbic,1'* Christopher J. Fennell,2 and Ken A. Dill2 1 Faculty of Chemistry and Chemical Technology, University of Ljubljana, Askerceva 5, SI-1000, Sloveni 2 Laufer Center for Physical and Quantitative Biology, Stony Brook University, 5252 SUNY, Stony Brook, NY 11794-5252, USA * Corresponding author: E-mail: tomaz.urbic@fkkt.uni-lj.si Received: 01-02-2012 Dedicated to Prof. Dr. Gorazd Vesnaver on the occasion of his 70h birthday Abstract The solvation free energies of ions are often computed using continuum theories, like the Born model. The Born model has the disadvantages that to fit experimental data, ionic radii are taken as adjustable parameters and you need to know the dielectric constant. We present here a more microscopic treatment for computing the free energies of ion solvation in water. Like the Born model, it gives an expression that is simple and can be computed quickly, but unlike the Born model, it uses true ionic radii and does not require inputting a dielectric constant. We show that the present model gives predictions for the free energies of transfer of alkali and halide ions into water that are in excellent agreement with recent experimentally derived estimates. Keywords: Free energy, ionic solvation, mean field theory 1. Introduction The free energies of ion solvation in water have been widely studied theoretically and experimentally.1-9 A simple classical way to compute ion-solvation free energies is using the Born model which gives the excess solvation free energy of a spherical ion of radius R and charge q as (1) By treating the solvent as a dielectric continuum with dielectric constant e0, the Born model has the advantage that it provides a closed-form expression, giving both physical insight and rapid computations. The drawbacks of the Born model are that: (1) to agree with experimental data, the ion radii are taken as adjustable parameters,510-12 and (2) in microscopically heterogeneous environments, the dielectric constant is rarely known, and is another free parameter. Related to these problems are the need, when developing atomically detailed forcefields, to adjust ion parameters, because refining the Lennard-Jones diameters is often a key step in getting computed ion-solvation pro- perties to match experiments.13-17 To avoid optimizing radii parameters, modeling efforts have introduced microscopic response functions in terms of fields of dipoles, but these often require several free parameters in order to reproduce experimental aqueous solvation properties.18-21 We develop here a model for computing the free energies of ion solvation in water. Our treatment involves simple assumptions and has the advantage that it does not require adjusting ion radii or inputting a dielectric constant. Rather than assuming the solvent is a dielectric continuum, here we treat the solvent by summing over the interactions of discrete water dipoles (see Fig. 1). a) N Figure 1. An illustration of averaged water dipoles in different solvation shells for (a) a cation and (b) an anion. 2. Theory In our approach we sum up the Coulombic interactions: (i) of the ionic charge interacting with dipole representations (cd) of water molecules in the surrounding solvent, and (ii) of all pairs of dipole-dipole interactions (dd) among these representations of water molecules. For each water site, the charge-dipole electrostatic energy Ucd is (2) where |i is the dipole moment of a water molecule at a given position in space, averaged over orientations, and E is the electric field acting on the water dipole at that position from the ionic charge. Fig. 1 shows our assumption that water molecules are arranged in spherically symmetric shells around the ion. Because we are assuming each water is already properly configurationally averaged, the per water charge-dipole component of the electrostatic Gibbs free energy of charging is (3) Here we assumed that the averaged induced solvent dipole is linearly dependent upon the charge of the ion, which is good assumption for all ions larger than lithium. We also approximated that we can split contributions to total free energy to contributions between pairs of particles and we neglected all three and more body interactions. To evaluate these free energy contributions, we need the configurationally averaged water dipole as a function of the field, For this purpose, we use a universal curve that we previously obtained by computer simulations of charged spherical solutes placed in the TIP3P model of water.22 Fig. 2 shows the configurationally averaged dipole moment for a TIP3P water molecule at the surface of negatively and positively charged spheres. Note that the curve is not symmetric with respect to the sign of the elec- -6e+10 -4e+10 -2e+10 0 2e+10 4e+10 6e+10 Electric Field (N/C) Figure 2. Water's average dipole charge separation vs. applied field E, determined from previous computer simulations of charged spheres in TIP3P water.22 Water's dipolar asymmetry is reflected in the different shapes of this curve for positive and negative fields. a m Midpoint of Dipole Charges (location of the point dipole) Figure 3. Location of averaged water dipoles near cations and anions. tric field. This is because water's dipole is asymmetrical with respect to water's center of mass (see Fig. 3). Water's dipolar asymmetry has long been known,10 and is usually handled in implicit-solvent models by treating the ionic radii as adjustable parameters and scaling the anion radii differently than the cation radii. Next, we compute the water dipole-dipole interactions using where Udd = 4K£o (4) (5) where | and are the averaged dipole vectors of the water molecules of interest, r is the separation between the two dipoles and u is the unit vector between the two interacting dipoles. Two water dipoles in a shell around an ion will repel each other, so this energy is positive. Having now obtained expressions for how each water molecule interacts with the ion and witho ther waters, we now sum over all the solvent-water molecules. We divide the space around the ion into shells as shown in Fig. 2. As indicated in the figure, the averaged water dipole moment will diminish with distance from the ion. It can be shown that dipoles do not interact across from one shell to another because each shell is closed, so the total electrostatic free energy of solvation is simply a sum over all M solvation shells, (6) Every shell contribution has two parts: (1) the free energy of the ion interacting with all the dipoles in the shell and (2) the free energy of each averaged water-water dipole interaction, summed over all pairs of dipoles in the shell, AG,- = NiAGcM + Ni(Ni~l)(AGdd(i))- (7) = 0 where Ni is average number of water molecules in i-th solvation shell. Note that the dd contribution is represented as an average because interactions between dipole pairs in a shell will come from a double sum over dipoles with different separations and relative orientations, while the cd term is identical for each water in the shell. Ni can be approximately calculated as the ratio between shell surface area and water molecule cross-section, (8) where ri is the distance from the center of the ion to the center of the i-th solvation shell and rw is the radius of a water molecule. Water is a soft core molecule and in order to approximate the hard core radius we use the Barker-Henderson thermodynamic perturbation theory23 2rw= Pil J o (9) Here, aw = 3:15061 A for TIP3P water (the water model used to derive the relation in Fig. 2) and uLJ is the Lennard-Jones component of the TIP3P water-water potential.24 The radius of the first solvation shell is calculated as the Pauling radius25 of an ion (rp) plus the radius of a water molecule (10) To obtain the radii r( of more distant shells, we assume that waters in close proximity to the ions form highly packed structures, so the distance between two solvation shells is equal to the height of a tetrahedron formed by water molecules 2rwy/6 n = n-1 + ■ 3 (11) consistent with computer simulations of the TIP3P mo-del.24,26 We calculate the dipole-dipole free energy contribution as an average of the free energy contribution over the relative positions of a second water molecule in the same shell without overlapping with the first water molecule positioned on the z-axis (AG,w(i)> = —J A GdddSl. (12) This integral can be calculated analytically. Given that all water molecules in a given shell have the same average dipole moment, ^ = p2 = we get (AG,/¿(f)} = 1 i + cos e. 4tt%Wîrf VI-cos0i (13) where is the closest angle the two water molecules can have without overlapping, (14) In this calculation we neglect the water-water correlation in the same solvation shell. Fig. 3 shows how we account for water's dipolar asymmetry. First, we calculate the distance between positive and negative charges in the averaged dipole from the curves in Fig. 2.22 For cations, we locate the negative charge on the oxygen atom at the center of the water molecule with the positive charge placed this distance further from the ion. The point dipole is located at the midpoint between these charges. For anions, we place the positive dipolar charge on the hydrogen atom with the negative charge at the calculated distance further from the ion. The point dipole was again placed at the midpoint of these charges. The explanation for this is that water is always pointing one hydrogen atom toward anion while second is freely rotatable around OH axis as shown on Fig. 4. All contributions pointing away from OH axis are being canceled when doing the average with respect to orientation. Figure 4. Orientational averaging with respect to OH axis gives us averaged dipole used in our calculations for anion. Finally, we note that the solvation of an ion is not purely electrostatic. To get the total free energy of solva-tion can be approximately calculated as sum of the free energy of transfer of anon-polar solute and the total electrostatic term, to get where we take AG in kcal/mol from:27 AGnp = 0.00542 *X4 + 0.92. SA is the surface area in A2 of a sphere of radius r1. (15) 3. Results and Discussion We tested our model by calculating the solvation free energies of the alkali and halide ions. For these calculations, our cutoff is M = 3 solvation shells, because beyond that point the electrostatic energies become smaller than kBT, where kB is Boltzmann's constant and T is temperature. Table 1 shows our computed electrostatic and non-polar contributions to free energy of solvation, as well as their resulting total. We compare our computed ion-solvation free energies to recent experimentally derived values.9 The agreement is remarkably good given the simplicity of the model. Fig. 5 shows that the model captures the trends of the solvation free energies with ion size, and captures the asymmetry between cations and anions. The anions are slightly over-favored relative to these experimentally derived quantities, though they are nearly identical to those derived by Marcus.2 Table 1. The Pauling radii in A and the calculated electrostatic, non-polar, and total solvation free energies in kcal/mol for the alkali and halide ion series alongside values derived from experimental thermodynamic data9 Ion rp -AGe -AG np -AG AGexpt Li+ 0.60 124.4 -1.2 123.2 128.4 Na+ 0.95 107.8 -1.3 106.5 103.2 K+ 1.33 90.8 -1.5 89.3 86.0 Rb+ 1.48 84.5 -1.5 83.0 80.6 Cs+ 1.65 77.9 -1.6 76.3 75.1 F- 1.36 113.4 -1.5 111.9 104.4 Cl- 1.81 85.5 -1.7 83.8 74.5 Br- 1.98 76.6 -1.7 74.9 68.3 I- 2.16 68.0 -1.8 66.2 59.9 One type of insight that we can get from the model, which is not obtainable from Born-like models, is the relative contribution to the solvation free energy from each water shell. Table 2 shows that for cations about 70% of the electrostatic contribution comes from first-shell waters, while for anions, about 85% is due to first-shell waters. This is consistent with explicit-simulation results of the hydration of solutes with strong electric fields.28 It also shows that while the second- and third- shell contributions are much smaller, they should not be neglected if one hopes to obtain accurate predictions of ion solvation free energies. Table 2. Percent contributions of first, second, and third shells of waters to the total electrostatic free energy of solvation Ion 1st shell % 2nd shell % 3rd shell % Li+ 70 22 8 Na+ 70 21 9 K+ 69 22 9 Rb+ 69 22 9 Cs+ 68 22 10 F- 83 15 2 Cl- 85 14 1 Br- 86 13 1 I- 86 13 1 Figure 5. The calculated and experimental9 AG values versus rp for the alkali (blue) and halide (red) ions. 4. Conclusions We presented here a more microscopic treatment for computing the free energies of ion solvation in water. The theory gives an expression that is simple and can be computed quickly, but unlike the Born model, it uses true ionic radii and does not require inputting a dielectric constant. We show that the present model gives predictions for the free energies of transfer of alkali and halide ions into water that are in excellent agreement with recent experimentally derived estimates. 5. Acknowledgments We appreciate the support of the Slovenian Research Agency (P1 0103-0201) and NIH Grant GM063592. 6. References 1. Webb, T. J. J. Am. Chem. Soc. 1926, 48, 2589-2603. 2. Marcus, Y. Biophys. Chem. 1994, 51, 111-127. 3. Hummer, G.; Pratt, L. R.; Garcia, A. E. J. Phys. Chem. 1996, 100, 1206-1215. 4. Tissandier, M. D.; Cowen, K. A.; Feng, W. Y.; Gundlach, E.; Cohen, M. H.; Earhart, A. D.; Coe, J. V. J. Phys. Chem. A 1998, 102, 7787-7794. 5. Schmid, R.; Miah, A. M.; Sapunov, V. N. Phys. Chem. Chem. Phys. 2000, 2, 97-102. 6. Grossfield, A.; Ren, P.; Ponder, J. W. J. Am. Chem. Soc. 2003, 125, 15671-15682. 7. Asthagiri, D.; Pratt, L. R.; Ashbaugh, H. S. J. Chem. Phys. 2003, 119, 2702-2708. 8. Rajamani, S.; Ghosh, T.; Garde, S. J. Chem. Phys. 2004, 120, 4457-4466. 9. Kelly, C. P.; Cramer, C. J.; Truhlar, D. G. J. Phys. Chem. B 2006, 110, 16066-16081. 10. Latimer, W. M.; Pitzer, K. S.; Slansky, C. M. J. Chem. Phys. 1939, 7, 108-111. 11. Rashin, A. A.; Honig, B. J. Phys. Chem. 1985, 89, 55885593. 12. Babu, C. S.; Lim, C. J. Phys. Chem. B 1999, 103, 79587968. 13. Âqvist, J. J. Phys. Chem. 1990, 94, 8021-8024. 14. Jensen, K. P.; Jorgensen, W. L. J. Chem. Theory Comput. 2006, 2, 1499-1509. 15. Joung, I. S.; T. E. Cheatham, I. J. Phys. Chem. B 2008, 112, 9020-9041. 16. Fyta, M.; Kalcher, I.; Dzubiella, J.; Vrbka, L.; Netz, R. R. J. Chem. Phys. 2010, 132, 024911. 17. Reif, M. M.; Hünenberger, P. H. J. Chem. Phys. 2011, 134, 144104. 18. Wertheim, M. S. J. Chem. Phys. 1971, 55, 4291-4298. 19. Florian, J.; Warshel, A. J. Phys. Chem. B 1999, 103, 10282-10288. 20. Jha, A. K.; Freed, K. F. J. Chem. Phys. 2008, 128, 034501. 21. Koehl, P.; Orland, H.; Delarue, M. J. Phys. Chem. B 2009, 113, 5694-5697. 22. Fennell, C. J.; Kehoe, C. W.; Dill, K. A. Proc. Natl. Acad. Sci. USA 2011, 108, 3234-3239. 23. Hansen, J. P.; McDonald, I. R. Theory of Simple Liquids; Academic Press: London, 1986. 24. Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem.Phys. 1983, 79, 926-935. 25. Pauling, L. J. Am. Chem. Soc. 1927, 49, 765-790. 26. Mark, P.; Nilsson, L. J. Phys. Chem. A 2001, 105, 99549960. 27. Rizzo, R. C.; Aynechi, T.; Case, D. A.; Kuntz, I. D. J. Chem. Theory Comput. 2006, 2, 128-139. 28. Mobley, D. L.; Barber II, A. E.; Fennell, C. J.; Dill, K. A. J. Phys. Chem. B 2008, 112, 2405-2414. Povzetek Solvatacijsko Gibbsovo prosto energijo velikokrat ocenimo s teorijami, kot je Bornov model. Ti modeli imajo pomanjkljivost, da moramo uporabiti prilagodljive parametre za ujemanje rezultatov z eksperimentom. Kot vhodni podatek pa potrebujemo tudi dielektrično konstanto topila. V prispevku bomo predstavili izračun Gibbsove proste energije na mikroskopskem nivoju. Izračun je preprost in računsko nezahteven. Uporabimo pa prave ionske radije in ne potrebujemo dielektrične konstante. Model zelo dobro napove proste energije transferja alkilnih in halidnih ionov v vodo.