Molecular Dynamics Simulations of the Structure and

Aug 4, 2009 - ... Andrew J. Masters. The Journal of Physical Chemistry B 2016 120 (23), 5183-5193 .... Published online 4 August 2009. Published in pr...
3 downloads 0 Views 3MB Size
11662

J. Phys. Chem. B 2009, 113, 11662–11671

Molecular Dynamics Simulations of the Structure and Thermodynamics of Carrier-Assisted Uranyl Ion Extraction Manori Jayasinghe† and Thomas L. Beck*,†,‡ Department of Chemistry and Department of Physics, UniVersity of Cincinnati, Cincinnati, Ohio 45221-0172 ReceiVed: April 15, 2009; ReVised Manuscript ReceiVed: June 8, 2009

We present molecular dynamics simulations of interfaces relevant to the selective chemical extraction of uranyl ions from aqueous solution. These molecular-level simulations model ion transfer in the PUREX process and in synthetic, selective membranes. We first present simulations of water/oil interfaces modified by incorporation of tributyl phosphate (TBP) into the oil phase (hexane). A range of concentrations is examined, from a single TBP molecule to values close to those utilized in the PUREX process. The TBP molecules exhibit strong interfacial activity, and the interface broadens relative to the water/oil case with increasing TBP concentrations. Additional structural features, including radial distribution functions and orientational distributions, are examined to elucidate the molecular ordering at the interface; the interface structure changes substantially with increasing TBP concentration. Finally, free-energy profiles are computed for (1) a single TBP molecule and a single uranyl nitrate complex [UO2(NO3)2] across the water/oil interface and (2) a UO2(NO3)2 · TBP2 complex across both water/oil and water/(oil+TBP) interfaces. The UO2(NO3)2 complex is strongly repelled from the water/oil interface, while the UO2(NO3)2 · TBP2 complex exhibits interfacial activity that decreases with increasing TBP concentration. The UO2(NO3)2 · TBP2 complex displays a net free-energy driving force for partitioning into the oil phase that increases with increasing TBP concentration. Introduction Understanding the behavior of uranyl ions in solution is crucial for the development of efficient nuclear waste management tools.1 Experimental results and molecular-level theories for the chemical speciation, the stability of the resulting molecular complexes, the structure and binding modes of active ligands, and the driving forces for partitioning at interfaces will aid the advance of extraction technologies. The well-known PUREX process2-4 is used for the remediation of contaminated water resulting from spent nuclear fuel. This process extracts uranium from a heavily concentrated nitric acid solution to an organic phase that consists of 30% (by volume) tributyl phosphate (TBP) in kerosene. The uranium ions are removed as uranyl nitrate complexes bound to TBP carrier molecules. In addition, a recently developed membrane separation technique5,6 utilizes a 30-50% TBP in oil (dodecane) phase to modify a gold-coated alumina membrane surface. It was observed in that work that the thin hydrophobic selectivity layer transported 100% of the uranium as neutral complexes followed by UO2(CH3COO)3- transport by diffusion across a 60 µm alumina membrane. Uranium exists in aqueous solution as the linear uranyl ion UO22+. This ion possesses very stable uranium-oxygen double bonds, leaving the oxygens largely unreactive. Due to the large atomic number of uranium, relativistic quantum effects are important.7-9 Uranium hybrid orbitals, influenced by relativistic effects, lead to strong chemical bonds to the two oxygens.10,11 The uranium atom in the uranyl complex possesses a substantial positive partial charge,12 and thus active electronegative ligands such as NO3-, sPdO, sCdO, and water bind to the metal in * To whom correspondence should be addressed. E-mail: thomas.beck@ uc.edu. † Department of Chemistry. ‡ Department of Physics.

the equatorial plane, leaving the UdO bond largely unaltered. Simulations,13,14 calculations,15 and experiments16,17 suggest that, in water, UO22+ forms a stable equatorial structure for the first hydration shell that contains five water molecules. An equatorial second shell consists of 10 water molecules coordinated via hydrogen bonds to the five first-shell water molecules. These structures are strongly solvated in water, and extraction thus requires a chelating agent in order to overcome the binding to water.18 Although spectroscopic studies19 and simulations20-22 suggest that the uranyl speciation at high concentrations of nitric acid (3-6 M in the PUREX process) consists of mono-, di-, and trinitrate complexes, extraction studies report the most probable extracted complex as the neutral UO2(NO3)2 · TBP2.22,23 The crystal structure of UO2(NO3)2 · (H2O)224,25 reveals that the two nitrates bind in a bidentate mode in the equatorial plane of the uranyl ion. Both entropic and enthalphic factors contribute to this binding motif in solution.20 In an aqueous environment, the particular coordination mode adopted in such complexes can have a subtle but important effect on their solvation, and ultimately on their extraction properties. In a study of UO2(NO3)2 · (H2O)2 with ab initio simulation, Bu¨hl et al.20 found that nitrate groups are quite flexible and well solvated in water, and the bonding to the uranyl ion is weakened upon solvation. They also noted substantial out-of-plane motions of the nitrate ligands in water. Finally, they observed a conformation with one nitrate ligand bound in a monodentate fashion; this conformer corresponds to a local minimum on the free-energy landscape in water, with a 1.4 kcal/mol transition energy between the bidentate and monodentate conformations. In early classical free-energy simulations, Guilbaud and Wipff26 proposed a new set of force field parameters for UO22+ that accounts for the difference in hydration free energies between UO22+ and Sr2+, as well as the coordination number of 5 in water. They further showed that, with this set of

10.1021/jp903470n CCC: $40.75  2009 American Chemical Society Published on Web 08/04/2009

Carrier-Assisted Uranyl Ion Extraction parameters, NO3- is loosely bound to UO22+ in water and dissociates at low concentration levels. Beudaert et al.27 presented initial simulations of low concentration TBP mixtures at the water/chloroform interface and found TBP surface activity with amphiphilic orientations of the TBP molecules. Using the force fields of ref 26, Baaden et al.28 modeled uranyl ion extraction systems by examining UO2(NO3)2 at a water/ organic(chloroform) interface with mixing/demixing molecular dynamics (MD) simulations. They showed that UO2(NO3)2, with or without constraints to a pseudo-D2h symmetrical structure, is repelled by the water/organic interface and migrates to the bulk aqueous phase. They also observed, for the first time, spontaneous complexation of the uranyl nitrate complex with one or two TBP molecules at the interface. Baaden et al.22 described the behavior of UO2(NO3)2 · TBP2 at water/(organic+TBP) interfaces for two TBP concentration levels in a more extensive mixing/demixing MD study. In this work, both chloroform and supercritical carbon dioxide (SC-CO2) were used for the organic phase. The mixed state was prepared by scaling down the electrostatic interactions by a factor of 100 and heating the system to 700 K. For a lower TBP concentration of 0.4 mol L-1 (78%/22% oil/TBP by volume for the uniform mixture), they observed a rapid separation of water and organic phase with the formation of an interface onto which most TBP molecules adsorbed. A spontaneous formation of 1:1 UO2(NO3)2 · TBP and 1:2 UO2(NO3)2 · TBP2 complexes near the interface was also observed. For a higher TBP concentration of 0.8 mol L-1 (64%/36% oil/TBP by volume prior to mixing), they observed relaxation to a microemulsion phase with neat water “bubbles” surrounded by a mixed organic/ TBP phase. In addition, they found that all of the uranyl nitrates formed 1:2 complexes with TBP at the liquid boundaries. Other simulation studies by Wipff and co-workers have extensively examined the role of nitric acid on the interfacial structure and extraction properties of the interfaces.29-32 In addition, that group has previously explored free-energy profiles for the extraction of the Cs+ ion by calixarene ligands.33 These simulations have led to significant insights into the molecular mechanisms of phase separation and metal ion extraction from aqueous solution. Besides the quantum and classical simulations discussed above, extensive quantum chemical research has been devoted to understanding the molecular interactions of UO22+ and its complexes in aqueous solution; thorough reviews are presented in refs 8 and 12. The studies have utilized methods ranging from the relativistic discrete variational Dirac-Slater approach,7 Hartree-Fock calculations,34 and density functional theory (DFT)12,35,36 to higher levels of electron correlation.12,37 As an example, an upper bound to the aqueous free energy of solvation of the uranyl cation was calculated with density functional theory and MP2 calculations by Gutowski and Dixon.12 Their value of -410 ( 5 kcal/mol is consistent with the best experimental value of -421 ( 15 kcal/mol. Other electronic structure calculations have addressed the possibility of alternative ligand binding motifs to uranium.38 Metal ion separations require recognition to achieve selectivity.5 A molecular-level investigation of the selectivity layer is important to understand the driving forces for extraction and hence to design membrane structures with specific ionic recognition properties. In this paper, we report molecular dynamics (MD) studies of the structure of the selectivity layer and the free-energy driving forces for the complexed uranyl ion. The selectivity layer is modeled as a liquid-liquid interface system39-45 consisting of water/oil and water/(oil+TBP) interfaces. Hexane is used to model the oil phase in our simulations

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11663 in order to mimic the organic phase employed in the PUREX process and in the selective membranes of ref 5. We have performed extensive MD simulations that employ the umbrella sampling technique46 to produce free-energy profiles for UO2(NO3)2 and UO2(NO3)2 · TBP2 across water/oil and water/ (oil+TBP) interfaces in order to better understand the solute partitioning. Further, we have investigated the structure of water/ (oil+TBP) interfaces in the absence of solutes to examine the molecular ordering which drives the extraction. The paper is organized as follows. In the next section, we describe the computational methodology employed in our simulations. We then discuss in detail the results for the structural properties and thermodynamic driving forces for these extraction interfaces. Finally, we present our conclusions and discussion of future research directions. Computational Methods Separate equilibrated aqueous and oil (or oil+TBP) systems were first prepared, and then joined at a liquid-liquid interface, followed by further equilibration runs. Even at the higher TBP concentrations, we observed the formation of a stable, although broad, water/organic interface over 10 ns time scales. Similar to the work of Baaden et al.,22 we did not attempt to model the highly acidic and high ionic strength environment of the aqueous solution present in the PUREX process; we note that the feed solution pH of 4 used in the modified membrane study of Jarvinen and co-workers5 is much higher than that used in the PUREX process, and thus our work is more relevant to that situation (although a high concentration LiNO3 solution was used in that study). Future work will examine the role of concentrated acid29-32 and high ionic strength on free-energy profiles for solute transfer. We note that the 50%/50% volume fraction for the oil/TBP mixture in ref 5 suggests that a stable interface impermeable to ions exists even at these high TBP concentrations, since pH levels and ion concentrations were maintained across the interface. Force Field Parameters. The potential energy of the system is described by a sum of bond, angle, and dihedral deformation energies and pairwise additive 1-6-12 (electrostatic + van der Waals) forces between nonbonded atoms for both intra- and intermolecular interactions. The parameters for UO22+ are the same as those in ref 26, while those for NO3- are the same as those in ref 22. The atomic charge distribution of TBP was taken from ref 22. All of the other parameters used in our simulations are from the CHARMM27 force field.47 Water was represented with the TIP3P model,48 while united-atom hexane was used to model the oil phase.49,50 Hexane was used for the oil phase to enhance diffusion rates (and thus shorten the equlibration times) relative to the dodecane fluid used in the experiments of ref 5. All bonds were constrained using the SHAKE algorithm. In concentrated nitrate solutions, UO22+ is typically bound to two nitrates to form a stable and neutral complex;5,22 we used that neutral stoichiometry in our computations. In our model, the secondary structures used for UO2(NO3)2 and the UO2(NO3)2 · TBP2 complex were constrained with distance harmonic constraints. In the course of the MD simulations, the harmonic constraints were used to keep the two nitrates close to the U atom in the equatorial plane. TBP is a strong complexing agent for UO2(NO3)2 in concentrated solutions and forms UO2(NO3)2 · TBP2 complexes.5,22 We have used distance harmonic constraints between the U atom and O atom of the sPdO bond in TBP in order to bind the two TBP molecules in the equatorial plane. Molecular Dynamics Setup. In building the complete simulation box of the water/oil system, we first created boxes

11664

J. Phys. Chem. B, Vol. 113, No. 34, 2009

of equilibrated water and equilibrated oil and then filled the z < 0 region with water and the z > 0 region with oil. A vacuum region was generated on the distant sides of the water and oil phases. Then, we relaxed and equilibrated the resulting box containing the two joined liquids. The resulting system contained 874 water molecules and 143 hexane molecules (5482 atoms) in a box of dimensions 30 Å × 30 Å × 100 Å in the x, y, and z directions. The dimensions of the water and oil layers are roughly 30 Å × 30 Å × 30 Å, respectively. Thus, the vacuum region is roughly 40 Å in the z direction. Consequently, the simulation box was large enough in the z direction to avoid significant interactions between the free interfaces. Next, we built the water/(oil+TBP) interface system starting with a box of oil and TBP molecules randomly placed at an initial concentration of 70%/30% by volume. The system was equilibrated for 1 ns. An equilibrated water box was then added to the equilibrated (oil+TBP) box followed by further equilibration. A vacuum region was inserted as in the water/oil system described above. The water/(oil+TBP) system consisted of 874 waters, 35 TBPs, and 105 hexane molecules in a 35 Å × 35 Å × 120 Å dimensional box in the x, y, and z directions. Periodic boundary conditions were applied in three dimensions. As discussed in the results section below, all of the TBP molecules in the above system moved to the water/(oil+TBP) interface within 1 ns, leaving no TBP molecules in the bulk oil phase. We thus subsequently prepared a system which generated a stable bulk TBP concentration in the oil phase. In order to generate two equivalent fluid interfaces, the vacuum region was removed from the simulations. Separate water and (oil+TBP) systems were prepared and joined at two interfaces under periodic boundary conditions. Enough TBP molecules were added such that, after extensive equilibration, an average bulk density of about 76%/24% oil/TBP by volume was observed. The system consisted of 1403 water molecules, 560 hexane molecules, and 112 TBP molecules. The resulting system was equilibrated for 10 ns to allow for potential long-time relaxation of the density distributions. The simulations of this system were conducted in the NPT ensemble (constant pressure of 1 atm). Thus, the system volume fluctuated so as to maintain the correct pressure. The Verlet propagator was used to integrate the trajectories in time with a 1 fs time step. Our calculations were performed in the NVT ensemble for the systems with a vacuum layer and the NPT ensemble for the system with two liquid-liquid interfaces. All of the simulations were carried out using the CHARMM package.47 The Lennard-Jones interactions were calculated with a unit-based cutoff for hexane and an oxygenbased cutoff for water at 15 Å. The long-ranged electrostatic interactions were calculated using the Ewald summation method. A proper treatment of long-ranged electrostatic interactions has been shown to be important for aqueous interfacial systems.51 Free-Energy Calculations. The umbrella sampling technique46 was used to enhance the sampling employed to generate the free-energy profiles for solute transport across the interfaces. The free-energy surface A(ξ) along the reaction coordinate ξ (the distance, perpendicular to the interface, between the center of mass of the solute and the center of mass of the water slab) was obtained from a series of simulations. The solute was moved along the z direction in sequentially overlapping windows separated by 2 Å. In each window, the solute was simulated in the presence of an artificial biasing window potential, k(ξ ξmin)2. The chosen range of the reaction coordinate was roughly -20 Å < z < 20 Å along the z axis. The solute was free to move in the x-y plane while being confined by the window

Jayasinghe and Beck potential along the z direction. A force constant of k ) 0.8 kcal/ mol · Å2 was chosen for all windows. This force constant allowed the solute to move across the 2 Å window several times in 1 ns; the biased probability distribution was accumulated for each window. The total simulation time to compute a complete free-energy profile was 28 ns. The weighted histogram analysis method (WHAM)52,53 was used to generate the unbiased probability distribution for the whole range of ξ. The free-energy profile was calculated as follows:

A(ξ) ) -kT ln F(ξ)

(1)

where k is Boltzmann’s constant, T is the temperature, and F(ξ) is the unbiased probability distribution generated by the WHAM process. Radial Distribution Functions. Various radial distribution functions (rdf’s) were computed to examine bulk and interfacial ordering in the 70%/30% (initial concentration) oil+TBP system discussed above. For the OO, OH, and HH rdf’s of water, we considered all water molecules for which z > -2.5 Å as interfacial water and all water molecules for which -15 Å < z < -10 Å as bulk water. For a given water atom in the interfacial region, we computed the number of atoms that are in a spherical shell of thickness dr ) 1 Å and radius r (ranging from 1 to 10 Å). We did this for all interfacial water molecules regardless of their z values. We then divided by the total number of interfacial atoms and by the number of atoms that would have been found inside a spherical shell of thickness dr and radius r in bulk water. Since the TBP molecules are highly surface active (below), we utilized a different geometry to compute interfacial water-TBP and TBP-TBP rdf’s. For these rdf’s, we used cylindrical shells with axes oriented parallel to the interface normal. All TBP molecules (since they reside in the interfacial region) and only interfacial water molecules were considered to examine the interactions at the interface. For a given atom, we counted the number of interfacial atoms in a cylindrical shell at r and with width dr ) 1 Å centered on the distinguished atom. Then, the rdf was averaged and normalized to yield a value of 1 at large distances. Results and Discussion We begin by discussing results concerning structural properties of the water/(oil+TBP) interface in the absence of the uranyl nitrate solutes. As discussed above, we first prepared an interfacial system whose oil phase consisted of a 70%/30% volume ratio of hexane and TBP with the TBP molecules randomly distributed throughout the hexane slab. Within 1 ns, all of the TBP molecules adsorbed to the liquid-liquid interface, exhibiting the strong interfacial activity of TBP. In this section, we will refer to this first system as the low density system. (For comparison, we also examined a single TBP molecule at the oil/water interface to explore the very low TBP concentration limit.) Later, we prepared a bulk system that, after extensive equilibration, displayed a 76%/24% oil/TBP volume ratio in the interior of the oil phase, with two equivalent interfaces in contact with water. We will refer to this system as the bulk system below. On the basis of the peak heights in the number density distributions near the interfaces for the low density and bulk systems (below), it appears that the low density simulations were close to the interface saturation level. Most of our structural analysis of the interfacial ordering was performed for the low density system. Following these structural studies, we examine free-energy profiles for the motion of uranyl nitrate and uranyl

Carrier-Assisted Uranyl Ion Extraction

Figure 1. (a) Density profiles for the low density water/(oil+TBP) interface. All TBP molecules segregated to the water/oil interface region. The water phase is on the left, and the oil phase is on the right. (b) Density profiles for the bulk water/(oil+TBP) interface, where the oil/ TBP bulk phase is roughly 76%/24% oil/TBP by volume. Here, the water phase has been replicated on the right side to show that two equivalent interfaces were created. The middle region is the oil+TBP mixture. (c) Density profiles for the bulk water/(oil+TBP) interfacial system, with a UO2(NO3)2 · TBP2 complex located at z ) 12 Å just inside the oil phase.

nitrate complexed with TBP molecules across the water/oil or water/(oil+TBP) interfaces to gain insights into the driving forces for the extraction. Structure of the Water/(Oil+TBP) Interface. We first present the density profiles for the low density and bulk water/ (oil+TBP) interfacial systems. Figure 1 presents the number density profiles for water, hexane carbons, TBP carbons, and TBP P atoms through the interfacial region for the (a) low density system, (b) bulk systems, and (c) bulk system with UO2(NO3)2 · TBP2 complex located at z ) 12 Å. It is apparent

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11665 from these profiles that TBP is highly interface active, consistent with earlier simulation studies.22,29-31 For the low density case of Figure 1a, all of the TBP molecules moved to the interfacial region during the 1 ns equilibration period and remained there for the rest of the simulations. The total simulation time for Figure 1a was 5 ns, at which point the profiles were wellconverged. The adsorbing TBP molecules significantly broaden the interface compared with the water/oil case. The TBP carbon distribution is quite broad, with some penetration into the water phase. The P atoms of the TBP molecules are relatively narrowly distributed at the interface. The bulk system displayed in Figure 1b also exhibits strong TBP interfacial activity at the two fluid interfaces. The interfacial peak heights for the TBP carbon density at the two interfaces are close in magnitude, suggesting good equilibration in the modeled system; the total simulation time for Figure 1b was 10 ns. We note that the TBP carbon peak maximum for the bulk system is slightly lower in magnitude compared with the low density case, reflecting some collective interaction of the interfacial and bulk TBP molecules. The TBP carbon density profile is not entirely flat in the middle of the oil phase but exhibits some structure. This structure persisted after 10 ns of simulation, so it appears to be stable. The size of a single TBP molecule is on the order of 10 Å, so it is clear that there is substantial solvation of TBP in the bulk-like environment in the middle of the oil phase. The average density of the TBP away from the interface regions is roughly 76%/24% by volume, close to the typical values of 70%/30% used in the PUREX process. Although the addition of TBP substantially broadens the liquid-liquid interfacial region, we observe stable, welldefined interfaces over 38 ns time scales (including the freeenergy calculations). Also, we do not observe any TBP layering on a 10 Å length scale (at the interface) in the low density system, as was observed for a lower density case in ref 30. Figure 1c presents the density profile for a 1 ns simulation of the bulk system with a UO2(NO3)2 · TBP2 complex located at z ) 12 Å. This figure is for reference below when we discuss the free-energy profiles for transfer of the complex across the water/(oil + TBP) interface. Two things are worth noting in this density profile. First, the oscillations in the hexane and TBP profiles in the interior of the bulk phase show that, over short times, the densities fluctuate significantly within the oil+TBP phase. Second, there exists a slight shoulder on the water density profile at the left interface (0 < z < 10 Å). That shoulder is due to an increased water density leading up to the complex that is located at z ) 12 Å on the oil phase side of the interface. This profile suggests that, to some extent, waters are dragged, at least over short distances, along with the complex as it passes from water into the oil+TBP phase. This effect, however, does not appreciably alter the water density profile. As noted above, the density profiles are relatively well-converged after 10 ns of simulation time. Molecular Ordering at the Liquid-Liquid Interface. We next discuss the molecular ordering at the water/(oil+TBP) interface through analysis of various distribution functions. For these results, we focus on the low density case shown in Figure 1a. Figure 2 displays the radial distribution functions for TBP at the interface, suggesting that the TBP molecules are structured in this region. In Figure 3a, we show the probability distribution for the angle θ between the PdO bond vector and the interface normal along the z direction. Although this is a broad distribution, the probability is peaked around 90°. Figure 3b displays the probability distribution for the angle ψ between the interface normal and the intermolecular direction vector between a

11666

J. Phys. Chem. B, Vol. 113, No. 34, 2009

Figure 2. TBP-TBP P-P, P-O, and O-O rdf’s for the low density system. These rdf’s suggest that the TBP interfacial molecules are well structured and interact strongly in the interfacial region.

phosphorus atom on one TBP and the PdO oxygen on other TBP molecules. This distribution is also peaked around 90°. These two distributions imply that the PdO bonds of the TBP molecules form a layered structure which is largely parallel to the interface plane. Further, the probability distribution of the angle φ between the intermolecular PdO bond-vector pairs (Figure 3c) exhibits a correlation peaked at 0° (two pairs are considered only if the P-P distance is less than the first valley of the P-P rdf at 7.5 Å, Figure 2). This figure indicates that the PdO dipoles of TBP molecules in the interfacial region are oriented parallel to each other. We believe that the origin of this PdO orientation is the intermolecular PdO dipole interaction. Our analysis of the TBP phase structure thus implies that the PdO dipoles form chain structures which lie in a slab region parallel to the interface. We also performed the same analysis for the interfacial TBP molecules in the bulk system, and found distributions similar to those presented here for the low density case; the distribution

Jayasinghe and Beck

Figure 3. (a) Probability distribution of the angle θ between the PdO vector and the interface normal along the z coordinate. At z ) 0, the PdO vector is on average perpendicular to the interface normal. This indicates that the PdO vectors are parallel to the interface plane. (b) Probability distribution of the angle ψ between the interface normal and the intermolecular PsO vector. The distribution is peaked near 90°, suggesting that the PdO groups lie in a relatively thin layer. (c) Probability distribution of the angle φ between the intermolecular PdO dipole vectors. The distribution is peaked at zero angle, indicating that the intermolecular PdO vectors are parallel to each other. This suggests chain-like structures for the PdO groups.

corresponding to Figure 3a is somewhat broader but centered around the same angle. We next examine the structural properties of bulk and interfacial water, and the water interactions with interfacial TBP molecules. The changes in the water environment can be observed in changes in the OO, OH, and HH rdf’s for water at the interface. Figure 4a shows a comparison of radial distribution functions of bulk and interface water. The local water structure at the interface, in terms of the shape of the rdf’s, is not significantly changed. The magnitude of the interfacial rdf’s is reduced, however, indicating some of the water density around a distinguished water is reduced, as expected at the interface. The water at the interface is also interacting with the polar groups of the TBP molecule. This is confirmed in Figure 4b, which displays the water-TBP rdf’s. These rdf’s suggest that

Carrier-Assisted Uranyl Ion Extraction

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11667

Figure 5. rdf for C(hexane)-C(TBP butyl groups) showing a regular hydrocarbon chain-chain hydrophobic interaction in the TBP-oil contact region.

Figure 4. (a) Water-water rdf’s for the low density system. The reduction of the height of the first peak for the interfacial waters suggests that the density of coordination shell water molecules is significantly reduced compared with bulk water. These interactions are replaced with water-TBP interactions. (b) Water-TBP rdf’s for the molecules in the interfacial region. The first letters in the labels of the rdf’s are for water, and the second letters are for TBP. The TBP O refers to the oxygen from the PdO group. These rdf’s show that the water-TBP interactions are significant.

there exists strong water-TBP interactions in the interfacial region. The water-TBP rdf’s also imply that individual waters bind tightly to specific TBP molecules, but there is a depletion in the water density at larger distances (5-10 Å). This indicates that the water solvation around TBP involves local clustering due to a complicated competition between water/TBP, TBP/ TBP, and TBP/oil interactions. In terms of the interactions of the hydrophobic parts of the TBP molecules and the oil phase, the rdf for the interactions of the carbons of TBP and the hexane carbons is shown in Figure 5. This rdf is typical for hydrocarbon liquids and shows that there is no special ordering of the hydrophobic components of the system. In order to explore the very low TBP concentration limit, we also performed simulations of a single TBP molecule at the water/oil interface. The PdO vector orientational probability

Figure 6. (a) Probability distribution for the angle between the PdO vector and the interface normal as a function of z for a single TBP molecule at the water/oil interface. This distribution indicates that TBP at very low concentration exhibits a typical amphiphilic orientation with the PdO group pointing into the water phase. (b) Relative free-energy profile for TBP across the water/oil interface for a single TBP molecule exhibiting a strong interface activity.

distribution (relative to the surface normal) as a function of the z coordinate is shown in Figure 6a. The peak at zero angle indicates an amphiphilic orientation of TBP at low concentration, in contrast to the orientations parallel to the interface at the higher TBP concentrations (Figure 3a). It is thus apparent that the nature of the molecular ordering at the interface changes qualitatively with increasing TBP concentration. We also computed the free-energy profile for passage of a single TBP molecule from the water phase to the oil phase (Figure 6b). A deep but broad interfacial well of 6 kcal/mol magnitude demonstrates the TBP interfacial activity. The net drop across the interface is roughly 1.5 kcal/mol. Solute Free-Energy Profiles and Orientations at the Liquid-Liquid Interface. The free-energy profile for UO2(NO3)2 across the water/oil interface and the free-energy profiles for UO2(NO3)2 · TBP2 across the water/oil and water/ (oil+TBP) interfaces were computed. The partition coefficients for these solutes in thermodynamic equilibrium can be calculated as follows:

11668

J. Phys. Chem. B, Vol. 113, No. 34, 2009 oil Kwater ) exp(-∆µex /kT)

Jayasinghe and Beck

(2)

where ∆µex is the free-energy change for transferring the solutes from water to the organic solvents. That free-energy change can be obtained as the difference between the plateau values on the water and oil sides of the free-energy profiles. On the basis of our analysis (below), the relative partition coefficient of UO2(NO3)2 from water to oil is almost zero and that of UO2(NO3)2 · TBP2 is 440. The partition coefficient of UO2(NO3)2TBP2 from water to the oil+TBP phase is increased to 2900 for the low density TBP system, and 35 800 for the bulk system. These partition coefficient values can be expected to have appreciable errors due to limitations of the force field representation, the lack of salt and/or acid in the feed solution, and sampling errors; we report them here to give some indication of changes in the extraction efficiency. The extractability of UO2(NO3)2 is thus enhanced considerably by complexing with two TBP molecules. It is further enhanced by increasing the TBP concentration levels in the organic phase. Figure 7 displays the free-energy profile for UO2(NO3)2 across the water/hexane interface, indicating a 20 kcal/mol free-energy barrier. The uranyl nitrate complex is thus strongly repelled from the water/oil interface, and the partition coefficient is effectively 0. UO2(NO3)2 is a large solute (diameter ≈ 10 Å) with a dipole moment of 2 D in vacuum when minimized with the CHARMM force field using the harmonic constraints employed in the MD simulations. During the MD simulations (with the harmonic constraints to maintain the integrity of the complex), we observed average dipole moments for UO2(NO3)2 of 16.0 D in water and 7.0 D in oil. Thus, the solvent plays a significant role in favoring polarized conformations of UO2(NO3)2, especially in water. From analysis of rdf’s between uranyl nitrate and inspection of snapshots of the system, we also noticed that two waters are always bound in the equatorial plane of the uranyl nitrate complex in the water phase. To further address the solvation effects on the dipole moment of the complex, quantum chemical calculations were performed at the DFT-B3LYP level with the Gaussian 03 code.54 The Stuttgart small core RECPs (core potentials) and the corresponding Stuttgart orbital basis sets for the U atom, and the 6-31G basis set for the N, O, and H atoms, were used for these calculations. The crystal structure24 was taken as the initial configuration, and the structure was optimized in vacuum while maintaining the D2h symmetry of the crystal. Due to the symmetry constraint, the resulting dipole of the complex was 0. Likewise, the dipole was 0 when the polarizable continuum model (PCM) for surrounding water solvent was included under the same symmetry constraints. One of the nitrate groups was then rotated about the symmetry axis of the OdUdO vector to determine the impact of the distortion on the total energy. In vacuum, distortion by 10° led to an energy increase of 0.6 kcal/ mol, while, with the continuum water model, the energy decreased by 4.2 kcal/mol. This qualitatively supports the notion that the condensed phase water stabilizes dipole fluctuations. Due to the artificial harmonic constraints in the simulations, and the limitations of the classical force field, we are not confident of the magnitudes of our computed average dipoles; it is clear, however, that the stabilization of the dipole fluctuations (especially in water) is involved in the large free-energy barrier. In addition, we observed that the first water coordination shell around UO2(NO3)2 is dragged to the oil phase with the complex, and a stable water finger connects the bulk water and the water shell around the complex via a hydrogen-bonding network

Figure 7. Relative free-energy profile for the uranyl nitrate (UO2(NO3)2) complex across the water/oil interface. The water phase is on the left side of the figure.

(Figure 8a). This water finger appears to limit the mobility of the complex, leading to potentially poor sampling in the oil phase. We then computed the free-energy profiles for the UO2(NO3)2 · TBP2 complex across the water/oil and water/(oil+TBP) (both low density and bulk systems) interfaces (Figure 9). The complex exhibits strong interfacial activity at the water/oil interface with a free-energy drop of 10 kcal/mol passing from bulk water to the interface. The net drop of about 4 kcal/mol from water to oil implies a net driving force favoring partitioning of the complex into the oil phase even with no TBP present. The free-energy profile for the low density oil+TBP system exhibits a broader interfacial region and a shallower free-energy well. (Here, we define the magnitude of the well as the height of the barrier to pass from the free-energy minimum into the oil phase.) The net driving force for partitioning increases to roughly 5 kcal/mol, indicating better solvation in the oil+TBP phase relative to pure oil. The trend continues with the bulk oil+TBP system; the interfacial free-energy well is shallower still, while the net free-energy drop into the oil+TBP phase is about 6.3 kcal/mol. We can thus conclude that increasing TBP concentration leads to a shallower (and broader) interfacial freeenergy well for the UO2(NO3)2 · TBP2 complex, but the net driving force for partitioning increases steadily with increasing TBP concentration in the oil phase. Figures 1c and 8b show that, as the complex is dragged into the oil phase from water, some waters remain in contact with the complex; similar behavior has been observed previously for related complexes at water/(chloroform+TBP) interfaces.33 Since the solubility of TBP in bulk water is about 1.5 × 10-3 M, and the uranium concentration in the aqueous feed solution in ref 5 is on the mM level, some complexation of uranyl nitrate by TBP may occur in the bulk phase prior to adsorption at the interface. Still, large fractions of the complexation events are expected to occur in the interfacial region, as has been discussed previously.22 Finally, we examined the orientations of the UO2(NO3)2 · TBP2 complex at the water/oil and low density water/(oil+TBP) systems. The probability distributions of the angle between the OdUdO vector and the interface normal for the two interface cases are shown in parts a and b of Figure 10, respectively. These results show that, at the water/oil interface, the uranyl ion tends to be oriented perpendicular to the interface plane, while, at the water/(oil+TBP) interface, the uranyl ion is oriented parallel to the interface. Again, we see that the addition of TBP qualitatively changes the molecular ordering in the interfacial region. Conclusions This paper has presented molecular dynamics simulations of the interface structure and the thermodynamics of solute transfer

Carrier-Assisted Uranyl Ion Extraction

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11669

Figure 8. (a) Snapshot of the uranyl nitrate complex in the oil phase near the water/oil interface. The water phase is on the left. Waters within a cutoff radius have been highlighted in the figure to emphasize the finger connecting the bulk water phase to the uranyl nitrate complex. The complex was dragged during the free-energy simulations from the water phase to the oil phase. (b) Shapshot of the UO2(NO3)2 · TBP2 complex in the oil+TBP phase near the water/(oil+TBP) interface. A water protrusion is apparent connecting the water phase to the complex embedded in the oil+TBP side of the interface. The protrusion shows up as a small shoulder on the water density profile in Figure 1c near the interface.

Figure 9. Free-energy profiles for the UO2(NO3)2 · TBP2 complex across the water/oil (dark line), water/(oil+TBP) low density system (dashed line), and the water/(oil+TBP) bulk system (dotted line). The interfacial well decreases in magnitude and broadens with increasing TBP concentration, and the driving force for partitioning into the oil phase increases with increasing TBP concentration. The water phase is on the left in the figure.

across liquid-liquid interfaces relevant to nuclear waste extraction. In particular, we examined water/oil and water/(oil+TBP)

interfaces for a range of TBP concentrations in the oil phase. TBP molecules are highly surface active at the water/oil interface. At very low TBP concentrations, the TBP molecules at the interface assume an amphiphilic orientation with the PdO group pointed into the aqueous phase. At higher concentrations, the interface broadens and the TBP molecules form chain structures parallel to the interface plane. There are significant interactions of the TBP molecules with water at the interface. The neutral uranyl nitrate complex UO2(NO3)2 is strongly repelled from the water/oil interface (as observed in ref 31 at a water/SC-CO2 interface), and the partition coefficient to move into the oil phase is negligibly small. When bound to TBP molecules in the predominant UO2(NO3)2 · TBP2 form, however, the complex is surface active at the water/oil interface with a significant free-energy well. There is also a net free-energy driving force for the complex to partition into the oil phase. With increasing TBP concentration in the oil phase, the interfacial free-energy well decreases in magnitude and broadens, and the net free-energy driving force increases for partitioning of the complex into the oil phase.

11670

J. Phys. Chem. B, Vol. 113, No. 34, 2009

Figure 10. (a) Orientational probability distribution for the angle between the OdUdO vector (in the UO2(NO3)2 · TBP2 complex) and the interface normal as a function of z for the water/oil interface. The OdUdO vector is oriented perpendicular to the plane of the interface on average. (b) The orientational probability distribution for the angle between the OdUdO vector (in the UO2(NO3)2 · TBP2 complex) and the interface normal as a function of z for the water/(oil+TBP) interface system (low density case). The OdUdO vector is oriented parallel to the interface plane on average.

Our results are broadly consistent with previous simulations of uranyl nitrate in related systems.22,29-32 The simulations by Wipff and co-workers have focused on water/chloroform and water/SC-CO2 interfaces, while the simulations reported here are for water/alkane interfaces as utilized in the PUREX process and in the selective membranes developed in ref 5. All of these organic phases exhibit well-defined interfaces when in contact with water and without modifiers. Our simulations agree with those from the Wipff group in several respects. First, for low TBP concentrations in the oil phase, all of the simulations indicate surface activity of the TBP molecules and an amphiphilic orientation with the phosphoryl sPdO group pointing on average into the aqueous phase. Second, increasing TBP concentration leads to a broadening of the interfacial region and increased disorder in the orientations of the TBP molecules at the interface (yet new structure develops, see below). Third, there exists strong interactions between water and the TBP molecules at the interface. Fourth, increasing TBP concentration leads to increased extraction of the complexed uranyl nitrate to the oil phase. Fifth, solvent molecules can be dragged (at least for short distances) along with the solute complexes when passing from the water phase to the oil phase. The simulations presented in this study also lead to some new insights into the uranyl nitrate extraction interfaces. First, the free-energy profiles computed here help to elucidate the thermodynamic driving forces for extraction. The interfacial activities of TBP and UO2(NO3)2 · TBP2 are apparent as freeenergy minima in the interfacial region. Since the free-energy barriers for passage from the interfacial minimum into the oil phase are several times kT, this helps to explain how the complexes may become trapped at the interface for rather long times during simulation studies. The increased extraction with increasing TBP concentration is apparent from the free-energy differences between the two bulk phases. Second, another feature of interest is the long-time (38 ns) stability of the interfaces simulated here. While hexane, chloroform, and SC-CO2 provide

Jayasinghe and Beck similar organic environments in several respects, there are also likely to be some differences. Comparing hexane and chloroform, the chloroform molecule has a dipole moment of 1.04 D, and the dielectric constant of the liquid is 4.8; the dielectric constant of liquid hexane is 1.9. In addition, the solubility of chloroform in water is 600 times larger than the solubility of hexane. Thus, stronger interactions between water and chloroform can be expected relative to water/hexane. Also, the SC-CO2 solvent incorporates more water than does the chloroform solvent,22 indicating that SC-CO2 interacts even more strongly with water due to the supercritical as opposed to liquid environments. These factors may help to rationalize the stable interfaces seen in the present study (even for high concentration TBP systems), as compared to the microemulsion behavior predicted in demixing simulations for high concentrations of TBP in chloroform.22 The experiments of ref 5, performed on water/alkane interfaces and at modest pH values, support the above discussion; no equilibration of pH or ion concentrations across the organic layer was observed, even at very high TBP concentrations (50% by volume). If the oil phase in that work rather existed as a microemulsion, charge transfer across the membrane would be expected. Thus, our results are consistent with those experiments. Third, by analysis of various radial distribution functions and orientational probability distributions, a structural picture of TBP at the interface emerges which differs somewhat from previous studies. As discussed above, a well-defined amphiphilic orientation is observed at low TBP concentrations. At high TBP concentrations, however, we observe a broad TBP layer that is somewhat disordered but also exhibits structure as evidenced by chain-like features oriented parallel to the interface. These chain-like structures are apparently driven by TBP dipole-dipole interactions. Fourth, the orientations of the UO2(NO3)2 · TBP2 complex at the interface change substantially with increasing TBP concentration in the oil phase. The present work is an initial step toward a better understanding of the free-energy driving forces for selective ion extraction of nuclear waste materials. The current simulations have not incorporated realistic acid and/or salt concentrations in the aqueous feed solution, and examining the impact of those low pH and/or high ionic strengths should help illuminate the influence of those factors on the extraction process. In particular, salting out of the UO2(NO3)2 · TBP2 complex from the aqueous solution can be expected to further drive the extraction process.31,55 Extensive simulations by Wipff and co-workers29-32 have begun to elucidate the important role of acid on the interface structure and TBP complexation. Extending those studies to computation of thermodynamic quantities will help to further quantify the complex driving forces occurring at these interfaces. The role of large, surface-active, hydrophobic anions in extraction22 also deserves further investigation. Such specificion or Hofmeister effects have received a great deal of recent attention.56 Finally, the present simulations have focused on the structure and thermodynamics for uranyl ion extraction only and have not directly addressed the problem of ion selectivity. This is a more complex issue that requires an explicit handling of chemical details, likely with quantum chemical methods.12,38 Several groups have initiated ab initio simulation studies of the uranyl ion and its complexes in aqueous solution,13,14,21 which is a first step in that direction. A recently developed quasichemical theory of solutions57-59 provides another avenue which should be pursued. That theory partitions the solvation free energy into inner-shell and outer-shell components. The parti-

Carrier-Assisted Uranyl Ion Extraction tioning is helpful both for a fundamental understanding of the solvation process and for the development of efficient numerical tools to compute the free energies. Future work will focus on applying the quasi-chemical theory to these challenging problems to try to develop a better fundamental understanding of selective ion extraction. Acknowledgment. We acknowledge the support of the National Science Foundation (NSF) through the Membrane Applied Science And Technology Center (MAST) at the Universities of Cincinnati and Colorado, and grant CHE0709560. We would like to thank Gordon Jarvinen and William Krantz for introducing us to this research area, and Lawrence Pratt, Steve Clarson, William Heineman, F. James Boerio, Joe Caruso, and Thomas Ridgway for many informative discussions. We acknowledge Alan Grossfield for permitting us to use his WHAM package for the free-energy calculations. Finally, we acknowldege a grant of computer time at the Ohio Supercomputer Center, where some of these calculations were performed. Note Added in Proof. Following submission of this paper, we became aware of a recent paper by X. Ye, S. Cui, V. de Almeida, and B. Khomami which presents closely related simulations of water/(oil+TBP) interfaces and uranyl ion extraction [J. Phys. Chem. B 2009, 113 (29), 9852-9862 (DOI: 10.1021/jp810796m)]. References and Notes (1) Department of Energy Five Year Plan FY 2007-FY 2011, Vol. 2, 2006. http://www.em.doe.gov/PDFs/170016EM_FYP_Final_3-6-06.pdf. (2) Choppin, G. R.; Khankhasayev, M. K. Chemical Separation Technologies and Related Methods of Nuclear Waste Management: Applications, Problems, and Research Needs, 1st ed.; Springer: New York, 1999; p 320. (3) Choppin, G. R.; Khankhasayev, M. K.; Plendl, H. S. Chemical Separations in Nuclear Waste Management: The State of the Art and a Look to the Future; Battelle Press: Columbus, OH, 2002; p 96. (4) Mills, A. L.; Logan, W. R. SolVent Extraction Chemistry; NorthHolland: Amsterdam, The Netherlands, 1967; p 322. (5) McCleskey, T. M.; Ehler, D. S.; Young, J. S.; Pesiri, D. R.; Jarvinen, G. D.; Sauer, N. N. J. Membr. Sci. 2002, 210, 273–278. (6) Jarvinen, G. D. Personal communication. (7) Hirata, M.; Monjyushiro, H.; Sekine, R.; Onoe, J.; Nakamatsu, H.; Mukoyama, T.; Adachi, H.; Takeuchi, K. J. Electron Spectrosc. Relat. Phenom. 1997, 83, 59–64. (8) Vallet, V.; Macak, P.; Wahlgren, U.; Grenthe, I. Theor. Chem. Acc. 2006, 115, 145–160. (9) Sonnenberg, J. L.; Hay, P. J.; Martin, R. L.; Bursten, B. E. Inorg. Chem. 2005, 44, 2255–2262. (10) Boncella, J. M. Nature 2008, 451, 250–252. (11) Arnold, P. L.; Patel, D.; Wilson, C.; Love, J. B. Nature 2008, 451, 315–7. (12) Gutowski, K. E.; Dixon, D. A. J. Phys. Chem. A 2006, 110, 8840– 56. (13) Bu¨hl, M.; Diss, R.; Wipff, G. J. Am. Chem. Soc. 2005, 127, 13506– 7. (14) Nichols, P.; Bylaska, E. J.; Schenter, G. K.; de Jong, W. J. Chem. Phys. 2008, 128, 4507. (15) Farkas, I.; Banyai, I.; Szabo, Z.; Wahlgren, U.; Grenthe, I. Inorg. Chem. 2000, 39, 799–805. (16) Aaberg, M.; Ferri, D.; Glaser, J.; Grenthe, I. Inorg. Chem. 1983, 22, 3986–3989. (17) Soderholm, L.; Skanthakumar, S.; Neuefeind, J. Anal. Bioanal. Chem. 2005, 383, 48–55. (18) Rozen, A. M. At. Energ. 1957, 2, 545–559.

J. Phys. Chem. B, Vol. 113, No. 34, 2009 11671 (19) Fratiello, A.; Kubo, V.; Lee, R. E.; Schuster, R. E. J. Phys. Chem. 1970, 74, 3726–3730. (20) Bu¨hl, M.; Kabrede, H.; Diss, R.; Wipff, G. J. Am. Chem. Soc. 2006, 128, 6357–68. (21) Bu¨hl, M.; Diss, R.; Wipff, G. Inorg. Chem. 2007, 46, 5196–206. (22) Baaden, M.; Schurhammer, R.; Wipff, G. J. Phys. Chem. B 2002, 106, 434–441. (23) Spence, R. Proc. R. Soc. London, Ser. A 1957, 243, 1–14. (24) Taylor, J. C.; Mueller, M. H. Acta Crystallogr. 1965, 19, 536– 543. (25) Mueller, M. H.; Dalley, N. K.; Simonsen, S. H. Inorg. Chem. 1971, 10, 323–328. (26) Guilbaud, P.; Wipff, G. J. Mol. Struct. 1996, 366, 55–63. (27) Beudaert, P.; Lamare, V.; Dozol, J.; Troxler, L.; Wipff, G. SolVent Extr. Ion Exch. 1998, 16, 597–618. (28) Baaden, M.; Berny, F.; Muzet, N.; Troxler, L.; Wipff, G. In Calixarenes for Separations; Lumetta, G. J., Rogers, R. D., Gopalan, A. S., Eds.; ACS Symposium Series 757; American Chemical Society: Washington, DC, 2000; p 71. (29) Baaden, M.; Berny, F.; Wipff, G. J. Mol. Liq. 2001, 90, 1–9. (30) Baaden, M.; Burgard, M.; Wipff, G. J. Phys. Chem. B 2001, 105, 11131–11141. (31) Schurhammer, R.; Wipff, G. In Separations and Processes using Supercritical Carbon Dioxide; Gopalan, A. S., Wai, C., Jacobs, H., Eds.; ACS Symposium Series 860; American Chemical Society: Washington, DC, 2003; pp 223-244. (32) Schurhammer, R.; Wipff, G. J. Phys. Chem. A 2005, 109, 5208– 5216. (33) Lauterbach, M.; Engler, E.; Muzet, N.; Troxler, L.; Wipff, G. J. Phys. Chem. B 1998, 102, 245–256. (34) Tsushima, S.; Suzuki, A. THEOCHEM 2000, 529, 21–25. (35) Spencer, S.; Gagliardi, L.; Handy, N. C.; Ioannou, A. G.; Skylaris, C. K.; Willetts, A.; Simper, A. M. J. Phys. Chem. A 1999, 103, 1831– 1837. (36) Tsushima, S.; Yang, T.; Suzuki, A. Chem. Phys. Lett. 2001, 334, 365–373. (37) Jackson, V.; Craciun, R.; Dixon, D.; Peterson, K.; de Jong, W. J. Phys. Chem. A 2008, 112, 4095–4099. (38) Hayton, T.; Boncella, J.; Scott, B.; Batista, E. J. Am. Chem. Soc. 2006, 128, 12622–12623. (39) Benjamin, I. Chem. ReV. 2006, 106, 1212–33. (40) Benjamin, I. Chem. ReV. 1996, 96, 1449–1476. (41) Benjamin, I.; Pohorille, A. J. Chem. Phys. 1993, 98, 236–42. (42) Pohorille, A.; Wilson, M. A. J. Chem. Phys. 1996, 104, 3760–73. (43) Wilson, M. A.; Pohorille, A. J. Am. Chem. Soc. 1994, 116, 1490– 501. (44) Pratt, L. R.; Pohorille, A. Chem. ReV. 2002, 102, 2671–92. (45) Nicolas, J. P.; de Souza, N. R. J. Chem. Phys. 2004, 120, 2464–9. (46) Torrie, G. M.; Valleau, J. P. Chem. Phys. Lett. 1974, 28, 578–581. (47) MacKerell, A.; et al. J. Phys. Chem. B 1998, 102, 3586–3616. (48) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926–935. (49) Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. J. Am. Chem. Soc. 1984, 106, 6638–6646. (50) Jorgensen, W. L.; Tirado-Rives, J. J. Am. Chem. Soc. 1988, 110, 1657–1666. (51) Feller, S.; Pastor, R.; Rojnuckarin, A.; Bogusz, S.; Brooks, B. J. Phys. Chem. 1996, 100, 17011–17020. (52) Kumar, S.; Rosenberg, J. M.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A. J. Comput. Chem. 1992, 13, 1011–1021. (53) Roux, B. Comput. Phys. Commun. 1995, 91, 275–282. (54) Frisch, M. J.; et al. Gaussian 03, revision B.05; Gaussian, Inc.: Wallingford, CT, 2004. (55) Ghosh, T.; Kalra, A.; Garde, S. J. Phys. Chem. B 2005, 109, 642– 651. (56) Kunz, W.; Nostro, P. L.; Ninham, B. Curr. Opin. Colloid Interface Sci. 2004, 9, 1–18. (57) Beck, T. L.; Paulaitis, M. E.; Pratt, L. R. The Potential Distribution Theorem and Models of Molecular Solutions; Cambridge University Press: New York, 2006. (58) Rogers, D.; Beck, T. J. Chem. Phys. 2008, 129, 134505. (59) Shah, J.; Asthagiri, D.; Pratt, L.; Paulaitis, M. J. Chem. Phys. 2007, 127, 144508.

JP903470N