A Density Functional Theory + U Study of Oxygen Vacancy Formation

7322

J. Phys. Chem. C 2009, 113, 7322–7328

A Density Functional Theory + U Study of Oxygen Vacancy Formation at the (110), (100), (101), and (001) Surfaces of Rutile TiO2 Benjamin J. Morgan* and Graeme W. Watson* School of Chemistry, Trinity College Dublin, Dublin 2, Ireland ReceiVed: December 21, 2008; ReVised Manuscript ReceiVed: February 6, 2009

Oxygen vacancy formation at the (110), (100), (001), and (101) surfaces of rutile TiO2 has been investigated using density functional theory with an on-site correction for strongly correlated systems (DFT + U). In agreement with experimental data, the reduced (110) surface shows an occupied defect state 0.7 eV below the bottom of the conduction band. The other reduced surfaces also show defect states in the band gap, with the defect state energies being strongly dependent on the choice of surface, and following the trend expected from crystal field arguments. For all the reduced surfaces, the excess charge associated with the defect states is primarily localized on two Ti sites neighboring the vacancy, formally reducing these to TiIII. For the (101) and (001) surfaces these Ti sites are geometrically inequivalent, and the corresponding gap states are separated in energy. Vacancy formation energies vary as ∆Evac(100) > ∆Evac(110) > ∆Evac(001) > ∆Evac(101). The variation in vacancy formation energy and gap state energies suggests potential differences between the surfaces in catalytic behavior for adsorbed reactant molecules. Introduction TiO2 is one of the most widely studied transition metal oxides1,2 and is often considered to be a model system because of the large number of published studies and the ease with which well-characterized sample surfaces can be prepared. TiO2 has potential applications in heterogeneous catalysis,3-5 gas sensor devices,6 and photocatalysis,7-10 where surface defects are active sites for the adsorption and dissociation of molecules.11-21 The most common polymorph is rutile, and for this phase the (110) surface is the lowest in energy2 and as such has been the subject of numerous studies. Because of its minimal energy the (110) surface is preferentially expressed in macroscopic crystal morphologies and has been assumed to dominate the surface chemistry. Minority surfaces however might display different catalytic behaviors;22 either increased yields or rates or improved selectivity, and may contribute significantly to the chemical behavior of samples. Where specific minority surfaces display desirable chemical behavior, it may be possible to develop methods to synthesize kinetically controlled samples that express these target surfaces in nonequilibrium proportions.23 Yang et al. recently reported one such technique, where the use of hydrofluoric acid as a morphology controlling agent in the synthesis of anatase TiO2 crystals resulted in the preferential expression of the normally disfavored (001) surface.24 The possibility of such directed synthesis of materials to achieve targeted surface chemistries means it is of interest to develop an understanding of the variation in redox chemistry across the range of surfaces. Bulk rutile TiO2 has a band gap of ∼3.0 eV,25 and defectfree surfaces have the same electronic structure as the bulk, with no surface states found in the band gap.26 Formation of surface oxygen vacancies is easily achieved through ion bombardment or controlled thermal annealing and quenching.17,27-30 The ease with which oxygen vacancies can be formed is thought to contribute to the catalytic behavior of TiO2, since this creates * To whom correspondence should be addressed. E-mail: [email protected] (B.J.M.); [email protected] (G.W.W.).

reactive sites into which molecules can adsorb. Oxygen vacancy formation at the (110) surface produces gap states that are seen in photoelectron spectra 2.3 eV above the valence band (∼0.75 eV below the conduction band).27,28 These states are also seen in energy loss spectroscopy as a peak related to Ti 3d-Ti 3d transitions in the region of 1-2 eV28,30-33 and in X-ray photoelectron spectroscopy as a shoulder on the lower binding energy side of the Ti 2p core feature,30 which disappears on subsequent exposure to O2.31 The defect-induced electronic states are believed to localize on neighboring Ti atoms, formally reducing these to TiIII,27,28 and recent theoretical studies of the (110) O vacancy support this interpretation,34-36 although these differ in their predictions of exact charge distribution. Data on the position of defect states on other surfaces is sparse. Photoelectron spectroscopy data have shown defect states for the (001) surface 1.0 eV below the conduction band,27 although the authors also reported low-energy electron diffraction data that suggested this surface had undergone a reconstruction, most likely to a facetted (110) structure. TiO2 is an example of a “strongly correlated” material having highly localized Ti 3d states, which are formally unoccupied in the stoichiometric system and are the predominant component of the conduction band. Electrons occupying these states experience strong mutual Coulombic repulsion, leading to highly correlated motion. Standard generalized-gradient approximation (GGA) and local-density approximation (LDA) functionals are parametrized to reproduce the interaction of a uniform electron gas and poorly describe the electronic energy in regions of space with strong local interactions. The exchange term used in standard GGA and LDA functionals is approximate and does not exactly cancel the self-term in the Coulombic sum. This results in the so-called self-interaction error, whereby electrons artificially repel themselves, causing partially occupied states to be excessively delocalized.37,38 For strongly correlated materials this error is acute, with standard density-functional theory (DFT) functionals often unable to achieve even qualitative agreement with experimental data,35,39-41 predicting metallic behavior for numerous systems which are experimentally

10.1021/jp811288n CCC: $40.75  2009 American Chemical Society Published on Web 04/01/2009

DFT Study of Oxygen Vacancy Formation observed to be insulating. A number of studies using standard DFT functionals have demonstrated this failing for O deficient TiO2.42-47 Two approaches that have been successfully applied to DFT calculations of the O-deficient TiO2 (110) surface to recover the experimentally observed features of the electronic structure are hybrid DFT,34 where a proportion of Hartree-Fock exchange is mixed with a standard DFT functional, and DFT + U,35,36 where the on-site Coulomb interaction is explicitly described within a Hubbard model.48 In this paper we report calculations of O vacancy formation at the (110), (100), (101), and (001) surfaces of rutile TiO2 and employ the DFT + U method to allow the strongly correlated nature of the Ti 3d orbitals occupied on reduction to be described. Localized defects are predicted for all four surfaces, and the energies and geometries of the associated electronic states are found to be strongly dependent on the surface termination. Methodology All calculations were performed using the DFT code VASP,49,50 in which valence electron states are described within a planewave basis set, with a cutoff of 500 eV. The PerdewBurke-Ernzerhof51 (PBE) exchange-correlation functional was supplemented by the rotationally invariant “+U” description of Dudarev et al.52 The projector augmented-wave method (PAW)53,54 was employed to treat valence-core interactions, with cores of [Ar] for Ti and [He] for O. Real space projection was employed for the PAW functions, with projection operators optimized to give an accuracy of ∼1 × 10-4 eV atom-1, and sufficient G vectors were included in the summation for Fourier transforms between real and reciprocal space that wrap-around errors were avoided. A U(Tid) value of 4.2 eV was applied to the Ti d states. This value has previously been shown to provide a description of O vacancies at the (110) surface that is in good agreement with the available spectroscopic data.28,30,35 A similar value of U ) 4.5 ( 0.5 eV was obtained by Bocquet et al.55 by optimizing the 2p core-level X-ray photoelectron spectroscopy (XPS) spectrum from a configuration interaction cluster model to reproduce experimental data. Satoshi et al.56 used an effective value of U ) 4.36 eV to describe the Ti 3d orbitals in LaTiO3 and SrTiO3. Calzado et al.36 proposed an ab initio value of U ) 5.5 ( 0.5 eV from LDA-based embedded cluster calculations. More recently Mattioli et al. have modeled oxygen vacancies in bulk rutile TiO2 with a U(Tid) value of 3.4 eV,57 calculated using the linear response approach of Cococcioni et al.38,58 It should be noted, however, that the optimal value of U varies with the choice of DFT functional and the details of the + U implementation, and caution should be taken in comparing values between different studies. The total energy and equilibrium volume of the bulk TiO2 rutile unit cell were computed with a Monkhorst-Pack grid of 4 × 4 × 4 k points. Zeropressure equilibrium lattice parameters were obtained by fitting a series of volume-energy data to the Murnaghan equation of state, giving values of a ) 4.669 Å and c ) 2.970 Å, compared with experimental parameters of a ) 4.593 Å and c ) 2.958 Å. The equilibrium bulk geometry was used to construct surfaces with slab geometries. These slabs were separated from their periodic images normal to the surface by a minimum vacuum gap of 15 Å, to give a pseudo-2D periodic system. The supercell expansions and lattice parameters used to construct each surface slab are listed in 1. All the surface slabs contained sufficient atoms to give a thickness of >14 Å, and Monkhorst-Pack

J. Phys. Chem. C, Vol. 113, No. 17, 2009 7323 TABLE 1: Slab Dimensions and Number of Atoms for the Four Stoichiometric Surface Systems surface

expansion

u/Å

V/Å

area/Å2

natoms

(110) (100) (101) (001)

4×2 2×5 2×2 3×3

11.68 9.34 11.07 14.01

13.21 14.85 9.34 14.01

154.29 138.70 103.39 196.28

240 330 144 297

k-point sampling of 2 × 2 × 1 was used, where the third reciprocal vector is aligned with the surface normal. Initial geometric relaxations were performed for each of the stoichiometric surfaces, in which all the atoms were allowed to move freely, and the minimization was deemed to be converged when the forces on each ion were smaller than 0.01 eV Å-1. After the stoichiometric surfaces had been relaxed, a single oxygen atom was removed from both surfaces of each slab, ensuring that an artificial dipole was not created perpendicular to the slab. These reduced slabs were then relaxed further, using the same convergence criteria. To calculate defect formation energies with respect to the stoichiometric surface, a single k-point (Γ-point) GGA calculation was performed on an isolated O2 molecule, using the same energy cutoff and convergence criterion as for the slab calculations, in a 15 Å side cubic cell. This gives a reference energy of -9.85 eV for O2(g). The PBE functional overbinds O2 by 1.4 eV per molecule, with reference to experiment,59 but no correction for this has been made in the values quoted in this work. While such a correction would alter the absolute formation energies it has no effect on the relative difference in formation energies between surfaces. Results The rutile structure is constructed from TiO6 octahedra, which form edge-sharing chains along the [001] direction and are corner sharing in the (001) plane. These octahedra are tetragonally distorted with axial Ti-O bonds longer than equatorial bonds. Different surface orientations determine how these octahedra are truncated at each of the four surfaces (Figure 1), resulting in variation in local geometry and coordination number for the surface atoms. For the low-energy (110) surface, Figure 1a, half of the surface Ti atoms are five-coordinated. The remaining six-coordinated Ti atoms are arranged in [001] rows with bridging surface oxygen atoms connecting them and completing their octahedral coordination. The (100) surface has only five-coordinate Ti atoms, Figure 1b, with two-coordinate O atoms bridging them. At the (110) surface the bridging Ti-O bonds are perpendicular to the surface plane. In contrast bridging oxygens at the (100) surface are bonded at an angle to the surface (61° in the relaxed stoichiometric structure). As is seen for the (100) surface, the (101) surface has 5-fold coordinated Ti atoms only, interconnected at the surface by two-coordinate bridging oxygen atoms. The five-coordinate surface Ti atoms at the (100) and (101) surfaces sit at the base center of TiO5 square pyramids, in both cases tilted from the surface normal. At the (100) surface these all have the same inclination, whereas the -O-Ti-O-Ti- chains that define the (101) surface are arranged in zigzag rows along the [101j] direction, giving the TiO5 square pyramids alternating tilt directions; Figure 1c. The truncated TiO6 octahedra that make up these (101) surface chains have alternately missing long and short Ti-O bonds, originally axial and equatorial local directions in the bulk; which makes alternating Ti surface sites along these chains enantiomeric. The (001) surface displays -O-Ti-O-Ti- rows running along the surface in the [110] direction, with alternating Ti atoms being six- and four-coordinate, Figure 1d. The six-coordinate surface

7324 J. Phys. Chem. C, Vol. 113, No. 17, 2009

Morgan and Watson

Figure 1. Coordination geometries at the (a) (110), (b) (100), (c), (101), and (d) (001) stoichiometric TiO2 rutile surfaces. Gray atoms are titanium, and red are oxygen. The blue polyhedra show the orientation of TiOn coordination environments.

Ti sites occupy the centers of TiO6 octahedra, whereas the removal of two oxygen atoms to give the four-coordinate surface sites results in these Ti atoms relaxing into the surface, to give distorted tetrahedral coordination. Every (001) surface O atom is bonded to one four-coordinate Ti and one six-coordinate Ti. Figure 2 shows slices through each of the four surfaces, normal to the surface plane. For each surface the oxygen considered for vacancy formation is indicated (Ov), along with the two neighboring Ti atoms. For the (101) and (001) surfaces, where the two Ti sites are inequivalent, these are distinguished as Tia, and Tib. The O atoms bonded to these neighboring Ti sites are also labeled. The Ti(a|b)-O distances for all the nearest neighbor O atoms are listed in Table 2, for the stoichiometric and reduced surfaces. Also listed are the coordination numbers for the Ti atoms neighboring the O vacancy, ncoord, and the mean bond length, jr. The calculated electronic densities of states (EDOS) for each of the four rutile surfaces are shown in Figure 3. All four EDOS show similar features, with distinct spin-polarized gap states appearing in the band gap. The EDOS for the reduced (110) surface shows a broad gap state (0.17 eV wide) centered 0.7 eV below the bottom of the conduction band, in agreement with the position of defect states observed experimentally.28,30 The charge density associated with this defect peak is shown in

Figure 2. Slices through the (a) (110), (b) (100), (c) (101), and (d) (001) rutile TiO2 surfaces. Each slice shows the oxygen removed for vacancy formation (Ov). The Ti atoms neighboring the Ov site are indicated, with Tia and Tib denoting inequivalent Ti sites for the (101) and (001) surfaces. The oxygen atoms coordinated to these Ti sites are also labeled as Oa,b,....

Figure 4a and is strongly localized on the two Ti atoms neighboring the O vacancy site, formally giving TiIII. The gap state for the reduced (100) surface is 1.4 eV below the bottom of the conduction band and is 0.12 eV wide. The narrowness of this feature compared to that for the (110) surface is due to the smaller overlap, and correspondingly weaker interaction between the two occupied Ti 3d orbitals. The Ti-Ti distance is only slightly larger for the (100) surface at 3.32 Å compared with 3.22 Å for the (110) surface. The relative orientation of the pairs of Ti orbitals, however, differs noticeably since this

DFT Study of Oxygen Vacancy Formation

J. Phys. Chem. C, Vol. 113, No. 17, 2009 7325

TABLE 2: Ti-O Bond Lengths for the Pairs of Ti sites (Tia and Tib) Neighboring the Oxygen Removed upon Vacancy Formation (Ov)a surface

site

ncoord

Ov

O2

O3

O4

O5

O6

jr

(110)stoic (110)red (100)stoic (100)red (101)stoic (101)red (101)stoic (101)red (001)stoic (001)red (001)stoic (001)red

Tia Tia Tia Tia Tia Tia Tib Tib Tia Tia Tib Tib

6 5 5 4 5 4 5 4 4 3 6 5

1.88Ov

1.88Oa 1.87Oa 1.89Oa 1.92Oa 1.92Oa 1.93Oa 1.88Ob 1.96Ob 1.87Oa 1.84Ob 1.95Oc 2.01 Of

2.06Ob 2.04Ob 1.94Ob 2.02Ob 1.92Oh 1.96Oc 1.92 Of 1.99Od 1.87Oa 1.95Oa 2.01Od 2.03Od

2.06Ob 2.04Ob 2.05Oc 1.98Oc 2.08Oc 2.03Og 2.08Oh 2.02 Of 1.91Ob 1.95Oa 2.01Od 2.03Od

2.09Oc 1.85Oc 2.05Od 2.00Od 2.08Og 2.08Od 2.08Od 2.04Oh

2.09Od 1.98Od

2.14Oe 2.05Oc

2.14 Of 2.22Oe

2.01 1.95 1.96 1.98 1.97 2.00 1.97 2.01 1.89 1.91 2.03 2.07

1.89Ov 1.88Ov 1.92Ov 1.91Ov 1.95Ov

a Ti coordination numbers, ncoord, and mean bond length, jr, for the stoichiometric and reduced (110), (100), (101), and (001) TiO2 rutile surfaces. The labile oxygen is denoted Ov, and all subscripts indicate sites in the local coordination environment as labeled in Figure 2.

Figure 3. Electronic DOS for the (a) (110), (b) (100), (c) (101), and (d) (001) surfaces. For the (101) and (001) surfaces, where there are two distinct gap states, the blue and yellow shades correspond to the same colored partial charge densities in Figure 4. The vertical dashed line indicates the top of the valence band, although the valence band in these EDOS appears to extend above zero eV due to a Gaussian smearing of 0.1 eV.

is determined by the local geometry of the TiIII sites. The (100) Ti 3d orbitals occupied by the excess electronic charge are titled 30° toward each other, whereas the (110) Ti 3d orbitals are tilted 50° toward each other to give greater overlap between orbital lobes (parts a and b of Figure 4). At the (101) surface the two-coordinate oxygen atoms occupy asymmetric sites: for a single oxygen vacancy the two neighboring Ti atoms are geometrically inequivalent (the O-Ti distances on the stoichiometric surface are 1.89 and 1.91 Å). Removal of a single surface oxygen to create a vacancy is accompanied by

a geometric relaxation that increases this asymmetry between the two neighboring Ti sites. Because of this, two distinct TiIII states are found for the reduced surface, both just above the valence band maximum, with a separation of 0.36 eV. The charge densities associated with these two peaks are shown in Figure 4c and are colored yellow and blue in correspondence with the peaks in Figure 3. The defect-peak-conduction-band gaps are 2.09 and 1.73 eV. For the (001) surface the two Ti sites neighboring the surface oxygen atoms are already inequivalent in the stoichiometric surface; one is initially six-coordinate, and the other is initially four-coordinate. Formation of oxygen vacancies again produces defect states, which have associated charge densities strongly localized on the Ti sites neighboring the vacancy site; Figure 4d. The three-coordinate surface site accommodates the charge associated with the lower energy defect state (yellow), and the bulklike five-coordinate Ti center accommodates the higher energy defect state (blue). It might be expected that a fourcoordinate stoichiometric surface Ti site would better stabilize the excess charge density introduced on reduction, leading to the charge density being associated with the three-coordinate site neighboring the vacancy, and a nearby stoichiometric fourcoordinate site. To examine this possibility a calculation was performed in which the bonds around a nearby four-coordinate site were lengthened in the starting structure, in an approximation of the expected local distortion were the excess charge to localize here. In addition the local magnetic moments were defined at the start of the calculation to encourage charge localization at this alternative site. This calculation however relaxed to give the same solution as described above, with the excess charge associated with the oxygen vacancy strongly localized on the two Ti sites neighboring the vacancy site. The lower peak lies at the same energy as the top of the valence band, with the higher energy peak split off by 0.40 eV, to give a 1.70 eV offset from the conduction band minimum. Despite being close in energy the valence band edge and defect band states can be distinguished by examining the distribution of the corresponding charge densities. The calculated surface energies for the stoichiometric surfaces with DFT(PBE) + U are listed in Table 3. The observed trend is Esurf(110) < Esurf(100) < Esurf(101) < Esurf(001). The same trend is obtained with pure PBE and is as seen in previous DFT calculations using different functionals.60-62 The order of increasing surface energies mirrors the variation in defect peak position offset from the conduction band minimum described above for the reduced surfaces.

7326 J. Phys. Chem. C, Vol. 113, No. 17, 2009

Morgan and Watson

Vacancy formation energies have been calculated according to ×
OO× + 2TiTi f VO•• + 2TiTi +

1 O 2 2(g)

(1)

where Kro¨ger-Vink notation has been used. Here the lattice site occupied is indicated by the subscript, and the superscript shows the electric charge of the species relative to that of the site it occupies in the stoichiometric crystal: × indicates no change to the charge, a single prime denotes a negative charge relative to the site in the stoichiometric material, a single dot denotes a positive charge, while a double dot indicates two positive charges relative to the charge of the appropriate standard site. These vacancy formation energies are listed in Table 4 and follow the trend ∆Evac(100) > ∆Evac(110) > ∆Evac(001) > ∆Evac(101). Previous calculations on the (110) surface have made predictions of the vacancy formation energy in the range of 2.6-4.5 eV,63 in broad agreement with the value of 3.66 eV presented here. It might be expected that oxygen vacancies form more readily on more unstable surfaces, with the vacancy formation energy following the same trend as the surface energies. This is not the case, however, with only weak correlation predicted. The vacancy formation energy is dependent on a number of interrelated factors, and identifying reasons for the predicted variation in O lability is not trivial. This energy has contributions from the binding energy of O in the stoichiometric surface, the geometric relaxation engendered when this O is removed, and the change in electronic energy that accompanies O vacancy formation as two electrons are transferred from the O to the neighboring Ti sites and the electronic structure responds to the change in electrostatic potential. Calculation of the O vacancy formation energy for a 3 × 2 × 4 bulk rutile cell (48 molecular units), with a 3 × 3 × 3 k-point mesh, and using the same calculation parameters gives 4.42 eV per vacancy. This suggests vacancy formation is more favorable at all the considered surfaces than in the bulk. This result is in general agreement with previous calculations that have compared surface and subsurface O vacancy formation. Quantitative differences in vacancy formation energies should not be compared directly; however, since these previous calculations neglected to correct for the self-interaction error for the occupied Ti 3d states and therefore have unphysical electronic structures.64,65 Discussion The electronic densities of states shown in Figure 3 show that the offset of the defect peaks from the conduction band decreases in the series (001) > (101) > (100) > (110). By consideration of the nearest neighbor bond lengths for the Ti sites which accommodate the excess charge on creation of the oxygen vacancies, this trend can be understood from crystal field considerations. The presence of O atoms in the coordination environment pushes up the energy of the occupied Ti 3d state associated with the reduced TiIII sites. In general, the smaller the Ti-O distance, the greater this effect, and we observe a correlation between the mean TiIII-O distance and the gap state-conduction band offset. The exception to this general trend is the three-coordinate TiIII site on the reduced (001) surface, which has a corresponding defect state at a lower energy than the occupied 3d states found on the other three surfaces, despite having a small mean Ti-O bond length. It is possible however that the low coordination number diminishes the net crystal field effect of the close O atoms. Crystal field considerations also

Figure 4. Partial charge density plots corresponding to the energies of the defect peaks in 3 for the reduced (a) (110), (b) (100), (c) (101), and (d) (001) surfaces. The isosurface levels are 0.05 e Å-3. For the (101) and (001) surfaces the occupied TiIII states are geometrically inequivalent and are colored yellow or blue in correspondence with the colors used to differentiate the defect peaks in the EDOS (Figure 3).

TABLE 3: Calculated Stoichiometric Surface Energies in J m-2 Using DFT(PBE) and DFT(PBE) + U for the (110), (100), (101), and (001) Surfacesa PBE + U PBE LDA60 LDA61 PBE61 PBE061 PW9162 (110) (100) (101) (001) a

0.83 0.92 1.18 1.54

0.58 0.67 0.99 1.31

0.89 1.12 1.40 1.65

0.89 1.20

0.42 0.69

0.55 0.83

1.88

1.39

1.59

0.50 0.69 1.03 1.25

Values from previous calculations using a range of functionals are included for comparison.60-62

DFT Study of Oxygen Vacancy Formation

J. Phys. Chem. C, Vol. 113, No. 17, 2009 7327

TABLE 4: Oxygen Vacancy Formation Energies in eV for the (110), (100), (101), and (001) Surfaces Calculated Using GGA + U surface (110) (100) (101) (001)

∆Ef(vac) 3.66 3.73 3.12 3.21

explain the splitting seen for the two TiIII sites on the reduced (100) surface. These sites have the same mean Ti-O bond length and are four-coordinate with approximate seesaw geometry. The low coordination however means the coordination environment has partially distorted toward tetrahedral geometry. For octahedral coordination crystal field effects place the t2g states below the eg states, and it is thus a t2g state which is preferentially occupied when TiIV is reduced, whereas the energetic order is reversed for tetrahedral coordination. For a four-coordinate TiIII state where two coordinated anions are removed, small distortions toward tetrahedral coordination increase the energy of the t2g states relative to the eg states. The larger the distortion the greater the shift in energy, with a large enough distortion expected to result in an eg state being occupied instead. The degree of distortion toward tetrahedral coordination can be estimated by measuring the O-Ti-O angle for the axial Ti-O bonds. This angle for the TiIII a site is 155°, and for the TibIII site is 169°, and the defect state associated with the former lies 0.36 eV above that of the latter. It might be expected that for the Ti sites that are reduced the change in electrostatic potential associated with the formal reduction from TiIV to TiIII will produce an increase in the local Ti-O bondlengths. This is not uniformly the case, however, and some of the bond lengths decrease upon vacancy formation. At the (110) surface the two Ti sites neighboring the vacancy are equivalent. These undergo a shortening in mean bond length, almost entirely due to the Ti sites moving away from the O vacancy site leading to a large decrease from 2.09 to 1.98 Å in the Ti-O bondlengths to the remaining bridging O sites (and the sites beneath these). There is little change in the mean bond length for the vacancy-neighboring Ti sites on the (100) and (101) surfaces, with an increase of 0.02 Å in both cases. For the (001) surface vacancy formation, the Tib mean bond length increases slightly from 2.03 to 2.07 Å, but this incorporates changes of opposite sign in the Ti-Oe and Ti-Of bonds which are equivalent in the stoichiometric system with lengths of 2.14 Å, one of which lengthens to 2.22 Å with the latter shortening to 2.01 Å, indicating a large distortion away from the octahedral lattice site with the Tib site being pulled into the bulk. The Tia site also undergoes a slight increase in mean bond length, again which incorporates large bond length changes in opposite directions. The application of DFT + U to the Ti d states is successful at recovering the localized polaronic nature of the occupied defect states, where the occupied and unoccupied states are well described by combinations of atomic-like Ti 3d functions. Since the valence band in TiO2 has partial Ti 3d character35 a secondary effect is to rehybridize and reduce the overlap with the Ti atomic-like functions that define the + U occupancies, resulting in an increased band gap. Whereas this demonstrates an improvement compared with pure DFT, the bandgap is still underestimated compared with experimental values, unless unphysically high values of U are used which exaggerate this effect (for a discussion of the effect of U on O vacancy-induced defect states in V2O5 when U is fitted to reproduce the

experimental bandgap see reference 66). Indeed, it would be surprising if such a simple correction was able to correct all the errors resulting from standard DFT functionals, especially those concerning the description of more diffuse electronic states. Since the occupied defect states and unoccupied conduction band states are both predominantly Ti 3d in character, it is expected that the distances between the defect state peaks and the bottom of the conduction band are more physically meaningful than the position of the defect states relative to the top of the valence band, for which substantial errors remain due to the use of DFT. The deviation between the experimentally observed valence-band-conduction-band gap and the DFT + U results is therefore expected to be mirrored by the valenceband-defect-state distance, with these defect states lying within the band gap for all four surfaces in experimental systems. In the absence of well-characterized experimental data for the minority TiO2 rutile surfaces higher-level theoretical methods such as GW calculations or DFT supplemented with the Gutzwiller variational approach67 or dynamical mean-field theory68 may indicate the position of these localized defect states more accurately. Concluding Remarks Density functional theory calculations have been performed on the (110), (100), (101), and (001) surfaces of rutile TiO2. For each surface stoichiometric surface energies have been calculated, and oxygen vacancy formation has been investigated. Surface energies increase as Esurf(110) < Esurf(100) < Esurf(101) < Esurf(001), while vacancy formation energies vary as ∆Evac(100) > ∆Evac(110) > ∆Evac(001) > ∆Evac(101), demonstrating only weak correlation between surface stability and oxygen lability. Oxygen vacancy formation requires the removal of a two-coordinate surface oxygen, and for all four surfaces the excess charge left at the surface localizes strongly at the two neighboring Ti sites, formally reducing these to TiIII and with the local geometry being strongly distorted. These localized states introduce new bands between 0.7 and 2.1 eV below the conduction band minimum. The precise position of these Ti 3d defect states varies strongly with the choice of surface, and follows the same trend as the stoichiometric surface energies, although there is no indication of a single parameter linking these phenomena. Since O-deficient surfaces possess an electron excess, these defect states constitute frontier orbitals, and are likely to be involved in redox reactions involving adsorbed molecules. Relative geometries of these reactive TiIII sites also vary between surfaces. The (110) and (100) reduced surfaces have pairs of symmetric occupied Ti orbitals, with variation in the angle these project from the surface, and the (101) and (001) reduced surfaces have asymmetric occupied Ti orbitals. These electronic and geometric variations suggest that O vacancy sites at different TiO2 surfaces may show strong variation in redox behavior, with possible corresponding modification of reactivity for experimental samples where nonequilibrium surface distributions are expressed. Acknowledgment. We acknowledge support from Science Foundation Ireland through the research frontiers programme (Grant No. 04/BR/C0216) and the principal investigate program (Grant No. 06/IN.1/I92). We also acknowledge the HEA and NDP for the PRTLI programmes IITAC and e-INIS and Trinity Centre for High Performance Computing for access to the TCHPC computational facilities. References and Notes (1) Diebold, U. Surf. Sci. Rep. 2003, 48, 53–229.

7328 J. Phys. Chem. C, Vol. 113, No. 17, 2009 (2) Henrich, V. E.; Cox, P. A. The Surface Science of Metal Oxides; Cambridge Univeristy Press: Cambridge, 1996. (3) Vannice, M. A.; Sudhakar, C. J. Phys. Chem. 1984, 88, 2429– 2432. (4) Raupp, G. B.; Dumesic, J. A. J. Phys. Chem. 1985, 89, 5240– 5246. (5) Linsebigler, A.; Lu, G. Q.; Yates, J. T. J. Chem. Phys. 1995, 103, 9438–9443. (6) Kim, I.-D.; Rothschild, A.; Lee, B. H.; Kim, D. Y.; Jo, S. M.; Tuller, H. L. Nano. Lett. 2006, 6, 2009–2013. (7) Fujishima, A.; Honda, K. Nature 1972, 238, 37&. (8) Linsebigler, A. L.; Lu, G.; Yates, J. T., Jr. Chem. ReV. 1995, 95, 735–758. (9) Fujishima, A.; Rao, T. N.; Tryk, D. A. J. Photochem. Photobiol. C 2000, 1, 1–21. (10) Thompson, T. L.; Yates, J. T., Jr. Chem. ReV. 2006, 106, 4428– 4453. (11) Henrich, V. E.; Dresselhaus, G.; Zeiger, H. J. Sol. Stat. Comm. 1977, 24, 623–626. (12) Smith, P. B.; Bernasek, S. L. Surf. Sci. 1987, 188, 241–254. (13) Kurtz, R. L.; Stockbauer, R.; Madey, T. E. Surf. Sci. 1989, 218, 178–200. (14) Smith, K. E.; Henrich, V. E. Surf. Sci. 1989, 217, 445–458. (15) Roman, E.; Desegovia, J. L. Surf. Sci. 1991, 251, 742–746. (16) Roman, E.; Desegovia, J. L.; Kurtz, R. L.; Stockbauer, R.; Madey, T. E. Surf. Sci. 1992, 273, 40–46. (17) Lu, G.; Linsebigler, A.; Yates, J. T., Jr. J. Phys. Chem. 1994, 98, 11733–11738. (18) Yamashita, H.; Ichihashi, Y.; Zhang, S. G.; Matsumara, Y.; Souma, Y.; Tatsumi, T.; Anpo, M. Appl. Surf. Sci. 1997, 121/122, 305–309. (19) Sorescu, D. C.; Rusu, C. N.; Yates, J. T., Jr. J. Phys. Chem. B 2000, 104, 4408–4417. (20) Rodriguez, J. A.; Jirsak, T.; Liu, G.; Hrbek, J.; Dvorak, J.; Maiti, A. J. Am. Chem. Soc. 2001, 123, 9597–9605. (21) Abad, J.; Bo¨hme, O.; Roma´n, E. Surf. Sci. 2004, 549, 134–142. (22) Ohno, T.; Sarukawa, K.; Matsumura, M. New J. Chem. 2002, 26, 1167–1170. (23) Tian, N.; Zhou, Z.-Y.; Sun, S.-G.; Ding, Y.; Wang, Z. L. Science 2007, 316, 732–735. (24) Yang, H. G.; Sun, C. H.; Qiao, S. Z.; Zou, J.; Liu, G.; Smith, S. C.; Cheng, H. M.; Lu, G. Q. Nature 2008, 453, 638–641. (25) Cronemeyer, D. C. Phys. ReV. 1952, 87, 876–886. (26) Henrich, V. E.; Kurtz, R. L. Phys. ReV. B 1981, 23, 6280–6287. (27) Tait, R. H.; Kasowski, R. V. Phys. ReV. B 1979, 20, 5178–5191. (28) Henrich, V. E.; Dresselhaus, G.; Zeiger, H. J. Phys. ReV. Lett. 1976, 36, 1335–1339. (29) Egdell, R. G.; Eriksen, S.; Flavell, W. R. Sol. Stat. Comm. 1986, 60, 835–838. (30) Go¨pel, W.; Anderson, A.; Frankel, D.; Jaehnig, M.; Phillips, K.; Scha¨fer, J. A.; Rocker, G. Surf. Sci. 1984, 139, 333–346. (31) Epling, W. S.; Peden, C. H. F.; Henderson, M. A.; Diebold, U. Surf. Sci. 1998, 412/413, 333–343. (32) Rocker, G.; Schaefer, J. A.; Go¨pel, W. Phys. ReV. B 1984, 30, 3704– 3708. (33) Eriksen, S.; Egdell, R. G. Surf. Sci. 1987, 180, 263–278. (34) Di Valentin, C.; Pacchioni, G.; Selloni, A. Phys. ReV. Lett. 2006, 97, 166803. (35) Morgan, B. J.; Watson, G. W. Surf. Sci. 2007, 601, 5034–5041.

Morgan and Watson (36) Calzado, C. J.; Herna´ndez, N. C.; Sanz, J. F. Phys. ReV. B 2008, 77, 045118. (37) Cohen, A. J.; Mori-Sa´nchez, P.; Yang, W. Science 2008, 321, 792– 794. (38) Cococcioni, M.; de Gironcoli, S. Phys. ReV. B 2005, 71, 035105. (39) Morgan, B. J.; Scanlon, D. O.; Watson, G. W. In press. (40) Scanlon, D. O.; Walsh, A.; Morgan, B. J.; Nolan, M.; Fearon, J.; Watson, G. W. J. Phys. Chem. C 2007, 111, 7971–7979. (41) Nolan, M.; Parker, S. C.; Watson, G. W. Surf. Sci. 2005, 595, 223– 232. (42) Ramamoorthy, M.; King-Smith, R. D.; Vanderbilt, D. Phys. ReV. B 1994, 49, 7709–7715. (43) Lindan, P. J. D.; Harrison, N. M.; Gillan, M. J.; White, J. A. Phys. ReV. B 1997, 55, 15919–15927. (44) Paxton, A. T.; Thieˆn-Nga, L. Phys. ReV. B 1998, 57, 1579–1584. (45) Rasmussen, M. D.; Molina, L. M.; Hammer, B. J. Chem. Phys. 2004, 120, 988–997. (46) Wang, Y.; Pillay, D.; Hwang, G. S. Phys. ReV. B 2004, 70, 193410. (47) Von Oertzen, G. U.; Gerson, A. R. Int. J. Quantum Chem. 2006, 106, 2054–2064. (48) Anisimov, I. V.; Zaanen, J.; Andersen, O. K. Phys. ReV. B 1991, 44, 943–954. (49) Kresse, G.; Furthmu¨ller, J. Phys. ReV. B 1996, 54, 11169–11186. (50) Kresse, G.; Furthmuller, J. Comput. Mater. Sci. 1996, 6, 15–50. (51) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865. (52) Dudarev, S. L.; Botton, G. A.; Savrasov, S. Y.; Humphreys, C. J.; Sutton, A. P. Phys. ReV. B 1998, 57, 1505. (53) Kresse, G.; Joubert, D. Phys. ReV. B 1999, 59, 1758–1775. (54) Blo¨chl, P. E. Phys. ReV. B 1994, 50, 17953. (55) Bocquet, A. E.; Mizokawa, T.; Morikawa, K.; Fujimori, A.; Barman, S. R.; Maiti, K.; Sarma, D. D.; Tokura, Y.; Onoda, M. Phys. ReV. B 1996, 53, 1161–1170. (56) Okamoto, S.; Millis, A. J.; Spaldin, N. A. Phys. ReV. Lett. 2006, 97, 056802. (57) Mattioli, G.; Filippone, F.; Alippi, P.; Amore Boapasta, A. Phys. ReV. B 2008, 78, 241201(R). (58) Kulik, H. J.; Cococcioni, M.; Scherlis, D. A.; Marzarzi, N. Phys. ReV. Lett. 2006, 97, 103001. (59) CRC Handbook of Chemistry and Physics, 79th ed.; Lide, D. R., Ed.; CRC Press, 1998. (60) Ramamoorthy, M.; Vanderbilt, D.; King-Smith, R. D. Phys. ReV. B 1994, 49, 16721–16727. (61) Labat, F.; Baranek, P.; Adamo, C. J. Chem. Theor. Comp. 2008, 4, 341–352. (62) Perron, H.; Domain, C.; Roques, J.; Drot, R.; Simoni, E.; Catalette, H. Theor. Chem. Acc. 2007, 117, 565–574. (63) Ganduglia-Pirovano, M. V.; Hofmann, A.; Sauer, J. Surf. Sci. Rep. 2007, 62, 219–270. (64) Oviedo, J.; San Miguel, M. A.; Sanz, J. F. J. Chem. Phys. 2004, 121, 7427–7433. (65) Pabisiak, T.; Kiejna, A. Sol. Stat. Comm. 2007, 144, 324–328. (66) Scanlon, D. O.; Walsh, A.; Morgan, B. J.; Watson, G. W. J. Phys. Chem. C 2008, 112, 9903–9911. (67) Deng, X. Y.; Wang, L.; Dai, X.; Fang, Z. Phys. ReV. B2009, 79, 075114. (68) Amadon, B.; Lecherman, F.; Georges, A.; Jollet, F.; Wehling, T. O.; Lichtenstein, A. I. Phys. ReV. B 2008, 77, 205112.

JP811288N