Ion Pair Potentials of Mean Force in the Interfacial Environment.

Fig. 1 summarizes how the effective interaction between oppositely charged ions is affected near an interface. Specifically, Fig. 1B shows the potential of mean force (PMF), W(r|z), between a M+−M− ion pair in bulk water and at six different z locations, sampling the region from 1 nm into the liquid side (z=1) to 1 nm on the vapor side of the interface (z=−1). The W(r|z) profiles in Fig. 1B have been set to zero at r=1.2 nm for reference. In bulk water, the M+−M− PMF shows features typical of solvent-mediated interactions (19, 20), with a contact minimum at r=0.39 nm, a solvent-separated minimum at r=0.55 nm, and a desolvation barrier between the two minima. W(r|z) becomes negligible in magnitude for ion separations larger than 1 nm. The effective interaction, W(r|z), between these small model ions, M+ and M−, at z=1 nm is already similar to that in bulk water. However, as the ion pair is moved toward the interface, i.e., z<1, W(r|z) becomes highly sensitive to the location of the ion pair relative to the interface. There is a monotonic decrease in the contact minimum, suggesting stronger water-mediated attraction between the M+−M− ion pair in the interfacial region. This is consistent with previous molecular dynamics and continuum electrostatic calculations that have shown an increased stability of a NaCl ion pair at the water–1,2-dichloroethane liquid–liquid interface (21). As the ion pair is drawn further toward the vapor side (z≤0), the effective ion–ion attraction is further strengthened, and W(r|z) displays long-range correlations spanning the box length. Such long-ranged correlations are not an artifact of electrostatic calculations or the simulation methods used. By performing electrostatic calculations using the particle mesh Ewald (PME) as well as the slab Ewald summation (22), by repeating calculations for a larger box size, and by using neutral hydrophilic solutes instead of ions in independent simulations, we show later in this paper (and in SI Text) that these long-ranged correlations arise from capillary fluctuations and interfacial deformation.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 1.

(A) The schematic shows an oppositely charged M+−M− ion pair separated by r and situated in a plane at a distance z from the interface. The z axis is normal to the plane of the average vapor–liquid interface (z=0), with z>0 being the liquid side. (B) PMF between M+ and M− ions calculated in various z planes. We set W(r|z)=0 in each z plane at r = 1.2 nm for reference. (C) The 3D W(r,z) PMF landscape for the M+−M− ion pair in which an ion pair separated by 1.2 nm and located at z=1 is used as a reference state; i.e., W(r=1.2,z=1)=0. This landscape highlights the work needed to move an ion pair from bulk water to the interfacial region.

The W(r|z) profiles in Fig. 1B are referenced to a convenient yet arbitrary choice of r=1.2 nm. We have also obtained the full 3D PMF landscape, W(r,z), where only a single state on the liquid side of the interface is used as a reference; i.e., W(r=1.2,z=1)=0 (Fig. 1C). This landscape highlights that work is required to move a M+−M− ion pair from water to the interface and to the vapor side, as reflected in the increase of W(r,z) along z for most r values. Once this work is expended, the effective ion–ion attraction is stronger than that in the bulk. For certain salts this strong attraction at the interface may even drive ion pairs from bulk water to the interface (23). Such an enhancement of ion density at the interface is thought to underlie the “Ray-Jones” effect (24, 25), especially at low bulk concentrations. The enhancement of ion density near the interface coupled with ion-induced quenching of capillary waves may lead to layering of ions in the direction perpendicular to the interface. Recent X-ray reflectivity experiments have shown layering of erbium chloride at the air–water interface (26).

What are the molecular origins of enhanced ion–ion attractions at the interface? A simple dielectric description of the interface would suggest a lower dielectric constant there and, correspondingly, lower screening of ion–ion interactions compared with that in bulk water. The observed stronger attraction between the M+−M− ions is qualitatively consistent with such an expectation. However, for like-charged ions, the dielectric description would suggest a stronger repulsion near the vapor–liquid interface, becoming more pronounced as the ions are drawn to the vapor side. Such enhanced repulsion is not observed in our simulations. Fig. 2 shows the W(r|z) profiles for like-charged M+−M+ and M−−M− ion pairs at various z locations. In bulk water, as expected, the direct ion–ion electrostatic repulsion dominates the short-range part (r<0.5 nm). At intermediate separations, 0.5<r<1nm, W(r|z) profiles for both ion pairs show small oscillations, reflecting the effects of solvent structuring in ion hydration shells. The first minimum in W(r|z) profiles at r≈0.55 nm is not a “contact minimum” but corresponds to a “water-bridged” configuration, in which a few hydration water molecules shared by both ions mediate their interaction. Structural details of such solvent-bridged minimum configurations have been reported previously for like-charged halide ion pairs (27, 28). PMFs for both ion pairs become negligible for r>1nm in bulk water.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 2.

(A and B) PMFs between like charges (A) M⁺ – M⁺ and (B) M⁻ – M⁻ at various z locations in the interfacial region. For reference, we set the PMF profile in each z plane to zero at r=1.2nm. The effective repulsion between like-charged ions is not enhanced. In fact, ion–ion interactions can become attractive on the vapor side as indicated by PMFs for z<0.

In contrast to the expectations of continuum dielectric models, the like-charged ion–ion PMFs display an overall strengthening of their effective attraction near the interface and especially on the vapor side. The stronger preference for water-bridged configurations is reflected in the deep minimum at r≈0.55 nm, in the M−−M− PMF at z=−1. Interestingly, some oscillations in the PMF are retained even at z=−1, suggesting that ions are partially hydrated on the vapor side of the interface.

Ion Solvation Thermodynamics.

To understand the above counterintuitive results, we first focus on single-ion solvation near the interface. We define ΔG(z) as the free energy or the reversible work required to move a solute from the liquid side of the interface (z=1) to a given z location. Fig. 3 shows ΔG(z) for M+ and M− ions and a neutral M0 solute. Water-to-vapor transfer of the neutral/hydrophobic (i.e., M0) solute is favorable, and its thermodynamics are well understood (16, 29⇓⇓–32). Indeed, ΔG(z) at z=−1 is negative and is favored by a positive entropy change, and its value matches the water-to-vapor transfer free energy calculated using Widom’s insertion method (33). This is expected as the M0 solute is completely dehydrated at z=−1 (Fig. 3 D and F). The ΔG(z) profile displays a minimum slightly to the vapor side of the interface, reflecting favorable energetic interactions of the M0 solute with the aqueous phase (3, 32).

• Download figure
• Open in new tab
• Download powerpoint

Fig. 3.

(A–C) Free energy (ΔG), enthalpy (ΔH), and entropy (−TΔS) of moving a solute from bulk water to a given z location for (A) M0 solute and (B) M+ and (C) M− ions. (D and E) Simulation snapshots show that the M0 solute is fully dehydrated (D), whereas the M+ ion remains partially hydrated (E), deforming the interface when pulled to z=−1. (F) The average number, 〈N〉, of water oxygen centers in the first hydration shell of the ions and of the neutral solute at various z locations. The hydration shell radius is defined by the first minimum in the solute–water radial distribution function in bulk water.

In contrast, ΔG(z) values for both M+ and M− ions increase monotonically to small positive values as the ions are moved from liquid to the interface. As the ions are drawn further into the vapor phase, ΔG(z) increases dramatically, reaching ≈+90 kJ/mol for both M+ and M− ions at z=−1. Details of the ΔG(z) profile depend on the size and the charge density of the ions. Even though there is no minimum in ΔG(z) profiles at the interface for M+ and M− ions studied here, a minimum does exist for ions having lower charge densities as shown by Otten et al. (7).

Despite the large amount of work required to pull the ions away from the liquid phase, as reflected in the large positive ΔG(z) at z=−1, its magnitude is only 20–25% of the liquid-to-vacuum transfer free energy of these ions (≈+380 kJ/mol for the M+ ion and ≈+450 kJ/mol for the M− ion). This observation suggests that the ions are not completely dehydrated at z=−1. Indeed, Fig. 3 E and F shows that the ions are able to maintain most of their hydration water molecules by deforming the vapor–liquid interface. The M− and M+ ions lose only about 2 and 4 molecules, of 8 and 10 hydration shell water molecules in the bulk, respectively. That the M− ion retains a larger fraction of its hydration water is consistent with the more favorable hydration of negative ions compared with positive ions of the same size (34⇓⇓–37). Interfacial deformation following the movement of ions across soft aqueous interfaces has been observed in previous simulation studies (3, 38).

Both interface deformation and partial dehydration of ions contribute to ΔG(z) for z<0. Below we estimate these contributions, using macroscopic interfacial thermodynamics and bulk coordination statistics to illustrate consistency with the physical picture presented above. Pulling an ion to a given z location on the vapor side increases the interfacial area by ΔA(z), and the corresponding free energy change is given by ΔGdef(z)≈γΔA(z), where γ is the water vapor–liquid surface tension. We calculate A(z) by applying the coarse-graining procedure of Willard and Chandler (39) to instantaneous configurations and estimate ΔA(z) as ΔA(z)=〈A(z)〉−〈A(z=1 )〉, where 〈〉 represents an ensemble average with the ion fixed at a given z location. Table 1 shows that for −0.5<z<0, the interface is slightly deformed, and ΔGdef(z) accounts for ∼50–60% of the overall free energy change, ΔG(z), for both the M+ and M− ions. As the ions are pulled further into the vapor phase, the dominant role of interface deformation becomes clear, with ΔGdef(z) contributing as much as 70–80% of ΔG(z).

The contribution to ΔG(z) from partial dehydration of ions can be estimated from the knowledge of bulk coordination statistics of ions (40). If P(Nsh) is the probability of observing Nsh water molecules in an ion’s hydration shell, then the reversible work required to remove δN water molecules from the shell is given approximately by −kBT⁡ln[P(Nsh∗−δN)/Pbulk(Nsh∗)], where Nsh∗ is the most probable number of hydration water molecules. At z=−1, δN = 2 and 4 for M− and M+ ions, respectively. For these δN values, using P(Nsh) for ions in bulk water, we estimate the contribution to ΔG(z) from partial dehydration to be ∼10 kJ/mol for both ions. Together, our estimates of contributions from interface deformation and partial dehydration account for ∼90% of ΔG(z) at z=−1.

Resolving ΔG(z) into entropy and enthalpy provides further insights into the interplay of different factors—ion–water vs. water–water interactions, water structure, capillary pinning, and interface deformation—in governing the ion solvation thermodynamics near the interface (Fig. 3). The entropy contribution (−TΔS) for moving the ions from bulk water (z=1) to the interface (z=0) is unfavorable (positive) for both M+ and M− ions. Otten et al. showed that this unfavorable entropy arises from pinning of the capillary fluctuations by ions at the interface (7). In contrast, the enthalpy contribution (ΔH) is favorable (negative) and displays a minimum near the interface, resulting from the balance of ion–water and water–water interactions in bulk and at the interface (7). As the ions are moved from z=0 to z=−1, ΔH increases sharply by more than 100 kJ/mol, and −TΔS decreases by ∼40–50 kJ/mol. Not only the signs of changes in ΔH (unfavorable) and −TΔS (favorable) but also their relative magnitudes are consistent with estimates from the temperature dependence of water vapor–liquid surface tension (17, 41). Creating interfacial area is enthalpically unfavorable and entropically favorable, with enthalpy dominating entropy by a factor of approximately 2. Our results together with those of Otten et al. (7) suggest that a favorable minimum in ΔH and an unfavorable maximum in −TΔS may be general features of ion solvation at aqueous interfaces. Whether ΔG(z) displays a minimum or not depends on their relative magnitudes, which in turn are influenced by ion size and charge density.

Thermodynamics of Ion–Ion Interactions.

The W(r|z) profiles in Figs. 1 and 2 characterize the distance-dependent interaction between the ions in a given z plane. For simplicity, we define ΔWa(z) as a measure of their association strength in that plane,ΔWa(z)=W(ra|z)−W(rd|z),[1]where ra and rd are the ion–ion separations in the associated and dissociated states, respectively, with more negative values of ΔWa(z) indicating stronger association. We use rd=1.2 nm for all of the ion pairs. For the M+−M− pair, we use ra=0.39 nm, corresponding to the contact minimum in the PMF. For the like-charged pairs, we use the solvent-bridged minimum as the associated state, with ra = 0.55 nm for M+−M+ and ra = 0.53 nm for M−−M− pairs, respectively.

Fig. 4A shows that the M+−M− ion pair association is favorable in water [ΔWa(z) < 0 at z = 1]. Both enthalpy (ΔHa(z)) and entropy (−TΔSa(z)) contributions to ΔWa(z) are negative (favorable) and have been shown to arise from direct ion–ion interactions and the release of restricted hydration shell water molecules upon association, respectively (42, 43). In contrast, the strength of association of like-charged ions into a solvent-bridged state is small (ΔWa(z)∼0) in bulk water, resulting from an approximate enthalpy–entropy compensation. This is consistent with the results of Zangi, who showed that for halide ion pairs, the enthalpy of solvent-bridged configurations is favorable and entropy is unfavorable (27). For ions with higher charge densities, the favorable enthalpic contribution can dominate, making the solvent-bridged state of the like-charged ions the most favorable configuration (27).

• Download figure
• Open in new tab
• Download powerpoint

Fig. 4.

(A–C) Enthalpy (ΔHa) and entropy (−TΔSa) contributions to the free energy of association (ΔWa, as defined by Eq. 1) for (A) M+−M−, (B) M+−M+, and (C) M−−M− ion pairs at various z locations. Inset in A focuses on the z locations that are on the liquid side of the interface.

Ion association at the interface is more favorable than in the bulk for all ion pairs as indicated by more negative values of ΔWa(z=0). Fig. 4 shows that ΔHa(z) and −TΔSa(z) values are correspondingly amplified at the interface. The molecular details underlying this amplification are likely subtly different for the like- and unlike-charged ion pairs. For example, the association of a M+−M− ion pair into the contact state is entropically favorable (Fig. 4A, Inset), consistent with the release of entropically restricted hydration water and reduced capillary pinning in the associated state. Whereas the association of M+−M+ and M−−M− pairs into solvent-bridged configurations remains entropically unfavorable (Fig. 4 B and C). Nevertheless, bulk-like trends with amplification are consistent with the well-hydrated state of ions at the interface.

All ion pairs display a strong tendency to associate on the vapor side. This association is dominated by interface deformation as reflected in the resolution of ΔWa(z) into ΔHa(z) and −TΔSa(z). Two separated ions deform the interface and create a larger interfacial area than does an associated pair. The decrease in interfacial area upon association is favored by enthalpy and opposed by entropy. Further, as shown for single-ion solvation in this region, the relative contributions of enthalpy and entropy are similar to those obtained from the temperature dependence of surface tension of water.

The Generality of the Mechanism Underlying Interface-Mediated Attraction Between Solutes on the Vapor Side.

That interface deformation governs the effective interactions on the vapor side of the interface is not unique to ions, but is observed for solutes having strong attractive interactions with water. We illustrate this by calculating the PMF between two neutral M0-like solutes, having ∼10 times the attraction with water compared with the M0 solute (σsolute−water=0.38 nm and ϵsolute−water=7.434 kJ/mol). To eliminate contributions from direct solute–solute attraction to the PMF, our test solutes interact with each other through the purely repulsive Weeks–Chandler–Andersen potential (44).

The W(r|z) profile for these solutes in bulk water shows oscillations with a well-defined solvent-separated minimum stabilized by the highly favorable solute–water interactions at r≈0.80 nm and a direct contact feature but not a contact minimum at r≈0.45 nm (Fig. 5A). The effective solute–solute attraction becomes stronger near the interface. At z=−1, the interface-mediated sticky solute–solute interaction is characterized by a broad ∼60-kJ/mol deep minimum, even though there is no direct solute–solute attraction. The interface-dominated association is further evident in the W(r|z) profiles at the vapor–liquid interface of an ethanol–water solution (Fig. 5B). Ethanol molecules are interfacially active, decrease the surface tension of water from 72 mN/m to 31 mN/m (45), and correspondingly reduce the driving force for solute association.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 5.

(A and B) The PMFs between two M0-like neutral solutes in (A) water and (B) 40% by mass ethanol solution at various z locations highlight the generality of the interface-dominated solute association on the vapor side of the interface. The solutes have strong attractive Lennard–Jones interactions with solvent molecules, but have repulsive Weeks–Chandler–Andersen (44) interactions with each other.

Implications for Peptides.

Many biomolecules and polyelectrolytes contain ionizable groups on their surfaces and are interfacially active (46). Naturally, these molecules will carry with them their ionizable groups to the interface. Our principal result—that solvent-mediated attraction between ions is enhanced and is sensitive to their location relative to the interface—has implications on biomolecular structure and function. To demonstrate this in a model system, we performed simulations of two versions of a Gly–(Leu)₅–Gly peptide—an end-neutral version, ⁰GL₅G⁰, where the N and C termini are capped with CH₃CO^– and CH₃NH₂^– groups, respectively, and an end-charged version, ⁺GL₅G⁻, where the termini are free amino (positive) and carboxyl (negative) groups, respectively. The AMBER-99SB (47) force field was used to represent the peptide. Fig. 6 shows snapshots from simulations of both peptides relative to the instantaneous aqueous interface. The hydrophobic leucine groups prefer to point to the vapor phase and make both the peptides interfacially active. The charged groups of the ⁺GL₅G⁻ peptide are present in the subsurface layer, thus pulling the peptide slightly to the liquid side of the interface.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 6.

The effect of charging the peptide ends on its structure at the interface. (A) Time dependence of the peptide C_α–C_α end-to-end distance, r1−7, over a 20-ns simulation trajectory for the end-neutral and end-charged versions. The representative interfacial conformations of ⁰GL₅G⁰ (α-helical) and ⁺GL₅G⁻ peptides (hairpin turn) obtained using a greedy clustering algorithm (48) are also shown. (B) The PMF along the peptide C_α–C_α end-to-end distance calculated using umbrella sampling simulations. (C) Simulation snapshots of the two versions of the peptide at the water vapor–liquid interface. The peptide backbone (cartoon), the leucine groups (sticks), C and N termini of the peptide (spacefill), and water molecules (oxygen, blue; hydrogen, gray) are shown. The Willard–Chandler instantaneous interface (39) is marked by a yellow mesh.

We performed 20-ns–long simulations of both peptides in bulk water and at the vapor–liquid interface and analyzed their structures, using a greedy clustering algorithm (48). To quantify the effects of charging the ends of the peptide, we calculated the PMF along the C_α–C_α end-to-end distance (r1−7), using umbrella sampling simulations. The end-neutral version, ⁰GL₅G⁰, folds into helical conformations at the interface, with average end-to-end distance of roughly 0.85 nm (Fig. 6A), at which the PMF also shows a clear minimum (Fig. 6B). Fig. 6 highlights the dramatic effect of charging the peptide ends. The ends of ⁺GL₅G⁻ peptide form a stable contact, with the hydrophobic leucine groups pointing toward the vapor phase, leading to a hairpin-turn–like structure at the interface. Correspondingly, the end-to-end PMF shows a deep minimum at ∼0.5 nm, as expected from the enhancement of effective ion–ion attraction in the vicinity of the interface shown in Fig. 1. An alternative perspective, focused on the effect of bringing the end-charged ⁺GL₅G⁻ peptide from bulk water to the interface is presented in the SI Text. The stronger water-mediated attraction between the oppositely charged ends of the peptide near the interface is clear in the differences between bulk and interfacial PMFs.