Molecular Dynamics Simulations of Uranyl and Uranyl Carbonate

Article pubs.acs.org/est

Molecular Dynamics Simulations of Uranyl and Uranyl Carbonate Adsorption at Aluminosilicate Surfaces Sebastien Kerisit* and Chongxuan Liu Physical Sciences Division, Pacific Northwest National Laboratory, Richland, Washington 99352, United States S Supporting Information *

ABSTRACT: Adsorption at mineral surfaces is a critical factor controlling the mobility of uranium(VI) in aqueous environments. Therefore, molecular dynamics (MD) simulations were performed to investigate uranyl(VI) adsorption onto two neutral aluminosilicate surfaces, namely, the orthoclase (001) surface and the octahedral aluminum sheet of the kaolinite (001) surface. Although uranyl preferentially adsorbs as a bidentate inner-sphere complex on both surfaces, the free energy of adsorption on the orthoclase surface (−15 kcal mol−1) is significantly more favorable than that at the kaolinite surface (−3 kcal mol−1), which is attributed to differences in surface functional groups and the ability of the orthoclase surface to release a surface potassium ion upon uranyl adsorption. The structures of the adsorbed complexes compare favorably with X-ray absorption spectroscopy results. Simulations of the adsorption of uranyl complexes with up to three carbonate ligands revealed that uranyl complexes coordinated to up to two carbonate ions are stable on the orthoclase surface whereas uranyl carbonate surface complexes are unfavored at the kaolinite surface. Combining the MD-derived equilibrium adsorption constants for orthoclase with aqueous equilibrium constants for uranyl carbonate species indicates the presence of adsorbed uranium complexes with one or two carbonates under alkaline conditions, in support of current uranium(VI) surface complexation models.

■

classical force fields (also termed potential models) allow the exploration of time and length scales sufficient for extracting not only structural information on uranyl surface complexes but also on the relative free energies of various adsorbed and dissolved uranyl species. The accuracy of such molecular-scale models is dependent on the potential parameters chosen to carry out the simulations. In previous work,27−30 we validated and applied potential parameters for modeling of water structure at feldspar surfaces, diffusion of uranyl carbonate species in aqueous solution, and diffusion and adsorption of uranyl and uranyl carbonate species in feldspar fractures. In particular, in a recent publication,31 we modified the potential parameters of Guilbaud and Wipff32,33 for simulating the uranyl ion in water to improve the agreement of the model with the experimental hydration free energy and water exchange rate. This improved potential model is expected to lead to more accurate predictions of the structure and thermodynamics of uranyl adsorption. In this work, we applied this improved potential model to the adsorption of uranyl and uranyl carbonate species at the kaolinite (001) and orthoclase (001) surfaces. Orthoclase

INTRODUCTION Uranium is a major contaminant at nuclear facilities across the U.S. Department of Energy (DOE) complex as a result of its use as a nuclear fuel for weapons production.1 In its oxidized form, U(VI), uranium is soluble and relatively mobile as the uranyl ion and therefore is of great concern for the environment, particularly in groundwater systems. Adsorption at mineral surfaces is a key factor affecting the mobility of U(VI) in the environment. Consequently, the ability to predict the fate of uranium in contaminated subsurface environments is dependent on the development of reactive transport models that include an accurate description of U(VI) surface complexation reactions at mineral surfaces. This has led to a sustained experimental effort2−11 to determine the nature and structure of U(VI) surface complexes using spectroscopic methods such as extended X-ray absorption fine structure spectroscopy (EXAFS) and time-resolved laser-induced fluorescence spectroscopy (TRLFS). However, interpretation of the experimental data can be complicated, for example, in cases where the absorber atom in EXAFS measurements is present in a mixture of different surface species. Therefore, this experimental effort has benefited from the application of computational modeling techniques, both classical12−18 and quantum-mechanical,19−26 to the investigation of uranyl adsorption on various mineral surfaces. In particular, molecular dynamics (MD) simulations based on © 2014 American Chemical Society

Received: Revised: Accepted: Published: 3899

December 3, 2013 February 18, 2014 March 1, 2014 March 3, 2014 dx.doi.org/10.1021/es405387c | Environ. Sci. Technol. 2014, 48, 3899−3907

Environmental Science & Technology

Article

these simulations, atoms are represented as point-charge particles that interact via 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 intramolecular angle bending. The short-range potential parameters and ionic charges used in this study are those of the SPC/E model43 for water, CLAYFF44 for orthoclase and kaolinite, Pavese et al.45 for the carbonate ion, and our modified version31 of the Guilbaud and Wipff model32,33 for uranyl. All of the parameters are reported in Tables SI−SVII in the Supporting Information and are the same as those employed in our previous publications on the orthoclase−water−uranyl−carbonate system27−30,40 with the exception of two sets of changes described in detail in the Supporting Information. The association of uranyl and carbonate ions was simulated in bulk water using simulation cells containing one uranyl carbonate species and 507−510 water molecules, depending on the size of the species of interest. Pre-equilibrated simulation cells were available from previous work.29 Uranyl adsorption was investigated (1) at the kaolinite (001)−water interface and (2) at the orthoclase (001)−water interface. Kaolinite is a layered aluminosilicate consisting of two-sheet layers bound by hydrogen bonds. Each layer consists of a tetrahedral silicon sheet and an octahedral aluminum sheet linked through bridging oxygen atoms. The (001) surface was considered because it is the cleavage plane of kaolinite and the dominant surface of kaolinite crystals. Density functional theory (DFT) calculations20 indicated that adsorption on the octahedral sheet is energetically favored over adsorption on the tetrahedral sheet. Accordingly, we focused solely on adsorption on the octahedral sheet. The perfect kaolinite (001) surface is polar.46 To cancel the surface dipole in our MD simulations, we adopted a configuration wherein half of the kaolinite layers were inverted (Figure 1a) so as to form a mirror plane parallel to the plane of the surface in the middle of the kaolinite slab. The octahedral aluminum sheet is exposed on both sides of the slab as a result. A change of basis was also employed to convert the simulation supercell to an orthogonal cell, as was done by Vasconcelos et al.47 The kaolinite slab was four layers thick and had surface dimensions of 25.75 Å × 26.80 Å. Orthoclase is the potassium end member of alkali feldspar, a framework aluminosilicate. Corner-sharing AlO4 and SiO4 tetrahedra form a three-dimensional lattice within which potassium ions fills a cavity to provide charge compensation for the substitution of Al for Si. The orthoclase (001) surface slab was generated as previously27 and was thus terminated by a full layer of hydroxyl groups. The surface dimensions were 25.69 Å × 25.93 Å. To create the interfaces, the direction normal to the surface was extended to create a gap of ∼25 Å between opposing surfaces that was filled with 500 or 556 water molecules for kaolinite or orthoclase, respectively (Figure 1). To equilibrate the interfaces, a 250 ps MD simulation was first run with the mineral slab fixed, and this was followed by a 1 ns MD simulation without any constraints before uranyl and uranyl carbonate species were introduced for the potential of mean force (PMF) calculations. One, two, three, or four water molecules were removed when UO 2 2 + , UO 2 CO 3 , UO2(CO3)22−, or UO2(CO3)34−, respectively, was inserted. The plane of the surface was used as the zero of distance and

(potassium feldspar, KAlSi3O8) was selected because feldspar has been shown to be an important mineral component of contaminated sediments at the DOE Hanford Site that preferentially sequesters U(VI) within its microsized fractures.34,35 Kaolinite [1:1 clay, Al2Si2O5(OH)4] is a common clay component of soils and is also one of the simplest clay minerals, with one octahedral aluminum sheet and one tetrahedral silicon sheet linked via hydrogen bonds. Electronic structure calculations20 have shown that adsorption on the octahedral aluminum sheet is energetically preferred over adsorption on the tetrahedral silicon sheet; therefore, the work described herein focused on adsorption on the octahedral aluminum sheet. Additionally, these two minerals exhibit two very different types of aluminosilicate surfaces: the orthoclase (001) surface shows only singly coordinated hydroxyl groups bonded to tetrahedral framework cations, whereas the kaolinite (001) surface shows only doubly coordinated hydroxyl groups bonded to octahedral aluminum atoms (Figure 1). Such a difference allows for an evaluation of the dependence of the thermodynamics of uranyl adsorption on surface structure.

Figure 1. Snapshots of the (a) kaolinite (001)−water and (b) orthoclase (001)−water interfaces. Si atoms are shown in yellow, Al atoms in pink, K atoms in green, O atoms in red, and H atoms in white.

An important factor that controls uranium adsorption on mineral surfaces is the presence of aqueous ligands such as carbonate. Indeed, U(VI) can form stable, soluble complexes with carbonate ions in soil and groundwater environments, which can lead to inhibition of uranyl adsorption. However, uranyl carbonate ternary surface complexes have been detected on mineral surfaces, mainly in EXAFS measurements,36−38,10 and are commonly used in surface complexation models.39−41 Therefore, one objective of the present work was to determine the effects of carbonate complexation on uranyl surface speciation and on the thermodynamics of uranyl adsorption on the two model mineral surfaces.

■

COMPUTATIONAL METHODS All of the calculations performed in this work were carried out with the molecular dynamics computer code DL_POLY.42 In 3900

dx.doi.org/10.1021/es405387c | Environ. Sci. Technol. 2014, 48, 3899−3907


Article

was defined as the plane that passes through the topmost octahedral Al sites or the topmost tetrahedral Si/Al sites for kaolinite or orthoclase, respectively. Surveys were performed to identify whether multiple innersphere adsorption sites/configurations were present at the two surfaces for uranyl and uranyl carbonate species, as described in detail in the Supporting Information. For the kaolinite (001) surface, only one configuration was identified in each case, whereas two configurations were found for adsorption of UO22+ and UO2CO3 at the orthoclase (001) surface. For these cases, the PMF simulations presented here considered only the energetically most stable configuration. All of the simulations were carried out at 298.15 K and zero applied pressure. The mineral−water interface simulations were carried out in the NVT ensemble (constant number of particles, constant volume, constant temperature), whereas the simulations of the uranyl carbonate species in bulk water were performed in the NPT ensemble (constant number of particles, constant pressure, constant temperature). The temperature and pressure were kept constant using the Nosé−Hoover thermostat48 and the Hoover barostat,49 respectively. The electrostatic interactions 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. For charged systems, the net cell charge was neutralized by a uniform background charge density. The Verlet leapfrog integration algorithm was used to integrate the equations of motion with a time step of 1 fs. The geometry of the water molecules was held fixed using the SHAKE algorithm.51

UO2 2 + + CO32 − → UO2 CO3

(1)

UO2 CO3 + CO32 − → UO2 (CO3)2 2 −

(2)

UO2 (CO3)2 2 − + CO32 − → UO2 (CO3)3 4 −

(3)

The equilibrium constant for association of species A and B, Ka, was calculated using the approach of Chialvo et al.,54 in which Ka in the infinite-dilution limit is defined as K a = 4π

∫r

rmax

∞ gAB (r )r 2 d r

min

(4)

where g∞ AB(r) is the radial distribution function in the infinite dilution limit, r is the distance between the centers of mass of A and B, and rmin and rmax define the distance range for which A and B are considered to be associated. g∞ AB(r) is calculated from the PMF WAB(r) according to ∞ WAB(r ) = −kT ln gAB (r )

(5)

To obtain each PMF, a series of simulations was performed in which r was constrained between 2.2 and 12.5 Å. The interval between data points was 0.1 Å, requiring a total of 104 runs for each species pair. At each distance, the mean force was obtained from a 1 ns simulation after a 100 ps equilibration period. The PMF was then obtained by integration of the mean force: WAB(r ) = WAB(r0) −

∫r

0

r

F (r ) d r

(6)

where r0 is the maximum separation distance considered in the simulations (i.e., 12.5 Å). Again, following Chialvo et al.,54 WAB(r0) was set to the continuum limiting value: q q WAB(r0) = A B εr0 (7)

■

RESULTS AND DISCUSSION This section is divided into four subsections describing the following studies: (1) the association constants of uranyl carbonate species in bulk water were calculated and compared to experimental data to evaluate the ability of the potential model to reproduce the energetics of the uranyl−carbonate− water system; (2) free energy calculations of uranyl adsorption on the (001) surfaces of kaolinite and orthoclase were performed to compare and contrast the thermodynamics of uranyl adsorption on two different types of aluminosilicate surfaces; (3) the adsorption of uranyl carbonate species on the two surfaces was studied to determine the effect of the presence of carbonate ions on the uranyl surface speciation; and (4) the results of the free energy calculations were used to estimate the U(VI) surface speciation as a function of pH. Uranyl−Carbonate Association in Aqueous Solution. As a reference, the hydration free energies of the isolated uranyl and carbonate ions were first calculated and compared to experimental data. The hydration free energies were obtained using the thermodynamic integration technique, and details of the calculation procedure can be found in a previous publication.31 Previous work on uranyl hydration31 showed that our model yielded good agreement with the experimental hydration free energy (−371 kcal mol−1, compared to −369 ± 15 kcal mol−1 based on the data of Gibson et al.52). When applied to the carbonate ion, the same approach yielded a calculated hydration free energy of −320 kcal mol−1, which is consistent with the experimental hydration free energy of −314 kcal mol−1 derived by Marcus.53 Next, the association constants for the reactions shown in eqs 1−3 below were calculated from PMF simulations:

where ε is the SPC/E dielectric constant under the same conditions (70.755) and qA and qB are the net charges on A and B, respectively. The PMFs are shown in Figure 2. The three reactions show similar PMFs with a very shallow minimum at approximately 7.1−7.4 Å, which corresponds to the solvent-separated ion pair, and a second, deeper minimum at 4.7−4.9 Å, which corresponds to the solvent-shared ion pair, wherein water

Figure 2. Potentials of mean force for the association of a carbonate ion with UO2(CO3)n2−2n, where n = 0, 1, or 2. 3901



Article

ion. Every other surface framework cation is bonded to a hydroxyl. The surface also exposes bridging oxygens. Water molecules in the first orthoclase hydration layer are localized in primarily three sites: above a surface potassium ion, near one surface hydroxyl, and between two surface hydroxyls. As a result, they have slightly different heights, which leads to a broad first peak. Details of the orthoclase (001)−water interface can be found elsewhere.27 The kaolinite surface is more atomically smooth than the orthoclase surface, as it is terminated only by doubly coordinated hydroxyls with all of the surface hydroxyl oxygens found at the same height. Water molecules in the first kaolinite hydration layer donate and/or accept hydrogen bonds with surface hydroxyls as well as between themselves and therefore are found in a narrower range of heights above the surface than for orthoclase. Because only the octahedral aluminum sheet, which has the gibbsite structure, was exposed to water in our simulations, a comparison of the calculated water density profile can be made with the Car−Parrinello MD simulations of the gibbsite− water interface reported by Lectez et al.,26 and these simulations also show a maximum density in the first hydration layer of approximately 4 times that of the bulk, in agreement with our calculations. In Figure 3, the coordination numbers of uranium with water oxygens (OW), hydroxyl oxygens (OH), and, for orthoclase, bridging oxygens (OB) are shown as functions of distance to correlate the PMF minima and maxima with changes in the uranyl coordination environment. The coordination numbers were obtained by integrating the radial distribution function (RDF) between uranium and each oxygen type up to its first minimum. As uranyl approaches the surface, the first minimum is found at 5.0 Å or 4.7 Å above the surface for kaolinite or orthoclase, respectively. At this distance, uranium is still coordinated to five water molecules, and therefore, this minimum corresponds to outer-sphere adsorption. For both orthoclase and kaolinite, there is then an energy barrier associated with the exchange of a water molecule for a surface hydroxyl in the uranyl coordination sphere, which leads to adsorption as a monodentate inner-sphere complex. This free energy barrier is substantially greater for adsorption on kaolinite (6.7 kcal mol−1) than on orthoclase (0.7 kcal mol−1). However, because the reaction coordinate is the distance along the normal to the surface in the PMF calculations, the free energy barriers for adsorption could be misrepresented if they also involve a significant component in the plane of the surface. As the uranyl ion approaches the surface further, it exchanges a second water molecule for a bridging oxygen on orthoclase (Figure 4a) or a second hydroxyl oxygen on kaolinite (Figure 4d). Adsorption on the orthoclase surface leads to much deeper free energy minima than adsorption on the kaolinite surface for both the monodentate (−7.3 vs 0.9 kcal mol−1) and bidentate (−15.2 vs −3.2 kcal mol−1) complexes. This can be explained by two factors: (1) surface potassium ions are observed to be released from the surface upon uranyl adsorption, which causes the development of a local negative charge on the surface, and (2) uranium directly binds to a singly coordinated hydroxyl and a bridging oxygen when adsorbing on the orthoclase surface, whereas it binds only to doubly coordinated hydroxyls on the kaolinite surface. The uranyl free energy profile for the orthoclase (001) surface is similar to that previously derived 30 except that the bidentate inner-sphere complex is now more stable than the

molecules directly coordinated to uranium via their oxygen atoms donate hydrogen bonds to the carbonate oxygen atoms. There are then two contact ion pair configurations, the first one at 3.4−3.6 Å and the second one at 2.9−3.0 Å, which correspond to mono- and bidentate configurations, respectively. The ion-pairing nomenclature used by Chialvo et al.54 and Marcus and Hefter56 is used here. The resulting association constants (log Ka) are 8.2, 6.7, and 4.7 for the reactions in eqs 1−3, respectively. The corresponding experimental values are 9.9, 6.7, and 5.2,57 respectively. The results are encouraging given that association constants in aqueous solutions are notoriously difficult to calculate accurately with potential-based MD simulations.58 Uranyl Adsorption on the Kaolinite and Orthoclase (001) Surfaces. The PMFs of a uranyl ion adsorbing on the (001) surfaces of orthoclase and kaolinite were determined from a series of MD simulations wherein the center of mass of the uranyl ion was constrained at specified distances r above the plane of the surface; r was varied from either 1.6 Å for orthoclase or 2.3 Å for kaolinite to 13.3 Å with 0.1 Å intervals. Each simulation was run for 250 ps after a 100 ps equilibration period. The resulting PMFs are shown in Figure 3 together

Figure 3. (top) Potentials of mean force for the adsorption of UO22+ on the kaolinite (001) and orthoclase (001) surfaces. Also shown are the water density profiles at the two surfaces. (middle, bottom) Coordination numbers of uranium with water oxygens (OW), hydroxyl oxygens (OH), and bridging oxygens (OB) for the orthoclase (middle) and kaolinite (bottom) surfaces.

with the water density profiles at the interfaces with the two minerals. Both water density profiles show a distinct first hydration layer. Integration of this first peak yields a water density of 8.6 molecules nm−2 for both orthoclase and kaolinite; however, the width and magnitude of the two peaks differ significantly. This difference is related to differences in the topography of the two mineral surfaces. At the orthoclase (001) surface, framework cations form 10-membered rings that create surface cavities, each of which is occupied by a potassium 3902



Article

Figure 4. Snapshots of adsorbed uranyl and uranyl carbonate inner-sphere complexes: (a) bidentate UO22+, (b) bidentate UO2CO3, and (c) bidentate UO2(CO3)22− complexes at the orthoclase (001) surface; (d) bidentate UO22+, (e) monodentate UO2CO3, and (f) monodentate UO2(CO3) complexes at the kaolinite (001) surface. Si atoms are shown in yellow, Al atoms in pink, K atoms in green, U atoms in blue, C atoms in grey, O atoms in red, and H atoms in white.

monodentate complex because of the increase in the partial charge on the uranium atom and the fact that in the bidentate configuration the uranium ion is bound directly to a bridging oxygen, which has a more negative partial charge than the other types of oxygen atom (Table S1 in the Supporting Information). The calculated structures of the uranyl ion adsorbed at the inner-sphere position on the kaolinite and orthoclase surfaces can be compared to EXAFS data published in the literature. To allow for an accurate comparison, additional 5 ns MD simulations were carried out in which the inner-sphere uranyl ions were not constrained. RDFs between the uranium atom and water oxygens, bridging oxygens, hydroxyl oxygens, and surface aluminum and silicon atoms were calculated to determine average interatomic distances and coordination numbers. There are no EXAFS studies of uranyl adsorption on orthoclase (potassium feldspar) in the literature, but Walter et al.6 reported structural parameters for uranyl adsorption on albite (sodium feldspar) from EXAFS measurements performed in air at pH 5.8 to 6.4. Focusing on the three pristine samples reported by Walter et al., as they more closely match the conditions employed in our simulations, fits to the EXAFS spectra of the three samples yielded 2.2 ± 0.3, 1.7 ± 0.2, and 2.5 ± 0.3 axial oxygens at 1.769 ± 0.007, 1.791 ± 0.004, and 1.785 ± 0.006 Å, respectively, compared with two axial oxygens at 1.86 Å in the MD simulations. In addition, Walter et al. found 6 ± 1, 6 ± 1, and 5 ± 1 equatorial oxygens at 2.34 ± 0.02, 2.36 ± 0.01, and 2.35 ± 0.02 Å, respectively. They noted the presence of some disorder in the equatorial oxygen shell on the basis of Debye−Waller factors greater than those obtained for aqueous uranyl species; however, a fit with split equatorial oxygen shells failed, and therefore, they concluded that the distortion within the equatorial oxygen shell was lower than their resolution of 0.19 Å. The MD simulations are in agreement with their findings. Indeed, they predict three water oxygens at 2.36 Å, one surface bridging oxygen also at

2.36 Å, and one surface hydroxyl oxygen at 2.27 Å. Finally, for one sample, the fit to the EXAFS spectrum was improved by adding 0.9 ± 0.2 Si or Al atom at 3.09 ± 0.01 Å, in accord with the simulations, which show one Al atom at 3.05 Å. Thompson et al.,59 Reich et al.,60 and Křepelová et al.8 reported EXAFS measurements of uranyl adsorbed on kaolinite at pH values from 5.0 to 8.5 in both air and Ar. Two axial oxygens were found at 1.77−1.81 Å in all of the samples. Again, the simulations predicted a slightly longer distance to the axial oxygens (1.85 Å). Fits to the EXAFS spectra also yielded five equatorial oxygens at distances between 2.34 ± 0.02 and 2.41 ± 0.02 Å. The MD simulations predicted three bonds to water oxygens at 2.34 Å and two bonds to surface hydroxyl oxygens at 2.44 Å, giving an average bond distance of 2.38 Å, which is within the experimentally derived range. Because the kaolinite aluminum octahedral sheet adopts the gibbsite structure, it is also informative to look at published EXAFS data for uranyl adsorbed on the gibbsite surface. Studies by Arnold et al.,61 Hattori et al.,62 and Gückel et al.63 carried out at pH values from 5.5 to 9.7 in both air and N2 found a similar structure for the equatorial oxygen shell with 5.0 ± 0.2 to 6.2 oxygens at distances from 2.373 ± 0.003 to 2.42 Å. The EXAFS data of Reich et al.60 and Křepelová et al.8 for kaolinite also showed the presence of one Al/Si atom at 3.06 ± 0.02 to 3.10 ± 0.02 Å and a second Al/Si atom at 3.26 ± 0.02 to 3.30 ± 0.02 Å, whereas the EXAFS data for all three gibbsite studies show only 0.8 ± 0.1 to 1.9 Al atoms at 3.33 to 3.405 ± 0.006 Å. The MD simulations predicted the presence of one Al atom at 3.42 Å. The similarities between the longer distances in the kaolinite EXAFS data, the gibbsite EXAFS data, and the MD prediction suggest that the two distances to Al atoms indicate the presence of two different adsorption sites on kaolinite, and indeed, this finding is consistent with the conclusions of Rösch and coworkers20,64,25 from their DFT studies of uranyl adsorption on basal and edge surfaces of kaolinite. The similarity between the shorter distances in the kaolinite EXAFS data and the MD prediction for orthoclase suggests that this shorter distance 3903



Article

could be assigned to uranyl adsorbed to a singly coordinated hydroxyl on kaolinite edge surfaces. Kaolinite edge sites are numerous enough to adsorb a significant fraction of the U(VI). Indeed, the maximum sorption density in the experiments of Křepelová et al.8 and Reich et al.60 (0.2 μmol/m2) is much less than the edge surface site density (3.4−6.1 μmol/m2) obtained from combining the kaolinite surface site density derived by Heidmann et al.65 (11.3−12.2 sites/nm2) and the percentage of kaolinite edge surfaces determined by Bickmore et al.66 (18− 30%). Effect of Carbonate Complexation on Uranyl Adsorption. To investigate how the formation of uranyl carbonate species affects the thermodynamics of adsorption on the two surfaces, constrained MD simulations similar to those described in the previous subsection were carried out for uranyl associated with one to three carbonate ions. The resulting PMFs are shown in Figures 5 and 6 together with the

Figure 6. (top) Potentials of mean force for the adsorption of UO2(CO3)n2−2n, where n = 0, 1, or 2, at the kaolinite (001) surface. (middle, bottom) Coordination numbers of uranium with carbonate oxygens (OC), water oxygens (OW), and hydroxyl oxygens (OH) for (middle) UO2CO3 and (bottom) UO2(CO3)22−.

UO2(CO3)22− adsorbs on the orthoclase (001) surface (Figure 4c), the steric hindrance due to the addition of a second carbonate ion to the uranyl carbonate complex leads to a decrease in the stability of the mono- and bidentate innersphere complexes; however, the adsorption free energy remains significantly favorable and comparable to the free energy of association of a third carbonate in the aqueous phase (Figure 2), which indicates that ternary uranyl complexes with up to two carbonate ions should form at the orthoclase (001) surface. In contrast, at the kaolinite surface, the adsorption of UO2(CO3)22− as an inner-sphere complex is unfavorable because of the steric hindrance caused by the two carbonate ions and the fact that uranium is only fivefold-coordinated in the equatorial plane compared with sixfold-coordinated for UO2(CO3)22− adsorption on the orthoclase (001) surface. UO2(CO3)34− adsorption on the kaolinite (001) surface was not considered, as the adsorption of UO2(CO3)22− was already found to be significantly unfavorable. Effect of pH and Implications for Surface Complexation Modeling. Following the approach of Chialvo et al.54 described above (with the 4π and r2 factors removed from eq 4 to account for the planar geometry), the equilibrium adsorption constants (log Kads) for adsorption on the orthoclase (001) surface,

Figure 5. (top) Potentials of mean force for the adsorption of UO2(CO3)n2−2n, where n = 0, 1, 2, or 3, on the orthoclase (001) surface. (middle, bottom) Coordination numbers of uranium with carbonate oxygens (OC), water oxygens (OW), hydroxyl oxygens (OH), and bridging oxygens (OB) for (middle) UO2CO3 and (bottom) UO2(CO3)22−.

coordination numbers of uranium with water oxygens, hydroxyl oxygens, bridging oxygens, and carbonate oxygens (OC). Snapshots of the inner-sphere complexes are shown in Figure 4. As was calculated in our previous work,30 the PMF of UO2CO3 adsorbing on the orthoclase (001) surface is similar to that of UO22+, which is due to the fact that UO2CO3 adsorbs to the surface with its UO22+ moiety exposed to the surface (Figure 4b). For the kaolinite (001) surface, the PMF of UO2CO3 differs slightly from that of UO22+, with the outersphere and inner-sphere adsorption sites being slightly destabilized by the presence of the carbonate ion. When

UO2 (CO3)n 2 − 2n (aq) → UO2 (CO3)n 2 − 2n (ads)

(8)

were calculated to be 10.61, 10.45, and 6.99 for n = 0, 1, and 2, respectively (the equilibrium adsorption constant for n = 3 was omitted because the corresponding PMF did not show a free energy minimum). As above for the uranyl−carbonate 3904



Article

limitation of these calculations is that they consider a single surface of a single mineral; however, as noted in the Introduction, feldspar minerals were shown to be an important mineral component of contaminated sediments at the DOE Hanford Site having a high affinity for U(VI).34,35 The calculations also did not consider the adsorption of hydrolyzed uranyl complexes, but those complexes have negligible concentrations above pH 8 in aqueous solutions in equilibrium with atmospheric CO2.68 Lastly, changes in surface protonation as a function of pH were not considered in the MD simulations, but this effect could be mitigated by the fact that the uranyl and uranyl carbonate complexes considered here all bind to the surface via their uranium atom. Despite these limitations, the approach presented in this work is a promising method for elucidating the nature and proportion of surface species and for providing a link between the spectroscopic determination of the microscopic structure of surface complexes and the derivation of surface complexation models from macroscopic experiments.

association in aqueous solution, the equilibrium adsorption constants were calculated in the infinite dilution limit. Kaolinite is not considered in this subsection because uranyl carbonate adsorption was found in the previous subsection to be mostly unfavorable. To investigate the effect of pH on the adsorption reactions through its influence on the aqueous speciation of uranyl and uranyl carbonate species, the three equilibrium adsorption constants were used in conjunction with experimental aqueous equilibrium constants for uranyl and uranyl carbonate species57 in the geochemical computer program PHREEQC.67 Several simplifying assumptions were made for these calculations: (1) the surface species described by eq 8 for the orthoclase (001) surface were assumed to be the only possible surface species; (2) the calculations did not include any modeling of the electrical double layer; and (3) because the MD simulations did not consider the adsorption of polynuclear uranium complexes, only aqueous speciation reactions involving mononuclear uranium complexes were considered. The total U(VI) concentration was set to 0.1 μM, 0.05 M NaNO3 was used to fix the ionic strength, and the aqueous carbonate concentration was assumed to be in equilibrium with the atmospheric CO2(g) partial pressure.68 Figure 7 shows the

■

ASSOCIATED CONTENT

S Supporting Information *

List of all short-range potential parameters and ionic charges, details of the changes made to the reference potential parameters, validation of uranyl−surface interactions, and search for surface adsorption sites. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This research was supported by the U.S. Department of Energy (DOE), Office of Biological and Environmental Research (OBER), 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 OBER and located at PNNL. PNNL is operated for the DOE by Battelle Memorial Institute under Contract DE-AC05-76RL01830.

Figure 7. (top) Aqueous U(VI) speciation as a function of pH based on experimentally derived equilibrium constants [U(VI)tot = 0.1 μM; pCO2 = 10−3.5 atm; I = 0.05 M NaNO3]. The sum of all uranyl hydroxyl complexes (n = 1−4) is shown as a single trace. (bottom) Surface U(VI) speciation on the orthoclase (001) surface as a function of pH based on MD-derived equilibrium constants (normalized to the total number of adsorbed uranium species).

■

REFERENCES

(1) Riley, R. G.; Zachara, J. M. Chemical Contaminants on DOE Lands and Selection of Contaminant Mixtures for Subsurface Science Research; U.S. Department of Energy, Office of Energy Research: Washington, DC, 1992. (2) Chisholm-Brause, C. J.; Conradson, S. D.; Buscher, C. T.; Eller, P. G.; Morris, D. E. Speciation of uranyl sorbed at multiple binding sites on montmorillonite. Geochim. Cosmochim. Acta 1994, 58, 3625− 3631. (3) Gabriel, U.; Charlet, L.; Schläpfer, C. W.; Vial, J. C.; Brachmann, A.; Geipel, G. Uranyl surface speciation on silica particles studied by time-resolved laser-induced fluorescence spectroscopy. J. Colloid Interface Sci. 2001, 239, 358−368. (4) Elzinga, E. J.; Tait, C. D.; Reeder, R. J.; Rector, K. D.; Donohoe, R. J.; Morris, D. E. Spectroscopic investigation of U(VI) sorption at the calcite−water interface. Geochim. Cosmochim. Acta 2004, 68, 2437−2448.

uranium aqueous speciation and the fraction of each surface species (normalized to the total number of adsorbed uranium species). In this model, UO2CO3 is the dominant adsorbed species from pH 6 to 9, even though it is the dominant aqueous species over a much narrower pH range. Above pH 9, UO2(CO3)22− becomes the dominant adsorbed species, while aqueous uranium is entirely present as UO2(CO3)34−. These results are consistent in general terms with surface complexation models derived from macroscopic experiments on contaminated sediments from the DOE Hanford Site, which consider uranyl surface complexes involving one or two carbonate ions under alkaline pH conditions.39,69 One 3905



Article

(5) Catalano, J. G.; Trainor, T. P.; Eng, P. J.; Waychunas, G. A.; Brown, G. E., Jr. CTR diffraction and grazing-incidence EXAFS study of U(VI) adsorption onto α-Al2O3 and α-Fe2O3 (11̅02) surfaces. Geochim. Cosmochim. Acta 2005, 69, 3555−3572. (6) Walter, M.; Arnold, T.; Geipel, G.; Scheinost, A.; Bernhard, G. An EXAFS and TRLFS investigation on uranium(VI) sorption to pristine and leached albite surfaces. J. Colloid Interface Sci. 2005, 282, 293−305. (7) Křepelová, A.; Brendler, V.; Sachs, S.; Baumann, N.; Bernhard, G. U(VI)−kaolinite surface complexation in absence and presence of humic acid studied by TRLFS. Environ. Sci. Technol. 2007, 41, 6142− 6147. (8) Křepelová, A.; Reich, T.; Sachs, S.; Drebert, J.; Bernhard, G. Structural characterization of U(VI) surface complexes on kaolinite in the presence of humic acid using EXAFS spectroscopy. J. Colloid Interface Sci. 2008, 319, 40−47. (9) 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. (10) Marques Fernandes, M.; Baeyens, B.; Dähn, R.; Scheinost, A. C.; Bradbury, M. H. U(VI) soprtion on montmorillonite in the absence and presence of carbonate: A macroscopic and microscopic study. Geochim. Cosmochim. Acta 2012, 93, 262−277. (11) Wang, Z.; Lee, S.-W.; Catalano, J. G.; Lezama-Pacheco, J. S.; Bargar, J. R.; Tebo, B. M.; Giammar, D. E. Adsorption of uranium(VI) to manganese oxides: X-ray absorption spectroscopy and surface complexation modeling. Environ. Sci. Technol. 2013, 47, 850−858. (12) 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. (13) 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. (14) 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. (15) 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. (16) 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. (17) Boily, J.-F.; Rosso, K. M. Crystallographic controls on uranyl binding at the quartz/water interface. Phys. Chem. Chem. Phys. 2011, 13, 7845−7851. (18) 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. (19) Drot, R.; Roques, J.; Simoni, E. Molecular approach of the uranyl/mineral interfacial phenomena. C.R. Chim. 2007, 10, 1078− 1091. (20) 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. (21) Sherman, D. M.; Peacock, C. L.; Hubbard, C. G. Surface complexation of U(VI) on goethite (α-FeOOH). Geochim. Cosmochim. Acta 2008, 72, 298−310. (22) Kremleva, A.; Krüger, S.; Rösch, N. Quantum chemical modeling of uranyl adsorption on mineral surfaces. Radiochim. Acta 2010, 98, 635−646. (23) Glezakou, V.-A.; deJong, W. A. Cluster-models for uranyl(VI) adsorption on α-alumina. J. Phys. Chem. A 2011, 115, 1257−1263. (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) 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. (26) 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, No. 154705. (27) Kerisit, S.; Liu, C.; Ilton, E. S. Molecular dynamics simulations of the orthoclase (001)− and (010)−water interfaces. Geochim. Cosmochim. Acta 2008, 72, 1481−1497. (28) Kerisit, S.; Liu, C. Molecular simulations of water and ion diffusion in nanosized mineral fractures. Environ. Sci. Technol. 2009, 43, 777−782. (29) Kerisit, S.; Liu, C. Molecular simulation of the diffusion of uranyl carbonate species in aqueous solution. Geochim. Cosmochim. Acta 2010, 74, 4937−4952. (30) Kerisit, S.; Liu, C. Diffusion and adsorption of uranyl carbonate species in nanosized mineral fractures. Environ. Sci. Technol. 2012, 46, 1632−1640. (31) Kerisit, S.; Liu, C. Structure, kinetics, and thermodynamics of the aqueous uranyl(VI) cation. J. Phys. Chem. A 2013, 117, 6421− 6432. (32) 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. (33) 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]arene6− and CMPO complexes. J. Mol. Struct. (THEOCHEM) 1996, 366, 55−63. (34) 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. (35) 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. (36) Bargar, J. R.; Reitmeyer, R.; Davis, J. A. Spectroscopic confirmation of uranium(VI)−carbonato adsorption complexes on hematite. Environ. Sci. Technol. 1999, 33, 2481−2484. (37) Bargar, J. R.; Reitmeyer, R.; Lenhart, J. J.; Davis, J. A. Characterization of U(VI)−carbonato ternary complexes on hematite: EXAFS and electrophoretic mobility measurements. Geochim. Cosmochim. Acta 2000, 64, 2737−2749. (38) Catalano, J. G.; Brown, G. E., Jr. Uranyl adsorption onto montmorillonite: Evaluation of binding sites and carbonate complexation. Geochim. Cosmochim. Acta 2005, 69, 2995−3005. (39) Liu, C.; Zachara, J. M.; Qafoku, N. P.; Wang, Z. Scale-dependent desorption of uranium from contaminated subsurface sediments. Water Resour. Res. 2008, 44, No. W08413. (40) 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. (41) Stoliker, D. L.; Liu, C.; Kent, D. B.; Zachara, J. M. Characterizing particle-scale equilibrium adsorption and kinetics of uranium(VI) desorption from U-contaminated sediments. Water Resour. Res. 2013, 49, 1163−1177. (42) Smith, W.; Forester, T. R. DL_POLY: A Package of Molecular Simulation Routines; Daresbury Laboratory: Warrington, U.K., 1996. (43) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269−6271. (44) Cygan, R. T.; Liang, J.-J.; Kalinichev, A. G. Molecular models of hydroxide, oxyhydroxide, and clay phases and the development of a general force field. J. Phys. Chem. B 2004, 108, 1255−1266. (45) Pavese, A.; Catti, M.; Parker, S. C.; Wall, A. Modelling of the thermal dependence of structural and elastic properties of calcite, CaCO3. Phys. Chem. Miner. 1996, 23, 89−93. (46) Hu, X. L.; Michaelides, A. The kaolinite (001) polar basal plane. Surf. Sci. 2010, 604, 111−117. 3906



Article

(47) Vasconcelos, I. F.; Bunker, B. A.; Cygan, R. T. Molecular dynamics modeling of ion adsorption to the basal surfaces of kaolinite. J. Phys. Chem. C 2007, 111, 6753−6762. (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, N. P. T. 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) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. Numerical integration of Cartesian equations of motion of a system with constraints: Molecular dynamics of n-alkanes. J. Comput. Phys. 1977, 23, 327−341. (52) 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. (53) Marcus, Y. Thermodynamics of solvation on ions. Part 5. Gibbs free energy of hydration at 298.15 K. J. Chem. Soc., Faraday Trans. 1991, 87, 2995−2999. (54) Chialvo, A. A.; Cummings, P. T.; Cochran, H. D.; Simonson, J. M.; Mesmer, R. E. Na+−Cl− ion pair association in supercritical water. J. Chem. Phys. 1995, 103, 9379−9387. (55) Reddy, M. R.; Berkowitz, M. The dielectric constant of SPC/E water. Chem. Phys. Lett. 1989, 155, 173−176. (56) Marcus, Y.; Hefter, G. Ion pairing. Chem. Rev. 2006, 106, 4585− 4621. (57) Guillaumont, R.; Fanghanel, T.; Fuger, J.; Grenthe, I.; Neck, V.; Palmer, D. A.; Rand, M. H. Update on the Chemical Thermodynamics of Uranium, Neptunium, Plutonium, Americium and Technetium; Elsevier: Amsterdam, 2003; Vol. 5. (58) Larentzos, J. P.; Criscenti, L. J. A molecular dynamics study of alkaline earth metal−chloride complexation in aqueous solution. J. Phys. Chem. B 2008, 112, 14243−14250. (59) Thompson, H. A.; Parks, G. A.; Brown, G. E., Jr. Structure and composition of uraniumVI sorption complexes at the kaolinite−water interface. In Adsorption of Metals by Geomedia; Jenne, E. A., Ed.; Academic Press: San Diego, CA, 1998; pp 349−370. (60) Reich, T.; Reich, T. Y.; Amayri, S.; Drebert, J.; Banik, N. L.; Buda, R. A.; Kratz, J. V.; Trautmann, N. Application of XAFS spectroscopy to actinide environmental science. AIP Conf. Proc. 2007, 882, 179−183. (61) Arnold, T.; Scheinost, A. C.; Baumann, N.; Brendler, V. Surface speciation of uranyl (VI) on gibbsite: A combined spectroscopic approach. In Institute of Radiochemistry Annual Report 2006; Bernhard, G., Ed.; Forschungszentrum Dresden Rossendorf: Dresden, Germany, 2007; p 53. (62) Hattori, T.; Saito, T.; Ishida, K.; Scheinost, A. C.; Tsuneda, T.; Nagasaki, S.; Tanaka, S. The structure of monomeric and dimeric uranyl adsorption complexes on gibbsite: A combined DFT and EXAFS study. Geochim. Cosmochim. Acta 2009, 73, 5975−5988. (63) Gückel, K.; Rossberg, A.; Brendler, V.; Foerstendorf, H. Binary and ternary surface complexes of U(VI) on the gibbsite/water interface studied by vibrational and EXAFS spectroscopy. Chem. Geol. 2012, 326−327, 27−35. (64) Martorell, B.; Kremleva, A.; Krüger, S.; Rösch, N. Density functional model study of uranyl adsorption on the solvated (001) surface of kaolinite. J. Phys. Chem. C 2010, 114, 13287−13294. (65) Heidmann, I.; Christl, I.; Leu, C.; Kretzschmar, R. Competitive sorption of protons and metal cations onto kaolinite: Experiments and modeling. J. Colloid Interface Sci. 2005, 282, 270−282. (66) Bickmore, B. R.; Nagy, K. L.; Sandlin, P. E.; Crater, T. S. Quantifying surface areas of clays by atomic force microscopy. Am. Mineral. 2002, 87, 780−783. (67) Parkhurst, D. L.; Appelo, C. A. J. User’s Guide to PHREEQC (Version 2)A Computer Program for Speciation, Batch-Reaction, OneDimensional Transport, and Inverse Geochemical Calculations; Water-

Resources Investigations Report 99-4259; U.S. Geological Survey: Denver, CO, 1999. (68) Zachara, J. M.; Brown, C. F.; Christensen, J. N.; Davis, J. A.; Dresel, P. E.; Liu, C. X.; Kelly, S. D.; McKinley, J. P.; Serne, R. J.; Um, W. A Site-Wide Perspective on Uranium Geochemistry at the Hanford Site; PNNL-17031; Pacific Northwest National Laboratory: Richland, WA, 2007. (69) Stoliker, D. L.; Kent, D. B.; Zachara, J. M. Quantifying differences in the impact of variable chemistry on equilibrium uranium(VI) adsorption properties of aquifer sediments. Environ. Sci. Technol. 2011, 45, 8733−8740.

3907


Molecular Dynamics Simulations of Uranyl and Uranyl Carbonate

Recommend Documents