Infrared and Raman Diagnostic Modeling of Phosphate Adsorption on Ceria Nanoparticles

Cerium dioxide (CeO2; ceria) nanoparticles (CeNPs) are promising nanozymes that show a variety of biological activity. Effective nanozymes need to retain their activity in the face of surface speciation in biological environments, and characterizing surface speciation is therefore critical to understanding and controlling the therapeutic capabilities of CeNPs. In particular, adsorbed phosphates can impact the enzymatic activity exploited to convert phosphate prodrugs into therapeutics in vivo and also define the early stages of the phosphate-scavenging processes that lead to the transformation of active CeO2 into inactive CePO4. In this work, we utilize ab initio lattice-dynamics calculations to study the interaction of phosphates with the three major surfaces of ceria and to predict the infrared (IR) and Raman spectral signatures of adsorbed phosphate species. We find that phosphates adsorb strongly to CeO2 surfaces in a range of stable binding configurations, of which 5-fold coordinated P species in a trigonal bipyramidal coordination may represent a stable intermediate in the early stages of phosphate scavenging. We find that the phosphate species show characteristic spectral fingerprints between 500 and 1500 cm–1, whereas the bare CeO2 surfaces show no active modes above 600 cm–1, and the 5-fold coordinated P species in particular show potential diagnostic P–O stretching modes between 650 and 700 cm–1 in both IR and Raman spectra. This comprehensive exploration of different binding modes for phosphates on CeO2 and the set of reference spectra provides an important step toward the experimental characterization of phosphate speciation and, ultimately, control of its impact on the performance of ceria nanozymes.


INTRODUCTION
Cerium dioxide (CeO 2 ; ceria) nanoparticles (CeNPs) are promising nanozymes (i.e., nanomaterials with enzymemimetic activity) with a range of potential applications in healthcare technologies. 1−6 CeNPs show high nanozymatic activity to a range of biological processes, including superoxide dismutases, 7−9 catalases, 8,10,11 phosphatases, 12−14 and oxidase, 15−17 together with high thermal and chemical stability and low toxicity. Nanozymatic activity is dictated by the surface composition and morphology of the CeNPs, which cannot easily be controlled once the NPs are in vivo. 18 In particular, in the cellular environment, CeNPs are directly exposed to phosphates as strong electrolytes at millimolar concentrations 19 and as part of phosphate-bearing molecules including NAD + /NADH, DNA, RNA, ATP, ADP, and phosphorylated proteins. 12 This can cause issues with the retention of phosphatase activity, since adsorbed phosphates can build up on the surface of the CeNPs and sequester the Ce 4+ binding sites, which then cannot be reduced to hydrolytically active Ce 3+ sites, resulting in reduced activity. 19 It has also been shown that retaining phosphatase activity can enhance simultaneous oxidase activity, where the coupling of oxidation and hydrolysis reactions uses ATP hydrolysis as an energy source to power the complementary reactions. 20 Furthermore, retaining phosphatase activity is crucial to applications such as decomposing the organophosphates used in chemical weapons 21 and activating phosphate prodrugs in nanozymeprodrug therapies. 22 The inevitable interaction between CeO 2 and phosphates can lead to excessive phosphate adsorption and the conversion of active CeO 2 to inactive CePO 4 , 23 and this "scavenging" or "poisoning" could be detrimental to the safe in vivo application of CeNP-based nanozymes. In particular, while these materials are known to have different activities toward the scavenging of reactive oxygen species, the underlying mechanisms of action are still unclear. 24−27 Broadly, phosphorus poisoning is also an issue with CeO 2 -based catalysts as phosphorus causes irreversible deactivation of three-way catalytic converters 28 and selective catalytic reduction of nitrogen oxides. 29 The incorporation of P-containing species into CeO 2 can also convert the materials into monazite, which leads to the deterioration of the oxygen storage capacity. 30 Although significant research effort has been directed toward enhancing phosphatase activity, there has been little investigation into establishing the particle morphology and surface speciation in vivo. The phosphatase activity may be related to the Lewis acidity of the cerium ions 12,31,32 or to the reducibility of the different surfaces of CeO 2 , 13,33−35 both of which also influence phosphate adsorption. Both also depend on the facets exposed by the NPs, with phosphate adsorption strengths in the order {111} < {110} < {100}, but with phosphatase activity showing the opposite ordering, such that octahedral CeNPs (dominant {111}) are the most active and cuboidal particles (dominant {100}) the least. 35 In order to predict and control the therapeutic activity of CeNPs, understanding the interaction of different surfaces with phosphate species is therefore crucial. In this work, we employ density-functional theory (DFT) calculations coupled with lattice-dynamics modeling to investigate the surface speciation of the major {100}, {110}, and {111} surfaces of CeO 2 in the presence of phosphates. We predict the vibrational fingerprints of the surface-bound phosphate species by simulating IR and Raman spectra and demonstrate that these spectroscopic techniques can be a sensitive probe for studying the surface speciation and, therefore, establishing its impact on nanozymatic activity.

METHODOLOGY
Calculations were performed using DFT as implemented using the Vienna ab initio simulation package (VASP) code. 36−38 All calculations were carried out using the Perdew−Burke− Ernzerhof (PBE) generalized gradient approximation (GGA) exchange−correlation functional. 39,40 A plane wave cutoff energy of 500 eV was used, which gives total energies converged to 4 meV per CeO 2 formula unit compared to a 20% larger cutoff. The ion cores were modeled using projector augmented-wave (PAW) pseudopotentials, 41,42 with frozen cores of [Kr] for Ce and [He] for O and P. The electron localization in the Ce 4f orbitals was accounted for using the DFT + U correction of Dudarev et al. 43 with U eff = 5 eV. The use of PBE + U is standard for CeO 2 , 26,44−48 and the next "tier" of functionals, i.e., hybrids, would be impractically expensive for these calculations.
The bulk conventional unit cell containing 4 CeO 2 formula units was minimized at constant pressure with electronic total energy and ionic force convergence criteria of 1 × 10 −5 eV and 1 × 10 −3 eV Å −1 , respectively. A k-point grid with 5 × 5 × 5 subdivisions was used to sample the Brillouin zone. The minimized CeO 2 bulk structure retains the Fm3̅ m space group with a lattice constant of 5.498 Å, which is a slight overestimation compared to the experimental value of 5.411 Å. 49 Slab models of the {100}, {110}, and {111} surfaces were generated from the minimized bulk CeO 2 unit cell using the METADISE code. 50 The {100} and {110} surface models correspond to √2 × √2 and √2 × 1 expansions, respectively, each with 7 surface layers and 28 CeO 2 formula units, and the {111} model corresponds to a 1 × 1 expansion with 5 surface layers and 20 CeO 2 formula units. A surface dipole on the {100} was avoided by shifting half the oxygen atoms from the top to the bottom of the slab. This is a common procedure in simulations of such surfaces. 46 We note that the oxygen surface layer on the {100} surface is quite flexible, with many configurations, as discussed extensively in the literature. 48,51−55 We also note that the surface structure and, therefore, parts of the vibrational spectrum, may be sensitive to the NP size. The thickness of our slab models is 12−18 Å, and XRD data has shown that the diffraction lines of CeNPs become similar to bulk crystals for particle sizes of around 40−50 nm. 49 However, calculations on thicker slabs would be prohibitively expensive and are beyond the scope of the current study.
The surface slab models were optimized at constant volume, with the top and bottom layers allowed to relax, using electronic total energy and ionic force convergence criteria of 1 × 10 −5 eV and 1 × 10 −2 eV Å −1 . A k-point grid with 2 × 2 × 1 subdivisions was used to sample the Brillouin zones. The surface energies of the {100}, {110}, and {111} surfaces were calculated to be 1.45, 1.07, and 0.70 J m −2 , respectively, in agreement with the literature. 45 The surface stability follows the order of {111} > {110} > {100}, with the {111} surface being the most stable.
As discussed below, for each surface, we attempted to prepare several models with phosphate species in different binding configurations. In some cases, we found that phosphoric acid sometimes dissociated following adsorption (i.e., its hydrogen atoms transferred to surface oxygen atoms to form hydroxyl groups), and the degree of dissociation observed in each of the models we prepared is listed in Table S2. Structure images were generated using VESTA. 56 The infrared (IR) and Raman spectra of the bulk and surface models were simulated using the procedure outlined in ref 57. This uses harmonic lattice-dynamics calculations at constant volume. Γ-Point phonon calculations were performed on each structure using the finite-displacement approach implemented in the Phonopy code, 58 with a displacement step of 5 × 10 −3 Å. Accurate single-point force calculations were performed using VASP with a tight electronic total-energy convergence criterion of 1 × 10 −9 eV. All of the configurations examined in this study did not show any imaginary modes, indicating them to be dynamically stable. The IR activities were determined using the Born effective-charge tensors computed using the densityfunctional perturbation theory (DFPT) 59 routines in VASP. The Raman activities were computed from the change in the macroscopic dielectric tensor when the atoms in the structure are displaced by ±5 × 10 −3 Å along each normal-mode vector, with the tensors computed using DFPT. The method used to compute the Raman spectra follows the approaches of Skelton et al. 57 and Schilling et al. 47 To assign spectral bands, animations of selected phonon modes were generated using Phonopy and visualized using the VMD software. 60

Structure and Energetics of Adsorbed Phosphates on CeO 2 Surfaces.
For each of the three major surfaces of CeO 2 , we attempted to construct models with adsorbed phosphate in several different binding configurations. We refer to the configurations using the notation {hkl}-xOP− yO surf , where {hkl} is the Miller index of the surface, x is the oxygen coordination of the P species, either 4 or 5 for tetrahedral or trigonal bipyramidal, and y is the adsorption Although we mainly focus on spectroscopy, we also discuss the energetics for completeness. The adsorption energies and optimized structures of the phosphate configurations examined in this work are presented in Figures 1 and 2, respectively. The values of the adsorption energies are presented in Table S1, and the extent of H 3 PO 4 dissociation in each of the models is listed in Table S2. The strength of the interactions between the adsorbed phosphates and CeO 2 surfaces is given by the adsorption energy (eq S1). Although, for a given surface, the adsorption of phosphate is generally stronger than the adsorption of water, which has binding energies of around −0.5, −1.2, and −1.6 eV on the {111}, {110}, and {100} surfaces, respectively, 45 substitution effects would need to be considered in the liquid phase. 61 3.1.1. 5-Fold Coordinated P Species (5OP). When 5-fold coordinated, phosphorus is directly bonded to surface oxygen (*O surf ) in a trigonal bipyramidal configuration (Figure 2a surfaces has also been reported and is related to the ease of forming oxygen vacancies at the surface, i.e., of pulling oxygen atoms out of the surface, but without creating a vacancy. 52 Based on this, it is possible that 5OP could be an intermediate in a surface scavenging mechanism (we note that this species is stable on all surfaces, as indicated by the negative binding energies and the absence of imaginary phonon modes, and hence would be classified as an intermediate).
As of yet, 5OP species have not been implicated experimentally, with only 4fold coordinated P species known to strongly adsorb to CeNPs. 48 This is consistent with the 4OP species being predicted to have larger binding energies than those of the 5OP species for a given surface ( Figure 1). We therefore suggest that, if involved in the scavenging mechanism, 5OP species may form directly on the surfaces.

4-Fold Coordinated P Species (4OP)
. The adsorption of phosphate as a 4-fold coordinated species (i.e., P with a tetrahedral oxygen coordination, 4OP) on the three surfaces of CeO 2 is found to be more stable than adsorption as a 5OP species (Figure 1). This indicates that direct bonding between the P atoms and surface O atoms is energetically less favorable than adsorption via the phosphoryl O atoms and surface Ce ions.
We attempted to investigate three different adsorption configurations of 4OP, viz., monodentate, bidentate, and tridentate (Figure 2d Figure 3, with the bands assigned in Table 1 (the relative intensities are provided in Table S4).
The IR spectrum of bulk CeO 2 (Figure 3a (Figure 3c). The {110} surface also has a feature at 282 cm −1 , which was assigned to ρCe−O surf and ρCe−O bulk modes, the former of which corresponds to surface oxygen atoms moving out of the surface plane.
The marked differences in the IR spectra of the three surfaces may be due to their different natures. The {110} is a Tasker type 1 surface, with O and Ce ions arranged in the same charge-neutral layer, whereas the {111} and {100} surfaces are Tasker types 2 and 3, respectively, and are both terminated by O ions. 64 The rocking mode in the {110} surface further highlights the ease with which O ions can be removed from the flat {110} surface, 46 and is consistent with the larger displacement of surface O ions when 5OP species are adsorbed (Figure 2b), and also with the facile predisposition of this surface to facet into {111}/{110} surfaces observed experimentally. 52 The spectra of all of the bare surfaces show strong stretching signals corresponding to the longitudinal Ce−O stretching modes in both the surface and bulk layers. In addition, the spectra of the {111} and {100} surface models show transverse     Tables 2  and S6. The spectra display characteristic fingerprints of the phosphate species between 500 and 1500 cm −1 (Figure 4). 44 This region overlaps with the IR spectrum of phosphoric acid (H 3 PO 4 ) in aqueous solution, where there are a number of features between 800 and 1200 cm −1 (ν s P−OH p , ν as P−OH p , δP−OH p , δP−(OH p ) 2 , and νP�O p ; Table 2 and Figure  S1). 65  Vibrations are labeled according to the scheme outlined in the text. The frequency ranges indicate where a given mode is seen across all binding modes on all three surfaces. Vibrational frequencies for molecular phosphate species in the gas phase and in solution are also given for comparison. The spectra also display peaks above 3000 cm −1 (Figure S2), corresponding to the stretching and wagging of the hydroxyl groups (νO−H surf /νO−H p and ωO−H p /ωO−H surf ) and the predicted frequencies of these vibrations are in good agreement with the literature. 65 We note, however, that the majority of the existing reports on the IR spectra of phosphate species, both from experimental 65,67 and computational studies, 65 are for isolated or aqueous phosphates and are therefore not directly comparable to our simulations, but can nonetheless still provide a useful point of comparison. The largest difference we see between the adsorbed and free phosphates is that the ρO p −H p signals from the absorbed phosphates appear at higher frequencies, viz. 500−1000 vs <400 cm −1 ( Table 2). 65 We attribute this to the interaction of the phosphate hydroxyl groups with surface Ce atoms, and we similarly predict the rocking modes of the surface hydroxyl groups (ρO surf −H surf ) to occur above 500 cm −1 due to their being more constrained.
The stretching of the P�O p bond (νP�O p ) has been reported experimentally to occur at 1225−1230 cm −1 in dimethyl methylphosphonate on the {111} CeO 2 surface 68 and at 1165−1178 cm −1 in aqueous H 3 PO 4 . 65,67 We predict this vibration to occur at 1250 cm −1 in an isolated H 3 PO 4 molecule, which compares well to the 1254 cm −1 predicted in a previous computational study. 44 For the 5OP structures, the direct interaction between phosphorus and a surface oxygen atom causes an elongation of the P�O p bond from 1.48 Å in H 3 PO 4 to 1.55−1.59 Å in the 5OP structures (Table  S3), which leads to a suppression of the band intensity in the IR spectra of the surface-bound 5OP species as reported in a previous study of phosphate adsorption to the CeO 2 {111} surface. 44 The IR spectra of the stable monodentate, bidentate, and tridentate 4OP adsorptions on the three surfaces of CeO 2 are shown in Figure 4d−j. In the following, we discuss the spectra with respect to the denticity of the binding.
In the spectrum of tridentate {111}-4OP−3O surf , the band at 821 cm −1 is a complex mode involving the phosphate and hydroxyl groups (νP−OH p , ρO surf −H surf ), and is consistent with previous PW91 + U calculations that only observed νP− O p for phosphates interacting with the {111} surface. A similar feature is observed in the {110}-4OP−3O surf spectrum at 833 cm −1 , but not in the {100}-4OP−3O surf spectrum. We also do not see any δP−OH p features above 900 cm −1 , in contrast to the feature at 939 cm −1 predicted for phosphate adsorbed to the {111} surface in previous work. 44 We predict ρO p −H p features at 948 and 965 cm −1 in the {111}-4OP−3O surf and {100}-4OP−3O surf spectra, but no corresponding feature in the {110}-4OP−3O surf spectrum.
The simulated IR spectra of the bidentate 4OP species all show ν s P−OH p bands in the region of 796−830 cm −1 and ν as P−OH p features between 853 and 927 cm −1 . Although the {111}-4OP−2O surf and {111}-4OP−3O surf spectra show similar vibrational signatures to those predicted in previous PW91 + U calculations, 44 we predict a greater variety of stretching P−O p modes between 800 and 1300 cm −1 , and we also do not predict any P−O p bending modes in this region. This is supported by experimental 65,67 and computational modeling 65 studies on phosphoric acid in aqueous solution.
As noted above, {100}-4OP−1O surf is the only stable monodentate adsorption configuration (Figure 4j). In this spectrum, the ν as P−O p asymmetric stretch produces a highintensity peak at 1152 cm −1 , consistent with the ν as P−O p between 1120 and 1200 cm −1 in the spectrum of aqueous phosphoric acid. 67 From our simulations, it is clear that the P−OH p stretch is a common fingerprint for the phosphate species. The band frequencies in phosphates adsorbed to the CeO 2 surfaces are ∼60 cm −1 lower in energy than the measured values of 820− 890 and 905−1388 cm −1 for the ν s P−OH p and ν as P−OH p of aqueous phosphoric acid. 65,67 We thus infer that adsorption of phosphate to ceria surfaces via a P−O bond weakens the P− OH bonds and thus shifts the corresponding stretching modes to lower frequencies.
3.2.3. Raman Spectra of Bulk and Bare Stoichiometric Surfaces. The simulated Raman spectrum of bulk CeO 2 shows a signal at 434 cm −1 (Figure 5a) from the ν s Ce−O bulk symmetric stretch, which is in agreement with the PBE + U value of 437 cm −147 but 30 cm −1 lower than the measured value of 464 cm −1 . 69 The simulated Raman spectra of the bare {100}, {110}, and {111} surfaces of CeO 2 are presented in Figure 5b−d and the main features are assigned in Table 3 Figure 6). The region below 600 cm −1 is dominated by intense signals from vibrations of the surfaces themselves (cf. Figure 5). The Raman active features related to hydroxyl groups (OH p and OH surf ) dominate the spectrum above 2000 cm −1 and are more intense than the signals from 500 to 1500 cm −1 , although the difference in intensity is not as pronounced as in the spectra of isolated phosphates modeled using B3LYP. 65 We attribute this to the interaction with the CeO 2 surfaces restricting the atomic motion in the νP−O p and νP−OH p stretches and, consequently, affecting the polarizability.
As in the IR spectra, the ρO p −H p phosphate hydroxyl rocking modes are shifted to higher frequencies of 916−1079 cm −1 compared to 158−376 cm −1 in isolated phosphates (cf. Table 2), 65 and are found in the same region as the ρO surf − H surf modes. This provides further evidence that the constraints imposed by the interactions between the phosphate species and ceria surfaces have a significant effect on the vibrational frequencies.
The intensities of the ν as P−OH p asymmetric stretches have been reported to be up to 20× and >100× less intense than ν s P−OH p symmetric stretches in aqueous phosphoric acid, based on calculations and experiments, respectively. 65 Our simulated spectrum of isolated H 3 PO 4 predicts a 7× difference in the Raman intensity ( Figure S1). In contrast to the {110}-5OP−1O surf and {100}-5OP−2O surf spectra, the spectrum of the {111}-5OP−1O surf configuration does not show a ν as P− OH p . However, as the ν s P−OH p for {111}-5OP−1O surf is relatively weak, it may be that ν as P−OH simply has too low an intensity to feature in the simulated spectrum.
In experiments, treating agglomerates of 3−4 nm stoichiometric CeO 2 NPs with 1 M NaH 2 PO 4 results in the appearance of a weak Raman signal at 965 cm −1 , which was assigned to symmetric stretching of the adsorbed phosphate ions. 70 We predict the Raman-active symmetric stretching modes of all our configurations to lie between 750 and 1077 cm −1 , and that these overlap with the asymmetric stretching modes between 904 and 1160 cm −1 .
Notably, the spectra of the 4OP configurations on the {100} surface show stronger features between 500−1500 cm −1 than those on the {110} and {111} surfaces. This could be due to the polar nature of the {100} surface, where surface oxygen atoms can rearrange to accommodate adsorbates. The higher intensity of the Raman signals of the {100} 4OP configurations

CONCLUSIONS
In summary, this study demonstrates that IR and Raman spectroscopy are potential diagnostic tools to determine the presence and binding modes of phosphate species on the surface of CeO 2 nanozymes and, in particular, highlights the role of modeling in predicting and assigning spectroscopic bands.
Our calculations show that phosphates can interact with ceria in both 5-fold coordinated trigonal bipyramidal and 4fold coordinated tetrahedral configurations. The 5-fold coordinated species are stable but higher in energy than the 4-fold coordinated species and may represent high-energy intermediates during the removal of O atoms from the surface in the phosphate scavenging of ceria nanozymes. This could be investigated using molecular-dynamics simulations, but such simulations are beyond the scope of this work. Among the 4fold coordinated configurations, binding modes that maximize the interaction between the phosphate and the surface are consistently more stable for a given surface, with tridentate configurations being the most energetically favorable.
The simulated spectra show features related to both the surface and the bulk. Notably, the {110} surface has a unique IR signal at ∼280 cm −1 , corresponding to the rocking of the surface Ce−O bonds where the oxygen atoms move out of the surface plane and also has strong Raman stretching signals around 180 cm −1 that are not seen in the spectra of the other surfaces. These signals are indicative of the ease with which surface oxygen atoms can be removed from the flat {110} surface and its facile predisposition to facet into {111}/{110} surfaces. 52 The simulated IR spectra of the adsorbed phosphates show characteristic fingerprints between 500 and 1500 cm −1 . In this region, the trigonal bipyramidal configurations show a distinct νP−*O surf signal between 650 and 700 cm −1 due to the interaction between the phosphorus and surface oxygen atoms, and this signal can potentially be used to identify such species experimentally. In general, we also find that while diagnostic features tend to be associated with strong signals in the IR spectra, they are often weaker in Raman spectra, which may make it challenging to identify and distinguish different phosphate species using Raman measurements.
Overall, this study provides a "stepping-stone" toward using IR and Raman spectroscopy as a diagnostic tool for characterizing the surface speciation of CeO 2 nanozymes, and future work should focus on more complex and realistic models, e.g., taking into account surface coverage and morphology and/or implicit or explicit solvent effects. ■ ASSOCIATED CONTENT
Energetics of phosphate adsorption, dissociation reactions of adsorbed phosphate species, structures of adsorbed phosphates, and vibrational fingerprints of isolated phosphate, bulk CeO 2 , and adsorbed phosphates on CeO 2 surfaces (PDF) Figure 6. Simulated Raman spectra of phosphate species adsorbed onto the {111}, {110}, and {100} stoichiometric surfaces of CeO 2 with the different binding modes shown by the inset structures. The spectra are normalized relative to each other such that the highest absolute intensity across all the spectra is set to unity.