Structure, Kinetics, and Thermodynamics of the ... - ACS Publications

Article pubs.acs.org/JPCA

Structure, Kinetics, and Thermodynamics of the Aqueous Uranyl(VI) Cation Sebastien Kerisit* and Chongxuan Liu Physical Sciences Division, Pacific Northwest National Laboratory, Richland, Washington 99352, United States ABSTRACT: In this work, molecular simulation techniques were employed to gain insight into the structural, kinetic, and thermodynamic properties of the uranyl(VI) cation (UO22+) in aqueous solution. The simulations made use of an atomistic potential model (force field) derived in this work and based on the model of Guilbaud and Wipff [J. Mol. Struct. (THEOCHEM) 1996, 366, 55−63]. Reactive flux and thermodynamic integration calculations show that the derived potential model yields predictions for the water exchange rate and free energy of hydration, respectively, that are in agreement with experimental data. The water binding energies, hydration shell structure, and selfdiffusion coefficient were also calculated and analyzed. Finally, a combination of metadynamics and transition path sampling simulations was employed to probe the mechanisms of water exchange reactions in the first hydration shell of the uranyl ion. These atomistic simulations indicate, based on two-dimensional free energy surfaces, that water exchanges follow an associative interchange mechanism. The nature and structure of the water exchange transition states were also determined. The improved potential model is expected to lead to more accurate predictions of uranyl adsorption energies at mineral surfaces using potential-based molecular dynamics simulations.

■

INTRODUCTION Uranium is one of the principal groundwater contaminants at nuclear facilities because of its use as nuclear fuel for weapons production. The ability to predict the fate of uranium in contaminated subsurface environments, in the context of remediation at nuclear facilities, relies heavily on the development of a thorough fundamental understanding of uranium transport and uranium reactivity with mineral surfaces. As a result, the speciation of uranium in contaminated sediments, particularly the surface complexation of the uranyl(VI) cation (its oxidized and soluble form) at relevant mineral surfaces, has attracted significant attention.1−10 Computational modeling can yield insights into uranyl complexation at mineral surfaces, and consequently, both classical11−19 and quantum mechanical20−27 simulation techniques have been applied to investigate uranyl adsorption at various mineral surfaces. Molecular dynamics simulations based on classical force fields constitute an appealing approach thanks to their ability to model large systems over long time scales. However, the accuracy of such simulations is dependent on the choice and derivation of the force field parameters. The force field parameters derived by Guilbaud and Wipff (GW)28,29 for simulating the uranyl ion and its interactions with water are the most commonly used in molecular dynamics simulations of uranyl adsorption at mineral−water interfaces.11,13−19 In addition, the GW model has been used on many occasions and in a wide variety of studies including those on uranyl liquid−liquid extraction,30−32 uranyl interactions with proteins33 and uptake by membranes,34,35 and uranyl complexation by carbonate36,37 and nitrate28,29,38 anions in aqueous solutions. In a previous publication,36 we reported that the water exchange rate in the first hydration shell of the uranyl ion, © 2013 American Chemical Society

as calculated with the GW model, is more than 2 orders of magnitude higher than the experimental value determined from 17 O nuclear magnetic resonance (NMR) measurements by Szabó et al.39 and Farkas et al.40 This difference between the calculated and measured water exchange rates indicates that the uranium−water interactions in the GW model are not sufficiently strong. Although the GW parameters were derived to reproduce the experimental difference in hydration free energies between UO22+ and Sr2+,29 the experimental value used in that study for the uranyl hydration free energy has since been shown to be too low.41−46 Therefore, the use of a low hydration free energy of uranyl in the derivation of the original model explains the erroneous high water exchange rate of the GW model compared to the NMR measurements. In this work, we aimed to (1) derive a more accurate model of the uranyl ion in aqueous solution by improving the agreement with the hydration free energy and the water exchange rate with respect to the GW model and (2) use the improved model to determine the structural, kinetic, and thermodynamic properties of the aqueous uranyl(VI) cation. In particular, the mechanism for water exchange in the first hydration shell of the uranyl ion was identified using a combination of rare-event techniques. The derivation of a new model with improved uranyl−water interactions will enable more accurate simulations of the energetics of adsorption of uranyl at mineral surfaces and will thus facilitate the development of uranyl surface complexation models. Received: May 8, 2013 Revised: June 26, 2013 Published: July 1, 2013 6421

dx.doi.org/10.1021/jp404594p | J. Phys. Chem. A 2013, 117, 6421−6432

The Journal of Physical Chemistry A

■

Article

COMPUTATIONAL METHODS All calculations presented in this work were carried out with the molecular dynamics computer code DL_POLY.47 Several techniques were employed to characterize the properties of the uranyl ion in aqueous solution. We have limited ourselves to describing only the basics of each technique necessary for introducing the simulation parameters; more details can be found in the references cited in this section. Standard Simulations. Standard molecular dynamics simulations were carried out to determine the equilibrium structural properties of the hydrated uranyl ion, as well as its diffusion properties. The simulations were performed in the NPT ensemble (constant number of particles, constant pressure, and constant temperature) at 300 K and zero applied pressure. The temperature and pressure were kept constant by the Nosé−Hoover thermostat48 and the Hoover barostat,49 respectively, with relaxation time parameters of 0.5 ps. The electrostatic forces were calculated by means of the Ewald summation method.50 A 9-Å cutoff was used for the short-range interactions and the real part of the Ewald sum. The Ewald sum parameters were chosen to achieve a relative error in the electrostatic energy of at most 10−7. The Verlet leapfrog integration algorithm was used to integrate the equations of motion with a time step of 0.001 ps. The simulation cells contained one uranyl ion and 511 water molecules, which corresponds to box lengths of approximately 25 Å. The geometry of the water molecules was held fixed using the SHAKE algorithm. The simulations were run for 10 ns after an equilibration period of 100 ps (equilibrated systems were available from a previous study36). The diffusion coefficient of the uranyl ion, D, was determined from the mean square displacement of the uranium atom. To compute D, a configuration was recorded every 0.2 ps, and a correlation time of 20 ps was used, following the approach used in a previous study.36 Reactive Flux. The reactive flux approach (ref 51 and references therein) was used to calculate the dynamics of water exchange in the first hydration shell of the uranyl ion. In this approach, the reaction coordinate is defined as the distance r between the centers of mass of the uranyl ion and the outgoing water molecule. The exchange rate constant, k, is written as k = κk

TST

W (r ) =

where κ is the transmission coefficient and k state rate constant, given by kBT r *2 e−βW (r*) = 2πμ ∫ r * r 2e−βW (r) dr 0

r

⎛r⎞ ⟨−FC(r )⟩ dr + 2kBT ln⎜ ⎟ ⎝ r0 ⎠

(4)

The transmission coefficient was determined from the plateau value of the normalized reactive flux, which can be computed as k(t ) =

⟨r(0) ̇ θ[r(t ) − r *]⟩c ⟨r(0) ⟩c ̇ θ[r(0)] ̇

(5)

where θ[x] is the Heaviside function, which is 1 if x is larger than 0 and 0 otherwise, and ṙ(0) is the initial uranyl−water velocity along the reaction coordinate. The subscript c means that the initial configurations were generated in the constrained reaction coordinate ensemble. To compute the transmission coefficient, a 5-ns simulation was performed in which the reaction coordinate was constrained to r* and a configuration was collected every 2 ps. Each of the 2500 configurations thus generated was then run both backward and forward in time for 2 ps, and the value of κ was determined by averaging k(t) over the last 0.5 ps. Thermodynamic Integration. Hydration free energies were calculated using the thermodynamic integration approach. The calculations were carried out in two stages, whereby the solvation cavity is grown by varying the short-range potential parameters in the first stage and the ion is slowly charged in the second stage (with the short-range parameters set to their full values). For the first stage, the integration was discretized into 12 intervals, whereby the dependence of the Hamiltonian H on the mixing parameter λ was defined as

Hλ = f (λ)H

(6)

where f (λ) = 1 − (1 − λ)k

(7)

and k equals 6. The following 12 values of λ were used following the approach described by Horinek et al.:53 0.00922, 0.04794, 0.11505, 0.20634, 0.31608, 0.43738, 0.56262, 0.68392, 0.79366, 0.88495, 0.95206, and 0.99078. For the second stage, the integration was discretized into 11 intervals, and the dependence of the Hamiltonian on λ was set to be linear

is the transition

kBT e−βWeff (r*) 2πμ ∫ r * e−βWeff (r) dr 0

f (λ ) = λ

(2)

(8)

where λ was varied from 0.0 to 1.0 in 0.1 intervals. Two sets of corrections were made to the calculated hydration free energies to enable comparison with experimental data.54−57 The first correction accounted for the effects of the finite size and periodicity of the system on the water polarization. The second correction accounted for the difference between the experimental and theoretical reference states of the gas phase. Finally, we note that there has been a significant debate in the literature regarding the most adequate summation method for determining the electrostatic potential at the ion center (refs 55−57 and references therein), whereby both the so-called Psummation scheme (based on the charges of all solvent particles within a sphere) and M-summation scheme (based on

where μ is the ion−water reduced mass; β is 1/kBT; r* is the transition state distance, which is taken to be the value of r for which W(r) reaches a maximum between the reactant (outgoing water molecule in the first hydration shell of the uranyl ion) and product (outgoing water molecule in the second hydration shell of the uranyl ion) states; and Weff(r), the centrifugally averaged effective potential, is defined as ⎛r ⎞ Weff (r ) = W (r ) − 2kBT ln⎜ ⎟ ⎝ r* ⎠

∫r

0

(1) TST

kTST =

in which the reaction coordinate was constrained between 2.2 and 3.2 Å in 0.025-Å intervals. Following the approach of Carter et al.,52 the mean constraint force, FC, from each run was then integrated to obtain the PMF, whereby the energy minimum of the first hydration shell, at r0, was set to the zero of energy

(3)

where W(r) is the potential of mean force (PMF). The PMF was obtained from a series of 1-ns constrained MD simulations 6422



Article

reaction coordinate, and only the reactant (when the outgoing water molecule is in the first hydration shell) and product (when the outgoing water molecule is in the second hydration shell) states need to be defined. Following the approach used in a previous study,68 the reactant and product states were defined based on the distance between the uranium ion and the oxygen atom of the outgoing water molecule, whereby the reactant state extends to the distance at which g(r) of the uranium− water oxygen radial distribution function falls below 1/3 and the product state begins where g(r) increases back up past 1/3. A Monte Carlo technique was used to sample the transition path ensemble and construct the path average by generating random walks in the space of trajectories following two basic trial moves called “shooting” and “shifting”. In a shooting move, a point along the existing trajectory was selected at random, next the atomic momenta were modified, and finally a new path was generated forward and backward from that point. The new momenta were drawn from a Gaussian distribution with a variance of 0.2 times the average kinetic energy per atom in the simulation and were then rescaled so that the temperature remains unchanged. In a shifting move, a segment was deleted from either end of the existing trajectory, and a new segment of the same length was grown from the opposite end. This had the effect of shifting the pathway forward or backward in time. The probability for shifting moves was set to 0.2, with a maximum temporal shift of 0.5 ps. The trial path was accepted or rejected depending on whether it connected the reactant and product states, and the process was repeated. The average of five TPS simulations, started with different configurations, was used to construct the path average. Initial trajectories were generated using configurations from a PMF run at the transition state distance whereby trajectories with random initial velocities were run forward and backward in time until a trajectory that linked the reactant and product states was found. In each TPS run, 10000 paths were attempted, for a total of 50000 paths. Each path was 2.4 ps long. Potential Models. In the simulations reported in this work, atoms are represented as point-charge particles that interact through long-range Coulombic forces and short-range interactions. The latter are described by parametrized functions and represent the repulsion between electron-charge clouds, van der Waals attraction forces, and many-body terms such as bond bending. In this work, the pairwise interaction energy, Vij, takes the form

the charges of all solvent molecules within a sphere) have been advocated. The free energies reported and discussed in this work are based on the P-summation scheme. As it is our intention to apply this approach in the future to the alkalineearth uranyl carbonate species studied in previous works,18,36 this choice will facilitate comparison with the hydration free energies calculated by Raiteri et al.58 for alkaline-earth and carbonate ions, who made use of the P-summation scheme. Metadynamics. The metadynamics technique introduced by Laio, Parrinello, and co-workers59−62 was used to explore the free energy landscape of water exchange around the uranyl ion. The metadynamics technique requires the definition of one or several collective variables (CVs), which are functions of the system coordinates and can distinguish between the reactant and product states. In the continuous metadynamics approach used here, a nonequilibrium molecular dynamics simulation is run during which the dynamics of the collective variables is driven by a history-dependent bias potential that consists of a sum of repulsive Gaussians added at time intervals τ and centered on the position of the collective variables at that time. Therefore, the bias potential, Vbias, at position s and time t takes the form Nτ ≤ t

∑

Vbias(s , t ) =

ti = 0, τ ,2τ ,..., Nτ

⎛ |s − s(t )|2 ⎞ i ⎟ H exp⎜ − 2ω 2 ⎠ ⎝

(9)

where H and ω are the height and width of the Gaussians, respectively. This history-dependent potential accumulates in energy wells and eventually allows the system to cross energy barriers and escape to nearby local minima. The free energy surface can be obtained from −Vbias(s,t) if the MD simulation is run long enough so that the bias potential fills all minima. Two collective variables were used in the metadynamics simulations carried out in this work. The first CV is the distance between the uranium ion and the oxygen atom of the outgoing water molecule (U−Ow‑out). The second CV is the coordination number of the uranium ion (CNU), calculated as N

CNU =

∑ i=1

12

() 1−( ) 1−

ri r0

ri r0

20

(10)

where N is the number of the spectator water molecules (i.e., not considering the outgoing water molecule), ri is the distance between the uranium ion and the corresponding water oxygen, and r0 was set to 2.905 Å for the GW model and 2.805 Å for the two new models. Gaussians were added every τ = 0.5 ps during a 10-ns MD simulation with a height of 0.115 kcal mol−1 and a width of 0.1 (for both collective variables). In addition, to prevent the system from sampling configurations that are not of interest (e.g., configurations in which the outgoing water molecule is found beyond the second hydration shell), a repulsive wall was added that acts as a harmonic restraining potential on the U−Ow‑out CV, with force constant 100 eV Å−1, when the value of this CV exceeded 6.0 Å. Transition Path Sampling. The transition path sampling (TPS) approach developed by Chandler and co-workers63−67 was employed to collect reactive trajectories for the water exchange reaction around the uranyl ion. The TPS approach consists of carrying out importance sampling of an ensemble of reactive trajectories (referred to as the transition path ensemble). This approach does not require the selection of a

⎛ r * ⎞12 ⎛ r * ⎞6 1 qiqj ij ij Vij = + εij⎜⎜ ⎟⎟ − 2εij⎜⎜ ⎟⎟ 4πε0 rij r r ⎝ ij ⎠ ⎝ ij ⎠

(11)

where ε0 is the permittivity of vacuum, rij is the interatomic distance between atoms i and j, qi is the charge carried by i εij =

εiεj

(12)

rij* = (ri* + r *j )/2

(13)

and

Additionally, a harmonic bond potential and a harmonic angle potential are also included to describe the U−O and O−U−O bond and angle, respectively, in the uranyl ion. The short-range potential parameters and ionic charges used in this study to model water are those of the SPC/E model.69 Two new sets of parameters for modeling the uranyl ion and its interactions with water are introduced in this work and compared with the model 6423



Article

of Guilbaud and Wipff.29 All of the parameters are reported in Table 1.

Gas-Phase Hydrated Uranyl Clusters. This work focuses on the properties of the uranyl ion in aqueous solutions; however, the structure and energetics of small gas-phase hydrated uranyl clusters allow for the new parameters to be evaluated against more accurate electronic structure calculations with varying levels of theory. Therefore, the aim of this subsection is to determine whether the new potential parameters can accurately describe the interactions between uranyl and water in the first hydration shell, by comparing the structural and energetic differences between coordination numbers 4, 5, and 6, and, more generally, whether they represent an improvement over the GW model. The structure and energies of [UO2(H2O)n]2+ gas-phase clusters, with n = 4−6, were compared with published electronic structure calculations. Tables 2 and 3 report the uranium−water oxygens distances (U−OW) and water binding energies, respectively, obtained with the three models and from five quantum mechanical studies.40,41,45,70,71 This list of quantum mechanical studies is not meant to be exhaustive, as a large number of such studies have been published to date; instead, the selected studies illustrate a range of levels of theory that have been applied to compute the structure and energies of these clusters. Table 2 shows that, for all three coordination numbers, the GW model yields U−OW distances that are in the upper range of those determined from electronic structure calculations, whereas model 2 is generally in the lower end of the range and model 1 is slightly below the lowest value. However, all three models reproduce fairly well the changes in U−OW in going from four to five and from five to six waters of hydration. The water binding energies from the five quantum mechanical (QM) studies (Table 3) span a wide range of energies, with those obtained by Farkas et al.40 from Møller− Plesset perturbation theory to the second order (MP2) [using Hartree−Fock- (HF-) optimized structures] as the lowest energies and the coupled-cluster (singles and doubles, CCSD) calculations of Cao and Balasubramanian71 (corrected for basisset superposition error) as the highest energies. The water binding energies obtained with models 1 and 2 generally fall in the middle of that range, whereas the GW energies are significantly lower. The three models give reasonable agreement with the energy difference between the tetra- and pentahydrates reported in the QM studies. Interestingly, the MP2 calculations40 and the DFT calculations with the Becke and Perdew (BP)41 and Becke−Lee−Yang−Parr (BLYP)70 exchange-correlation potentials predict a positive energy difference between the penta- and hexahydrates, whereas the CCSD calculations71 predict a negative difference. The water binding energies from the three models are in agreement with the trend of the CCSD calculations. Notably, model 2 gives very good agreement with the binding energy differences obtained from the CCSD calculations.

Table 1. Model Parameters Used in This Worka parameter

GW

model 1

model 2

qU (e) qO (e) εU (kcal mol−1) εO (kcal mol−1) r*U (Å) r*O (Å)

2.500 −0.250 0.400 0.200 1.58 1.75

3.250 −0.625 0.120 0.200 1.60 1.75

3.500 −0.750 0.300 0.200 1.58 1.75

The harmonic bond potential [VU−O = kU−O/2(rU−O − r0)2, where kU−O = 43.364 eV Å−2 and r0 = 1.8 Å] and harmonic angle potential [VO−U−O = kO−U−O/2(θO−U−O − θ0)2, where kO−U−O = 13.009 eV rad−2 and θ0 = 180°] derived by Guilbaud and Wipff (GW)28,29 were used for all three models. The SPC/E parameters are as follows: qO = −0.8476 e, qH = 0.4238 e, εO = 0.155 kcal mol−1, and rO* = 1.777 Å, with the O−H bond distances and H−O−H bond angle fixed at 1.0 Å and 109.47°, respectively. a

■

RESULTS AND DISCUSSION Derivation of New Model Parameters. As discussed previously, the objective of the new derivation was to strengthen the uranyl−water interactions in an attempt to obtain a better representation of the water exchange rate and hydration free energy. One way to strengthen the uranyl−water interactions is to increase the uranium partial charge, qU; therefore, two series of 1-ns MD simulations were run in which qU was set to 3.25 and 3.50 e and the values of εU and rU* were varied from 0.12 to 0.30 kcal mol−1 and from 1.50 to 1.74 Å, respectively, while the values of εO and r*O from the GW model were kept constant. To evaluate the quality of any new set of parameters, the calculated hydration energies and position of the first maximum of the uranium−water oxygen radial distribution function were compared with experimental data. We note that, to avoid performing a more computationally expensive thermodynamic integration for each pair of parameters and, thus, to accelerate the derivation process, the average potential energy was used instead of the free energy for comparison with the experimental hydration free energies. For a few selected cases, we observed that the ratio of hydration energy to hydration free energy was fairly constant, which indicates that this approach is appropriate for identifying trends. Thermodynamic integration calculations were performed with selected parameters as described later in this article. For each series, the pair of parameters (i.e., εU and rU*) that yielded the best overall agreement with the two evaluation properties while retaining a coordination number of 5.0 was selected for further investigation, as reported in Table 1. The models in which qU was 3.25 and 3.50 e are denoted as models 1 and 2, respectively.

Table 2. Uranium−Water Oxygen Distances (Å) of [UO2(H2O)n]2+ Gas-Phase Clusters (n = 4−6) n=4 n=5 n=6 Δ(5−4) Δ(6−5) a

GW

model 1

model 2

HFa

LDAb

BLYPc

B3LYPd

CCSDe

2.49 2.53 2.63 0.04 0.10

2.33 2.38 2.50 0.05 0.12

2.39 2.43 2.52 0.04 0.09

2.52 2.57 2.70 0.05 0.13

2.36 2.42 2.53 0.06 0.11

2.48 2.55 2.69 0.07 0.14

2.44 2.50 − 0.06 −

2.47 2.50−2.52 2.51−2.69 ∼0.04 −

Farkas et al.40 bMoskaleva et al.41 cSpencer et al.70 dGutowski and Dixon.45 eCao and Balasubramanian.71 6424



Article

Table 3. Water Binding Energies (kcal mol−1) of [UO2(H2O)n]2+ Gas-Phase Clusters (n = 4−6) n=4 n=5 n=6 Δ(5−4) Δ(6−5)

GW

model 1

model 2

MP2a

BPb

BLYPc

B3LYPd

CCSDe

−180.4 −211.9 −226.2 −31.5 −14.3

−234.7 −270.8 −281.7 −36.1 −10.9

−234.1 −273.6 −290.1 −39.5 −16.5

−214.0 −243.8 −241.4 −29.8 2.4

−237 −263 −258 −26.0 5.0

−233.1 −256.6 −252.6 −23.5 4.0

−223.9 −248.4 − −24.5 −

−259.7 −298.4 −316.4 −38.7 −17.9

a

Farkas et al.,40 from HF-optimized structures. bMoskaleva et al.,41 from LDA-optimized structures. cSpencer et al.70 dGutowski and Dixon.45 eCao and Balasubramanian.71

Overall, models 1 and 2 provide significant improvement over the GW model with respect to the magnitude of the water binding energies. This is achieved by increasing the magnitude of the partial charges on the uranium ion, which leads to shorter uranium−water distances. For model 1, this causes the U−OW distances to fall slightly below the range determined in quantum mechanical calculations. However, the distance differences as well as the energy differences, when compared with the CCSD calculations, are consistent with the QM calculations. Hydration Shell Structures. Standard 10-ns MD simulations were carried out with models 1 and 2, as well as with the GW model. Figure 1 compares three radial distribution functions (RDFs), namely, those for uranium−water oxygen (U−OW), uranyl oxygen−water hydrogen (OU−HW), and uranyl oxygen−water oxygen (OU−OW), with those obtained from the ab initio molecular dynamics simulation (AIMD) of Nichols et al.72 [specifically, Car−Parrinello DFT simulations with the gradient-corrected Perdew−Burke−Ernzerhof (PBE) exchange-correlation functional]. Additionally, the values of the position, magnitude, and coordination number of the first peak of each RDF are reported in Table 4. The coordination numbers were calculated by integrating each RDF up to its first minimum. The position of the first maximum in the U−OW RDF of model 2 (2.40 Å) agrees exactly with that obtained from the AIMD simulation, whereas model 1 and the GW model under(2.36 Å) and over- (2.49 Å) estimate this distance, respectively. Other AIMD simulations by Bühl et al.73,74 resulted in mean U−OW distances of 2.47 ± 0.09 and 2.48 ± 0.10 Å from AIMD, whereas the quantum mechanics/molecular mechanics MD simulation of Frick et al.75 and the MD simulation of Hagberg et al.76 (with a force field derived from multiconfigurational QM calculations) yielded values of 2.49 and 2.40 Å, respectively. Experimental mean U−OW distances from X-ray absorption fine structure spectroscopy (XAFS), X-ray diffraction (XRD), and X-ray scattering (XS) of aqueous uranyl(VI) range from 2.40 to 2.42 Å.77−84 Although all three models exhibit a coordination number of exactly 5.0 for the first hydration shell as obtained by Nichols et al.,72 the magnitude of the first peak is greater than that seen in the AIMD simulation, reflecting a narrower range of distances in the first hydration shell. The AIMD simulation also appears to show a slight asymmetry in the first peak, whereas the bonding environment is more symmetrical in the MD simulations. Note that the AIMD simulation of Nichols et al.72 was for a smaller simulation cell (122 H2O) and shorter simulated time (9 ps) than the MD simulations carried out in this work. Experimental coordination numbers obtained from XAFS vary from 4.5 to 5.3.78−81,83 Because the uncertainties in the coordination numbers reported in these studies range approximately from 0.3 to 0.8, the coordination numbers of

Figure 1. (Top) Uranium−water oxygen (U−OW), (middle) uranyl oxygen−water hydrogen (OU−HW), and (bottom) uranyl oxygen− water oxygen (OU−OW) radial distribution functions (RDFs) obtained with the original GW model and the two models derived in this work. Also shown for comparison are the RDFs from the AIMD simulation of Nichols et al.72

exactly 5.0 obtained with the three models agree with the available XAFS data within uncertainties. However, it should be noted that recent high-energy XS data provide evidence for a small proportion of uranyl ions in 4-fold coordination.82,84 Based on the shape and position of the second RDF peak, models 1 and 2 have very similar second hydration shells that are in better agreement with the AIMD simulation than the GW model. The OU−HW RDF from the AIMD simulation 6425



Article

uranyl and water diffusion coefficients. The uranyl diffusion coefficient of Marx and Bischoff,85 which was determined to be the most adequate for comparison with the MD simulations for reasons discussed in detail in the same previous publication,36 yields a D0/DH2O ratio of 0.330. Although the relative diffusion coefficients of models 1 and 2 decrease slightly, they remain within the uncertainty of the experimental value, thus showing that the new parametrization remains consistent with the experimental data. Free Energies of Hydration and Water Exchange Rates. As discussed previously, the main deficiency of the GW model was shown to be its low hydration free energy and high water exchange rate. Table 5 lists the free energies of hydration

Table 4. Structural and Diffusion Properties of the Hydrated Uranyl Ion Obtained with the Three Models properties

GW

model 1

model 2

AIMDa

dU−OW (Å)

2.49

2.36

2.40

2.40

g(r) first maximum CN dOU−HW (Å)

11.9

15.9

15.9

7.8

5.0 −

5.0 2.08

5.0 1.90

5.0 1.96

−

0.1

0.4

0.2

− 3.08

0.3 2.98

0.7 3.01

0.3 2.96

3.6

4.5

4.3

3.5

12.5 0.94 ± 0.03 0.33 ± 0.04

14.7 0.85 ± 0.04 0.30 ± 0.03

14.5 0.88 ± 0.03 0.31 ± 0.03

14.8 − −

g(r) first maximum CN dOU−OW (Å) g(r) first maximum CN D0 (10−9 m2 s−1) D0/DH2O a

Table 5. Hydration Free Energies and Water Exchange Properties of the Hydrated Uranyl Ion Obtained with the Three Models

Nichols et al.72

shows a small peak centered at 1.96 Å that corresponds to a small number of water molecules forming hydrogen bonds with the uranyl oxygens. This peak is not present in the GW simulation (and there is therefore no entries for OU−HW for the GW model in Table 4), but because the uranyl oxygen charge increases in model 1, this peak begins to appear and becomes more defined as the uranyl oxygen charge is increased further in model 2. Integration of this first peak shows that model 1 reproduces best the extent of hydrogen bonding observed in the AIMD simulation. Farther away, models 1 and 2 again provide a better agreement with the AIMD simulation than the GW model. In particular, the subsequent two peaks, which correspond to the two hydrogen atoms of the water molecules in the hydration shell of the uranyl oxygens, are more separated and in closer agreement with the AIMD simulation in the model 1 and model 2 simulations than in the GW simulation. For the OU−OW RDF, although models 1 and 2 show a first peak with a greater magnitude than the GW and AIMD simulations, the coordination number and peak position are improved in models 1 and 2 over the GW model with respect to the AIMD simulation (Table 4). Overall, models 1 and 2 show improvement over the GW model for reproducing the structure of the uranyl ion’s hydration shells. The first hydration shell of the uranium ion is similar for models 1 and 2, although with a mean U−OW distance for model 2 that is more consistent with XAFS data and the AIMD simulation of Nichols et al.,72 whereas the uranyl oxygen hydration shell is more accurately described by model 1. Self-Diffusion Coefficients. The size-independent diffusion coefficient, D0, decreases for models 1 and 2 compared to the GW model (Table 4). Because of the tighter hydration shell in models 1 and 2 compared to the GW model, one would expect the diffusion coefficient of the new models to increase with respect to that calculated for the GW model based on the Stokes−Einstein equation; however, the difference in size is only slight, and the calculated difference in diffusion coefficient is more likely due to the stronger interactions of the uranyl ion with water in models 1 and 2 as a result of the increase in magnitude of the partial charges of the uranium and uranyl oxygens. As noted in a previous publication,36 the sizeindependent diffusion coefficient of the SPC/E water model is greater than that determined experimentally. Therefore, comparison with experimental data is based on the ratio of the

property

GW

model 1

model 2

−1 ΔGabs hyd (kcal mol ) ‡ −1 ΔG (kcal mol ) (PMF) kTST (s−1) κ k (s−1) ΔG‡ (kcal mol−1) (metadynamics/ TPS)

−304 4.6 3.9 × 109 0.188 7.4 × 108 4.4

−371 8.4 7.8 × 106 0.194 1.5 × 106 9.1

−377 5.7 7.8 × 108 0.149 1.2 × 108 6.4

calculated for the three models using thermodynamic integration. Based on the discussions and references contained in Gibson et al.42 and Rai et al.,46 the absolute free energy of hydration of the uranyl ion was calculated as abs abs ΔG hyd [UO2 2 +] = ΔHhyd [UO2 2 +] − T ΔS hyd[UO2 2 +]

(14)

where the absolute hydration enthalpy is defined as abs conv ΔHhyd [UO2 2 +] = ΔHhyd [UO2 2 +] + 2ΔHfabs[H +(aq)]

= ΔHfconv[UO2 2 +(aq)] − ΔHf [UO2 2 +(g)] abs + 2(ΔHf [H+(g)] + ΔHhyd [H+])

(15)

where conversion from conventional to absolute energy scale is needed to allow direct comparison with the results of the MD simulations. Using the values for ΔHf[UO22+(g)] (289 ± 5 kcal mol−1) and ΔHconv [UO22+(aq)] (−249 ± 1 kcal mol−1) calculated by f 86 Marcus based on the lattice enthalpy of UO2F2, together with the literature values of 36786 and −261 ± 242 kcal mol−1 for + ΔHf[H+(g)] and ΔHabs hyd[H ], respectively, Marcus reported a 2+ value of −325 ± 5 kcal mol−1 for ΔHabs hyd[UO2 ]. In the same 2+ publication, Marcus also calculated ΔShyd[UO2 ] to be −79 ± 2+ 1 cal K−1 mol−1, which yields ΔGabs hyd[UO2 ] = −302 ± 6 kcal −1 mol . The free energy of hydration obtained from thermodynamic integration with the GW model is in good agreement with this value. This finding is expected because the free energy from Marcus’ work was used in the derivation of the GW model.29 However, using the more recent value of ΔHf[UO22+(g)] determined by Gibson et al.42 (364 ± 15 kcal mol−1), together with the literature values of ΔH fconv [UO 2 2+ (aq)] and ΔHf[H+(g)] used by Gibson et al.42 (−244 ± 1 and 366 kcal 2+ mol−1, respectively), the magnitude of ΔHabs hyd[UO2 ] increases −1 42 2+ to −398 ± 15 kcal mol , which translates to ΔGabs hyd[UO2 ] 6426



Article

Figure 2. (Left) Potential of mean force and (right) transmission coefficient functions obtained from the reactive flux calculations of the original GW model and the two models derived in this work.

= −369 ± 15 kcal mol−1 (using a more appropriate value of ΔShyd[UO22+] = −96 ± 1 cal K−1 mol−1 advocated by Gibson et al.42 and Rai et al.46). The results of the thermodynamic integration with models 1 and 2 are in very good agreement with the latest value. It should also be noted that an alternative + −1 value of ΔHabs hyd[H ] (−275 kcal mol ) was reported recently 87 by Tissandier et al. Using this alternative value, 42 2+ ΔGabs becomes hyd[UO2 ] based on the work of Gibson et al. −1 −397 ± 15 kcal mol . Therefore, the hydration free energies obtained with models 1 and 2 fall within the range of values + obtained from the two reported values of ΔHabs hyd[H ]. Several research groups have calculated the uranyl hydration free energy from electronic structure calculations of uranyl−water clusters embedded in a water continuum model. Moskaleva et al.41 calculated a value of −420 kcal mol−1 using DFT with the Becke and Perdew exchange-correlation functional, which they reported should increase to −370 kcal mol−1 when considering zero-point energy and thermal corrections. Shamov and Schreckenbach43,44 obtained values of −383.4 and −413.5 kcal mol−1 from PBE and B3LYP calculations (the two sets of calculations also differed in their approach to the treatment of relativistic effects). Finally, Gutowski and Dixon45 predicted free energies ranging from −406 to −418 kcal mol−1 from DFT/B3LYP and MP2 calculations with a varying number of water molecules in the second hydration shell. These calculations support the more recent derivation of 2+ 2+ ΔGabs hyd[UO2 ] obtained from the value of ΔHf[UO2 (g)] 42 determined by Gibson et al. The water exchange rate was calculated for each of the three models using the reactive flux approach described in the Computational Methods section. The PMFs thus obtained are shown in Figure 2, together with the results for the transmission coefficient calculations. The water exchange properties are summarized in Table 5. The low free energy of hydration of the GW model translates into a water exchange rate (k = 7.4 × 108 s−1) higher than that obtained experimentally by Grenthe and co-workers39,40 using 17O NMR measurements (k = 1.4 × 106 and 1.30 ± 0.05 × 106 s−1). Consequently, the GW model also underestimates the activation free energy for water exchange (ΔG‡ = 4.6 kcal mol−1) with respect to that determined by Farkas et al.40 (ΔG‡ = 9.1 ± 0.5 kcal mol−1). Table 5 shows that model 1 (ΔG‡ = 8.4 kcal mol−1 and k = 1.5 × 106 s−1) provides a significant improvement over the GW model and good agreement with the experimental exchange rate and activation free energy.

Although model 2 has a free energy of hydration similar to that of model 1, it shows a lower free energy of activation and, consequently, a much higher exchange rate. This can be explained by considering the water exchange mechanisms discussed in the next subsection. Water Exchange Mechanism. A classification of ligand substitution mechanisms was first introduced by Langford and Gray.88 This classification consists of three main substitution mechanisms: associative (A), dissociative (D), and interchange. When applied to the water exchange reaction around the pentacoordinated uranyl ion, the A and D mechanisms proceed through hexa- and tetra-coordinated intermediates, respectively, whereas the interchange mechanism involves exchanges for which no intermediate is formed. The interchange mechanism can be further divided into associative interchange (Ia), dissociative interchange (Id), and concerted interchange (I) depending on whether bond formation, bond breaking, or neither dominates, respectively. Computational studies published to date that have investigated water exchange mechanisms around the uranyl ion can be divided into two main families. The first family consists of ab initio cluster models of the uranyl ion with one or two hydration shells either in the gas phase or in solution as simulated by a continuum model. The advantage of ab initio cluster models is that they can be treated with highly accurate electronic structure methods and, indeed, a wide variety of such methods has been applied to this subject, including DFT and post-HF methods.89−95 The consensus that has emerged from this body of work is that, given a sufficient description of the solvent beyond the first hydration shell, the A mechanism shows a lower activation energy than the D mechanism and the activation parameters of the A mechanism are in closer agreement with the experimental data of Farkas et al.40 than those of the D mechanism. Notably, Vallet et al.89 found the I and A mechanisms to have similar activation energies, and Rotzinger90 concluded, from his calculations of the lifetime of the associative intermediates, that an associative interchange mechanism could be possible (based on calculations carried out with one of the two hydration treatments considered in that work). In the second family of models, the solvent is treated explicitly using a number of water molecules in a periodic simulation cell at a finite temperature.73,96,97 Specifically, Bühl and co-workers73,96 used constrained Car−Parrinello molecular dynamics simulations to isolate the dissociative and associative 6427



Article

interchange mechanisms. They found the associative interchange mechanism to show a lower free energy of activation (6.7 versus 10.8 kcal mol−1). Atta-Fynn et al.97 carried out metadynamics AIMD simulations with the coordination number of the uranium ion with all water oxygens as the collective variable and concluded that the water exchange most likely followed a dissociative mechanism. Although this result might seem at odds with the calculations of Bühl and coworkers, the choice of collective variable in Atta-Fynn et al. did not allow them to consider interchange mechanisms, which prevents a direct comparison with the findings of Bühl and coworkers. The free energy surfaces obtained from the metadynamics simulations carried out in this work with the three models are shown in Figure 3. As described in the Computational Methods section, the collective variables used to describe these surfaces are the coordination number of the uranium ion with spectator water molecules (i.e., excluding the outgoing water molecule) (CNU) and the distance between the uranium ion and the oxygen atom of the outgoing water molecule (U−Ow‑out). All three models show two minima: a first minimum that corresponds to the reactant configuration (CNU is approximately 4 and U−Ow‑out is short) and a second minimum that corresponds to configurations in which the outgoing water molecule is in the second hydration shell and a new water molecule has entered the first hydration shell (CNU increases to approximately 5 and U−Ow‑out is long). The metadynamics simulations do not show any intermediate minima, which means that the reactions proceed through an interchange mechanism. In addition, because CNU increases before U− Ow‑out increases significantly, the simulations indicate an associative interchange mechanism, which is consistent with the majority of the quantum mechanical calculations discussed herein. In Figure 3, the average paths from the transition path sampling calculations are superimposed on the metadynamics free energy surfaces. The two approaches yield consistent results. The minima and maxima in free energy along these paths were used to determine the activation free energies predicted by the metadynamics simulations. As shown in Table 5, the activation free energies thus calculated are consistent with those obtained from the reactive flux simulations. The combination of metadynamics and transition path sampling calculations allows for determining the transition state configurations predicted by the three models. Based on the metadynamics free energy surfaces, configurations from the transition path sampling calculations are considered as transition states when 4.5 ≤ CNU ≤ 4.7 and 2.75 Å ≤ U− Ow‑out ≤ 2.95 Å for model 1 and when 4.55 ≤ CNU ≤ 4.75 and 2.85 Å ≤ U−Ow‑out ≤ 3.05 Å for model 2 and the GW model. From this pool of transition state configurations, the distribution of the distance between uranium and the oxygen of the incoming water molecule was calculated and is shown in Figure 4a for each of the three models. Figure 4a shows that the incoming water molecule is much closer to the uranium ion in the transition state configurations obtained with model 2 than in those obtained with model 1. Although water exchanges follow an interchange mechanism for both models, the associative character is more pronounced for model 2. Figure 4b,c shows examples of transition states obtained with the two models. The distance distribution obtained with the original GW model is intermediate between those of the two new models.

Figure 3. Free energy surfaces (in kcal mol−1) obtained from the metadynamics carried out with (top) the original GW model, (middle) model 1, and (bottom) model 2. Coordination numbers are between the uranium ion and the oxygen ions of the spectator water molecules, and distances are between the uranium ion and the oxygen ion of the outgoing water molecule. Also shown are the corresponding average paths calculated with the transition path sampling approach.

This analysis helps explain the variation in activation free energies reported in Table 5 and discussed in this work. As the uranium and uranyl oxygen partial charges are increased from the GW model to model 1, the uranyl−water interactions strengthen, resulting in an increase in activation free energy for water exchange and a decrease in the water exchange rate. However, as the uranium and uranyl oxygen partial charges are 6428



Article

Figure 4. (a) Probability distribution of the distance between uranium and the oxygen of the incoming water molecule for the original GW model and the two models derived in this work. (b,c) Examples of transition state configurations (only the water molecules in the first hydration shell are shown for clarity) obtained with models (b) 1 and (c) 2. Uranium is shown in blue, oxygen in red, and hydrogen in white. The oxygens of the incoming and outgoing water molecules are shown in yellow and green, respectively.

Notes

increased further from model 1 to model 2, the attractive uranium−water oxygen and uranyl oxygen−water hydrogen interactions in the transition state configurations also increase further, stabilizing the transition state and thus reducing the activation free energy for water exchange, which explains the increase in water exchange rate for model 2 with respect to model 1.

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors acknowledge Dr. Patrick Nichols for providing the radial distribution functions from his ab initio molecular dynamics simulation of the uranyl(VI) cation in water. This research was supported by the U.S. Department of Energy (DOE) Biological and Environmental Research (BER) Division through the Subsurface Biogeochemistry Research (SBR) Program of the Science Focus Area (SFA) at Pacific Northwest National Laboratory (PNNL). The computer simulations were performed in part using the Molecular Science Computing (MSC) facilities in the William R. Wiley Environmental Molecular Sciences Laboratory (EMSL), a national scientific user facility sponsored by the DOE’s Office of Biological and Environmental Research (OBER) and located at Pacific Northwest National Laboratory (PNNL). PNNL is operated for the DOE by Battelle Memorial Institute under Contract DE-AC05-76RL01830.

■

CONCLUSIONS The potential parameters of the uranyl model of Guilbaud and Wipff were modified to strengthen the uranyl−water interactions following the observation that this model underestimates the experimental water exchange rate and hydration free energy. Two set of potential parameters were derived by increasing the uranium ion partial charge from 2.50 e in the GW model to 3.25 e (model 1) or 3.50 e (model 2) and readjusting the uranium ε and r* parameters. Both models showed improvement over the GW model with respect to the strength of the uranyl−water interactions, as indicated by the water binding energies of uranyl hydrate gas-phase clusters and the hydration free energies from thermodynamic integration calculations. Model 1 showed slightly shorter uranium−water oxygen distances (approximately −2%) than found experimentally from XAFS, XRD, and XS measurements and from AIMD; however, model 1 reproduces well the hydration shell structure of the uranyl oxygens and, importantly, yields good agreement with the experimental water exchange rate and hydration free energy. Model 2 gives uranium−water oxygen distances that are more consistent with experimental and quantum mechanical data and also yields a hydration free energy within the uncertainty of the experiment data, but the extent of hydrogen bonding between water molecules and uranyl oxygens is likely overestimated, which causes a decrease in the activation free energy for water exchange. Overall, model 1 performs best and will be used in our future work aimed at investigating the energetics of uranyl surface complexation species at mineral surfaces relevant to contaminated subsurface sediments.

■

■

REFERENCES

(1) Zheng, H.; Singh, A.; Basak, S.; Ulrich, K.-U.; Sahu, M.; Biswas, P.; Catalano, J. G.; Giammar, D. E. Nanoscale Size Effects on Uranium(VI) Adsorption to Hematite. Environ. Sci. Technol. 2009, 43, 1373−1378. (2) Rossberg, A.; Ulrich, K.-U.; Weiss, S.; Tsushima, S.; Hiemstra, T.; Scheinost, A. C. Identification of Uranyl Surface Complexes on Ferrihydrite: Advanced EXAFS Data Analysis and CD-MUSIC Modeling. Environ. Sci. Technol. 2009, 43, 1400−1406. (3) Catalano, J. G.; Heald, S. M.; Zachara, J. M.; Brown, G. E., Jr. Spectroscopic and Diffraction Study of Uranium Speciation in Contaminated Vadose Zone Sediments from the Hanford Site, Washington State. Environ. Sci. Technol. 2004, 38, 2822−2828. (4) Froideval, A.; Del Nero, M.; Gaillard, C.; Barillon, R.; Rossini, I.; Hazemann, J. L. Uranyl Sorption Species at Low Coverage on AlHydroxide: TRFLS and XAFS Studies. Geochim. Cosmochim. Acta 2006, 70, 5270−5284. (5) Malin, J. N.; Holland, J. G.; Saslow, S. A.; Geiger, F. M. U(VI) Adsorption and Speciation at the Acidic Silica/Water Interface Studied by Resonant and Nonresonant Second Harmonic Generation. J. Phys. Chem. C 2011, 115, 13353−13360. (6) Bargar, J. R.; Reitmeyer, R.; Lenhart, J. J.; Davis, J. A. Characterization of U(VI)-Carbonato Ternary Complexes on

AUTHOR INFORMATION

Corresponding Author

*Tel.: 509-371-6382. E-mail: [email protected]. 6429



Article

Hematite: EXAFS and Electrophoretic Mobility Measurements. Geochim. Cosmochim. Acta 2000, 64, 2737−2749. (7) Chardon, E.; Bosback, D.; Bryan, N. D.; Lyon, I. C.; Marquardt, C.; Römer, J.; Schild, D.; Vaughan, D. J.; Wincott, P. L.; Wogelius, R. A.; Livens, F. R. Reactions of the Feldspar Surface with Metal Ions: Sorption of Pb(II), U(VI) and Np(V), and Surface Analytical Studies of Reaction with Pb(II) and U(VI). Geochim. Cosmochim. Acta 2008, 72, 288−297. (8) Liu, C.; Shang, J.; Kerisit, S.; Zachara, J. M.; Zhu, W. ScaleDependent Rates of Uranyl Surface Complexation Reaction in Sediments. Geochim. Cosmochim. Acta 2013, 105, 326−341. (9) McKinley, J. P.; Zachara, J. M.; Liu, C. X.; Heald, S. C.; Prenitzer, B. I.; Kempshall, B. W. Microscale Controls on the Fate of Contaminant Uranium in the Vadose Zone, Hanford Site, Washington. Geochim. Cosmochim. Acta 2006, 70, 1873−1887. (10) Liu, C.; Zachara, J. M.; Qafoku, O.; McKinley, J. P.; Heald, S. M.; Wang, Z. Dissolution of Uranyl Microprecipitates in Subsurface Sediments at Hanford Site, USA. Geochim. Cosmochim. Acta 2004, 68, 4519−4537. (11) Boily, J.-F.; Rosso, K. M. Crystallographic Controls on Uranyl Binding at the Quartz/Water Interface. Phys. Chem. Chem. Phys. 2011, 13, 7845−7851. (12) Steele, H. M.; Wright, K.; Hillier, I. H. Modelling the Adsorption of Uranyl on the Surface of Goethite. Geochim. Cosmochim. Acta 2002, 66, 1305−1310. (13) Greathouse, J. A.; Cygan, R. T. Molecular Dynamics Simulation of Uranyl(VI) Adsorption Equilibria onto an External Montmorillonite Surface. Phys. Chem. Chem. Phys. 2005, 7, 3580−3586. (14) Greathouse, J. A.; Cygan, R. T. Water Structure and Aqueous Uranyl(VI) Adsorption Equilibria onto External Surfaces of Beidellite, Montmorillonite, and Pyrophyllite: Results from Molecular Simulations. Environ. Sci. Technol. 2006, 40, 3865−3871. (15) Greathouse, J. A.; O’Brien, R. J.; Bemis, G.; Pabalan, R. T. Molecular Dynamics Study of Aqueous Uranyl Interactions with Quartz (010). J. Phys. Chem. B 2002, 106, 1646−1655. (16) Zaidan, O. F.; Greathouse, J. A.; Pabalan, R. T. Monte Carlo and Molecular Dynamics Simulation of Uranyl Adsorption on Montmorillonite Clay. Clays Clay Miner. 2003, 51, 372−381. (17) Doudou, S.; Vaughan, D. J.; Livens, F. R.; Burton, N. A. Atomistic Simulations of Calcium Uranyl(VI) Carbonate Adsorption on Calcite and Stepped-Calcite Surfaces. Environ. Sci. Technol. 2012, 46, 7587−7594. (18) Kerisit, S.; Liu, C. Diffusion and Adsorption of Uranyl Carbonate Species in Nanosized Mineral Fractures. Environ. Sci. Technol. 2012, 46, 1632−1640. (19) Greathouse, J. A.; Stellalevinsohn, H. R.; Denecke, M. A.; Bauer, A.; Pabalan, R. T. Uranyl Surface Complexes in a Mixed-Charge Montmorillonite: Monte Carlo Computer Simulation and Polarized XAFS Results. Clays Clay Miner. 2005, 53, 278−286. (20) Kremleva, A.; Krüger, S.; Rösch, N. Quantum Chemical Modeling of Uranyl Adsorption on Mineral Surfaces. Radiochim. Acta 2010, 98, 635−646. (21) Kremleva, A.; Krüger, S.; Rösch, N. Density Functional Model Studies of Uranyl Adsorption on (001) Surfaces of Kaolinite. Langmuir 2008, 24, 9515−9524. (22) Kremleva, A.; Martorell, B.; Krüger, S.; Rösch, N. Uranyl Adsorption on Solvated Edge Surfaces of Pyrophyllite: A DFT Model Study. Phys. Chem. Chem. Phys. 2012, 14, 5815−5823. (23) Sherman, D. M.; Peacock, C. L.; Hubbard, C. G. Surface Complexation of U(VI) on Goethite (α-FeOOH). Geochim. Cosmochim. Acta 2008, 72, 298−310. (24) Skomurski, F. N.; Ilton, E. S.; Engelhard, M. H.; Arey, B. W.; Rosso, K. M. Heterogeneous Reduction of U6+ by Structural Fe2+ from Theory and Experiment. Geochim. Cosmochim. Acta 2011, 75, 7277− 7290. (25) Glezakou, V.-A.; deJong, W. A. Cluster Models for Uranyl(VI) Adsorption on α-Alumina. J. Phys. Chem. A 2011, 115, 1257−1263.

(26) Drot, R.; Roques, J.; Simoni, E. Molecular Approach of the Uranyl/Mineral Interfacial Phenomena. C. R. Chim. 2007, 10, 1078− 1091. (27) Lectez, S.; Roques, J.; Salanne, M.; Simoni, E. Car−Parrinello Molecular Dynamics Study of the Uranyl Behaviour at the Gibbsite/ Water Interface. J. Chem. Phys. 2012, 137, 154705. (28) Guilbaud, P.; Wipff, G. Hydration of UO22+ Cation and Its NO3− and 18-Crown-6 Adducts Studied by Molecular-Dynamics Simulations. J. Phys. Chem. 1993, 97, 5685−5692. (29) Guilbaud, P.; Wipff, G. Force Field Representation of the UO22+ Cation from Free Energy MD Simulations in Water. Tests on Its 18Crown-6 and NO3− Adducts, and on Its Calix[6]arene(6-) and CMPO Complexes. J. Mol. Struct. (THEOCHEM) 1996, 366, 55−63. (30) Jayasinghe, M.; Beck, T. L. Molecular Dynamics Simulations of the Structure and Thermodynamics of Carrier-Assisted Uranyl Ion Extraction. J. Phys. Chem. B 2009, 113, 11662−11671. (31) Ye, X.; Cui, S.; de Almeida, V.; Khomami, B. Interfacial Complex Formation in Uranyl Extraction by Tributyl Phosphate in Dodecane Diluent: A Molecular Dynamics Study. J. Phys. Chem. B 2009, 113, 9852−9862. (32) Bühl, M.; Wipff, G. Insights into Uranyl Chemistry from Molecular Dynamics Simulations. ChemPhysChem 2011, 12, 3095− 3105. (33) Pible, O.; Guilbaud, P.; Pellequer, J.-L.; Vidaud, C.; Quéméneur, E. Structural Insights into Protein−Uranyl Interaction: Towards an in Silico Detection Method. Biochimie 2006, 88, 1631−1638. (34) Lins, R. D.; Vorpagel, E. R.; Guglielmi, M.; Straatsma, T. P. Computer Simulation of Uranyl Uptake by the Rough Lipopolysaccharide Membrane of Pseudomonas Aeruginosa. Biomacromolecules 2008, 9, 29−35. (35) Lin, Y.-W.; Liao, L.-F. Probing Interactions between Uranyl Ions and Lipid Membrane by Molecular Dynamics Simulation. Comput. Theor. Chem. 2011, 976, 130−134. (36) Kerisit, S.; Liu, C. X. Molecular Simulation of the Diffusion of Uranyl Carbonate Species in Aqueous Solution. Geochim. Cosmochim. Acta 2010, 74, 4937−4952. (37) Doudou, S.; Arumugam, K.; Vaughan, D. J.; Livens, F. R.; Burton, N. A. Investigation of Ligand Exchange Reactions in Aqueous Uranyl Carbonate Complexes Using Computational Approaches. Phys. Chem. Chem. Phys. 2011, 13, 11402−11411. (38) Ye, X.; Smith, R. B.; Cui, S.; de Almeida, V.; Khomami, B. Influence of Nitric Acid on Uranyl Nitrate Association in Aqueous Solutions: A Molecular Dynamics Simulation Study. Solvent Extr. Ion Exch. 2010, 28, 1−18. (39) Szabó, Z.; Glaser, J.; Grenthe, I. Kinetics of Ligand Exchange Reactions for Uranyl(2+) Fluoride Complexes in Aqueous Solution. Inorg. Chem. 1996, 35, 2036−2044. (40) Farkas, I.; Bányai, I.; Szabó, Z.; Grenthe, I. Rates and Mechanisms of Water Exchange of UO22+(aq) and UO2(oxalate)F(H2O)2−: A Variable-Temperature 17O and 19F NMR Study. Inorg. Chem. 2000, 39, 799−805. (41) Moskaleva, L. V.; Krüger, S.; Spörl, A.; Rösch, N. Role of Solvation in the Reduction of the Uranyl Dication by Water: A Density Functional Study. Inorg. Chem. 2004, 43, 4080−4090. (42) Gibson, J. K.; Haire, R. G.; Santos, M.; Marçalo, J.; Pires de Matos, A. Oxidation Studies of Dipositive Actinide Ions, An2+ (An = Th, U, Np, Pu, Am) in the Gas Phase: Synthesis and Characterization of the Isolated Uranyl, Neptunyl, and Plutonyl Ions UO22+(g), NpO22+(g), and PuO22+(g). J. Phys. Chem. A 2005, 109, 2768−2781. (43) Shamov, G. A.; Schreckenbach, G. Density Functional Studies of Actinyl Aquo Complexes Studied Using Small-Core Effective Core Potentials and a Scalar Four-Component Relativistic Method. J. Phys. Chem. A 2005, 109, 10961−10974. (44) Shamov, G. A.; Schreckenbach, G. Density Functional Studies of Actinyl Aquo Complexes Studied Using Small-Core Effective Core Potentials and a Scalar Four-Component Relativistic Method. J. Phys. Chem. A 2006, 110, 12072. 6430



Article

(45) Gutowski, K. E.; Dixon, D. A. Predicting the Energy of Water Exchange Reaction and Free Energy of Solvation for the Uranyl Ion in Aqueous Solution. J. Phys. Chem. A 2006, 110, 8840−8856. (46) Rai, N.; Tiwari, S. P.; Maginn, E. J. Force Field Development for Actinyl Ions via Quantum Mechanical Calculations: An Approach to Account for Many Body Solvation Effects. J. Phys. Chem. B 2012, 116, 10885−10897. (47) Smith, W.; Forester, T. R. DL_POLY_2.0: A General Purpose Parallel Molecular Dynamics Simulation Package. J. Mol. Graphics 1996, 14, 136−141. (48) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695−1697. (49) Melchionna, S.; Ciccotti, G.; Holian, B. L. Hoover NPT Dynamics for Systems Varying in Shape and Size. Mol. Phys. 1993, 78, 533−544. (50) Ewald, P. P. Die Berechnung optischer und elektrostatischer Gitterpotentiale. Ann. Phys. 1921, 64, 253−287. (51) Rey, R.; Hynes, J. T. Hydration Shell Exchange Kinetics: An MD Study for Na+(aq). J. Phys. Chem. 1996, 100, 5611−5615. (52) Carter, E. A.; Ciccotti, G.; Hynes, J. T.; Kapral, R. Constrained Reaction Coordinate Dynamics for the Simulation of Rare Events. Chem. Phys. Lett. 1989, 156, 472−476. (53) Horinek, D.; Mamatkulov, S. I.; Netz, R. R. Rational Design of Ion Force Fields Based on Thermodynamic Solvation Properties. J. Chem. Phys. 2009, 130, 124507. (54) Hummer, G.; Pratt, L. R.; Garcia, A. E. Ion Sizes and Finite-Size Corrections for Ionic-Solvation Free Energies. J. Chem. Phys. 1997, 107, 9275−9277. (55) Hummer, G.; Pratt, L. R.; Garcia, A. E.; Berne, B. J.; Rick, S. W. Electrostatic Potentials and Free Energies of Solvation of Polar and Charged Molecules. J. Phys. Chem. B 1997, 101, 3017−3020. (56) Kastenholz, M. A.; Hünenberg, P. H. Computation of Methodology-Independent Ionic Solvation Free Energies from Molecular Simulations. I. The Electrostatic Potential in Molecular Liquids. J. Chem. Phys. 2006, 124, 124106. (57) Kastenholz, M. A.; Hünenberg, P. H. Computation of Methodology-Independent Ionic Solvation Free Energies from Molecular Simulations. II. The Hydration Free Energy of the Sodium Cation. J. Chem. Phys. 2006, 124, 224501. (58) Raiteri, P.; Gale, J. D.; Quigley, D.; Rodger, P. M. Derivation of an Accurate Force Field for Simulating the Growth of Calcium Carbonate from Aqueous Solution: A New Model for the Calcite− Water Interface. J. Phys. Chem. C 2010, 114, 5997−6010. (59) Laio, A.; Parrinello, M. Escaping Free-Energy Minima. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 12552−12556. (60) Laio, A.; Rodriguez-Fortea, A.; Gervasio, F. L.; Ceccarelli, M.; Parrinello, M. Assessing the Accuracy of Metadynamics. J. Phys. Chem. B 2005, 109, 6714−6721. (61) Laio, A.; Parrinello, M. Computing Free Energies and Accelerating Rare Events with Metadynamics. Lect. Notes Phys. 2006, 703, 315−347. (62) Iannuzzi, M.; Laio, A.; Parrinello, M. Efficient Exploration of Reactive Potential Energy Surfaces Using Car−Parrinello Molecular Dynamics. Phys. Rev. Lett. 2003, 90, 238302. (63) Dellago, C.; Bolhuis, P. G.; Csajka, F. S.; Chandler, D. Transition Path Sampling and the Calculation of Rate Constants. J. Chem. Phys. 1998, 108, 1964−1977. (64) Bolhuis, P. G.; Chandler, D.; Dellago, C.; Geissler, P. L. Transition Path Sampling: Throwing Ropes over Rough Mountain Passes, in the Dark. Annu. Rev. Phys. Chem. 2002, 53, 291−318. (65) Dellago, C.; Bolhuis, P. G.; Chandler, D. Efficient Transition Path Sampling: Application to Lennard-Jones Cluster Rearrangements. J. Chem. Phys. 1998, 108, 9236−9245. (66) Bolhuis, P. G.; Dellago, C.; Chandler, D. Sampling Ensembles of Deterministic Transition Pathways. Faraday Discuss. 1998, 110, 421− 436. (67) Dellago, C.; Bolhuis, P. G.; Chandler, D. On the Calculation of Reaction Rate Constants in the Transition Path Ensemble. J. Chem. Phys. 1999, 110, 6617−6625.

(68) Kerisit, S.; Rosso, K. M. Transition Path Sampling of Water Exchange Rates and Mechanisms around Aqueous Ions. J. Chem. Phys. 2009, 131, 114512. (69) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The Missing Term in Effective Pair Potentials. J. Phys. Chem. 1987, 91, 6269−6271. (70) Spencer, S.; Gagliardi, L.; Handy, N. C.; Ioannou, A. G.; Skylaris, C.-K.; Willetts, A. Hydration of UO22+ and PuO22+. J. Phys. Chem. A 1999, 103, 1831−1837. (71) Cao, Z.; Balasubramanian, K. Theoretical Studies of UO2(H2O)n2+, NpO2(H2O)n+, and PuO2(H2O)n2+ Complexes (n = 4−6) in Aqueous Solution and Gas Phase. J. Chem. Phys. 2005, 123, 114309. (72) Nichols, P.; Bylaska, E. J.; Schenter, G. K.; de Jong, W. Equatorial and Apical Solvent Shells of the UO22+ Ion. J. Chem. Phys. 2008, 128, 124507. (73) Bühl, M.; Diss, R.; Wipff, G. Coordination Environment of Aqueous Uranyl(VI) Ion. J. Am. Chem. Soc. 2005, 127, 13506−13507. (74) Bühl, M.; Kabrede, H.; Diss, R.; Wipff, G. Effect of Hydration on Coordination Properties of Uranyl(VI) Complexes. A FirstPrinciples Molecular Dynamics Study. J. Am. Chem. Soc. 2006, 128, 6357−6368. (75) Frick, R. J.; Hofer, T. S.; Pribil, A. B.; Randolf, B. R.; Rode, B. M. Structure and Dynamics of the UO22+ Ion in Aqueous Solution: An Ab Initio QMCF MD Study. J. Phys. Chem. A 2009, 113, 12496−12503. (76) Hagberg, D.; Karlstrom, G.; Ross, B. O.; Gagliardi, L. The Coordination of Uranyl in Water: A Combined Quantum Chemical and Molecular Simulation Study. J. Am. Chem. Soc. 2005, 127, 14250− 14256. (77) Åberg, M.; Ferri, D.; Glaser, J.; Grenthe, I. Structure of the Hydrated Dioxouranium(VI) Ion in Aqueous Solution. An X-ray Diffraction and 1H NMR Study. Inorg. Chem. 1983, 22, 3986−3989. (78) Allen, P. G.; Bucher, J. J.; Shuh, D. K.; Edelstein, N. M.; Reich, T. Investigation of Aquo and Chloro Complexes of UO22+, NpO2+, Np4+, and Pu3+ by X-ray Absorption Fine Structure Spectroscopy. Inorg. Chem. 1997, 36, 4676−4683. (79) Thompson, H. A.; Brown, G. E., Jr.; Parks, G. A. XAFS Spectroscopic Study of Uranyl Coordination in Solids an Aqueous Solution. Am. Mineral. 1997, 82, 483−496. (80) Wahlgren, U.; Moll, H.; Grenthe, I.; Schimmelpfennig, B.; Maron, L.; Vallet, V.; Gropen, O. Structure of Uranium(VI) in Strong Alkaline Solutions. A Combined Theoretical and Experimental Investigation. J. Phys. Chem. A 1999, 103, 8257−8264. (81) Sémon, L.; Boehme, C.; Billard, I.; Hennig, C.; Lützenkirchen, K.; Reich, T.; Roßberg, A.; Rossini, I. W.; Wipff, G. Do Perchlorate and Triflate Anions Bind to the Uranyl Cation in an Acidic Aqueous Medium? A Combined EXAFS and Quantum Mechanical Investigation. ChemPhysChem 2001, 2, 591−598. (82) Neuefeind, J.; Soderholm, L.; Skanthakumar, S. Experimental Coordination Environment of Uranyl(VI) in Aqueous Solution. J. Phys. Chem. A 2004, 108, 2733−2739. (83) Hennig, C.; Tutschku, J.; Rossberg, A.; Bernhard, G.; Scheinost, A. C. Comparative EXAFS Investigation of Uranium(VI) and -(IV) Aquo Chloro Complexes in Solution Using a Newly Developed Spectroelectrochemical Cell. Inorg. Chem. 2005, 44, 6655−6661. (84) Soderholm, L.; Skanthakumar, S.; Neuefeind, J. Determination of Actinide Speciation in Solution Using High-Energy X-ray Scattering. Anal. Bioanal. Chem. 2005, 383, 48−55. (85) Marx, G.; Bischoff, H. Transport Processes of Actinides in Electrolyte Solutions. I. Determination of Ionic Mobilities of Uranium in Aqueous Solutions at 25 °C by the Radioisotope Method. J. Radioanal. Chem. 1976, 30, 567−581. (86) Marcus, Y. Some Thermodynamic Data Concerning the Dioxouranium(VI) Ion and Its Compounds and Reactions. J. Inorg. Nucl. Chem. 1975, 37, 493−501. (87) Tissandier, M. D.; Cowen, K. A.; Feng, W. Y.; Gundlach, E.; Cohen, M. H.; Earhart, A. D.; Coe, J. V.; Tuttle, T. R. J. The Proton’s Absolute Aqueous Enthalpy and Gibbs Free Energy of Solvation from Cluster-Ion Solvation Data. J. Phys. Chem. A 1998, 102, 7787−7794. 6431



Article

(88) Langford, C. H.; Gray, H. B. Ligand Substitution Processes; W. A. Benjamin, Inc.: New York, 1965. (89) Vallet, V.; Wahlgren, U.; Schimmelpfennig, B.; Szabó, Z.; Grenthe, I. The Mechanism of Water Exchange in [UO2(H2O)5]2+ and [UO2(oxalate)2(H2O)]2−, as Studied by Quantum Chemical Methods. J. Am. Chem. Soc. 2001, 123, 11999−12008. (90) Rotzinger, F. P. The Water-Exchange Mechanism of the [UO2(OH2)5]2+ Ion Revisited: The Importance of a Proper Treatment of Electron Correlation. Chem.Eur. J. 2007, 13, 800−811. (91) Vallet, V.; Wahlgren, U.; Grenthe, I. Comment on “The WaterExchange Mechanism of the [UO2(OH2)5]2+ Ion Revisited: The Importance of a Proper Treatment of Electron Correlation”. Chem. Eur. J. 2007, 13, 10294−10297. Rotzinger, F. P. Chem.Eur. J. 2007, 13, 800. (92) Rotzinger, F. P. Reply to the comment on “The Water-Exchange Mechanism of the [UO2(OH2)5]2+ Ion Revisited: The Importance of a Proper Treatment of Electron Correlation”. Chem.Eur. J. 2007, 13, 10298−10302. Rotzinger, F. P. Chem.Eur. J. 2007, 13, 800. (93) Tsushima, S. Hydration and Water-Exchange Mechanism of the UO22+ Ion Revisted: The Validity of the “n + 1” Model. J. Phys. Chem. A 2007, 111, 3613−3617. (94) Wåhlin, P.; Danilo, C.; Vallet, V.; Réal, F.; Flament, J. P.; Wahlgren, U. An Investigation of the Accuracy of Different DFT Functionals on the Water Exchange Reaction in Hydrated Uranyl(VI) in the Ground State and First Excited State. J. Chem. Theory Comput. 2008, 4, 569−577. (95) Wåhlin, P.; Schimmelpfennig, B.; Wahlgren, U.; Grenthe, I.; Vallet, V. On the Combined Used of Discrete Solvent Models and Continuum Descriptions of Solvent Effects in Ligand Exchange Reactions: A Case Study of the Uranyl(VI) Aquo Ion. Theor. Chem. Acc. 2009, 124, 377−384. (96) Bühl, M.; Kabrede, H. Mechanism of Water Exchange in Aqueous Uranyl(VI) Ion. A Density Functional Molecular Dynamics Study. Inorg. Chem. 2006, 45, 3834−3836. (97) Atta-Fynn, R.; Bylaska, E. J.; de Jong, W. A. Free Energies and Mechanisms of Water Exchange around Uranyl from First Principles Molecular Dynamics. MRS Proc. 2012, 1383, mrsf11-1383-a07-06.

6432


Structure, Kinetics, and Thermodynamics of the ... - ACS Publications

Recommend Documents