The Use of the “+U” Correction in Describing Defect States at Metal Oxide Surfaces: Oxygen Vacancies on CeO2 and TiO2, and Li-doping of MgO

There are many examples of defects in strongly-correlated metal oxides for which density functional theory predicts electronic structures that qualitatively disagree with experimental data. This behaviour arises from the self-interaction error inherent to standard density functionals, and is demonstrated by both pand n-type systems where the defect state is a small polaron associated with host lattice atoms. An approximate correction is to describe the electron–electron interactions in the orbitals of interest within the DFT+U formalism. This gives improved descriptions for systems where the states of interest are well represented by atomic-like orbitals. The qualitative failure of standard DFT and corresponding improvement achieved with DFT+U is illustrated for cases where the defect state is primarily associated with localised cation f and d states (O vacancies in CeO2 and TiO2) and anion p states (Li-doped MgO). [DOI: 10.1380/ejssnt.2009.389]


I. INTRODUCTION
Numerous metal oxide surfaces are technologically useful as heterogeneous catalysts. In many cases surface defects are essential to these catalytic processes, either modifying the general chemical properties of the surface or being directly involved in the catalytic mechanism. Defects both intrinsic-due to non-stoichiometryand extrinsic-due to doping-are thought to be active sites for catalysis. To model the interaction of reactant molecules with defective surfaces, it is first necessary to correctly describe the relevant surface defects in isolation. Often the computational method of choice for studying metal oxide surface systems is density functional theory. Yet there are a number of important defect systems where this methodology fails to give even a qualitatively correct description of the geometric and electronic structure, predicting delocalised electron or hole states for systems that are observed experimentally to be well localised.
This failure of standard DFT functionals to even qualitatively describe the electronic structure of specific metal oxides is historically related to Mott insulators; the archetype being NiO [1]. These materials are expected by conventional band theory to be conducting, but are found experimentally to be insulating. This is frequently the case for so-called strongly correlated systems, which have localised valence orbitals; typically metal d or rareearth metal f states. This localisation entails a strong on-site Coulomb interaction that opens an energy gap between occupied and unoccupied states resulting in insulating behaviour. The DFT+U formalism arose from attempts to describe these materials with a modified Hamiltonian in which the on-site electron-electron interaction described by the standard DFT expressions is replaced with a "+U " term that explicitly describes this Coulomb interaction [1][2][3] where λ Iσ i are occupation numbers for an orthogonal set of localised orbitals on atom I with angular momentum σ [3,4]. This correction penalises partial occupation of the orbitals of interest and hence favours disproportionation to fully occupied or fully vacant orbitals.
To explain the utility of the DFT+U formalism it is necessary to examine how standard DFT functionals incorrectly describe the change in energy of a system when basis functions are partially occupied. For an isolated atom allowed to transfer charge to or from some external reservoir, the change in energy of the system when part of an electron is added or removed is linear in δE/δn elec (this is necessarily so since partial occupation is described quantum mechanically as a time average over states with integer electron numbers), and is discontinuous at integer occupation values [5]. This behaviour is reproduced by the true (unknown) DFT functional, yet standard DFT functionals (LDA/GGA) deviate from this behaviour, having continuous derivatives with respect to n elec , and artificially lowering the energy for partially occupied states [6,7]. This gives rise to the self-interaction error, which is exhibited as an electron repelling itself, and results in artificial excessive delocalisation for partially occupied states. In the atomic limit the deviation is quadratic in n elec , and adding a quadratic "+U " term can correct for this, restoring the discontinuity at integer electron numbers for the derivative of the energy functional, and opening a gap of size U between occupied and unoccupied orbitals [3].
For extended (crystalline) wavefunctions in solids the situation is complicated since fractional occupation of atomic orbital functions can arise from hybridisation (and in fact must for those crystal states where all sites are equivalent by symmetry). A good approximate correction can still be achieved, however, for those Kohn-Sham states that are well represented by atom-centered atomiclike orbitals. Defects in strongly correlated metal oxides are often characterised by small polaronic states, where the electronic state associated with the defect is localised at a single atomic center, and hence are expected to be particularly well described by a single atomic-like atom-centered orbital (i.e. there is minimal hybridization). DFT+U is therefore expected to be particularly suited to an improved description of these states.
To use the DFT+U formalism it is important to determine a suitable value of U . First principles prescriptions have been developed for this purpose [6,8]. Such ab initio values of U are geometry dependent however, and although this is expected to be a small effect in the bulk, it could be a large factor for systems with low symmetry, such as surfaces and those containing defects [9]; different U values may in theory be required for each inequivalent site, and would vary during large structural relaxations. This makes a self-consistent ab initio approach computationally prohibitive for all but the simplest systems. The degree to which a single ab initio +U term corrects the descriptive error with which one is concerned additionally depends on how well the Kohn-Sham orbitals which determine that property are represented by appropriate atomic orbitals. This is exhibited even for uniform systems as different U values being required to achieve agreement with experimental data for different system properties [10,11], and can be acute when contrasting properties dependent on localised and extended states in defective systems. For example, a larger value of U is often required to fit band gaps to experimental data than to correctly localise a polaronic defect state [12,13]. Thus even a "mean-field" ab initio U value derived from stoichiometric bulk calculations is not guaranteed to provide an appropriate description of the properties of interest. A pragmatic approach adopted here, since we are interested primarily in the description of the defect state, is to select U based on the relative position in the band gap of the defect state energy level, in the spirit of an optimised hybrid DFT functional [14]; see e.g. ref. [15].
In this paper we describe three defective metal oxide systems that are widely used in catalytic technologies, and which illustrate the applicability of the DFT+U methodology to obtaining descriptions that are in qualitative agreement with experimental data [16]. These examples are O deficient CeO 2 and TiO 2 ; which are n-type with the defect state expected to be localised on Ce 5f / Ti 3d cation orbitals [17][18][19][20]; and Li-doped MgO; a p-type system with a defect state experimentally observed as a localised O 2p state [21].
All calculations have been performed using the periodic plane-wave DFT code VASP [22,23], with the interaction between core and valence electrons implemented through the projector augmented wave (PAW) approach [24]. The Perdew-Burke-Ernzerhof (PBE) [25] GGA exchangecorrelation functional was employed. In each system the relaxed geometries and electronic structures obtained using GGA and GGA+U functionals were compared, with the GGA+U calculations constructed by supplementing the PBE functional with the Dudarev "+U " term [4] applied to the atomic orbitals expected to characterise the defect states. The calculation models had slab geometries, where two equivalent surfaces of the bulk system are exposed and separated by a large vacuum gap to give a pseudo-2D construction [26].  CeO 2 plays an important role in automobile catalytic converters, where heterogeneous catalysis enables the degradation of environmental pollutants such as CO and NO 2 [28]. These processes involve O vacancies at the CeO 2 surface: CO is oxidised by the stoichiometric surface, to produce an oxygen vacancy; NO 2 binds into vacancy sites, restoring the surface stoichiometry and departing as NO. Oxygen vacancy formation on CeO 2 is formally associated with the reduction of two cerium atoms neighbouring the vacancy site from Ce 4+ to Ce 3+ : Experimental data have demonstrated that partial reduction of ceria surfaces results in the appearance of a peak in the valence band photoemission spectrum [17][18][19], which has been attributed to the presence of Ce 3+ ions. For the (111) surface, this peak lies 1.5 eV above the valence band.
A second metal oxide with important catalytic properties is TiO 2 , with applications in photocatalytic degradation of atmospheric pollutants [29]. Again, oxygen vacancies are formed readily, and these are expected to perform as active sites in catalytic processes involving adsorbed molecules. The formation of oxygen vacancies is accompanied by the appearance of an occupied gap state ∼ 0.8 eV below the bottom of the conduction band [30,31]. EPR data have been used to assign this peak to localised Ti 3+ states [32,33], formally making the reduction although there is debate over the geometric position of these reduced sites [15,34,35]. Many theoretical studies have shown that oxygen vacancies for both CeO 2 [36][37][38], and TiO 2 [39][40][41][42] are predicted by standard DFT to give delocalised electronic structures in disagreement with the experimental data: the excess electronic charge associated with surface reduction is predicted to occupy states at the bottom of the conduction band, indicating delocalisation, rather than isolated gap states as seen experimentally. The upper panels of Figs. 1 and 2 show the calculated electronic densities of states (EDOSs) for oxygen deficient CeO 2 (111) and TiO 2 (110), which illustrate this failure of standard DFT. In both cases GGA does not predict the experimentally observed gap state, instead placing the excess charge associated with the defect at the bottom of the conduction band, giving a metallic system, as for previous DFT calculations [36,43]. The charge densities associated with these occupied conduction band states are shown in the upper panels of Figs. 3 and 4. For the reduced CeO 2 (111) surface this charge is delocalised over all the Ce sites, whereas for the reduced TiO 2 (110) surface the delocalised excess charge is concentrated at subsurface Ti sites arranged in [001] rows.
The EDOS and excess charge density obtained for the reduced CeO 2 (111) surface using GGA+U (U = 5.0 eV) are shown in the lower panels of Figs. 1 and 3. A gap opens between the occupied and unoccupied Ce 5f states leading to a distinct peak in the band gap, in agreement with experimental UPS data, and the excess charge is found to be strongly localised at the two Ce atoms neighbouring the vacancy, formally giving a reduction to Ce 3+ . A similar difference is seen for the reduced TiO 2 (110) surface when described by GGA+U (U = 4.2 eV) (Figs. 2 and 4). Again a distinct gap state peak is seen in the EDOS ∼ 0.7 eV below the conduction band in good agreement with experimental observations. This defect state corresponds to a localisation of charge at the two Ti atoms neighbouring the vacancy site, formally giving Ti 3+ . In both reduced systems the change of functional from GGA to GGA+U has consequences for the optimised structure. For the TiO 2 (110) surface these are quantitative: The Ti-Ti distance for the two Ti atoms neighbouring the defect site decreases from 3.39Å to 3.22Å [15]. GGA+U predicts strongly reduced Ti 3+ ions that are larger than those partially reduced in the GGA calculation, leading to a greater repulsive overlap between these Ti ions and the remaining neighbouring O atoms that increases Ti-O distances and pushes the Ti ions together. There is also a weaker electrostatic repulsion between the two ions than for the GGA Ti (4−δ)+ species. For CeO 2 a qualitative difference is observed. The O atoms that terminate the CeO 2 (111) surface occupy the apices of trigonal pyramids; being coordinated to three equivalent Ce atoms. In the GGA calculation the excess charge is distributed over the entire system, and these three Ce atoms are equally partially reduced to Ce (4−δ)+ so retain their equivalence. With the GGA+U calculation two are reduced to Ce 3+ with the third remaining Ce 4+ . This improved geometric and electronic description of the oxygen vacancy site in CeO 2 has subsequently been used as a starting point in examining the molecular adsorption of NO 2 [44,45]; as a step towards a full description of the  Both Mott insulators, and metal oxides that exhibit ntype small polarons (such as CeO 2−x and TiO 2−x ), are strongly correlated systems characterised by partially occupied cation states. It is therefore perhaps unsurprising that the DFT+U methodology developed to describe the former is effective in providing qualitatively improved descriptions of the latter. The DFT self-interaction error can also lead to unrealistic delocalisation in systems where it is the anion states that are partially occupied. In systems such as metal oxides with p-type defects, where these are also expected to be characterised by small polaronic states, DFT+U is again applicable to providing descriptions in agreement with experimental data.
Li-doped magnesium oxide catalyzes the coupling of methane with high (73%) selectivity [46]. The commonly accepted mechanism proceeds via hydrogen extraction from methane to form methyl radicals which combine in the gas phase. High-resolution electron energy loss spectroscopy of Li-doped MgO shows a feature at temperatures above 1100 K at 1.6 eV, which has been attributed to the formation of O • O holes at the surface [47]. These are the proposed reactive sites for the extraction of hydrogen from CH 4 . The presence of trapped holes has been confirmed by EPR [21,[48][49][50] and by reaction with NO [46,51], and through electron nuclear double resonance techniques, with hole formation proceeding as Decomposition of the magnetic hyperfine interaction between Li 7 and O • O indicates a Li-O bond length of ∼ 2.5Å, compared with an Mg-O distance of 2.11Å [21]. This increased bond length is also predicted by unrestricted Hartee-Fock cluster calculations, which show localisation of the hole on the O atom [52]. DFT does not predict a localised hole state, instead distributing the unpaired electron over all the O atoms neighbouring the Li dopant [53][54][55].
That DFT delocalises hole states that are expected from experimental data to be localised has been seen in other p-type systems, for example Al-doped SiO 2 [56,57]. MgO is highly ionic and therefore the O 2p orbitals are well localised and again this problem can be approached within the DFT+U methodology, by applying a suitable +U term to the O 2p states. Figure 5 shows the EDOS for GGA and GGA+U (U = 7.0 eV) calculations [54]. With GGA the hole state is at the top of the valence band. The spin density for the system shows the charge distribution associated with this state (Fig. 6: upper panel), which is delocalised over all the O atoms directly bonded to the Li dopant. All the Li-O bond lengths are equal. The GGA+U EDOS in the lower panel of Fig. 5 instead shows two new features. An empty hole state is ∼ 1.0 eV above the top of the valence band, and a corresponding occupied state with opposite spin-polarity lies below the bottom of the valence band; both suggesting a localised hole state. This is confirmed by the spin density for the GGA+U calculation ( Fig. 6: lower panel). The excess spin is localised on a single O atom neighbouring the Li dopant, and the distance between this O atom and the Li has increased, moving the Li dopant off the symmetric lattice site, as expected from ENDOR data [21].

IV. SUMMARY
There are numerous technologically useful surface systems where chemically active defects are characterised by small polaronic states. Standard DFT typically fails to even qualitatively predict electronic structures for these systems that are in agreement with available experimental data. This is a consequence of the self-interaction error inherent to standard DFT functionals, which favours excessive delocalisation of partially occupied states. For systems where the state of interest is well represented by atom centered atomic-like orbitals, such as the small polarons considered here, a good approximate correction is provided by the DFT+U methodology. This defines a new effective Hamiltonian that explicitly describes the Coulombic on-site interaction in the atomic orbitals of interest. Selecting a suitable value of U requires care, and the degree to which the states that determine the properties of interest are represented by the localised orbitals to which the +U correction is applied is relevant both to the physical justification of the DFT+U approach, and the magnitude of U necessary to achieve the required outcome. DFT+U may be an inappropriate approach to correcting specific deficiencies of DFT in cases where the self interaction error does not play a role in problem being considered, i.e. the degree of electronic localisation is not related to the properties with which one is concerned. In the approach described here U is optimised so that the experimentally observed position of the defect state relative to the valence and conduction bands is reproduced [15]. Doing so recovers the qualitative features seen in experimental UPS spectra (by design), and also has consequences for the geometry of the defect state. Often the localisation of the electron/hole state is accompanied by strong local geometric distortions; these often include a lowering of the symmetry of the defect site, as atoms neighbouring the site are differentially reduced or oxidised. Such charge-localisation-driven lowering of symmetry is seen here for O vacancies at the CeO 2 (111) surface, and for Li doping at the MgO (100) surface.