from Density Functional Theory - American Chemical Society

Dec 9, 2008 - theory (DFT) combined with the Marcus/Holstein theory of electron/polaron transfer. We previously applied a similar methodology to model...
5 downloads 0 Views 2MB Size

J. Phys. Chem. C 2009, 113, 346–358

Intrinsic Hole Migration Rates in TiO2 from Density Functional Theory N. Aaron Deskins* and Michel Dupuis Chemical & Materials Sciences DiVision, Pacific Northwest National Laboratory, Battelle BouleVard, K1-83, Richland, Washington 99354 ReceiVed: April 3, 2008; ReVised Manuscript ReceiVed: October 17, 2008

Formation and migration of hole polarons in bulk rutile and anatase TiO2 were modeled using density functional theory (DFT) combined with the Marcus/Holstein theory of electron/polaron transfer. We previously applied a similar methodology to model electron polarons and extended the approach to hole polarons. Holes were formed by removal of an O(2p) valence electron, and the quantum mechanical characterization of hole hopping (reorganization energy and electronic coupling) was carried out with the DFT+U method, a method that corrects for self-interaction errors and facilitates charge localization, combined with Hartree-Fock cluster calculations. Several elementary hole transfer processes along various directions were investigated for both rutile and anatase. Hopping along most directions in rutile and along one direction in anatase was found to be adiabatic in character, i.e., thermal processes coupled to phonons. The activation energies for these processes were found to be about twice as large as the activation energy for electron hopping in rutile [Deskins and Dupuis, Phys. ReV. B 2007, 75, 195212], in agreement with experiment that indicates electron diffusivity in TiO2 to be faster than hole diffusivity. Comparison of our calculated rutile hole mobility with experiment shows good agreement (theoretical value of 0.16 cm2/V · s at 1300 K). Lattice distortions around hole polarons were found to be larger than around electron polarons. Our results showed also that holes are thermodynamically more stable in the rutile phase, while electrons are more stable in the anatase phase. A hole trapping site with hemibond structure (i.e., a “molecular” polaron with the hole shared between two nonbonded oxygen sites) was also identified in anatase. We also characterized the formation of hole and electron polarons at the (110) surface. The energy barriers for hole and electron hopping at the surface are larger than in the bulk, by ∼0.07 and ∼0.11 eV, respectively. These studies form the basis for further development of models to describe polaron transport in TiO2 structures such as surfaces, interfaces, or nonperfect crystals. I. Introduction Titania, TiO2, is an important photoactive material that has been extensively studied since the discovery of its ability to split water1 in the 1970s. TiO2 is able to adsorb light to create excited electrons and use this energy gain for solar energy production2 or pollutant decomposition.3,4 Several factors limit the efficiency of such photoactive materials, including the trapping and recombination of charge carriers.5,6 Fundamental characterization of charge transport phenomena is therefore important for developing models and defining strategies to control and lower recombination rates. In this article, we applied density functional theory (DFT) to atomic level modeling of hole transport in TiO2. When titania adsorbs light to create excited conduction electrons, positively charged holes are also created. These conduction electrons and holes can migrate to surfaces and participate in surface reactions;4 they can become trapped at defect sites and also recombine for charge neutralization. The general description of these excited electrons and holes revolves around the polaron model, which has both experimental7-10 and theoretical11,12 evidence. These charge carriers become confined at localized sites and can then migrate through the material by hopping to adjacent sites. In TiO2 electrons localize at Ti4+ cations to form Ti3+, while holes form O- ions from O2- ions. There are several review articles available focusing on polaron * To whom correspondence [email protected].





formation in a variety of materials,13-17 including a review by Schirmer discussing hole polarons in oxide materials.18 A comprehensive understanding of electron/hole stability and mobility is also fundamental to understanding the photocatalytic activity of mixed-phase catalysts, such as the Degussa P25 catalyst. A mixture of the two titania phases has been shown to be more active for several reactions than the individual phases.19-22 One proposed explanation for the increased catalytic behavior of mixed-phase catalysts is that electrons selectively migrate from anatase to rutile, leading to an excess of holes on anatase.23 The excess holes react with surface adsorbed species and a net increase in reaction rate is seen compared to anatase alone. Another theory proposes the opposite: electrons selectively migrate to anatase with a similar reaction rate increase on rutile.24 Accurate polaron models at the atomic level may help address the issue of the stability of the polarons in the different phases. Other experimental work25 has shown that the diffusion of electrons is faster than the diffusion of holes in nanocrystalline TiO2. We previously modeled26 the formation of electron polarons in defect-free bulk rutile and anatase, including the transfer of the polarons along various directions. Our previous results represented intrinsic transfer rates for stoichiometric TiO2 while ignoring important “real” effects, such as defects, applied voltage, or photoexcitation. Other workers have modeled polarons in TiO2 clusters,27,28 but our work was the first step in an effort to characterize charge transport in bulk titania from first-principles. The results indicated that electron transport in rutile is thermal in character and the calculated activation energy

10.1021/jp802903c CCC: $40.75  2009 American Chemical Society Published on Web 12/09/2008

Intrinsic Hole Migration Rates in TiO2 of 0.09 eV was in good accord with the experimentally derived activation energy of 0.07 eV.9,15 In contrast we found electron transport in anatase to be nonadiabatic in character. We report here a similar study describing hole polaron transport in these materials. Direct experimental detection of holes has proven difficult due to the presence of defect sites that can trap holes, and theory is key to characterize nontrapped hole structures and diffusion. In this context we have already been able to use our DFT results to develop empirical potentials that describe well polaron structures and formation in the bulk as well as at surfaces and interfaces.29 We note that in prior work we successfully modeled electron and hole transport in chromia (R-Cr2O3),30 hematite (R-Fe2O3),31 and iron oxide32 using a similar methodology as the one contained herein and obtained excellent agreement with experiment for these materials, giving credence to our predictions for titania. For example our calculations reproduced very well the strong anisotropy in diffusivity in the basel plane compared to the c direction. There are other reports in the literature dealing with computational modeling of polarons using related theories and approaches, most often based on cluster or embedded cluster models and dealing with a variety of oxide materials.33-39 Computational issues often associated with embedded cluster models include the finite size of the clusters treated at the first principles level of theory and the less rigorous treatment of longrange effects. In contrast our approach employed in our earlier work and here is based on the periodic DFT+U methodology, which accounts consistently for long-range electrostatic effects while affording a treatment of electron correlation. We present in this work a computational characterization of hole polaron structures and intrinsic polaron hopping in both perfect bulk rutile and anatase. We also selected directions and sites to model electron and hole transfer at the (110) rutile surface. For the bulk and surface we applied the Marcus/Holstein theory40,41 of polaronic systems, calculated energy barriers and electronic couplings, two key parameters that govern polaron hopping, and extracted polaron hopping rates and polaron mobility. An overview of the formulation of Marcus/Holstein theory used herein to characterize polaron hopping is given in section II. Section III gives the details of our computational approach. Section IV contains our results which allow a comparison of hole hopping through rutile and anatase as well as a comparison of hole vs electron hopping. We also present in this section results for polaron hopping at the rutile (110) surface and the description of a hole trapping site in anatase. II. Model Background We described the hole polaron as an O- species. We explored the possibility of the hole binding at Ti to form Ti5+ but found this species to be unstable. These observations are consistent with the band structure of both rutile and anatase that shows the valence band to involve O(2p) states. A schematic representation of the electronic structure of the states involved in a basic polaron transfer process is depicted in Figure 1, which shows three relevant states. In equilibrium configuration qA, the missing electron is localized on one particular O site, and ΨA represents the state denoted {O--O2-}, while in configuration qB, the hole has moved to an adjacent site, and ΨB represents the state {O2--O-}. The crossing point or “transition” state between the initial and final states is designated qC. The labels qA, qB, and qC are short-hand notations representing the polaronic configurations. To describe the transfer process, we make use of Marcus/Holstein theory of which several reviews exist.13,40,41

J. Phys. Chem. C, Vol. 113, No. 1, 2009 347

Figure 1. Schematic diagram of h+ polaron transfer. In the initial state Ψa with structure qA, the hole is localized on the left O ion, while in the final state Ψb with structure qB, the hole is localized on the right O ion. At the transition state Ψc with structure qC, the hole is shared between the two O ions.

Figure 2. General features of Marcus/Emin/Holstein/Austin/Mott theory for symmetric polaron transfer. The potential energy surfaces of the initial state, ΨA, and final state, ΨB, are shown with equilibrium structures qA and qB. In our case the initial state is O--O2- and the final state O2--O-. The coincidence state or transition state between the two states is shown as qC. The reorganization energy, λ, corresponds to the energy of the final state ΨB at the geometric configuration qA. The energy barrier is shown as ∆G*. The adiabatic energy curves are shown as dashed lines, with the electronic coupling matrix element, VAB, given as half the energy difference between the two adiabatic states.

While Marcus theory provided the context for our work, several authors have put forth similar descriptions of charge transfer in solids. Emin and Holstein, as summarized by Austin and Mott,15 have given theoretical formulations for the description of polaron movement in solids. Levich, Dogonadze, and Kouznetsov have also given a picture similar to Marcus’ for electron transfer (see Dogonadze42 and the references therein). We note that several concepts in Emin/Holstein/Austion/Mott theory (EHAM) and in Marcus theory bear resemblance in spite of their different terminology. In the review by Austin and Mott it is notable that the diagram of the potential energy surfaces for the initial and final states that enter the EHAM polaron hopping theory is identical to Marcus diagram for electron transfer reproduced below (Figure 2). Both EHAM and Marcus consider the electronic states as residing in parabolic energy surfaces (see Figure 2 of the text by Austin and Mott). Electronic coupling J in EHAM’s theory is denoted VAB in Marcus’ theory, and the energy barrier is termed WH in EHAM’s theory rather than ∆G* in Marcus’ theory. Within these theories there are three key parameters: the reorganization energy λ, the diabatic energy barrier ∆G*, and the electronic coupling element VAB. These parameters and the general features of the theories are illustrated in Figure 2. Under

348 J. Phys. Chem. C, Vol. 113, No. 1, 2009

Deskins and Dupuis

these theories the initial and final states are considered to be individual states that vibrationally fluctuate about their equilibrium structures. The energy curves are assumed parabolic in shape. If an initial state has enough energy to surmount the energy barrier then the charge may be transferred to an adjacent site, and this gives rise to the final state. The energy barrier depends strongly on the interaction between the initial and final state wave functions. This interaction is quantified by the electronic coupling of the two states, VAB. In the absence of significant interaction between the initial and final states at qC (small VAB), the energy barrier is simply ∆G*. This is known as a nonadiabatic or diabatic transfer, and the main mechanism for electron transfer occurs via tunneling between the final and initial states. If there is strong coupling (large VAB), the energy barrier is lowered by the amount VAB. This is the adiabatic regime, and the hopping mechanism occurs via thermal, phononassisted hopping as the charge can easily move between initial and final states. In what follows we will refer in a generic way to the “energy barrier” of a charge-transfer process when we are not interested in emphasizing the nonadiabatic or adiabatic character of the process. This energy barrier corresponds to the height of the crossing point at qC in Figure 2. When needed, we will use the term “nonadiabatic energy barrier” for nonadiabatic processes and “activation energy” for adiabatic processes (akin to activated processes in chemical reactions). The other key parameter, λ, is the reorganization energy; it corresponds to the energy required to move the electron to the site of the final state while retaining the initial state geometry. In the diabatic limit with quadratic potentials, the reorganization energy was formally shown by Marcus to be ∼4∆G*. Whether the polaron transfer has adiabatic or diabatic character is approximately determined by the transmission probability κ. This transmission coefficient gives the probability of transfer from the initial state to the final state and the semiclassical treatment43 is given by


2P012 (1 + P012)


The probability of conversion to the final state per passage 0 , is given below for the selfthrough the intersection, P12 exchange process,44 based on the work of Landau45 and Zener46,47 2 P012 ) 1 - exp [-(VAB /hυn)(π3/λkBT)1⁄2]


In the above equation υn is the frequency for the phonon that assists in a given hopping process. It is specific to each of the elementary hopping processes along the various directions. The specific phonon modes of interest here are the breathing modes around the O and/or Ti atoms of the polaron sites in the initial and final states. When the two sites involved in the hole polaron transfer share bridging titanium atoms, the phonon or vibration mode that is relevant to the transfer process exhibits an asymmetric displacement of the bridging atoms. We capture this phonon mode in our linearized approximation of the polaron transfer coordinate specific to each process considered, as described later. An accurate definition of the specific hopping phonon modes could be obtained from a full determination of the force constants and the phonons for the system at hand. However because of the computational cost associated with the calculation of the phonon spectrum and more importantly because of the approximate nature of the criteria in eq 1, we chose instead to approximate the phonon mode frequency of interest with a single frequency that we use for all the cases

investigated here, i.e., the longitudinal optic mode phonon frequency for both rutile48 and anatase,49 corresponding to hυn values of 0.10 and 0.11 eV for rutile and anatase, respectively. We believe that this is a sufficiently good approximation for our purpose. In any case a value near 1 for κ indicates adiabatic transfer, while a small value for κ indicates a nonadiabatic process. In the adiabatic case the rate of polaron transfer is given as44

[ ]

-∆G*ad ket ) υn exp kBT


In eq 3 the strongest dependence of ket lies in the activation energy, given by

λ ∆G*ad ) - VAB 4


which further justifies our approximation of υn since rates will be dominated by the exponential term in eq 3. In contrast the diabatic transfer rate is given by40

ket )


1 2π (∆G0 + λ)2 |VAB|2 exp p 4λkBT √4πλkBT



For our case all O sites are equivalent, so ∆Go is zero. Beyond calculating the transfer rates, we can calculate the polaron mobility, µ, from the Einstein relation50


eD kBT


The diffusion coefficient, D, is calculated from51,52

D ) R2nket


R in the above equation is the distance between transfer sites, and n is the number of neighboring electronic accepting sites. In section III below, we give the details for the calculation of the key parameters of the theory, mainly the reorganization energy and the electronic coupling element. III. Computational Methodology We employed the plane-wave-based DFT code Vienna ab initio Software Package (VASP)53-56 for our calculations. We described the core electrons using the projector augmented-wave (PAW) method,57 and the valence electron configurations were Ti:3p64s23d2 and O:2s22p4. All calculations were performed allowing spin polarization. The exchange-correlation functional for this work was the functional of Perdew, Burke, and Ernzerhof (PBE).58,59 A finite temperature scheme with Gaussian smearing was applied to describe the electron occupancies. This standard method speeds up convergence, while introducing little error in the calculation due to the small degree of smearing. We tested our model for convergence and found it sufficient to describe the hole structures. We used supercells, 3 × 3 × 3 for rutile and 3 × 3 × 2 for anatase, for the calculations, as shown in Figure 3. Both had been determined to be adequate in our earlier work involving electron polarons.26 We also tested these supercell sizes against selected calculations using a 4 × 4 × 4 rutile cell. The geometries around the hole were nearly identical for the 3 × 3 × 3 and 4 × 4 × 4 cells. The calculated energy barriers differed by only 0.01 eV, further validating that cell size effects are small. Use of our supercells leads to cell lengths of 14.01 and 8.91 Å for the rutile cell and 11.46 and 19.42 Å for the anatase cell. These cells have 162 and 216 atoms, respectively. We also modeled a four layer (110) 3 × 2

Intrinsic Hole Migration Rates in TiO2

J. Phys. Chem. C, Vol. 113, No. 1, 2009 349

Figure 4. Charge density difference plot of the h+ polaron structure in rutile. The charge difference is taken between the cell with a missing electron and the cell without the missing electron. The density plot is shown for an isocharge of -0.19 e-/Å3.

Figure 3. Bulk structures of TiO2 used in this study. (a) 3 × 3 × 3 bulk rutile supercell. (b) 3 × 3 × 2 bulk anatase supercell. Ti atoms are depicted as gray spheres, while O atoms are red spheres.

surface with approximately 17 Å of vacuum space and a total of 144 atoms. Our calculations gave lattice parameters of 4.65 and 2.97 Å for rutile and 3.81 and 9.73 Å for anatase, all in good agreement with experiment.60 For our large bulk supercells we used 2 × 2 × 2 and 2 × 2 × 1 k-point meshes for rutile and anatase, and a 2 × 2 × 1 k-point mesh for the surface calculations. The method for creating the polaron structures involved removing an electron from bulk super cells. To this end we employed the DFT+U level of theory, and then allowed the system’s lattice to relax in response to this missing electron. The relaxed lattice effectively modeled the polaron structure with a localized hole. We experienced problems with standard DFT (GGA) in describing the polaron structures. The formal positive charge and spin created by removal of an electron and described by DFT becomes delocalized throughout the cell and not confined to any particular O ion. This is a well-known deficiency of DFT, a manifestation of the self-interaction of unpaired electrons in theories using functionals of the density in the generalized gradient approximation. To overcome the difficulty in localizing hole polarons, the Hartree-Fock method has been used by different researchers.61,62 Several approaches to correcting for electron self-interaction within DFT have been proposed, including the definition and inclusion of a self-interaction correction term, the ad hoc mixing

of a fraction of the exact Hartree-Fock exchange (a very popular approach in chemistry), or the inclusion of a Hubbardlike U correction term in the Kohn-Sham operator that reduces self-interaction errors. This is the DFT+U method of Liechtenstein, Anisimov, and Zaanen,63 and we used it herein. A similar method was used by Maxisch, Zhou, and Ceder35 to describe polaron transfer in LixFePO4 or Nolan and Watson64 to describe hole formation in Al-doped silica. Within this framework, the parameters U (Coulomb interaction) and J (exchange interaction) become the main determining factors for the magnitude of the correction. There is no universal value of U that makes DFT+U work equally well for all systems, and a common approach to finding the “optimal” parameters for a given system is to fix J and vary the U value. Analysis of properties of interest as a function of the effective U value (Ueff ) U - J) suggests the appropriate parameters. We chose a J value of 1 eV for our work. We detail our efforts to find an optimal Ueff under section IV. Verification of the polaron structure is necessary. We performed Bader’s atoms in molecules charge analysis65 using the program of Henkelman, Arnaldsson, and Jo´nsson.66 For example, the charges in the 3 × 3 × 3 rutile cell were found to be -0.74 e for the O polaron site, compared to a charge of -1.16 e for the nonpolaron O ions. In addition, a plot of the charge density difference (cdd) of the 3 × 3 × 3 rutile cell with and without the missing electron is given in Figure 4. As can be seen, the missing electron is taken from a 2p state on the O ion. The correct localization of the hole allows calculation of energy barriers and reorganization energies according to the following description. With the hole localized on initial and final sites that exhibit lattice distortions, we used a linear interpolation scheme to define a polaron transfer pathway

q(x) ) xqA + (1 - x)qB


qA, qB, and qC are vectors of the three-dimensional coordinates of all the atoms in the unit cell. In eq 8, x varies between 0 and 1, 0 e x e 1. x ) 1 corresponds to the final state; x ) 0 corresponds to the initial state; x ) 0.5 corresponds to the midpoint that happens to be the crossing point or “transition” state, owing to translation or symmetry equivalence of the final and initial states. This approach was shown67 to give a good representation of the polaron transfer coordinates in our earlier work on other metal oxides as it is consistent with the antiphase breathing phonon modes involved in the polaron transfer. The

350 J. Phys. Chem. C, Vol. 113, No. 1, 2009

Deskins and Dupuis

linearization of the coordinates proved acceptable in these earlier cases, and we expect it to work in the present case as well, owing to the small changes in O-Ti bond distances associated with the polaron distortion and hopping. As we mentioned earlier this linearized coordinate captures the physics of the specific phonon mode associated with the polaron transfer, the antiphase breathing polaronic expansion around the oxidized oxygen site and the contraction around the reduced oxygen site, those being associated with the changing negative charge on the oxygen atoms from O- to O2- and O2- to O-, respectively, for the initial and final states. The energy of the structure at qC is taken as the free energy of the crossing point or “transition state”, depending on the nonadiabatic or thermal regime of the polaron transfer. We assumed the enthalpic and entropic changes of the free energy to be negligible, which is typical for solid state calculations. We evaluated the reorganization energy as 4∆G*, assuming parabolic potential energy surfaces. That this is indeed the case will be shown later. We used the molecular code NWChem68,69 to calculate the VAB terms. NWChem contains a module based on the method of Farazdel et al.67 for calculating the electronic coupling. The formula used for this calculation is



HAB - SAB(HAA + HBB)/2 2 1 - SAB



with HAB ) 〈ΨA|H|ΨB〉 and SAB ) 〈ΨA|ΨB〉. We used cluster models for this work and extracted these clusters from the converged supercell periodic structures at the crossing-points, x ) 0.5, (the linearized pathway methodology was not applied to the clusters; only single-point calculations of the A and B states at the crossing point structures were carried out). In these clusters, we terminated the appropriate oxygen ions with hydrogen atoms along Ti-O bond lengths to maintain charge neutrality. The O-H distances were fixed at 0.96 Å when we performed these static calculations. We considered several cluster sizes, which we detail in section IV. These calculations were performed at the unrestricted Hartree-Fock (UHF) level with 6-311G basis sets68 for O and H and Ahlrichs’ VTZ basis set70 for Ti. We tested basis sets with more diffuse and polarization functions, but these larger basis sets had minimal effect on the VAB values. Hole localization for these cluster calculations was verified by Mulliken population analysis. IV. Results and Discussion A. Bulk Hole Calculations. 1. Energy Barriers from Periodic Calculations. We first set out to determine an appropriate U value for hole localization in TiO2. There is no universal U value that holds for every material or even for every property of interest. For example, studies detailing how different properties of CeO2 are affected by the choice of U are available.71,72 In our previous study26 involving electron polarons, we used the band gap and energy barriers as a guide to selecting a U value. In the current work, applying the DFT+U correction to the O 2p orbitals changed the band gap only slightly. The conduction bands are primarily composed of Ti d orbitals, and the effect of the DFT+U is to raise these bands relative to the valence bands, and thus applying DFT+U to the O 2p orbitals has very little effect on the band gap. Thus we used the calculated energy barriers as a function of Ueff, as detailed below, to choose a U value. An alternative approach could have been to choose a Ueff value based on location of the hole state in the band gap. Low Ueff values (6 eV) give gap states near the conduction band edge.

Figure 5. h+ transfer directions for (a) bulk rutile and (b) bulk anatase. The localized hole moves to an adjacent O atom in all cases.

While low Ueff values give more “physically” correct gap states, hole localization can be difficult, if not impossible (as we shall see below) at the initial and transition states. This is a problem inherent with the DFT+U method, and the U value must often be chosen to match a specific property. We have confidence in our choice of Ueff because, as shall be discussed further on, our rates (the property of interest in this work) show good agreement (vide infra) with experiment in the preference of hole and electron mobility rates (holes < electrons), as well as hole mobility rates and also because our previous work involving electron polarons showed good agreement with experiment. There are several possible hole transfers to nearest neighboring O atoms, as shown in Figure 5. For rutile we identified five unique transfers (labeled A-E), while for anatase there are six unique transfers (A-F). Our procedure involved localizing the hole at an O site, localizing it at one of the adjacent sites and then interpolating between the lattice distortions associated with the two holes, as discussed in section III, to get the crossing point or “transition” state (x ) 0.5 in eq 8). At times we refer to the interpolated path as the “reaction coordinate”. The results of one such calculation are displayed in Figure 6, which shows the energy along the reaction coordinate for transfer along direction B in rutile. The curves are parabolic with a nearly perfect quadratic polynomial fit, indicating that the application of Marcus/Holstein theory in the harmonic approximation is justified for this work. We calculated the energy barriers along the various directions in bulk rutile for two different Ueff values of 6 and 10 eV. Table 1 summarizes the results. Note that the energy barriers are the actual calculated barriers prior to applying the nonadiabatic/ adiabatic criterion and selecting the appropriate expression yielding activation energies or nonadiabatic energy barriers (eqs 3 or 5). The table indicates that the energy barriers generally increase with increasing Ueff. We attribute the anomalous results for transfer along direction B to incomplete charge localization for the Ueff value of 6 eV as we experienced difficulties in localizing the hole at the “transition” state in that particular case as indicated by the Bader charge (-1.10 e) of the O atom on

Intrinsic Hole Migration Rates in TiO2

J. Phys. Chem. C, Vol. 113, No. 1, 2009 351 TABLE 2: Calculated Reorganization Energies, Energy Barriers, and Electronic Coupling for Hole Transfer (Ueff ) 10 eV)


anatase Figure 6. h+ transfer energies and energy barrier along direction B in bulk rutile for Ueff ) 10 eV.

TABLE 1: Calculated Energy Barriers (eV) for Bulk Rutile along the Different Directions for Different Ueff Values direction

Ueff ) 6 eV

Ueff ) 10 eV


0.46 0.69 0.33 0.51 0.52

0.55 0.54 0.50 0.59 0.62

which we attempted to localize the hole. In contrast, O atoms with localized holes had charges near -0.75 e. At a Ueff value of 10 eV, we were able to achieve charge localization for the “transition” state along direction B, and the energy barrier for direction B is of similar magnitude to the energy barriers along the other directions. Figure 7 shows the energy barrier as a function of Ueff along the A direction in rutile. Again the energy barrier increases with increasing Ueff value, but there appear two distinct regimes in the figure. The slope is much larger below a Ueff value of 6 eV compared to above 6 eV. Below Ueff values of 3.5 eV localization of the hole proved impossible. At the larger Ueff values, the energy barrier increases little with increasing Ueff. On the basis of Table 1 and Figure 7, we used a Ueff value of 10 eV for the remainder of this work. We shall see that this value gives rates in good agreement with experiment. Once we selected the U value we calculated the energy barriers for all the selected directions for rutile and anatase. The energy barriers for hole transfer in rutile and anatase along with the reorganization energies, the electronic coupling values, and symmetry operations mapping initial and final states are given in Table 2. The anatase results are similar to the rutile results with calculated energy barriers near 0.5 eV. Barriers of this magnitude are typical for small polarons.13 Similar energy barriers for polaron transfer in different polymorphs of the same

Figure 7. Energy barriers (x ) 0.5) at different Ueff values. Results are along the [001] direction: direction A in Figure 5a.


Ereorg (eV)

Ebarrier (eV)

VAB (eV)

symmetry of transfer






2.16 2.00

0.54 0.50

0.29 0.34


2.35 2.49 2.37 2.04

0.59 0.62 0.59 0.51

0.34 0.03 0.00 0.00






1.68 2.35

0.42 0.59

0.25 0.06





translation/no symmetry mirror improper rotation/ mirror inversion mirror improper rotation improper rotation/ mirror translation/no symmetry inversion translation/no symmetry improper rotation

material were seen for electron transfer in TiO226 and Cr2O3.30 In our work involving electron polarons in TiO2, we attributed this observation to the similar local geometries and phonon structure of rutile and anatase. Rutile and anatase have similar Ti-O and Ti-Ti distances in the Ti octahedral geometries. The phonon frequencies are also similar, 811, 806, 500, and 458 cm-1 for rutile48 and 876, 755, 435, and 367 cm-1 for anatase.49 Thus, an argument linking energy barriers to lattice distortions suggests that the similarity in energy barriers between rutile and anatase is reasonable, both for electron and hole polarons, albeit differences exist for electron and hole polarons in a given material. In Figure 8 are displayed the structures (bond lengths and changes in positions relative to the neutral system) of the hole in the initial state for all the rutile and anatase processes considered in this work. The pictures illustrate that holes in rutile have significant structural asymmetry while holes in anatase have significant symmetric character. This observation has implications when defining the transfer pathways for the various directions as discussed below. The bond lengths in rutile lie between 1.82 and 2.17 Å, while in anatase they lie between 1.85 and 2.29 Å. Thus, both phases have similar bond lengths near the hole. The changes in the bond lengths upon hole formation relative to neutral bulk (with no holes) are shown in Figure 9. The bond length changes vary between -0.15 and +0.20 Å in rutile and between -0.12 and +0.29 Å in anatase. Distortions of similar magnitude (0.3 Å) have been observed for other point defects in TiO2 (O vacancies and Ti interstitials).73 These changes indicate that the hole distorts the crystal structure significantly. These distortions appear to propagate farther in rutile than in anatase, with bond distortions in rutile as large as +0.14 Å three bond distances away from the hole. While we can rationalize the lengthening of Ti-O bonds on the basis of electrostatic arguments, weakening the Ti-O bonds both for electron polarons in our earlier work and hole polarons in this work, we surmise that removing a valence electron from an oxygen atom strongly weakens the electrostatic interactions between O and Ti that are at the heart of the structure of stoichiometric TiO2, more so than adding an extra electron in an “empty” d atomic state of Ti. This results in lattice distortions associated with hole polarons that are larger and spread farther than those for electron polarons. We further modeled hole formation using a 4 × 4 × 4 rutile cell, and we obtained results

352 J. Phys. Chem. C, Vol. 113, No. 1, 2009

Deskins and Dupuis

Figure 8. Geometry near the hole in (a) rutile and (b) anatase. All numbers have units of Å. The H indicates the oxygen with the hole. Also shown are the structures for the unperturbed bulk with arrows indicating changes in atom locations upon hole formation in (c) rutile and (d) anatase.

consistent with the 3 × 3 × 3 cell results for the distortions near the hole or in the first shell region. Distortions in the second shell region and beyond, while smaller to the first shell, did not reach zero for both 3 × 3 × 3 and 4 × 4 × 4 cells as we would expect for small polarons. Presumably for large enough supercells the distortions should decay to zero far enough from the polaronic site. We further compared our work with our molecular dynamics simulations of hole polarons in TiO2.29 The classical potential parameters were derived by fitting the calculated energy barriers only, not the polaron lattice distortions. Despite only fitting to our energy barriers, the first and second shell distortions calculated from molecular dynamics calculations using a 5 × 5 × 6 cell were in good agreement with our current results using a 3 × 3 × 3 cell. In the 5 × 5 × 6 cell the distortions did reach nearly zero far from the hole. DFT calculations on 5 × 5 × 6 cells are currently too computationally demanding unfortunately. Other work on FeO74 showed that distortions are largest near the polaron and are minimized over 8 Å from the polaron. On the basis of the comparison of our DFT results and the MD simulations for smaller supercells and the MD results for larger supercells, we surmise that long-range distortions do not appear to affect hole structures and hopping energies significantly and that the first shell region dominates the structure/barriers. We would point out that large lattice distortions have also been reported for SiO2, a material with much less ionic character than TiO2 here.75 An important aspect of calculating the energy barriers is to ensure that appropriate symmetry operations map the initial and final hole sites. This is not an issue for anatase since the holes have a symmetric structure, but for rutile care must be taken. For example, the crystallographic group operator that relates the initial and final states for transfer A is a translation operator

(no symmetry). The initial and final holes differ simply by translation of the hole 2.97 Å along the [001] direction. In Table 2, we list the symmetry operations that map the initial state and the final state. This is used to define the final state structure (without calculating it explicitly) and to define the transfer pathway. If the wrong symmetry relationship between initial and final states is used, then the transfer pathway is incorrect. Care must be taken with respect to symmetry when modeling charge transfer. 2. Electron Coupling from Cluster Calculations. We extracted cluster models of the crossing-point structures (x ) 0.5) from the bulk periodic calculations described above and calculated the electronic coupling matrix elements, VAB. In our earlier work on electron polarons,26 we determined that clusters with four Ti atoms (Ti4) were sufficient to calculate the VAB values. We also tested the effect of cluster size on the hole VAB calculations. For transfer along direction A in rutile we obtained a VAB value of 0.00 eV for a Ti5 cluster and 0.02 eV for a Ti9 cluster. Also, for transfer along direction B in rutile we obtained a VAB value of 0.29 for a Ti4 cluster and 0.25 eV for a Ti6 cluster. These results suggest that very large clusters are not needed to obtain reasonable VAB values for the hole structures. The exception to this occurs when localization of the hole proved difficult. We found that localization of the hole was easier when we used symmetric clusters and for this we needed to use slightly larger clusters. This was the case for direction F in anatase. The clusters used for the VAB calculations are shown in Figures 10 and 11 along with the VAB values already gathered in Table 2. We obtained VAB values for rutile in the range of 0.00 eV (very small) and 0.34 eV. Two of the VAB values (directions A and E) are near 0 eV, while the other three values

Intrinsic Hole Migration Rates in TiO2

J. Phys. Chem. C, Vol. 113, No. 1, 2009 353

Figure 10. Clusters used for VAB calculations of h+ transfer in rutile. Grey spheres represent Ti atoms, red spheres O atoms, and white spheres H atoms. Figure 9. Geometry changes induced upon hole formation in (a) rutile and (b) anatase relative to the bulk structure without a hole. All numbers have units of angstroms. The H indicates the oxygen with the hole.

are near 0.3 eV. VAB values in anatase range between ∼0 and 0.48 eV. Most of the values in anatase are near ∼0 eV, and only two values (direction D and F) have any significant magnitude. These findings indicate that hole transfer processes in rutile have mostly adiabatic character and the processes in anatase have mostly nonadiabatic character. Similar behavior was seen for electron transfer; electron polaron transfer was found to be adiabatic in one rutile direction and nonadiabatic in the other direction, while electron polaron transfer was found to be nonadiabatic in anatase in all cases. The direction-specific hopping rates and mobilities were calculated using the appropriate formula determined on the basis of the transmission coefficient, all data are presented in Table 3 and discussed below. Electronic coupling elements are many-electrons quantities and their analysis is not straightforward although we have identified three factors that might affect the magnitude of VAB. They are the O-O distance, the degree of bridging Ti atoms, and the orientation of atomic states (orbitals) involved in the overlap. The localized hole state dominates the many-electron overlap between the initial and final states and therefore the electronic coupling between states A and B. In both rutile and anatase the directions (B in rutile and D in anatase) that have the smallest O-O distance have significant VAB values. Nearly all of the hole transfer directions have bridging Ti atoms and thus it appears that the degree of Ti bridging does not contribute to distinguishing the different VAB values. The effect of the orientation of the p orbitals, the hole states, however, can be seen for transfer along direction F in anatase. The O atoms involved in the hole transfer are nearly in the same plane (see

Figure 5) with a relatively short O-O distance of 2.74 Å. The angle for the O-Ti-O bond structure in this plane is 90°, so that the O p orbitals on the two O ions are nearly perpendicular or parallel to each other, resulting in a large value of the coupling term VAB. Figure 4 shows an example of such orbitals. This orientation for transfer along direction F allows for very favorable overlap between the p orbitals and leads to the very large VAB value (0.48 eV). This illustrates that when orbitals are aligned favorably a large electronic coupling results due to this overlap. In conclusion, the magnitude of the calculated electronic coupling along the various directions suggest that hopping along directions A and E in rutile is nonadiabatic in character, and adiabatic along directions B, C, and D. For anatase hopping in all directions except direction D is nonadiabatic in character. The adiabatic transfers should be the dominant modes of hole transport. The corresponding nonadiabatic or adiabatic energy barriers along with derived rates are given in Table 3 and discussed below. 3. Polaron Transfer Rates. We have calculated the transfer rates using eqs 1-7, as given in Table 3. A wide range of rates are seen in TiO2. The energy barriers given in Table 3 are now either nonadiabatic energy barriers for nonadiabatic processes (κ , 1) or activation energies for adiabatic processes (κ ≈ 1). In the nonadiabatic cases the energy barriers are the calculated barriers given in Table 2, while for the adiabatic processes the activation energies are derived from eq 4 and in general are significantly smaller than the actual calculated energy barriers. The barriers range from 0.16 to 0.62 eV in rutile and from 0.17 to 0.59 eV in anatase. Three of the five barriers in rutile are smaller than 0.25 eV and adiabatic in character, while the other two are larger than 0.50 eV and nonadiabatic in character. Only one barrier in anatase is small (0.17 eV) for the adiabatic process in the D direction, while the other four are larger than 0.50 eV,

354 J. Phys. Chem. C, Vol. 113, No. 1, 2009

Figure 11. Clusters used for VAB calculations of h+ transfer in anatase. Grey spheres represent Ti atoms, red spheres O atoms, and white spheres H atoms.

corresponding to nonadiabatic transfers in these directions. Thus hole transfers B, C, and D in rutile and D in anatase are adiabatic in character, while transfers A and E in rutile and A, B, C, and E in anatase are nonadiabatic. The rates, diffusion coefficients, and mobilities are typically several orders of magnitudes larger for the transfers with adiabatic character. Note that transfer along direction F in anatase represents a special type of hole polaron involving multiple O sites (termed “molecular polaron”), so we have not included the transfer rates in Table 3. A specific discussion of this hole polaron is given below. We have provided rates for two different temperatures (298 and 1300 K) in Table 3. Of note is that the nonadiabatic rates/ mobilities increase much more with temperature than the adiabatic rates and are only approximately 1 order of magnitude smaller than the adiabatic rates at 1300 K. Both adiabatic and nonadiabatic terms depend exponentially on E/kT, where E is the activation barrier or energy barrier (see eqs 3 and 5). The energy barriers are larger for nonadiabatic processes than the activation energies for adiabatic processes, which leads to larger slope of ket as a function of temperature. There is a preexponential term in eq 5 that depends on temperature, but the entire expression is dominated by the exponential term. These results suggest that, despite the large difference in activation energies, at high temperatures the nonadiabatic processes may become important. Charge diffusion is much more anisotropic (direction dependent) at low temperatures. Experimental measurements of charge mobility in TiO2 show a wide variety of values for both electron and hole mobilities (see Table 1 of Bak et al.76). There is no standard charge carrier mobility data because of the wide variety in TiO2 samples and measuring techniques. Defects in TiO2, ambient conditions (such

Deskins and Dupuis as O2 concentration), and temperature are all factors that affect charge mobility. Creation of n-type TiO2 is also easier than p-type TiO2, so most of the available charge mobility data is for electrons. Nonetheless, the Hall mobility, which occurs under an applied magnetic field, was measured to be 0.1 cm2/V · s at 1300 K for holes,77 while drift mobility, which occurs under applied electrical voltage, was measured to be 180 cm2/V · s at 1000 K also for holes.78 In contrast the recent results by Bak et al. give a hole mobility of 0.3 cm2/V · s at 1300 K.79 At 1300 K our largest intrinsic mobility in rutile was 0.16 cm2/V · s (direction C), which is in quite good agreement with the experimental data available considering the assumptions of our model. In summary, we identified in this work several directions (B, C, and D) of hole transfer in rutile that correspond to adiabatic processes and that are several orders of magnitude faster than along directions A and E for which hopping is of nonadiabatic character. For anatase one direction (D) corresponds to adiabatic hopping and is several orders of magnitude faster than nonadiabatic hole hopping along the A, B, C, and E directions in anatase. Finally we note that the adiabatic activation energy calculated for electron polaron hopping in rutile26 was in good accord with the experimentally derived activation energy (0.09 eV calculated and 0.07 eV experimentally,9,15 respectively). Electron hopping in anatase was found to be slow (in accord with the predicted nonadiabaticity of electron hopping in anatase). Here the smallest calculated activation energy for hole hopping in rutile or anatase is ∼0.16 eV, in accord with the experimental observation25 that hole transport in TiO2 is slower than electron transport. B. A “Molecular Polaron” Trapping Site in Anatase. A special finding was obtained in the case of hole transfer in the direction F in anatase. While a premise of this work is that holes are localized on one O atom only (albeit Marcus/Holstein theory does not require it), the calculations suggested the existence of delocalized hole trapping sites involving multiple O atoms in anatase. Similar two-site holes have been predicted for other oxides, such as KNbO3,80 Al-doped BaTiO3,81 or SiO2.33 For all the clusters under consideration in this work we used the Mulliken population analysis to verify charge localization. For the anatase F direction, we also confirmed that the initial and final wave functions calculated at the crossing point qC exhibited only one unpaired electron localized on the initial and final state O atoms respectively using a singular value decomposition to maximize the overlap between alpha and beta orbitals82 of the spin-polarized wave functions and determined the “best” unpaired hole orbital. The singular value decomposition and the Mulliken analysis both confirmed that the final and initial states had well-defined localized electrons and thus validated the calculation of the coupling element as discussed here. For the F direction in anatase, while the initial and final states were single O hole states, the electronic coupling between the two states was found to be very large (0.48 eV in our standard excised Ti7 cluster and 0.54 eV for a larger Ti10 cluster model). Accordingly, this large electronic coupling made the energy of the hole at the qC configuration barely larger if not lower than the energies of the initial and final states. (The fully calculated activation energy was +0.03 eV for the Ti7 cluster model and -0.03 eV for the larger Ti10 cluster model). The appropriate potential energy surface associated with the F direction is depicted in Figure 12, which shows that the lowest energy state for the F process corresponds to the hole being shared (trapped) between the two O sites originally thought to be involved in the F hopping process and adopting structure qC

Intrinsic Hole Migration Rates in TiO2

J. Phys. Chem. C, Vol. 113, No. 1, 2009 355

TABLE 3: Hole Transfer Parameters and Associated Rates 298 K phase/direction

∆G* (eV)


transfer mode

rutile/A B C D E anatase/A B C D E

0.55 0.25 0.16 0.25 0.62 0.59 0.51 0.52 0.17 0.59

0.00 1.00 1.00 1.00 0.30 0.00 0.00 0.04 1.00 0.69

nonadiabatic adiabatic adiabatic adiabatic nonadiabatic nonadiabatic nonadiabatic nonadiabatic adiabatic nonadiabatic


ket (s ) 3.75 × 102 1.41 × 109 4.17 × 1010 1.36 × 109 2.13 × 102 2.04 × 10-4 5.40 × 10-2 1.43 × 103 3.91 × 1010 4.63 × 103

as the low energy polaron structure. Figures 5 and 11 show the two O atoms involved in this F process. The hole polaron involved in the F process now appears to be better described as a delocalized hole trapping site involving at least two oxygen atoms. Indeed from any given O atom, there are in fact three other equivalent O atoms on to which the hole can become trapped suggesting that the hole may be delocalized even more extensively, among five O atoms. Attempts to converge a solution where a hole was shared between two O sites in a Ti7 cluster were unsuccessful. More work would be needed to characterize more precisely the possibility of this delocalized trapping site and the hole delocalized over two or more O atoms. In any case the large electronic coupling between the single O atom final and initial states as calculated here is responsible for the hole being delocalized over at least two O atoms, giving rise to the delocalized trapping sites (molecular polaron). As pointed out earlier the favorable orientation and overlap of the O(2p) states between the two O atoms is responsible for the large electronic coupling. The finding of the very large electronic coupling negates the description of the F process as involving hole localization on single O atoms, albeit the model is completely valid for all the other processes considered in this work. In a broader sense, the finding does not negate the Marcus/ Holstein two state model of hole polaron transfer. What it negates is a single O hole polaron model for the F process. A proper description of a hole polaron for the F process would simply require an initial and/or final state electronic wave function of the hole that would be delocalized over multiple O atoms (two or five or whatever the case may be), albeit such wave functions may be challenging to obtain computationally in practice. C. Comparison of Bulk Hole and Electron Polaron Results. We can now compare the current results for bulk holes with our previous results on electron polarons in bulk TiO2.26 The electron and hole polarons show different structural changes in rutile and anatase with regards to symmetry. In rutile the addition of an electron creates a symmetric polaron (symmetric lattice relaxation), while holes give rise to asymmetric polarons

Figure 12. Illustration of charge transfer when VAB is large. In this case the lowest energy state occurs at qC.

1300 K -1

D (cm /s)


µ (cm /V · s)

ket(s )

D (cm2/s)

µ (cm2/V · s)

6.60 × 10-13 9.22 × 10-7 1.32 × 10-4 4.18 × 10-6 3.45 × 10-13 3.25 × 10-19 1.03 × 10-16 4.14 × 10-12 4.77 × 10-5 1.34 × 10-11

2.57 × 10-11 3.59 × 10-5 5.13 × 10-3 1.63 × 10-4 1.34 × 10-11 1.27 × 10-17 4.01 × 10-15 1.61 × 10-10 1.85 × 10-3 5.22 × 10-10

2.86 × 109 2.60 × 1012 5.65 × 1012 2.58 × 1012 1.33 × 1010 5.17 × 103 1.15 × 105 4.42 × 109 5.91 × 1012 1.01 × 1011

5.03 × 10-6 1.70 × 10-3 1.78 × 10-2 7.91 × 10-3 2.15 × 10-5 8.23 × 10-12 2.19 × 10-10 1.28 × 10-5 7.15 × 10-3 2.92 × 10-4

4.49 × 10-5 1.51 × 10-2 1.59 × 10-1 7.06 × 10-2 1.92 × 10-4 7.35 × 10-11 1.96 × 10-9 1.15 × 10-4 6.39 × 10-2 2.61 × 10-3


(see Figure 8). In anatase we observe that both the electron polaronic structure and hole structure are symmetric. The differences in polaron relaxation symmetry between the two phases are related to their structures. Rutile has the smaller molar volume of the two phases. Both rutile and anatase have as basic building blocks TiO6 octahedrons. In rutile the octahedrons arrange so that octahedrons share two edges, while in anatase the octahedrons share four edges. Rutile also has TiO6 octahedrons that are less distorted than the octahedrons of anatase. Combined, these structural and crystallographic features may facilitate the stronger breaking of the symmetry in rutile compared to anatase upon removal of an electron. Next, we analyze the bond length changes near the two different polarons. The bond distortions near the polaron in rutile for the Ti3+ electron polarons and O- hole polarons both show increases in bond length relative to bond lengths in chargeneutral bulk. We rationalize this finding on the basis of electrostatic arguments. The increase in Ti-O distances for Ti3+ in electron polarons is consistent with the addition of an electron decreasing the electrostatic attraction between the Ti and O ions. Removing an electron on the O ion is also consistent with a decrease in the electrostatic attraction between the Ti and O ions. This weakening leads to a longer Ti-O bond length. In anatase we observe similar behavior in that the Ti-O bond lengths for both types of polarons increased. While both O- and Ti3+ species lead to increased bond lengths, the changes are larger in magnitude for hole polarons. The largest Ti-O bond length change associated with electron polarons in rutile is +0.10 Å, with other changes being in the range of +0.04 to -0.03 Å. The largest Ti-O bond length changes associated with a hole in rutile is +0.20 Å, with many bond lengths changes being near (0.1 Å. Similar results are seen in anatase; the largest bond length change in the hole polaron structure is +0.29 Å, compared to +0.15 Å in the electron polaron structure. Besides larger bond lengths in the hole structures, the distortions reach out farther. For the hole polarons the absolute changes in the third closest Ti-O bond lengths are still as high as 0.14 Å, whereas in the electron polaron the third closest Ti-O bond length changes are small, with absolute changes being less than 0.02 Å. These results suggest that hole polarons in TiO2 are “larger” than electron polarons. The activation energies also show different behavior for the two types of polarons. The nonadiabatic energy barriers for the Ti3+ polarons are all near 0.3 eV, while the nonadiabatic barriers for the O- polarons are all near 0.5-0.6 eV. One reason hole barriers could be larger is that the distortions are larger in hole systems and require more energy to distort to the crossing point. The larger barriers result in hole diffusion being slower than electron diffusion along the specific directions considered. This

356 J. Phys. Chem. C, Vol. 113, No. 1, 2009

Deskins and Dupuis

Figure 13. (110) rutile surface cell used for h+ transfer along bridging row oxygens and e- transfer along Ti5c rows. The transfer of holes and electrons from one surface atom to an adjacent atom, as reported in this work, is shown by the arrows. 25

finding is consistent with experiment, which shows that diffusion coefficients of electrons are larger than the diffusion coefficients of holes. There is also the question of polaron stability in the two phases of titania. The mixed-phase Degussa P25 catalyst is thought to be more reactive than single-phase catalysts because the presence of two phases leads to increased charge separation. There is likely a selective preference of the polarons to one of the two phases which leads to this increased charge separation and lowers the recombination rate.23,24 We calculated the relative stability of the polarons in the two phases using the following transfer reactions.

rutilehole + anataseneutral f rutileneutral + anatasehole (10) rutileelectron + anataseneutral f rutileneutral + anataseelectron

(11) We found that the energy for hole “transfer” from rutile to anatase according to eq 10 is +1.3 eV, while the energy for the electron “transfer” according to eq 11 is -0.71 eV. These results suggest that the holes prefer bulk rutile, while electrons prefer bulk anatase. This thermodynamic drive is a contributing factor to charge separation that would increase the reactivity of mixed-phase catalysts. We must acknowledge that these comparisons are for the perfect bulk phases, and defects or surface states may have a very significant effect (especially for nanosized particles) on the stability of the polarons in the different phases. D. Polarons at the (110) Rutile Surface. Surfaces play an important role in charge transport and trapping for particles as well as providing sites for charge transfer to adsorbate molecules. There is also evidence that polarons prefer to move toward surfaces.27,29 For these reasons we modeled the transfer of a hole on a rutile (110) surface. We used a 4-layer thick slab and a 2 × 3 surface cell, as shown in Figure 13. The hole was localized at a bridging oxygen (OB) and then transferred to an

Figure 14. Changes in geometry between the initial and transition state for polaron transfer along the rutile (110) surface. (a) h+ polaron; (b) electron polaron. The O atom with the hole polaron is indicated by H, while the Ti atom with the electron polaron is indicated by e-. All numbers are given in units of angstroms. Numbers in parentheses indicate the change in bond lengths for a polaron in the bulk structure.

adjacent OB. The surface OB atoms are coordinated to two Ti atoms, while bulk O atoms are coordinated to three Ti atoms. The transfer of a hole along the OB rows is equivalent to transfer along the [001] direction in the bulk or transfer along direction A. We calculated an energy barrier of 0.62 eV, compared to a value of 0.55 eV in the bulk. The surface atoms have fewer coordinating atoms contributing to the increase in activation energy (0.07 eV) for hole transfer on the surface as there are fewer neighbors onto which the polaron distortion is distributed. In Figure 14a we show the geometry changes that occur in taking the surface from the initial to transition state. Also shown are the bond length changes for hole transfer in the bulk. In general the geometrical distortions are larger for transfer on the surface compared to the bulk. Several bond length changes are near (0.1 Å on the surface while the bond length changes are typically smaller in the bulk. These structural changes in the surface as a result of the unsaturated coordination are consistent with the increase in energy barrier (0.07 eV) on the surface compared to the bulk. We also modeled the localization and transfer of an electron between surface five-coordinated Ti ions, as shown in Figure 13. The energy barrier for this transfer is 0.40 eV, compared to 0.29 eV between six-coordinated bulk Ti. Again, a higher energy barrier is observed at the surface, and we assign this to the difference between the surface and bulk distortions. We show in Figure 14b the bond distance changes that occur in taking the electron polaron from the initial to transition state. The bond distances in the plane of the surface are similar for both the surface and bulk polaron. The bond distances perpendicular to the surface show the greatest difference between the surface and bulk with bond distance changes being (0.10 Å for the surface and only -0.01 Å in the bulk structure. These larger structural changes are in line with the 0.11 eV difference in activation energy between the surface and bulk. In summary, these results illustrate that localization of charge at a surface is possible, opening the door for modeling the interaction of charged surfaces and adsorbates as well as further studies of charge transport at TiO2 surfaces.

Intrinsic Hole Migration Rates in TiO2 V. Conclusions By use of DFT+U we have localized holes at O ions within rutile and anatase TiO2 to form O- species. Transfer pathways were interpolated to obtain energy barriers for transfer of a hole to an adjacent site. We calculated several parameters within Marcus/Holstein theory to obtain direction-specific rates of charge transfer. We applied the U correction to the O 2p orbitals to ensure hole localization. Our results show that most of the calculated energy barriers for hole transfer in rutile and anatase fall between 0.5 and 0.6 eV. Calculated VAB values typically range between 0.0 and 0.3 eV. These results show that hole transfer in rutile is of adiabatic character along many directions, while hole transfer in anatase is adiabatic along only one direction. Thus hole polaron transfer in TiO2 is adiabatic (thermal, phonon coupled) for both rutile and anatase. The activation energy for hole polaron transfer in rutile and anatase is slightly higher than the activation energy for thermal electron polaron transfer in rutile, a finding consistent with experimental observation of electron polaron transfer being faster than hole polaron transfer. We recall also that the electron polaron transfer activation energy previously reported by us26 of 0.09 eV was in good accord with the experimentally derived activation energy of 0.07 eV. Our calculated hole mobility in rutile (0.16 cm2/ V · s) was also in good agreement with experimental hole mobility rates. Structural analysis shows that the lattice distortions (Ti-O distances) induced by the polaron are greater for holes compared to electrons. The magnitudes of the distortions are also greater in rutile. Our calculations show that the holes are thermodynamically more stable in rutile, while the electrons are more stable in anatase. We also calculated activation barriers for hole and electron transfer on the (110) rutile surface. Hole transfer along the surface has an activation energy 0.07 eV larger than the comparable bulk value, while electron transfer along the surface has an activation energy 0.11 eV larger than the bulk value. Finally a possible hole trapping site in anatase was identified as a “molecular polaron”, whereby a hole is shared between two or more O atoms as a result of a very large electronic coupling between these atoms. We have now characterized intrinsic hole and electron transport in defect-free bulk TiO2 using DFT, a common electronic structure method. There are other issues important for polaron transport and recombination, such as surface effects or particle size effects. We are currently using our DFT results to further develop empirical potentials to describe polaron transport so we can also model surfaces, interfaces, and other nonideal structures. This work will allow further characterization of charge transport in TiO2, which is important for designing and synthesizing efficient photocatalysts. Acknowledgment. We acknowledge Sebastien Kerisit and Kevin Rosso for stimulating discussions. Funding was provided by the Department of Energy, Office of Basic Energy Sciences. Computational resources were provided by the Molecular Science Computing Facility located at the Environmental Molecular Science Laboratory in Richland, WA and the National Energy Research Scientific Computing Center at Lawrence Berkeley National Laboratory. All work was performed at Pacific Northwest National Laboratory (PNNL). Battelle operates PNNL for the U.S. Department of Energy. References and Notes (1) Fujishima, A.; Honda, K. Nature 1972, 238, 37. (2) Kalyanasundaram, K.; Gratzel, M. Coord. Chem. ReV. 1998, 177, 347.

J. Phys. Chem. C, Vol. 113, No. 1, 2009 357 (3) Carp, O.; Huisman, C. L.; Reller, A. Prog. Solid State Chem. 2004, 32, 33. (4) Linsebigler, A. L.; Lu, G. Q.; Yates, J. T. Chem ReV. 1995, 95, 735. (5) Furube, A.; Asahi, T.; Masuhara, H.; Yamashita, H.; Anpo, M. J. Phys. Chem. B 1999, 103, 3120. (6) Frank, A. J.; Kopidakis, N.; van de Lagemaat, J. Coord. Chem. ReV. 2004, 248, 1165. (7) Yagi, E.; Hasiguti, R. R.; Aono, M. Phys. ReV. B 1996, 54, 7945. (8) Nowotny, J.; Radecka, M.; Rekas, M. J. Phys. Chem. Solids 1997, 58, 927. (9) Bogomolov, V. N.; Kudinov, E. K.; Firsov, Y. A. SoV. Phys. Solid State 1968, 9, 2502. (10) Bogomolov, V. N.; Kudinov, E. K.; Mirlin, D. N.; Firsov, Y. A. SoV. Phys. Solid State 1968, 9, 1630. (11) Schelling, P. K.; Halley, J. W. Phys. ReV. B 2000, 62, 3241. (12) Nelson, J.; Chandler, R. E. Coord. Chem. ReV. 2004, 248, 1181. (13) Shluger, A. L.; Stoneham, A. M. J. Phys.: Condens. Matter 1993, 5, 3049. (14) Holstein, T. Ann. Phys. 2000, 281, 725. (15) Austin, I. G.; Mott, N. F. AdV. Phys. 2001, 50, 757. (16) Honig, J. M. IBM J. Res. DeV. 1970, 14, 232. (17) Stoneham, A. M.; Gavartin, J.; Shluger, A. L.; Kimmel, A. V.; Ramo, D. M.; Ronnow, H. M.; Aeppli, G.; Renner, C. J. Phys.: Condens. Matter 2007, 19. (18) Schirmer, O. F. J. Phys.: Condens. Matter 2006, 18, R667. (19) Ohno, T.; Sarukawa, K.; Tokieda, K.; Matsumura, M. J. Catal. 2001, 203, 82. (20) Yan, M. C.; Chen, F.; Zhang, J. L.; Anpo, M. J. Phys. Chem. B 2005, 109, 8673. (21) Ohno, T.; Tokieda, K.; Higashida, S.; Matsumura, M. Appl. Catal. A 2003, 244, 383. (22) Wu, C. Y.; Yue, Y. H.; Deng, X. Y.; Hua, W. M.; Gao, Z. Catal. Today 2004, 863, 93–95. (23) Kawahara, T.; Konishi, Y.; Tada, H.; Tohge, N.; Nishii, J.; Ito, S. Angew. Chem., Int. Ed. 2002, 41, 2811. (24) Hurum, D. C.; Agrios, A. G.; Gray, K. A.; Rajh, T.; Thurnauer, M. C. J. Phys. Chem. B 2003, 107, 4545. (25) Duzhko, V.; Timoshenko, V. Y.; Koch, F.; Dittrich, T. Phys. ReV. B 2001, 64, 075204. (26) Deskins, N. A.; Dupuis, M. Phys. ReV. B 2007, 75, 195212. (27) Qu, Z. W.; Kroes, G. J. J. Phys. Chem. B 2006, 110, 8998. (28) Li, S. G.; Dixon, D. A. J. Phys. Chem. A 2008, 112, 6646. (29) Kerisit, S.; Deskins, N. A.; Rosso, K. M.; Dupuis, M. J. Phys. Chem. C 2008, 112, 7678. (30) Iordanova, N.; Dupuis, M.; Rosso, K. M. J. Chem. Phys. 2005, 123, 074710. (31) Iordanova, N.; Dupuis, M.; Rosso, K. M. J. Chem. Phys. 2005, 122, 144305. (32) Rosso, K. M.; Dupuis, M. J. Chem. Phys. 2004, 120, 7050. (33) Sicolo, S.; Palma, G.; Di Valentin, C.; Pacchioni, G. Phys. ReV. B 2007, 76, 075121. (34) Ramo, D. M.; Shluger, A. L.; Gavartin, J. L.; Bersuker, G. Phys. ReV. Lett. 2007, 99. (35) Maxisch, T.; Zhou, F.; Ceder, G. Phys. ReV. B 2006, 73, 104301. (36) Nolan, M.; Watson, G. W. Surf. Sci. 2005, 586, 25. (37) Mackrodt, W. C.; Middlemiss, D. S. J. Phys.: Condens. Matter 2004, 16, S2811. (38) Mackrodt, W. C.; Harrison, N. M.; Saunders, V. R.; Allan, N. L.; Towler, M. D. Chem. Phys. Lett. 1996, 250, 66. (39) McKenna, K. P.; Sushko, P. V.; Shluger, A. L. J. Am. Chem. Soc. 2007, 129, 8600. (40) Marcus, R. A.; Sutin, N. Biochim. Biophys. Acta 1985, 811, 265. (41) Marcus, R. A. ReV. Mod. Phys. 1993, 65, 599. (42) Dogonadze, R. R. Theory of Molecular Electrode Kinetics. In Reactions of Molecules at Electrodes; Hush, N. S., Ed.; Wiley-Interscience: London, 1971; pp 135. (43) Ulstrup, J. Charge Transfer Processes in Condensed Media. In Lecture Notes in Chemistry; Springer: New York, 1979; Vol. 10; pp 175. (44) Brunschwig, B. S.; Logan, J.; Newton, M. D.; Sutin, N. J. Am. Chem. Soc. 1980, 102, 5798. (45) Landau, L. Phys. Z. Sowjetunion 1932, 2, 24. (46) Zener, C. Proc. R. Soc. London, Ser. A 1932, 137, 696. (47) Zener, C. Proc. R. Soc. London, Ser. A 1933, 140, 660. (48) Porto, S. P. S.; Fleury, P. A.; Damen, T. C. Phys. ReV. 1967, 154, 522. (49) Gonzalez, R. J.; Zallen, R.; Berger, H. Phys. ReV. B 1997, 55, 7014. (50) Einstein, A. InVestigations on the Theory of the Brownian MoVement; Dover: New York, 1956. (51) Goodenough, J. B. Metallic Oxides. In Progress in Solid State Chemistry; Reiss, H. , Ed., 1971; Vol. 5; pp 145. (52) Heikes, R. R.; Johnston, W. D. J. Chem. Phys. 1957, 26, 582.

358 J. Phys. Chem. C, Vol. 113, No. 1, 2009 (53) Kresse, G.; Furthmuller, J. Comput. Mater. Sci. 1996, 6, 15. (54) Kresse, G.; Furthmuller, J. Phys. ReV. B 1996, 54, 11169. (55) Kresse, G.; Hafner, J. Phys. ReV. B 1993, 47, 558. (56) Kresse, G.; Hafner, J. Phys. ReV. B 1994, 49, 14251. (57) Kresse, G.; Joubert, D. Phys. ReV. B 1999, 59, 1758. (58) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865. (59) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1997, 78, 1396. (60) Howard, C. J.; Sabine, T. M.; Dickson, F. Acta Crystallogr. B 1991, 47, 462. (61) Pinto, H.; Stashans, A. Phys. ReV. B 2002, 65, 134304. (62) Gavartin, J. L.; Shluger, A. L. Phys. ReV. B 2001, 64, 245111. (63) Liechtenstein, A. I.; Anisimov, V. I.; Zaanen, J. Phys. ReV. B 1995, 52, R5467. (64) Nolan, M.; Parker, S. C.; Watson, G. W. Surf. Sci. 2006, 600, 175. (65) Bader, R. F. W. Acc. Chem. Res. 1985, 18, 9. (66) Henkelman, G.; Arnaldsson, A.; Jonsson, H. Comput. Mater. Sci. 2006, 36, 354. (67) Farazdel, A.; Dupuis, M.; Clementi, E.; Aviram, A. J. Am. Chem. Soc. 1990, 112, 4206. (68) Apra`, E.; Windus, T. L.; Straatsma, T. P.; Bylaska, E. J.; de Jong, W.; Hirata, S.; Valiev, M.; Hackler, M.; Pollack, L.; Kowalski, K.; Harrison, R.; Dupuis, M.; Smith, D. M. A.; Nieplocha, J.; V., T.; Krishnan, M.; Auer, A. A.; Brown, E.; Cisneros, G.; Fann, G.; Fruchtl, H.; Garza, J.; Hirao, K.; Kendall, R.; Nichols, J.; Tsemekhman, K.; Wolinski, K.; Anchell, J.; Bernholdt, D.; Borowski, P.; Clark, T.; Clerc, D.; Dachsel, H.; Deegan, M.; Dyall, K.; Elwood, D.; Glendening, E.; Gutowski, M.; Hess, A.; Jaffe, J.; Johnson, B.; Ju, J.; Kobayashi, R.; Kutteh, R.; Lin, Z.; Littlefield, R.; Long, X.; Meng, B.; Nakajima, T.; Niu, S.; Rosing, M.; Sandrone, G.; Stave, M.; Taylor, H.; Thomas, G.; van Lenthe, J.; Wong, A.; Zhang, Z. NWChem, A Computational Chemistry Package

Deskins and Dupuis for Parallel Computers, version 4.7; Pacific Northwest National Laboratory: Richland, WA 99352-0999, USA., 2005. (69) Kendall, R. A.; Apra, E.; Bernholdt, D. E.; Bylaska, E. J.; Dupuis, M.; Fann, G. I.; Harrison, R. J.; Ju, J. L.; Nichols, J. A.; Nieplocha, J.; Straatsma, T. P.; Windus, T. L.; Wong, A. T. Comput. Phys. Commun. 2000, 128, 260. (70) Schafer, A.; Horn, H.; Ahlrichs, R. J. Chem. Phys. 1992, 97, 2571. (71) Andersson, D. A.; Simak, S. I.; Johansson, B.; Abrikosov, I. A.; Skorodumova, N. V. Phys. ReV. B 2007, 75. (72) Loschen, C.; Carrasco, J.; Neyman, K. M.; Illas, F. Phys. ReV. B 2007, 75, 035115. (73) Cho, E.; Han, S.; Ahn, H. S.; Lee, K. R.; Kim, S. K.; Hwang, C. S. Phys. ReV. B 2006, 73. (74) Kerisit, S.; Rosso, K. M. J. Chem. Phys. 2005, 123. (75) Sulimov, V. B.; Sushko, P. V.; Edwards, A. H.; Shluger, A. L.; Stoneham, A. M. Phys. ReV. B 2002, 66. (76) Bak, T.; Nowotny, J.; Rekas, M.; Sorrell, C. C. J. Phys. Chem. Solids 2003, 64, 1069. (77) Odier, P.; Baumard, J. F.; Panis, D.; Anthony, A. M. J. Solid State Chem. 1975, 12, 324. (78) Yahia, J. Phys. ReV. 1963, 130, 1711. (79) Bak, T.; Nowotny, M. K.; Sheppard, L. R.; Nowotny, J. J. Phys. Chem. C 2008, 112, 12981. (80) Kotomin, E. A.; Eglitis, R. I.; Postnikov, A. V.; Borstel, G.; Christensen, N. E. Phys. ReV. B 1999, 60, 1. (81) Varnhorst, T.; Schirmer, O. F.; Krose, H.; Scharfschwerdt, R.; Kool, T. W. Phys. ReV. B 1996, 53, 116. (82) King, H. F.; Stanton, R. E.; Kim, H.; Wyatt, R. E.; Parr, R. G. J. Chem. Phys. 1967, 47, 1936.