Quantum mechanical calculation of aqueuous uranium complexes: carbonate, phosphate, organic and biomolecular species

Background Quantum mechanical calculations were performed on a variety of uranium species representing U(VI), U(V), U(IV), U-carbonates, U-phosphates, U-oxalates, U-catecholates, U-phosphodiesters, U-phosphorylated N-acetyl-glucosamine (NAG), and U-2-Keto-3-doxyoctanoate (KDO) with explicit solvation by H2O molecules. These models represent major U species in natural waters and complexes on bacterial surfaces. The model results are compared to observed EXAFS, IR, Raman and NMR spectra. Results Agreement between experiment and theory is acceptable in most cases, and the reasons for discrepancies are discussed. Calculated Gibbs free energies are used to constrain which configurations are most likely to be stable under circumneutral pH conditions. Reduction of U(VI) to U(IV) is examined for the U-carbonate and U-catechol complexes. Conclusion Results on the potential energy differences between U(V)- and U(IV)-carbonate complexes suggest that the cause of slower disproportionation in this system is electrostatic repulsion between UO2 [CO3]35- ions that must approach one another to form U(VI) and U(IV) rather than a change in thermodynamic stability. Calculations on U-catechol species are consistent with the observation that UO22+ can oxidize catechol and form quinone-like species. In addition, outer-sphere complexation is predicted to be the most stable for U-catechol interactions based on calculated energies and comparison to 13C NMR spectra. Outer-sphere complexes (i.e., ion pairs bridged by water molecules) are predicted to be comparable in Gibbs free energy to inner-sphere complexes for a model carboxylic acid. Complexation of uranyl to phosphorus-containing groups in extracellular polymeric substances is predicted to favor phosphonate groups, such as that found in phosphorylated NAG, rather than phosphodiesters, such as those in nucleic acids.


Background
The toxicity and radioactivity of U makes it a potentially hazardous element in the environment. In areas of high U concentrations, understanding U chemistry is imperative in order to predict its fate, transport, and risk. Uranium is capable of forming a wide variety of aqueous and surface complexes. Furthermore, redox reactions, mainly between U(VI) and U(IV), are common in subsurface environments (e.g., [1]).
Research has focused on the environmental chemistry of U with the goal of managing and remediating U-contaminated sites in the most effective manner ( [2][3][4] and references therein). Recent studies have probed the molecularlevel structures and processes that influence the overall behavior of U in the environment (e.g., [5]). Both analytical and theoretical studies have discussed complexation with numerous ligands [6][7][8][9] and the redox reactions between U(VI) and U(IV) (e.g., [10][11][12][13][14][15]). Computational chemistry is an important complement to experimental studies of U chemistry because this methodology can provide information that is not available via experiment, especially for transient species and those with short kinetic lifetimes. In order for the molecular modeling to be useful however, one must demonstrate that the computational methodology produces accurate results compared to known experimental data.
Before one can simulate structures, thermodynamics and kinetics with confidence, a computational methodology must be tested against observation. Environmental chemists are interested in U complexation and redox reactions, so this study focused on evaluating the ability of quantum mechanical calculations to reproduce experimental data on aqueous U complexes and redox chemistry. Specifically, models of aqueous U(VI), U(V), and U(IV) were generated and compared with experiment and previous calculations. Uranium complexes with inorganic (carbonate and phosphate), organic (oxalate and catechol), and biological (phosphodiester, phosphorylated glucosamine, and the 2-Keto-3-deoxyoctanoate) ligands were modeled and analyzed in light of previous observations. The model results are compared to interatomic distances from EXAFS, observed vibrational frequencies, and 13 C and 17 O NMR chemical shifts. Calculations on the observed oxidation of catechol by U(VI) are also presented.
The models were created including explicit H 2 O molecules around the complex to account for H-bonding by aqueous solutions. In this paper, a H-bond is considered to exist if the H---O distance is less than or equal to 2.0 Å and if the O-H---O angle is greater than 120°. These criteria are similar to those used by others (e.g., [26]) and are useful for identifying significant shifts in the calculated O-H stretching frequencies [27]. In general, initial models of solvation were created by positioning H 2  Previous work [28,29] has shown that including the H 2 O molecules in the primary solvation shell of UO 2 2+ is important for obtaining accurate structures, vibrational frequencies and energetics. This study (as in [29]) investigates the effects of adding a second solvation shell to the hydrated UO 2 2+ cation. The number of H 2 O molecules was chosen to be at least the minimal number necessary to form one H-bond to each of the possible H-bonding atoms in the U coordination sphere (e.g., 2 H 2 O molecules for each U-OH 2 group). In some cases (e.g., UO 2oxalate), an increasing number of H 2 O molecules were included in the model to assess the effects of explicit solvation on the predicted interatomic distances, vibrational frequencies, and NMR chemical shifts. Energy minimizations were generally carried out with the default criteria in Gaussian 03. When imaginary frequencies were calculated from an energy minimized structure, a re-optimization of the structure was performed with the "Opt = Tight" option until a structure with no imaginary frequencies was found. Although we have obtained potential energy minima, there is no guarantee that each configuration is the global minimum because the potential energy surface of these models will be complicated due to many possible Hbonding configurations. Any energy minimization scheme is unlikely to find the global potential energy minimum, so molecular dynamics simulations would be useful in the future to investigate configuration space at the temperature of interest and determine average configurations for these models.
Calculated results were compared to observed EXAFS, IR, Raman and NMR spectra. Model interatomic distances were directly compared to the values extracted via analysis of EXAFS data and vibrational spectra. Theoretical vibrational frequencies were compared to observed values for uranyl without a scaling factor applied because the appropriate value is not known for this computational methodology. For the vibrations of the ligands, a scaling factor of 0.96 was applied as determined by Wong [30] for B3LYP/ 6-31G(d) with the assumption that the p-functions added to the H atoms do not significantly affect the vibrational frequencies of the C-C and C-O bonds. This assumption is justified based on the work of Andersson and Uvdal [31] who showed scaling factors varying only between 0.96 and 0.97 for the 6-311G basis set with addition of polarization and diffuse functions.
Chemical shieldings were calculated with the gaugeincluding atomic orbital (GIAO) method [32]. The 13 C chemical shift values, δ 13 C, were calculated relative to the chemical shielding calculated for C in tetramethyl silane (TMS) using a model TMS structure energy minimized with the B3LYP/6-31G(d,p) method to be as similar to the uranium complex results as possible. The 17 O reference values were the average of O shielding values in a 19 H 2 O cluster which worked well in the past for predicting 17 O chemical shifts [33] and (CD 3 ) 2 CO which has been used experimentally as well [34]. The 13 C and 17 O chemical shieldings were calculated using the HF/6-31G(d,p) method [17] for both the TMS and uranium complexes. Tests were run using B3LYP/6-31G(d,p) and B3LYP/VTZ [35] for the aqueous UO 2 2+ and uranyl-oxalate complexes, but the agreement with observed δ 13 C values was not as good as for the HF/6-31G(d,p) method. However, the B3LYP/VTZ method resulted in significantly better agreement with the observed δ 17 O values. In general, the accuracy of calculated chemical shifts with different methods vary among elements (e.g., [36]), but a reason as to why this should be is a matter for further investigation and beyond the scope of this paper. There are also issues especially regarding the electric field gradient when including U in the NMR calculations due to the relativistic nature of the electrons in this atom [37,38].
Atomic charges were calculated with the natural bond orbital analysis (NBO) method of Glendening et al [39]. A modified version of natural population analysis was used for the explicitly solvated uranyl-catechol models that includes U 6d electrons in the valence space [40].
Starting structures were based on previous interpretations of experimental data where available [41][42][43][44][45][46][47]. The initial aqueous U models were UO 2 2+ (OH 2 ) 6 10 •20(H 2 O) (in both singlet and triplet states). In this notation, the first solvation sphere is denoted by "(OH 2 ) n " and the second solvation sphere is denoted by "•n(H 2 O)". Thus, two initial coordination numbers (CN) were tested for U(VI): the first with CN=8 and the second with CN=10. Comparisons of the calculated structures and potential energies generated from these initial configurations were used to test which coordination state should be favored. Starting the aqueous structure in a higher coordination state is a test of whether the lower coordination state observed (i.e., CN=7) can be reproduced or whether a local minimum would be found. Gibbs free energies of each model were calculated using the integrated equation formalism polarized continuum model (IEFPCM) of Cancès et al. [48]. The Gibbs free energies include the electrostatic energy interaction between the model and the continuum as well as the cavitation, dispersion and repulsion energy terms. The thermal and zero-point energy contributions to the Gibbs free energy as estimated in the gas-phase frequency calculations were added to the IEFPCM Gibbs free energies to correct for temperature.

Raman
A Nicolet Almega model dispersive Raman spectrometer from Thermo Scientific with a 785 nm laser was used to obtain spectra from 3446 to 111 cm -1 . The final spectra for each sample were the result of 128 scan averaging at high resolution (4.8-8.9 cm -1 ) with a 3.1 micron spot size. Liquid Raman spectra were collected with 100% laser intensity and were baseline corrected by subtracting a water spectrum from the background.

NMR
Solid-state 13 C CP/MAS NMR spectra of a freeze-dried 1:1 uranyl-catechol solution were obtained at 100.6 MHz with a 400 MHz Varian Inova spectrometer with an optimized tangential ramp of the transverse 13 C field over a contact time of 2 ms. The sample was spun at 8.0 kHz in a 5 mm Varian Jakobsen-style probe assembly. The 13 C chemical shifts are referenced to TMS via an external adamantane standard set to δ 13 C = 38.6 ppm. After re-dissolution of this sample in DMSO-d 6 solution-state 13 C and 1 H NMR spectra were obtained on a 600 MHz Varian Inova spectrometer fitted with a 10 mm broadband probe, using pulse delays of 1 and 10 s respectively.

U aqueous
The aqueous species of U(VI), U(V) and U(IV) were modeled to test the ability of the computational methodology to predict accurate interatomic distances as compared to those derived via EXAFS spectroscopy. The various bond lengths are listed in Table 1 [49,50]. In general, observed and calculated interatomic distances are in good agreement. The U(VI) aq model is discussed first followed by U(V) aq then U(IV) aq .

U(VI) aqueous
The experimental U=O and U-OH 2 distances of U(VI) aq are 1.76 ± 0.02 and 2.41 ± 0.02 Å, respectively, with a coordination number (CN) of 7 to 8 counting both the axial and equatorial O atoms [41,[51][52][53], and the calculated values are 1.79 and 2.44 Å, respectively, with a CN=7 (i.e., 2 axial doubly bonded O atoms and 5 H 2 O molecules). Interpretations of the experimental data have also been confirmed by high-level quantum calculations with continuum and explicit solvation by Gutowski and Dixon [29]. These model first solvation shell values are also similar to DFT calculations [54]  Gal et al. [56] observed IR peaks at 963, 253 and 160 cm -1 and Raman peaks at 874 and 198 cm -1 for this species associated with uranyl motion. Model frequencies involving U-OH 2 motion were predicted at 152 to 188 cm -1 (these modes were predicted to have minimal IR intensity; Fig. 1a), UO 2 2+ motion at 246 to 267 cm -1 , U=O stretches at 833 to 848 cm -1 (only the 848 cm -1 mode had significant IR activity), and U=O stretches at 949 cm -1 . A number of vibrations in the 150 to 300 cm -1 range have uranyl motion character, but the Raman activity of these modes is weak. The only strong Raman vibration is a triplet predicted at 834, 840 and 848 cm -1 (Fig. 1b). Peaks at 879 and 895 cm -1 have minor motion of the uranyl ion and are weak Raman scatterers. The present results are similar to previously calculated U=O stretching frequencies (≈820 and 890 cm -1 - [54]; 915 and 1010 cm -1 - [29]). This correspondence is considered a fair reproduction of the observed frequencies, and therefore the model is a reasonable representation of the solvation environment. Sources of error are the choice of basis set, electron correlation method, solvation and differences in the experimental and model reference state. Using the same structure and the B3LYP/VTZ method [35], the calculated O ax and O eq δ 17 O values are 1128 and -522 ppm, respectively; both of which are closer to observation. Thus, basis set and electron correlation have a more significant effect for the calculated δ 17 O, and this helps confirm the reasonable accuracy of the model structure.
Kubicki and Sykes [57] were able to reproduce δ 17 O values in aluminosilicate glasses to within a few ppm using the B3LYP/6-311+G(d,p) method, so the error here may be associated with relativistic effects due to bonding with the U atom. Solvation plays a significant role as evidenced by the large range of model δ 17 O values depending on details of the U-O and H-bonding. Lastly, the experimental reference states are the bulk liquid phases of water and (CD 3 ) 2 CO whereas the model reference states are a 19 H 2 O cluster and gas-phase (CD 3 ) 2 CO. The gas-phase to bulk liquid shift for 17 O chemical shielding can be as large as 36 ppm [58], so this discrepancy contributes to the overal errors mentioned above as well.
Two deprotonated U(VI) monomer species were also modeled. These models are useful for estimating the behavior of aqueous U(VI) as pH increases because UO 2 2+ undergoes hydrolysis above pH 5 [59,60]. These species are also the likely precursors of uranium oligomers such as (UO 2 ) 2 (OH) 2 2+ that form in solution and as such can be used in future studies of uranium polymerization reactions.
As expected, the U-OH bond lengths (Table 1) are intermediate between the U=O and U-OH 2 bond distances. These different bond lengths may not be distinguishable as distinct peaks in EXAFS spectra but would decrease the average U-O eq bonds by approximately 0.05 Å. Recently, Müller et al. [61] have observed an infrared peak at 922 cm -1 in water at pH 4 and assigned this peak to monomeric UO 2 UO 2 2+ in parallel rather than perpendicular as in the lower energy configuration significantly underestimates the frequency of the observed peak. However, the dominant IR peak in the UO 2 OH(H 2 O) 5 •14(H 2 O) is calculated to be at 909 cm -1closer to the observed peak. Since the assignment in Müller et al. [61] was based upon an interpolation of this mode between UO 2 2+ (aq) and U(OH) 4 2-(aq) , the actual assignment may be the UO 2 (OH) + species.
As is commonly observed, the calculated ΔG's for hydrolysis have significant discrepancies with experiment. Using a ΔG hyd for H + as +1098 kJ/mol [62] and the calculated G's in Table 1, the ΔG hydrolysis are equal to +130 and +282 kJ/ mol as compared with the experimental values of +29 and +59 kJ/mol [6]. Such a large discrepancy is large for this level of theory and may indicate that oligomers are being formed simultaneously with hydrolysis. Because we have not modeled the oligomerization reaction, if the oligomers have a lower free energy the observed ΔG hydrolysis would be less than calculated here consistent with the discrepancy above. Oligomerization is possible depending on the concentration of uranyl in solution [61].

U(V) aqueous
A literature search did not result in any published EXAFS studies of aqueous U(V). This is likely due to the instability of this species when uncomplexed. Renshaw et al. [44] have published interatomic distance for U(V)O 2 CO 3 complexes in solution and Chen et al. [63] have determined the crystal structure of K(UO)Si 2 O 6 in which U is pentavalent. Mizuouka and Ikeda [64] also found stable U(V) in dimethyl sulfoxide (DMSO) solution. Using these species as a guide, one expects U=O ax distances of 1.82 to 2.07 Å, respectively. The calculated changes in these parameters for aqueous U(V) ( Table 1) compared to U(VI) are consistent with this magnitude of elongation. One reason for the elongation of the U=O ax bonds in these calculations is that the H-bonding to the O atoms in these bonds has increased significantly over that found in the U(VI) aq model where only one H-bond of 1.92 Å formed. Hbonding to the O ax atoms increases the U=O ax bond length as electron density is shifted from the U=O ax bond to the H-bond. In contrast, the U(V) model had two Hbonds of 1.76 and 1.86 to one O ax atom and another of 1.58 Å to the second O ax atom.
The increase in H-bonding is a result of the increased electron density of the O ax atoms in the U(V) versus U(VI) species (the atomic charge as calculated with natural population analysis changes from -0.93 in U(VI) to -1.15 in U(V)). Although U(V) aq can still be represented as a UO 2 + ion, these changes in U=O ax bond length and H-bonding foreshadow the change to a aquo-U(IV) species where only H 2 O molecules are bound to the U atom. The U(V)-O eq bonds were found to be 2.47 Å in [44]. With U=O ax of 1.87 Å and U-O eq of 2.52 Å, the model results appear to be a reasonable prediction for this species and consistent with MP2 (customized basis set) results (U=O ax = 1.81 and U-O eq = 2.51) of Tsushima et al. [65] and the DFT (small-core ECP for U and 6-31G(d) for other atoms) results of Shamov and Schreckenbach [54] (U=O ax = 1.83 to 1.91 and U-O eq = 2.52 to 2.59).
Vibrational spectra are available for U(V) aq , in Best et al. [66] who observed a band attributed to aqueous U(V) at 910 cm -1 . No calculated vibrational modes associated with U(V)-O vibrations with significant IR intensity were found in these model calculations. Best et al. [66] noted that the presence of this band was highly sensitive to pH and that U(V) could be stabilized by dimerization with U(VI) species, but neither effect was modeled in these calculations. Hence, this assignment remains uncertain.
Our calculated IR frequencies are compared to U(V) in other phases to test whether the calculations are missing this vibration or whether the observed peak could have other origins. U=O vibrations are likely to be de-coupled from other bonds, so their frequencies should be similar (although not exactly the same) in the aqueous and crystalline phases. Chen et al. (2005) reported an IR frequency in a uranium silicate near 941 cm -1 that was assigned to a U-O stretching mode of UO 2 + . Mizuouka and Ikeda [64] observed an IR peak at 775 cm -1 that was assigned to asymmetric O=U=O stretching in U(V) in DMSO. A doublet with significant IR intensity exists in the model spectrum at 779 + 783 cm -1 . This mode does have a component O=U=O asymmetric stretch coupled to H 2 O motions which was also predicted in Shamov and Schreckenbach [54]. The closest IR peak with dominantly O=U=O asymmetric stretching motion in the theoretical spectrum is predicted to occur at 826 cm -1 . This vibration was predicted to occur at 786 cm -1 in a UO 2 + (OH 2 ) 5 •12(H 2 O) model by Shamov and Schreckenbach [54], and this model frequency was noted to change significantly with solvation. The fact that the observed modes in these various phases have close counterparts in the U(V) aq model calculations suggest that our assumption of de-coupling for these modes is justified; however, the model calculations fail to predict any peaks for U(V)aq near 910 or 949 cm -1 , so the extension of the U=O bond from the U(VI) to U(V) species may be exaggerated in the model results, thus lowering the bond force constant more than observed.
Although calculated frequencies generally are scaled for comparison with experiment to account for anharmonic and basis set effects, the correct scaling factor for uranyl vibrations and the method used here is not known, so no scaling was used. However, the close agreement of the observed and model results suggests that the model realistically represents the UO 2 + ion, especially considering that the aqueous UO 2 + model is being compared to a crystalline uranium silicate or U(V) complexes in DMSO.

U(IV) aqueous
The observed U-OH 2 bond lengths for aqueous U(IV)are 2.41 ± 0.02 Å [41,67], and there are 9 U-O(H 2 ) bonds averaging 2.47 Å from these computations (Table 1). (Note: This value was found for the triplet state. In the singlet state, U(IV) aq is predicted to have 8 U-O(H 2 ) bonds averaging 2.42 Å, but the energy is approximately +190 kJ/ mol higher in the singlet state.). A discrepancy on the CN exists, however, in that Moll et al. [67] concluded that U(IV) should be 10-fold coordinated.
In consideration of this interpretation, the model U(IV) aq species was generated with 10 H 2 O molecules initially bonded to the U(IV) cation. A minimum energy configuration for a CN of 10 could not be found. Attempts were made to find a stable U(IV) with 10 H 2 O molecules bonded using B3LYP/6-311+G(d,p) and with the LANL2DZ/6-31G(d,p) methods, but the O---O distances in these complexes were so short that a stable electron density could not be found. The resulting optimized structure for the B3LYP/6-311+G(d,p) U(IV) aq was similar to the B3LYP/6-31G(d,p) structure and experiment (Hennig, 2005), so we conclude that the predicted structure is not a strong function of basis set above 6-31G(d,p) (U-O(H 2 ) = 2.44 and 2.42 Å, respectively). The LANL2DZ/6-31G(d,p) structure resulted in average U-O(H 2 ) distances of 2.10 Å, so we conclude that the LANL2DZ is not as accurate as the ECP60MWB. In all three cases, the 9-fold coordinate U(IV) aq resulted in a stable energy minimum with no imaginary frequencies and the 10-fold coordinate U(IV) aq was unstable.
One would expect that these energy minimizations at 0 K would tend to overestimate the stability of higher coordinate species because thermal energy would favor the release of H 2 O from the coordination sphere which has a positive enthalpy [68]. Energy minimization for this model within a self-consistent reaction field could be performed to test whether long-range solvation could stabilize the higher coordination state.

UO 2 -Inorganic complexes
Uranyl complexes with a variety of inorganic anions (e.g., carbonate, phosphate, silicate). The uranyl-carbonate complexes are particularly important in the environment (e.g., [69,70]), so this study focused mainly on these Uinorganic complexes. In addition, phosphate or phosphoryl groups may be important in the chemistry of uranium with bacterial surfaces [71], so the uranyl-phosphate complex was also modeled here as a precursor to the organophosphate calculations.

UO 2 (CO 3 ) 2 2-
Uranyl dicarbonate is a common complex in typical groundwaters around pH 8 [70]. UO 2 (OH 2 ) 2 (CO 3 ) 2 2-•28(H 2 O) was modeled in two configurations with respect to the carbonate groups: trans and cis. Both optimized structures produced interatomic U-O (2.43 Å) and U---C (2.88 Å) distances in close agreement with experiment (Table 1). A similar computational methodology was used in Majumdar et al. [72] to calculate structures of the UO 2 (CO 3 ) complex, but these results predicted U-O bond lengths of ≈2.3 Å. Majumdar and Balasubramanian [73] used B3LYP, MP2 and coupled cluster doubles (CCD) methods with COSMO solvation [74] to model UO 2 (CO 3 ) 2 2-(with and without Li + and Na + ) and predicted U-O(C) bonds of 2.34 Å without Li + or Na + . This distance increased to 2.35-2.37 Å with the charge-balancing alkalis (No U---C distance was reported in this paper). The difference of approximately 0.1 Å between these models and observation is the result of neglecting the H-bonding due to explicit solvation of the carbonate group. Lack of H-bonding causes a stronger attraction of the [CO 3 ] 2group to the uranyl and will overestimate the complexation energy as well as underestimate the bond lengths. These differences highlight the need to include at least one solvation sphere when modeling aqueous anions to account for the strong H-bonding that occurs to these species [75].
The trans uranyl-dicarbonate species was predicted to be -9 kJ/mol lower in Gibbs free energy ( Table 1). The trans configuration should be more stable because it reduces the electrostatic repulsion between the negative carbonate groups, but the free energy difference is not large compared when one considers thermal and entropic effects. Because including diffuse functions in the basis set is relatively more important for anionic species, potential energy minimizations were also run with the B3LYP/6-31+G(d,p) method on the H, C and O atoms for these two configurations. In this case, the calculated ΔG = +23 kJ/ mol in favor of the cis configuration. In both cases, the calculated ΔG is within the range of error in the computational methodology employed, so we expect that both configurations could exist in solution.
A δ 13 C of approximately 165.8 ppm has been observed for UO 2 (CO 3 ) 3 4at pH 8 [78] and 168.2 ppm at pH 12 [79]. Differences in ionic strength between these two experiments could also contribute to this change in observed δ 13 C. Also, a value of 170.5 ppm was observed by Brucher et al. [80] at pH 9. The calculated value is 168.5 ppm, consistent with observation considering the accuracy of the calculations is a few ppm. The former is in somewhat closer agreement with the observed 215 ppm value, and the latter is off the scale in the published spectrum [81]. One can conclude that the computational methodology and solvation model are reasonably accurate for predicting δ 13 C values in this complex, but that solvation effects complicate prediction of the δ 17 O values.
Anionic complexes such as this can be more difficult to model due to the strong H-bonding that occurs to the O atoms of the anion [75]. The increased solubility of these U-carbonate complexes is due to this strong H-bonding to the carbonate groups. In our model, both the O eq and the O distal atoms form H-bonds to the surrounding H 2 O molecules. Furthermore, the increased electron density of this complex over the U(VI) aq species led to the formation of shorter H-bonds (≈1.86 Å) to the U-O ax than were predicted for the uncomplexed UO 2 2+ ion. These model values are slightly less accurate compared to experiment than those in de Jong et al. [6] based on small core Stuttgart relativistic effective core potentials and local density approx-  (Table 1; [42]). H-bonds to the O atoms of the UO 2 2+ group are relatively weak. Reduction of the U atom to U(V) causes a slight twisting of the carbonate ligands as proposed by Docrat et al. [43]. H-bonding to the uranyl O atoms becomes relatively stronger compared to the analogous U(VI) complex. (c) Addition of Ca 2+ ions to charge-balance this model aqueous species results in a configuration close to the observed crystal structure of the mineral liebigite [45].
imation DFT calculations with diffuse functions added to the basis sets, but the results in this study are still considered in reasonable agreement with observation.
The addition of Ca 2+ to aqueous solutions of uranyl-carbonates has a dramatic effect on increasing their stability and decreasing reactivity [49,[82][83][84]. To model this effect, two Ca 2+ ions were added to the UO 2 (CO 3 ) 3 (Fig. 2c). Two structures were energy minimized that varied in their Ca 2+ positions: one with both Ca 2+ ions approximately in the plane of the carbonate groups (Fig. 2c) and another with both Ca 2+ out of the equatorial plane. Both had similar structures with regard to the uranyl-carbonate interatomic distances (Table 1). One exception was the U---Ca distances where the complex with two Ca 2+ ions in-plane predicted distances of 4.05 Å, much closer to the experimental 3.95 Å compared to the 3.68 Å of the out-ofplane complex. The calculated stability of the in-plane structure is also -120 kJ/mol lower in Gibbs free energy than the latter (Table 1). This suggests a fairly strong regular structure for this ternary complex which will influence how it can react with mineral and bacterial surfaces [71,84].

5-
Docrat et al. [43] reported the structure of aqueous UO 2 [CO 3 ] 3 5-(i.e., the U(V) version of the uranyl triscarbonate species; Fig. 2b) based on EXAFS, and the interatomic distances are listed in Table 1. Starting with the model optimized structure of UO 2 [CO 3 ] 3 4-, an electron was added and energy minimization performed. The calculated interatomic distances are reasonably close to the experimental values (Table 1; [43,44,50]). The model results are in closer agreement with the values of Docrat et al. [44] and Ikeda et al. [50] who found shorter U=O ax bonds and longer U-O eq , U---C and U---O distal bonds than Renshaw et al. [44], but the discrepancies are still generally larger between observation and model than between the two experimental studies.
The observed δ 13 C NMR chemical shift for UO 2 [CO 3 ] 3 5is 106.7 ppm [79] compared to the 113.6 ppm average value calculated in our model which is in reasonable agreement for a paramagnetic complex considering that pseudo-contact shifts were not considered. These differences with experiment indicate that the U(V)-triscarbonato complex is not modeled as accurately as the corresponding U(VI) complex, but calculating NMR shifts can be problematic in U-containing systems where relativisitic effects may be significant [85], especially for an open shell system such as UO 2 [CO 3 ] 3  Similar changes in interatomic distances and H-bonding were predicted for the U(V)-carbonate complex as were calculated for the aqueous U(V) species compared to U(VI). This result is consistent with the suggestion of Tsushima et al. [65] and Ikeda et al. [50] that H-bonding to axial oxygens should increase as the U atom is reduced from U(VI) to U(V).
Aqueous U(V) rapidly disproportionates into U(IV) and U(VI), but the UO 2 [CO 3 ] 3 5complex can be stable for 2 h [43]. Recently, Ikeda et al. [50] have shown that the redox between U(VI)-and U(V)-triscarbonato species occurs "quasi-reversibly" and that U(V)-carbonate species can be stable for up to two weeks in a sealed cuvette. This change in electrochemical behavior with complexation has been suggested to be due to changes in the conformation of the uranyl-carbonate complexes with reduction from U(VI) to U(V) [43]. However, the similarity in measured interatomic distances via EXAFS [43,44] makes it difficult to explain the change in disproportionation kinetics. One possibility is that the carbonate group twists and the O eq atoms are displaced from the U-O eq plane defined in the UO 2 [CO 3 ] 3 4complex [43].
The model results indicate that such a twist does occur because the U-O-C-O torsions change from 163, 170, and 175° for the UO 2 (CO 3 ) 3 4complex to 160, 170, and 173°f or the UO 2 [CO 3 ] 3 5complex. This twist is not dramatic and the UO 2 (CO 3 ) 3 4complex does not have all six O eq atoms in the same plane to begin with. However, the static nature of these energy minimizations does not allow one to estimate the dynamic distortion that will occur in each type of complex at finite temperatures. Quantum mechanical molecular dynamics simulations of these complexes would be useful for examining this phenomenon.
The calculated Gibbs free energy differences between U(VI) aq /U(V) aq and UO 2 1). (Note that these are reduction half-reactions and do not include the energy of the electron, so these results do not imply that the reduced species are more stable overall. See [65] for methods and results on reduction energies between U(VI) and U(V).) The calculated reduction energies are not expected to be highly accurate at this level of theory and size of basis set, but the large relative change should be qualitatively correct. The direction of the change is unexpected because it would suggest that the U(V) aq species is relatively more stable with respect to re-oxidation and/or disproportionation compared to the carbonate complex. The lower calculated energy difference of the UO 2 (CO 3 ) 3 4-•28(H 2 O) complex compared to the uncomplexed U(VI) aq model relative to the U(V) version of each species is in direct contrast to the longer stability of the U(V)-carbonate complex. The observed metastability of the UO 2 (CO 3 ) 3 5-•28(H 2 O) complex is likely due to kinetic factors rather than thermodynamics. Surrounding the U(V) with carbonate anions may hinder the reaction between two U(V) atoms as electrostatic repulsion raises the activation energy barrier between anions. Modeling of the disproportionation reaction of complexed and uncomplexed U(V) species could help address this hypothesis.

Uranyl phosphate
The calculated structure of the uranyl-phosphate complex in Majumdar et al. [86] is close to the observed interatomic distances [46]. In this instance, the ligand has been protonated (i.e., HPO 4 2-) in the previous theoretical study, and the U=O ax and U-O(P) distances match experiment to within 0.01 Å. The U-O(H 2 ) distance is 0.08 Å larger than observed which may be the result of the model U existing in 7-fold coordination compared to 6-fold coordination extracted from the experimental data. One weak bond, such as the U-OH 2 bonds, may break at room temperature to form fewer, but shorter, bonds compared to the 0 K energy minimization structure obtained computationally. Another discrepancy is that the U---P distance is 0.15 Å longer in the model of Majumdar et al. [86] compared to the value obtained by Locock et al. [46].
The protonation state of the phosphate ligand will vary with pH (pK a 's of 2.  [46]. Two attempts at optimizing a monodentate UO 2 -HPO 4 complex were made using initial configurations that should favor the monodentate complex (i.e., nearly linear U-O-P angles and strong Hbonding to the phosphate O atoms not bonded to the U). However, both resulted in conversion to a bidentate configuration (Fig. 3a). The approximate potential energy dif-  ference between the mono-and bidentate configurations was greater than -1000 kJ/mol. Even considering the partial optimization of the monodentate structure and all the approximations and errors inherent in these model calculations, this difference is extremely large and suggests that the monodentate configuration is not stable for UO 2 -HPO 4 . The possibility exists for a monodenate UO 2 -H 2 PO 4 complex as discussed below.
The possibility of a 1:2 uranyl phosphate complex also exists. Two HPO 4 2ions could bond to uranyl to form a monodentate complex, so this structure was used as an initial configuration. As in the 1:1 complex, one of the HPO 4 groups re-oriented itself to form a bidentate configuration during the energy minimization (Fig. 3b). The second HPO 4 remained in a monodentate configuration throughout the energy minimization. The monodentate configuration had a shorter U-O(P) distance (2.22 Å) and a longer U---P distance (3.78 Å) than experimentally     4 group. Based on the results discussed above and the fact that the U-O(P) and U---P distances do not agree with observation, it is likely that this monodentate configuration is metastable.
The relative stability of this monodentate HPO 4 2ligand was checked by running an energy minimization of the same model with both HPO 4 2groups starting in a bidentate configuration. The model with both HPO 4 2groups bidentate was -101 kJ/mol lower in Gibbs free energy than the mixed monodentate/bidentate model (Table 1). Again, these results indicate that the monodentate configuration is less stable for HPO 4 2ligands than the bidentate configuration.
A stable monodentate configuration for the 1:1 UO 2 -H 2 PO 4 model complex was found. The U---P interatomic distance in this configuration was 3.64 Å in poor agreement with Locock et al. [46] (Table 1) but in excellent agreement with a value for U---P for uranyl adsorbed to phosphoryl groups on bacterial surfaces at low pH (3.64 ± 0.01 Å; [47]). The Kelly et al. [47] value was interpreted as a monodentate phosphoryl bonded to uranyl. The Gibbs free energy of the monodentate configuration was +9 kJ/mol higher than the same model system in a bidentate configuration (UO 2 H 2 PO 4 (bi); Table 1). The above results all indicate that the mechanism of ligand binding to uranyl ions is sensitive to pH and the number of ligands involved. These factors need to be considered as spectroscopic data collected under specific conditions is applied to surface complexation modeling and to field studies.

IR and Raman
Quantum mechanical calculations can be useful in helping to identify speciation and protonation changes in conjunction with IR and Raman spectroscopies as well. By predicting vibrational frequencies and modes, assignment of observed spectra to specific complexes can be done with greater confidence.
Allen et al. [78] observed a ν 1 (O=U=O) peak at 813 cm -1 in Raman spectra of aqueous UO 2 (CO 3 ) 3 4-; this vibrational mode is predicted to be at 871 cm -1 in these model calculations (average of two modes at 865 and 877 cm -1 ). Model frequencies may be approximately 4% too high due to neglecting anharmonicity [30]), so this agreement is reasonable. However, recent work has suggested that the overestimate of the vibrational frequencies for the uranyl ion is more likely due to basis set and electron correlation effects rather than anharmonicity [89]. Bargar et al. [42] measured ATR FTIR spectra of aqueous UO 2 (CO 3 ) 3 4and found peaks near 1360 and 1510 cm -1 associated with the carbonate vibrational modes (symmetric and asymmetric ν 3 modes, respectively). The model UO 2 (CO 3 ) 3 4-•28(H 2 O) used in this study resulted in frequencies of 1302, 1372, 1427 and 1545 cm -1 for these modes with the highest IR intensity in this range (Fig. 1c).
At first, the disagreement would seem to be significant, but when one considers that the observed bands have fullwidths at half-maximum of approximately 50 cm -1 , then the possibility can be considered that these observed bands are mixtures of the calculated bands. The IR intensity weighted-averages of the two sets of frequencies result in calculated values 1338 and 1503 cm -1 that are as close as can be expected to the observed frequencies for modeling aqueous phase complexes at this level of theory [90]. Without including H 2 O molecules of solvation around the complex, de Jong et al. [6] calculated sets of peaks at 1356 + 1383 and 1508 + 1533 cm -1 . Hence, either solvation is not critical for calculating the UO 2 (CO 3 ) 3 4frequencies or there are compensating effects in the DFT approach used by de Jong et al. (2005). Previous work has shown that strong H-bonding can dramatically affect carbonate vibrational frequencies [75], so the reason for the excellent agreement in this case with the gas-phase calculations should be investigated further by performing a systematic study of explicit and implicit solvation on this model.
Although neither IR nor Raman spectra are available for uranyl-phosphate aqueous solutions, vibrational frequencies calculated (Fig. 1d) [91]. Calculated peaks around 900 cm -1 were not observed in the Raman spectrum of this mineral nor did the calculations predict observed peaks between 952 and 976 cm -1 [91]. Considering the fact that a model aqueous complex is being compared to a mineral, the level of agreement is remarkable and suggests that the bonding environment of the uranyl-phosphate complex is being modeled accurately.

UO 2 -Organic
Complexation of UO 2 2+ with oxalate was modeled because oxalate is a simple organic ligand that forms strong complexes with cations. In addition, oxalate is a common organic acid in nature. Alliot et al. [92] studied the role that uranyl-oxalate complexation plays in the adsorption of uranium to alumina and found that making the charge of the complex neutral or negative would cause uranium to desorb due to higher solubility of the aqueous complex. Previous quantum mechanical calculations have been performed on uranyl-oxalate complexes by Vallet et al. [93] using Hartree-Fock gas-phase energy minimizations and single point energies based on MP2 calculations in the conductor-like polarized continuum model (CPCM) for solvent effects. Our structural results including explication solvation by H 2 O molecules will be compared to these earlier calculations.   Fig. 5c. The strongest IR adsorbing modes group into three regions in the 500 to 2000 cm -1 part of the spectrum shown (higher frequencies are neglected because they are dominated by motions of H 2 O): 700 to 900, 1200 to 1400, and 1600 to 1800 cm -1 . The peak at 1797 cm -1 is the C=O stretch. The group of peaks just below this sharp peak in Fig. 5a   mechanical calculations in Tsushima et al. [94] predicted a 1.69 Å U=O ax bond, but this underestimate is likely due to the fact that solvating H 2 O molecules were not included in this early study. The 1144 and 1270 cm -1 peaks compare favorably with the calculated peaks at 1184 and 1275 cm -1 (neglecting anharmonicity in the calculations tends to increase the predicted frequencies). The peak at 1450 cm -1 in the experimental spectrum has no equivalent in the model uranyl-oxalate complex and is due to uncomplexed oxalate. A peak at 1483 cm -1 has been observed for aqueous Na-oxalate by attenuated total reflectance FTIR [95].
IR spectra have been collected by Giesting et al. [96] on a crystal of (UO 2 ) 2 C 2 O 4 (OH) 2 (H 2 O) 2 but the uranyloxalate configuration is not similar to the aqueous model because the oxalate ligand is bonded to two uranyl groups. Thus, the vibrational spectra will be completely different from the case where two C-O(U) and two C=O bonds exist in the oxalate. A literature search revealed no other IR or Raman spectra of uranyl-oxalate complexes, so these calculations stand as predictions in the event this type of study is conducted in the future.
Values for δ 13 C of the uranyl-oxalate complexes are listed in Table 2 and the NMR spectra shown in Figure 6. The HF/6-31G(d,p) method was used to calculate chemical shifts relative to TMS and the experimental values were measured relative to trimethylsylyl propionate, but the difference in chemical shielding between TMS and TMP is only -0.08 ppm [97]. Hartree-Fock calculations have proven to be more accurate for reproducing δ 13 C values in other calculations of organic complexes [98], so this method was used instead of a DFT method such as B3LYP. (B3LYP was tested but results did not match experiment as closely for test cases.) The observed δ 13 (Table 2). Thus, the liquid-state NMR experiment appears to be observing the singly deprotonated oxalic acid species, a time-averaged value of the singly and doubly deprotonated oxalic acid, or the method is accurate to only 5 ppm. The experimental uranyl-oxalate δ 13 C (169 ppm) and calculated (167 ppm) values are in agreement within computational error for the [UO 2 (OH 2 ) 4 ] 2+ -Ox•30(H 2 O) model. Since the measurement was at pH 3, then the spectrum corresponds to H 2 Ox (pK a = 1.27) and HOx -(pK a = 4.28) in rapid exchange equilibrium (time-average). Having an observed chemical shift that is more shielded than calculated for HOxis consistent with this.
Solvation appears to be important for obtaining more accurate chemical shifts because the models with no or few H 2 O molecules are less accurate using the same basis set and method. Note that the model value is the average of two distinct values for the oxalate C atoms: 171 and 164 ppm. Reproduction of the observed vibrational frequencies and 13 C chemical shifts is evidence that the oxalate binds to UO 2 2+ in the bidentate fashion as previously described by Vallet et al. [93].

UO 2 -Catecholate
Uranyl has been shown to interact with catechol in aqueous solutions. Catechol is an organic ligand generated by plants [99] and an analog for common functional groups within humic substances [100]. Consequently, modeling uranyl-catechol complexation is important for environments with significant amounts of natural organic matter present. Additionally, Eng [101] has shown that U(VI) can oxidize catechol and form organic polymers. Although the U(VI) was subsequently re-oxidized in these aerobic experiments, one can presume that the U(VI) is reduced to U(V) or U(IV) in the process of oxidizing the catechol. Calculations on the structures and 13 C NMR chemical shifts of UO 2 2+ -catechol complexes are presented here and results on the uranyl-catechol redox reaction are discussed.

Solvated with 4 and 28 H 2 O
The catechol ligand does not deprotonate in aqueous solution with UO 2 2+ as evidenced by NMR spectroscopy (Table 2; Fig. 7). Signal from both phenol protons on the NMR spectra of oxalic acid and uranium-oxalate aqueous solutions Figure 6 NMR spectra of oxalic acid and uranium-oxalate aqueous solutions.

Oxalic acid
Uranyl-oxalate complex catechol are present as indicted by the approximately 1:2 ratio of intensity relative to the C-H (Fig. 7) which is suggestive of an outer-sphere uranyl-catechol complex. This interpretation is supported by observation of 13 C chemical shifts for freeze-dried uranyl-catechol solutions that are nearly identical (within 0.7 ppm) to solid catechol. Outer-sphere complexation leads to the formation of significantly weaker bonds between the catechol and uranyl compared to the oxalate ligand, and the U-O distances are approximately equal between the U-O(H 2 ) and U-O(C) bonds in the UO 2 2+ -catechol complexes ( Table 2). Notably, addition of 28 H 2 O molecules around the bare [UO 2 (OH 2 ) 4 ] 2+ -catechol complexes decreases both of the equatorial bond types by approximately 0.15 Å. This is opposite to the effect predicted for carbonate and oxalate ligands which lose electron density to surrounding H 2 O molecules, and we have no explanation for this bond shortening at this time. The weakest of the U-O(H 2 ) bonds is also broken with the addition of H 2 O's of solvation to lower the CN to 7 in the solvated case (Table 2).
Catechol has three pairs of equivalent carbon atoms with distinct δ 13 C values (Table 2). For uncomplexed catechol in DMSO solvent, these are found at 116.1 (C atoms next to the COH groups), 119.8 (for the C atoms opposite the COH groups), and 145.6 ppm (for the phenolic C atoms), respectively, which are very similar to values of 115.7, 121.7, and 142.7 ppm for the neat solid [102]. HF/6-31G(d,p) NMR calculations on the B3LYP/6-31G(d,p) optimized structures result in predicted values of 115.5, 121.1, and 144.3 ppm for each of these peaks (Table 2), so this methodology is expected to reproduce δ 13 C values to approximately 2 ppm. Calculated values for the innersphere [UO 2 (OH 2 ) 4 ] 2+ -H 2 Cat•28(H 2 O) complex are in disagreement with the observed 13 C chemical shifts by 4 to 6 ppm (Table 2). Thus, this model with direct bonding of the phenol O atoms to the uranyl cation is unlikely to exist in solution. The other inner-sphere models, [UO(OH)(OH 2 ) 4 ] 2+ -HCat•28(H 2 O) and [U(OH) 2 (OH 2 ) 4 ] 2+ -Cat•28(H 2 O), result in even larger discrepancies from experiment predicting six separate chemical shifts ranging from 113 to 342 ppm, so another type of uranyl-catechol complex must give rise to the observed 13 C NMR spectra that show δ 13 C values shifted by approximately 0.5 ppm from the catechol-only solution (Fig. 7).
In a search for the structure of the uranyl-catechol complex, the structure of a model uranyl-catechol outersphere complex was also calculated. The energy minimized structure is shown in Fig. 9d. Agreement between model and observed δ 13 C values is much improved over the inner-sphere and quinone models (Fig. 9e). Although significant errors remain for three of the δ 13 C values, the other three C atoms have calculated chemical shifts that are reasonably close to observation ( Table 2). Furthermore, the stability of this species is calculated to be approximately -80 kJ/mol more favorable in potential energy (Table 2) than the inner-sphere [UO 2 (OH 2 ) 4 ] 2+ -H 2 Cat•28(H 2 O) complex shown in Fig. 9a. Apparently, the uranyl bonds to the phenolic groups are weaker than the UO2 2+ -OH 2 bonds that are broken to form them.
One can conclude from the NMR and energy results that the [UO 2 (OH 2 ) 4 ] 2+ -H 2 Cat (OS) model is the closest rep-NMR spectra of uranyl-catechol complexes freeze-dried solid 13 Figure 9 (
resentation of those studied here of the actual uranyl-catechol interactions in low pH experimental solutions.
Modeling the structure of an outer-sphere complex is more challenging than inner-sphere complexes that are dominated by covalent bonding because the outer-sphere complex is likely to be much more flexible and experience a wider range of configurations. Thus, conformation searching for the uranyl-catechol outer-sphere complex would be useful.

U(VI) reduction and catechol oxidation
As mentioned above, oxidation of catechol has been observed in conjunction with uranyl-catechol complex formation [101]. To begin modeling this redox reaction, the [UO 2 (OH 2 ) 4 ] 2+ -H 2 Cat•28(H 2 O) complex that had been energy minimized in the singlet state (Fig. 9a) was also energy minimized in the triplet (  Table 2). The lower energy of the quintet state is consistent with the prediction above that the U(IV) aq species is in the triplet rather than singlet state.
The computed charge on the U atom changes from 3.42 to 3.29 to 3.19 using natural bond orbital analysis [103][104][105][106][107] as the multiplicity of the complex changed from singlet to triplet to quintuplet. The charge on the UO 2 2+ group was predicted to decrease from +1.10 to +0.65 to +0.31 (Table  3) reflecting a reduction of the uranyl cation (Note that the NPA charges on the uranium atom only change by 0.1 electron.). The atomic spin densities on U also increased from 0 to 1.27 to 2.10 in the singlet, triplet and quintet state, respectively. The remaining spin density in each case was located on the C and O atoms of the original catechol. This result is consistent with the observation of polymerization in these solutions [101] via a radical mechanism.
These calculations predict a thermodynamic driving force for the reduction reaction; but, because the U(VI)-catechol complex is stable enough to give rise to Raman and NMR spectra in low pH solutions, the activation energy barrier must be significant at room temperature as well consistent with the observation that oxidation/polymerization requires days to occur at room temperature [101].
Modeling the reaction pathway is beyond the scope of this paper, but the reaction mechanism is likely to involve a proton-coupled electron transfer because the stable species is [U(OH) 2 (OH 2 ) 4 ] 2+ -catecholate where the doublybonded O atoms in the uranyl cation have become protonated via transfer of H atoms from the phenol groups in the catechol (Fig. 9b). This pathway makes sense in light of the experimental evidence that a quinone-like species is produced when catechol is oxidized in this reaction [101].
Using the outer-sphere uranyl-catechol model (Fig. 9d) as the reactant and a U (IV) (OH) 2 (OH 2 ) 5 -quinone•28(H 2 O) (Fig. 9e) as the product, the overall reaction ΔG is estimated to be approximately -65 kJ/mol (Table 2); however, our model results suggest that the triplet and quintet states in Table 2 are at a lower Gibbs free energy than the outersphere U (IV) (OH) 2 (OH 2 ) 5 -quinone•28(H 2 O) complex.

UO 2 -Biological
Adsorption of uranium to bacterial surfaces is an important aqueous geochemical process (see [108]). This adsorption is generally thought to occur to functional groups located in the extracellular polymeric substances (EPS) located around many bacteria. Water may be present throughout EPS [109], so this section of the paper A modified l607 routine was used for the explicitly solvated models to include the U 6d electrons in the valence space for the natural population analysis [40]. The multiplicities (singlet, triplet and quintet) correspond to U formal oxidation states of +6, +5 and +4. The "U Charge" and "UO 2 2+ Charge" headings indicate the actual charge calculated using NBO. Note that the protonation state of the catechol (H 2 Catechol, HCatechol and Catechol) indicates a H + transfer to the uranyl group simultaneous with the electron transfer.
examines aqueous uranium complexation with EPS functional groups.
Previous work has concluded that the phosphoryl and carboxlate groups are the most important for binding uranyl [47,110,111]. Consequently, we have modeled two types of phosphoryl groups, one that mimics phosphoryl backbones in nucleic acids (OrgPO 4 ; see also [112]) and another that represents phosphoryl group in glucosamine (GlcNPO 4 ) which is thought to be an important component in EPS for metal adsorption [113]. The carboxylate group was modeled in the 2-Keto-3-deoxyoctanoate (KDO) component of EPS which has also been suggested as playing a role in uranyl adsorption to bacteria [114,115].

UO 2 -Phosphodiester
Previous workers have suggested that monodentate complexes of uranyl and bacterial EPS phosphoryl groups are the predominant mechanism of bonding (e.g., [47,71,111]). These calculations sought to compare the relative stabilities of monodentate and bidentate configurations. Model results in this study strongly indicate that the bidentate complex should be favored to the nucleic acid phosphoryl groups. Three initial configurations were subjected to energy minimizations using the same methodology described above for the inorganic and organic Ucomplexes. Two initial monodentate structures were used (e.g., Fig. 10a); in both cases, the energy minimization led to a bidentate structure (e.g., Fig. 10b) that was at least -520 kJ/mol lower in Gibbs free energy. Although relaxa-  [47]showing monodentate phosphoryl bonding to U on bacterial surfaces at low pH.
tion of atoms in the model other than those involved in forming the second bond between the phosphoryl group and the uranyl account for some of this lowering in energy, it is clear that the bidentate configuration is lower in Gibbs free energy by much more that the expected errors in the calculations. (Note: The energy minimizations beginning with three different starting configurations resulted in final structures with energies within ± 25 kJ/mol of the average. This order of magnitude is used as an estimate of the uncertainty in the Gibbs free energies resulting from these energy minimization procedures.) Although the calculated energetics clearly favor the bidentate configuration, the calculated U---P distance in these complexes is approximately 3.3 Å (Table 4) whereas the observed value is 3.64 to 3.68 Å [111,116]. This discrepancy could be due to the fact that the Kelly et al. [47,111] experiments were conducted at a pH of 1.7 to 4.8, so the phosphoryl groups were probably protonated which prevented bidentate bonding. The model structures of UO 2 -OrgPO 4 monodentate (Table 4) are also consistent with a shorter U-O eq bond of 2.30 to 2.37 Å measured by Koban et al. [117] using EXAFS in solutions of pH 3.5 to 5.5. Francis et al. [116] performed EXAFS on samples reacted at pH 5 and obtained similar results as in Kelly et al. [111], but the expected pK a for phosphoryl groups on bacterial surfaces is approximately 7.2 [118], so experimental conditions may not have reached a state where the phosphoryl was predominantly deprotonated.
In order to investigate the possibility that a protonated phosphoryl (i.e., mimicking low pH conditions) favors the monodentate uranyl complex, a uranyl-organophosphate complex (UO 2 -OrgHPO 4 ; Fig. 10c) was investigated. A stable monodentate configuration was obtained when the phosphoryl group was protonated which resulted in a U---P distance of 3.8 Å (  [116,119]. These complexes were not modeled in this study and should be included in future work.

UO 2 -GlcNPO 4
Another possible explanation of the discrepancies in interatomic distances between EXAFS and these calculations may be that uranyl binds to a phosphoryl group unlike the phosphodiester model. To examine this possibility, models of the glucosamine phosphoryl groups, GlcNPO 4 (Fig. 11) were also studied which represents a phosphonate compound. Phosphonates have been shown to complex U(VI) effectively [113,120]. The UO 2 -GlcNPO 4 model calculations also predicted that the bidentate configuration had a lower Gibbs free energy than the monodentate configuration but only by approximately -30 kJ/mol as compared with a difference of -520 kJ/mol for the phosphodiester (Table 4). In addition, a monodentate configuration was found in this case that represented a local potential energy minimum, so the calculated difference is less of an estimate in this case compared to the UO 2 -OrgPO 4 models (Table 4).
Furthermore, the U---P distance calculated for monodentate UO 2 -GlcNPO 4 is in excellent agreement with observation [47]. The U-O(P) bonds for the mondodentate and bidentate configurations are 2.23 and 2.43 Å, respectively. Thus, the U-O bond shortening observed by Koban et al. [117] is predicted, and the average of the two configurations results in a U-O(P) bond length close to that observed (i.e., 2.30 to 2.37 Å). Consequently, the conclusion that uranyl will bond to phosphoryl groups in a mondentate manner at low pH and convert to a bidentate complex at circumneutral pH is supported by these results because the model structures are similar to observation at lower pH, but theoretical potential energies of the bidentate configuration are significantly lower. These results also suggest that uranyl will favor glucosamine phosphoryl groups over phosphodiester groups. Micro-environmental conditions (such as lower dielectric constant or lower pH) within the EPS region could be responsible for stabilizing the monodentate over the bidentate configuration because these factors were not considered in the model calculations of this study.

UO 2 -KDO
As in the case of uranyl-phosphoryl binding, uranyl-carboxylate complexes may either form mono-or bidentate structures (Fig. 12). Previously, the bidentate configuration of uranyl with acetate has been modeled [6] and compared to experimental interatomic distances and frequencies. The relative potential energies of the mono-and bidentate UO 2 -KDO complexes and the model U---C distances are consistent with the interpretation that the bidentate configuration is favored (

Figure 11
(a) The monodentate UO 2 -GlcNPO 4 complex has a calcu-lated potential energy +70 kJ/mol higher than the corre-sponding (b) bidentate complex, but the former has interatomic distances in closer agreement with EXAFS data on uranyl adsorbed onto bacterial surfaces [47] Figure 11 (a) The monodentate UO 2 -GlcNPO 4 complex has a calculated potential energy +70 kJ/mol higher than the corresponding (b) bidentate complex, but the former has interatomic distances in closer agreement with EXAFS data on uranyl adsorbed onto bacterial surfaces [47]. These two facts suggest that protonation of the phosphoryl group at low pH stabilizes the monodentate configuration and that the bidentate configuration should be more stable at circumneutral pH.  [6] study and the UO 2 -KDO values calculated in this study). The main differences are for the U-O bonds to the carboxylate group and the U---C distance. However, these differences are mainly due to the fact that the uranyl coordination sphere was not completed in the 1:1 complex modeled in de Jong et al. [6]. Their reported values for the U-O and U---C distances increase toward the values reported here in the 1:3 model uranyl-acetate complex and in the EXAFS data on solids [121].
Although the monodentate complex is predicted to be of slightly higher Gibbs free energy, it may be difficult to rule out the existence of this complex based on EXAFS results. The model interatomic distances in the monodentate complex are similar to those predicted for the bidentate complex except for the U---C distance (i.e., U=O = 1.76, C-O = 1.26, U-O = 2.34, and U---C = 3.45 Å). The latter is longer by approximately 0.5 Å, but this type of interatomic distance can be difficult to measure with EXAFS in materials such as bacterial EPS.
A third possible configuration is the outer-sphere pair formed by the uranyl cation and KDO (Fig. 12c). The outer-sphere model results in a higher potential energy than both the mono-and bidentate configurations by +12 kJ/mol. This value is well within computational accuracy and suggests that an equilbrium exists between inner-and outer-sphere species. Previously, it has been common to assume that inner-sphere species are in a lower energy state than outer-sphere, but recent work by Catalano et al. [122] on arsenate adsorption to hematite has demonstrated that this is not necessarily the case. Consideration of these outer-sphere species must be included in bacterial complexation modeling.
Vibrational spectroscopy could be more diagnostic in this case because the splitting of the carboxylate vibrational modes changes depending on whether the bonding is mono-or bidentate [123]. For the bidentate complex (Fig.  13a), the calculated C-O s mode is approximately 1500 cm -1 (there are four modes involving symmetric C-O stretching coupled with CH 3 motion in the KDO; Fig. 13), and the C-O as mode is at 1612 cm -1 . The C-O s and C-O as frequencies observed via ATR FTIR for bidentate uranyl acetate in solution are 1467 and 1527 cm -1 , respectively. If the theoretical frequencies are scaled by 0.96 (the scale factor for organic molecules using B3LYP/6-31G(d) method; [124]), this results in errors of 27 and 20 cm -1 , respectively, for these two modes. Considering the difficulties in model aqueous uranyl-organic complexes, and the fact that the model results are for UO 2 -KDO rather than the experimental uranyl acetate, this size of error is satisfactory.
The symmetric and asymmetric U=O vibrations are predicted to be approximately 900 (two frequencies at 895 and 906 cm -1 ) and 983 cm -1 (Fig. 13a). Both modes are relatively strong in this model, but the U=O as vibration UO 2 -KDO (a) monodentate complex is most consistent with observed interatomic distances [110], but the (b) bidentate configuration is calculated to be a lower energy state using the methodology discussed in this paper Figure 12 UO 2 -KDO (a) monodentate complex is most consistent with observed interatomic distances [110], but the (b) bidentate configuration is calculated to be a lower energy state using the methodology discussed in this paper. (c) An outer-sphere configuration has a higher potential energy than either of the two inner-sphere model complexes.
has about twice as much IR intensity. The calculated frequencies are close to those predicted by de Jong et al.
(2005) for the 1:2 uranyl-acetate complex, but they represent a much smaller splitting between the U=O s and U=O as modes than the 1:1 complex in de Jong et al. [6]. However, the theoretical vibrational frequencies are consistent with the interpretation of de Jong et al. [6] that adding more H 2 O molecules of solvation would lower the calculated splitting.
Measured U=O vibrational frequencies for the 1:1 complex in solution are 861 and 954 cm -1 . [125]. Hence, the values calculated here are 3 to 5% higher than observed as expected, and the experimental and model splitting between the U=O s and U=O as modes are 93 and 83 cm -1 , respectively. One can conclude that the model is representing the vibrations of the uranyl-carboxylate complex fairly accurately, especially considering that the KDO molecule has been used here instead of acetate.
For the monodentate complex, the U=O s and U=O as modes were calculated to occur at 894 and 982 cm -1 (Fig.  13b); and hence, these vibrations may not be helpful in distinguishing between mono-and bidentate complexation because they do not change frequency between the two model complexes. Consequently, tracking changes in the U=O modes may not be useful in distinguishing complex structures. On the other hand, the C-O s and C-O as modes were computed at 1425 (three modes between 1420 and 1428 cm -1 ) and 1676 cm -1 (two modes at 1646 and 1686 cm -1 ). Estimating a scale factor or 0.96 based on the comparison to experiment above and assuming these calculations do as accurate a job on the monodentate complex vibrations, the C-O s and C-O as stretches should occur near 1415 and 1610 cm -1 . The calculated splitting between C-O s and C-O as in the monodentate case is almost 200 cm -1 -more than 100 cm -1 greater than the observed splitting of these bands for uranyl acetate [125]. Hence, vibrational spectroscopy can be used as a relatively quick and inexpensive method for determining uranyl carboxylate bonding configurations. The success of timeresolved laser-induced fluorescence spectroscopy (TRLFS; e.g., [113] and references therein) also suggests that future quantum mechanical calculations focusing on predicting fluorescence spectra would be worthwhile.

Conclusion
The aqueous speciation calculations demonstrate the accuracy of the computational methods employed for predicting these structures because the model results are consistent with available experimental data. The U(V) aq species maintains the UO 2 moiety, but the increased Hbonding to the axial O atoms is a precursor to the protonation of these atoms in the reduced U(IV) aq species. U(IV) aq is more stable in the triplet state over the singlet state by approximately 190 kJ/mol in these model calculations.
Uranyl-carbonate model complexes reproduce experimental EXAFS and NMR results provided explicit solvation (i.e., including H 2 O molecules) is included. The stability and structure of the ternary Ca-UO 2 -CO 3 aqueous complex is also predicted. Kinetic hindrance of U(V) disproportionation by the presence of carbonate was not consistent with the hypothesis that distortions of the U-O-C-O torsion angle stabilizes U(V); instead, the slower disproportionation rate is likely due to electrostatic repulsion between the highly negatively charged U(V)triscarbonate complexes (-5 e -). Model results are shown to be consistent with spectroscopic results on uranyl-organic complexes as well provided the first solvation shell around these complexes is included in the model. In particular, the NMR spectra collected in this study are consistent with an outer-sphere uranyl-catechol complex. The oxidation of catechol by UO 2 2+ (aq) was shown to occur through a H-radical mechanism as two phenolic H atoms are transferred in sequence to the axial O atoms of the UO 2 2+ . This results in a U(IV) aq and quinone. The intermediate quinone radical species can explain the observation of catechol oxidation and polymerization in the presence U(VI) in aqueous solutions [101].
For uranyl-cell surface complexation, uranyl is predicted to favor binding at phosphonate groups rather than phosphodiester groups. Although the inner-sphere bidentate configuration is predicted to have the lowest Gibbs free energy in these models, the differences between these configurations and outer-sphere associations is relatively small suggesting that a significant portion of the observed complexation could involve outer-sphere binding.