Hole Polaron Migration in Bulk Phases of TiO2 Using Hybrid Density Functional Theory

Understanding charge-carrier transport in semiconductors is vital to the improvement of material performance for various applications in optoelectronics and photochemistry. Here, we use hybrid density functional theory to model small hole polaron transport in the anatase, brookite, and TiO2-B phases of titanium dioxide and determine the rates of site-to-site hopping as well as thermal ionization into the valance band and retrapping. We find that the hole polaron mobility increases in the order TiO2-B < anatase < brookite and there are distinct differences in the character of hole polaron migration in each phase. As well as having fundamental interest, these results have implications for applications of TiO2 in photocatalysis and photoelectrochemistry, which we discuss.

ABSTRACT: Understanding charge-carrier transport in semiconductors is vital to the improvement of material performance for various applications in optoelectronics and photochemistry. Here, we use hybrid density functional theory to model small hole polaron transport in the anatase, brookite, and TiO 2 -B phases of titanium dioxide and determine the rates of site-to-site hopping as well as thermal ionization into the valance band and retrapping. We find that the hole polaron mobility increases in the order TiO 2 -B < anatase < brookite and there are distinct differences in the character of hole polaron migration in each phase. As well as having fundamental interest, these results have implications for applications of TiO 2 in photocatalysis and photoelectrochemistry, which we discuss.

INTRODUCTION
Titanium dioxide (TiO 2 ) is a non-toxic and earth-abundant semiconductor that finds a number of photocatalytic and renewable energy applications, for example, as an electron transport layer in solar cell devices, 1 a photocatalyst for removing harmful NO x species from the environment, and a water splitting catalyst for hydrogen production. 2 TiO 2 has a number of polymorphs with anatase, rutile, brookite, and TiO 2 -B being the best studied. A wide variety of nanostructures including nanocrystals of different shapes, nanorods, and nanosheets are also commonly explored, offering additional parameters for material optimization. Many of the applications above rely directly or indirectly on the efficient transport of charge carriers through the material. For example, for photochemical applications, electrons and holes that are photogenerated in the bulk need to diffuse efficiently to surfaces where they can interact with adsorbates to facilitate chemical reactions. 3,4 For efficient solar cell devices, electrons must diffuse effectively through nanoporous titania layers to the electrode interfaces. 5 The strong electron−phonon coupling in TiO 2 means that in many cases, photogenerated electrons and holes can self-trap to form small polarons (quasiparticles consisting of a charge carrier and an associated lattice deformation). 6,7 Since the mobility of small polarons is much lower than that of band-like charge carriers (or large polarons) such self-trapping can enhance rates of recombination (radiative and non-radiative), negatively affecting material performance for applications. 8−10 However, the dominant mechanisms of charge transport in TiO 2 are still not well understood. Even if small polarons can form, it does not necessarily mean the mechanism of long range carrier diffusion is purely polaronic (i.e., hopping of polarons between neighboring lattice sites). Instead, it may be more favorable for the small polaron to thermally ionize before propagating as a band-like carrier and subsequently retrapping. 11−13 It is challenging to probe such processes experimentally and a deep understanding is currently missing. Aside from its fundamental importance, atomistic insights into these effects could ultimately help in the development of strategies to control and improve charge-carrier mobility for applications.
In this article, we tackle this problem theoretically and investigate the mobility of holes in the three common phases of TiO 2 where hole polarons are predicted to be stable (anatase, brookite, and TiO 2 -B). This is a challenging problem since it requires approximations to the many-body problem that are sufficiently accurate to correctly describe the crystal and electronic structures as well as the stability of localized versus delocalized charge carriers. Here, we use a hybrid density functional theory (DFT) approach, which has been explicitly parameterized for accurate description of small polarons in TiO 2 using a known property of exact DFT. 14 Using this approach, single site small polaron hole localization is favorable only in anatase, brookite, and TiO 2 -B and so we focus on these three polymorphs for simulation of hole mobility. Hole trapping in rutile on the other hand is not predicted to be favorable. Multisite (or molecular) hole polarons are also predicted to be stable in TiO 2 -H and TiO 2 -R but are not considered here since they are less relevant experimentally and hopping of multisite polarons would be considerably more complex to model. Our key findings are that long-range diffusion of holes in brookite proceeds via thermal ionization and band-like conduction since the thermal ionization energy of the small hole polaron (0.12 eV) is small compared to the activation energies for polaron hopping. The thermal ionization energy for small hole polarons in anatase is larger (0.24 eV) and therefore at low temperatures, holes diffuse by polaronic hopping, with preferential diffusion in the [100] and [010] directions, but at room temperature the migration has a mixed character. TiO 2 -B has the largest ionization energy (0.53 eV), suggesting polaronic hole hopping to be more dominant. The anisotropy in polaron mobility and differences between the three phases have important implications for understanding and designing materials for applications, which we discuss.
The remainder of the article is organized in the following way. In Section 2, we discuss the importance of hole polarons TiO 2 with respect to previous experimental studies and applications. We also highlight the challenges associated with accurate modeling of small polarons and discuss previous theoretical studies on TiO 2 . In Section 3, we present the computational methods employed before presenting the results in Section 4. Finally, we present a detailed discussion and our conclusions.

BACKGROUND
TiO 2 is usually intrinsically n-type due to the presence of oxygen vacancies or titanium interstitials, which act as shallow donors. However, holes can be generated on photoexcitation and their subsequent dynamics is important for understanding a number of processes of fundamental and practical significance. For example, the migration of photoinduced holes through the lattice to surfaces is essential for photocatalysis. 15 A number of studies have highlighted facet-dependent reactivity in nanocrystal systems, which has been interpreted in terms of the different stabilities of hole polarons on different surfaces. 16,17 Although relatively unexplored, anisotropy in hole diffusion may also play a role. The dynamics of photoinduced holes are critical in nanostructured TiO 2 photoelectrodes for water splitting. 18,19 Complex features in temperature-dependent photoluminescence spectra can also be interpreted in terms of carrier dynamics, where evidence suggests that in general holes are less mobile than electrons and the determining factor for recombination. 20,21 The examples given above provide clear evidence of the importance of hole polaron migration in TiO 2 and need for deeper atomistic insights.
There have been numerous DFT investigations into the stability and electronic properties of hole polarons in TiO 2 in perfect and amorphous crystals as well as at extended defects. 8,9,14,22−34 Standard local or semi-local exchange correlation functionals fail to describe hole polarons correctly due to significant residual self-interaction error (SIE), which tends to favor delocalized electronic states. Therefore, most of the above studies have been performed using modifications to standard DFT such as DFT + U or hybrid functionals incorporating some proportion of Hartree−Fock (HF) exchange (hybrid DFT). Predicted polaron properties are sensitive to the parameterization of these methods (i.e., the U in DFT + U or the fraction of HF exchange), with some methods overlocalizing and others underlocalizing. In principle, the generalized Koopmans condition (a known property of the exact functional) provides a means to constrain the choice of parameters in the functional, ensuring a more complete cancellation of the SIE 35−38 It has been recently shown for model systems in 1D (for which the many-body problem can be solved exactly) that hybrid functionals parameterized in this way provide ionization energies and electron densities in excellent agreement with exact solutions. 39 This approach has been employed to predict the properties of stable electrons and hole polarons in all known phases of TiO 2 , where the use of large supercells (>500 atoms) also minimizes finite size effects as well as at surfaces, twin boundaries, and nanocrystals. 14,31,40−42 DFT calculations of polaron hopping in TiO 2 are less common. A few studies using the Marcus−Emin−Holstein− Austin−Mott (MEHAM) formalism 6,43−45 together with DFT + U have explored electron and hole polaron migration in anatase and rutile. 22 For holes in anatase, hopping is predicted to be non-adiabatic for most intersite pathways with activation energies in the range 0.17−0.59 eV. 26 However, the size of the SIE in these calculations is unknown (and could be significant given the large Hubbard U value used of 10 eV). Although hybrid functional methods have been used together with MEHAM theory to investigate electron transfer process in some materials, for example, MgO and HfO 2 , 46−48 it is yet to be applied to bulk phases of TiO 2 (including brookite or TiO 2 -B). Furthermore, previous studies have not compared activation energies for thermal ionization of the polaron against activation energies for hole hopping at the same level of theory which is needed in order to determine the migration mechanism.

METHODS
Hybrid DFT calculations using the generalized gradient approximation were carried out using the CP2K simulation package. 49 Exact HF exchange is mixed into the exchange− correlation functional to overcome the issue of the SIE that is well known in DFT. We use a truncated PBE0 hybrid-DFT exchange−correlation functional that includes long-range corrections to the interaction potential (PBE0-TR-LRC) with a global 1/r dependence. 49 This defines a range of separations in the electron integrals to implement the HF exact exchange, and standard PBE is used outside this defined range. The truncation radius (R c ) must be less than half the length of the smallest supercell lattice vector to ensure that there is no interaction between neighboring cells, and we set our radius to 6.0 Å, shown previously to give converged structural and electrical properties. 14 The percentage of HF exact exchange (α) to include in these calculations was stringently parameterized by satisfying Koopmans' condition with a 0.05 eV tolerance for electron and hole polarons in each of the bulk TiO 2 phases, providing a band gap within 3% of the experimental value. 14 This was achieved by first optimizing the bulk crystal structure and hole polaron geometry for different choices of α. Then, for each case, the charge transition energy associated with adding an electron vertically (without geometry optimization) is computed and compared with the corresponding highest occupied eigenvalue. To satisfy Koopmans' condition (a necessary requirement for an exact functional), these two energies should be equal, providing a means to determine optimal values for each phase. We refer the reader to ref 14 where full details on this approach are provided. The optimal α for each phase is: 11.5% for anatase, 10.5% for brookite, and 12% for TiO 2 -B. Triple ζ basis sets were used for both titanium and oxygen for accurate calculations 50,51 and the The Journal of Physical Chemistry C pubs.acs.org/JPCC Article Goedecker−Teter−Hutter (GTH) pseudopotentials for both species available within CP2K. 52−54 A multi-grid approach for mapping products of Gaussians onto a real-space integration grid is used in CP2K, where the wide and smooth Gaussian functions are mapped onto a coarser grid, and the electron density is mapped onto the finest grid. The plane wave energy cut-off, a reference grid which controls the Gaussian mapping onto the multi-grid, is set to 60 Ry. Five multi-grids are used, and the plane wave cut-off is sufficiently converged at 600 Ry for the finest level of the multi-grid. The cell vectors and bulk geometries were optimized using the Broyden−Fletcher− Goldfarb−Shanno method where the forces were converged to less than 8 × 10 −4 Ry/a 0 (0.02 eV/Å) and the electronic convergence was set to 1 × 10 −6 Ry per self-consistent field (SCF) cycle.
Our previous work has shown that small hole polarons can form in anatase, brookite, and TiO 2 -B, and the calculated bulk structures and corresponding hole polaron geometries are taken as starting points for this investigation, see Elmaslmane and McKenna for computational details. 14 The supercells we use are approximately cubic in shape and contain 600 (anatase), 576 (brookite), and 720 (TiO 2 -B) atoms with the length of lattice constants in the range 16−25 Å. Since only the hole is introduced to the system without an accompanying electron, finite size effects from the electrostatic interaction of the hole with periodic images are present. However, we have performed careful supercell scaling calculations in our previous work and shown that the large supercells we employ are sufficient to minimize finite size effects. 14 Furthermore, calculation of barriers to hole hopping are far less affected by finite size effects since the hole remains localized throughout transition. To model polaron hopping, we employ the MEHAM model, which is able to describe the charge-transfer processes in both the adiabatic and non-adiabatic limits and is well suited to describing small polaron hopping in TiO 2 . 6,43−45 There are many previous applications of the approach to model polaron  Here, we employ MEHAM with the common approximation that the hole hopping process can be described as a onedimensional process coupling to a single generalized phonon mode. We note that although we consider a single phonon mode, it includes displacements of all atoms and the wavefunction is not constrained in any way, so it is not the case that adjacent atoms are not involved (in general they all are). The reasonableness of the one-dimensional approximation is discussed in detail in the work by Alkauskas et al. 57 and references therein. By considering such hopping processes between all inequivalent sites, we can describe hopping via adjacent atoms by chaining together these individual transitions we discuss in Section 4.3.
For each of the TiO 2 phases, we consider the hopping of small polarons between different anion sites along various pathways as indicated in Figure 1. These pathways include all paths with inter-anion separations of less than 4 Å. With the pathways identified, the first step is to obtain the optimized geometries for the polaron localized on both the initial and final anion sites (denoted by the set of coordinates R i a and R i b , respectively, where i is an index labeling each atom). Atomic structures for hole polarons localized on different sites are obtained by introducing a precursor potential well for polaron trapping by displacing neighboring Ti atoms around a given O site by 0.3 Å followed by full geometry optimization. As shown recently in a separate study, this approach is highly efficient for obtaining ground state structures of polarons. 58 We then linearly interpolate between the initial and final geometries to define a one-dimensional reaction coordinate for polaron hopping The corresponding one-dimensional potential energy surface (PES) is then computed as E[R(t)]. One can calculate diabatic PESs (where the hole remains localized on one site or the other as t varies) by initializing calculations using the converged wavefunction from the adjacent point in the pathway (starting from t = 0 and increasing as well as starting at t = 1 and decreasing). The adiabatic PES (where the hole is delocalized across the sites) can be computed by performing the same calculation with random initialization of the wavefunction. From these PESs, we are able to compute or estimate the various parameters needed to evaluate the rate of polaron hopping using MEHAM theory, including the adiabatic activation energies, the tunneling matrix elements, and the effective one-dimensional phonon-mode frequency coupling to the charge-transfer process (details given in the following section). We note that it is important to have well controlled self-interaction errors to accurately calculate the PESs required. In the next section, we present results for various polaron hopping pathways in each of the bulk phases of titania as well detail the approach to determine the various parameters used to calculated the polaron migration rates. All structural images and spin density plots are visualized using the VESTA software. 59 Figure 1. In anatase, all anion sites are equivalent, with oxygen atoms coordinated to three titanium atoms in a trigonal planar geometry. The hole small polaron in anatase is predicted to localize primarily onto a single oxygen atom. In brookite, there are two inequivalent oxygen sites in the bulk structure. The more stable O 1 site and less stable O 2 (difference in energy 0.03 eV) are indicated by blue and green spheres, respectively in Figure 1b. Similar to anatase, both oxygen sites are coordinated to three titanium atoms in a trigonal planar geometry. Of the four inequivalent oxygen sites in the TiO 2 -B phase, small hole polarons are found to localize on only two of them. The more stable O 1 site and less stable O 2 (difference in energy 0.11 eV) are indicated by blue and green spheres, respectively, in Figure 1c. Site O 1 again adopts a threecoordinated trigonal planar geometry, while O 2 is coordinated to two titanium atoms in a linear configuration.
For all of the pathways identified in Figure 1, we perform linear interpolations to obtain the PESs for small polaron hopping (as described in Methods and eq 1). For each pathway, we calculate the total energy at intervals Δt = 0.1. To illustrate our approach, we show a concrete example in Figure 2 that corresponds to the hole hopping pathway (ii) in anatase (see Figure 1a). The spin density is also shown in the figure highlighting that at the transition state (t = 0.5), the hole polaron is symmetrically localized across the two sites, corresponding to an adiabatic solution. The adiabatic activation energy (E act ) can therefore be determined directly from this PES. Initially, we attempted to also compute the diabatic PES at the transition state (t = 0.5) by carefully initializing the calculation with a charge density corresponding to hole localization on one of the sites only. However, while we were able to obtain symmetry broken solutions, there was still a significant adiabatic character in the solution (i.e., partial delocalization over the two sites). As an alternative, we computed an estimate of the diabatic energy at the transition state by fitting a harmonic function to the points t < 0.4 where the spin density shows a diabatic-like solution (black curves in Figure 2a). The fitted diabatic PESs provide estimates for the diabatic energy E dia at t = 0.5, the reorganization energy λ at t = 1, and the effective phonon frequency (Ω) corresponding to this onedimensional electron transfer model where Q is the generalized configuration coordinate and M i are the atomic masses. The electronic coupling constant H ab is then estimated as the difference between the adiabatic and diabatic energies at the transition state. We note a more accurate approach would be to obtain the diabatic solutions at the transition state using constrained DFT and then diagonalizing to obtain the electronic coupling constant as we have done in previous studies. 46−48 However, the high computational cost of these calculations and large number of pathways make this prohibitive at the present time.
With the various parameters calculated or estimated above, we proceed to calculate the polaron migration rate (k et ) using the MEHAM theory as follows where the electronic transmission coefficient is P LZ is the Landau−Zener transition probability P 1 exp( 2 ) LZ = − − πγ (6) and γ is the adiabaticity parameter defined as  Table 1 shows the calculated parameters for pathway (ii) in anatase together with those of the other pathways in anatase and the other phases considered. The calculated adiabatic activation energy for pathway (ii) is 0.20 eV with H ab = 109 meV. The pathways for each phase in Figure 1 are labeled according to the corresponding activation energies from smallest to largest. Pathway (i) with the smallest activation energy (0.12 eV) and reorganization energy (0.93 eV) in anatase corresponds to hopping between neighboring oxygen sites separated by 2.80 Å. The calculated effective phonon frequencies (Ω) are also found to be consistent with typical infrared active modes measured experimentally in anatase (e.g., the E u (1) mode at 358 cm −161,62 ).
Turning to brookite, the situation becomes more complex since there are two inequivalent oxygen sites and more pathways to consider. In the cases of hopping between inequivalent oxygen sites, the PES for hole hopping becomes asymmetric. For example, Figure 3a shows the calculated PES for pathway (i) in brookite. The energy difference between the two end points (ΔG 0 ) must be taken into account when calculating the hole hopping transfer rate using MEHAM theory. 26,63 In the following, we consider hole hopping from the higher to the lower energy site. The activation energy for hole hopping in the For asymmetric paths, the average for λ and Ω is given along with the rate in both directions. See Figure 1 for definition of the pathways.
The Journal of Physical Chemistry C pubs.acs.org/JPCC Article reverse direction is then simply given by E act − ΔG 0 . The transition state for asymmetric pathways does not occur at t = 0.5; hence, we need to account for this in the determination of E act by evaluating the adiabatic energy at the transition state once the fits for the diabatic PESs have been performed. However, in all cases considered, the transition state is within 0.014 (brookite) and 0.023 (TiO 2 -B) of 0.5 and so the adiabatic energy at the mid-point is a very good approximation. In general, for asymmetric pathways, the diabatic PESs also have different reorganization energies and vibrational frequencies. For simplicity, here we simply take an average of these quantities to evaluate κ el , which should be a reasonable approximation since these values do not differ significantly. The brookite structure consists of alternating (001) oriented layers of O 1 and O 2 oxygen sites as shown in Figure 2b. Therefore, hole migration in the [100] and [010] directions involves hopping between equivalent oxygen sites, while migration in the [001] direction involves hopping between inequivalent sites. Table 1 shows the calculated hole hopping parameters for all pathways in brookite. The lowest energy pathway (i) for the hole migration in brookite is along the [001] direction between two different O sites with an activation energy of 0.08 eV (or 0.13 eV in the reverse direction). Direct hopping between the most stable O 1 sitespathways (iv) and (vi)are associated with higher activation energies of 0.15 and 0.24 eV. The calculated effective phonon frequencies are again in good agreement with experimental data. 64 The structure for the TiO 2 -B phase is much more complicated than either anatase or brookite since there are four inequivalent O sites; however, hole polarons localize on only two of them, with three-coordinated site O 1 more stable than the twocoordinated site O 2 by 0.11 eV (Figure 1c). For the other two sites, delocalized holes are more favorable than the localized holes for the value of α that satisfies Koopmans' condition. The calculated hole hopping parameters for all pathways in TiO 2 -B are shown in Table 1. The most favorable pathway (i) involves hopping from site O 2 to site O 1 with an activation energy of 0.04 eV (0.15 eV in the reverse direction). Direct hopping between O 1 sites is possible via only two pathways: pathway (vi) in the [010] direction with an activation energy of 0.39 eV and pathway (vii) with an activation energy of 0.41 eV. These activation energies are sufficiently high that long-range diffusion in TiO 2 -B is likely to proceed via the less stable O 2 sites (discussed further below). Finally, we note that all of the pathways in Table 1 are adiabatic (κ el ∼ 1) except for brookite pathway (viii) which as discussed below is in practice not accessible in any case.
4.2. Thermal Ionization of Hole Polarons in Anatase, Brookite, and TiO 2 -B. To compute the energy barriers to thermal ionization of polarons into the valence band, E b I , we compute the corresponding adiabatic PES by linearly interpolating between the polaron geometry and the bulk geometry, defining a corresponding configuration-coordinate, Q (where Q = 0 is the geometry of the ideal bulk crystal and Q = ΔQ is the equilibrium geometry of the polaron). The energy was evaluated between the range −ΔQ to 2ΔQ (to capture the local minima at Q = 0 and Q = ΔQ as well as the curvature on either side) with a spacing of 0.1ΔQ between points. From these PESs, we compute E b I (which is the barrier to cross from the minima at Q = ΔQ to the minima at Q = 0) and the barrier to self-trapping E b ST , which is the barrier for the reverse process (where one exists). Figure 4 shows the computed PES for anatase, where E b I = 0.24 eV and E b ST = 0.03 eV. The calculated values for the ionization and trapping barriers for the other phases and hole polaron sites are given in Table 2. For all phases of TiO 2 , the barrier to trapping is much lower than the barrier to ionization, and in the case of TiO 2 -B, there is no barrier at all for trapping. The ionization barrier is relatively large for TiO 2 -B, suggesting that the polarons are highly stable for these phases, while for brookite the ionization barriers are much lower.
The reported values for ionization and trapping have strong implications in p-type doping for semiconductors. For example, in TiO 2 -B, we can see that there is not only a relatively high   The Journal of Physical Chemistry C pubs.acs.org/JPCC Article barrier to ionization, but there is no barrier to trapping. There remains no barrier to trapping even when the spacing between points in configuration-coordinate space is reduced to 0.01ΔQ in the region of interest. This means that TiO 2 -B would be difficult to dope p-type. The opposite can be said for brookite as the ionization barrier is much lower, while the trapping barrier is higher, suggesting it could be accessible for p-type doping. 4.3. Long-Range Polaronic Diffusion in Anatase, Brookite, and TiO 2 -B. Comparison of the calculated barriers to thermal ionization of hole polarons into the valence band against the calculated barriers for hole hopping allows one to make immediate conclusions on the mechanism of hole transport in the three phases of TiO 2 considered. First, as noted above, the barrier to thermally ionize hole polarons in brookite is very small (0.12 eV) and almost all hole hopping pathways have larger activation energies. The barrier to selftrapping is also in comparison large (0.05 eV). This suggests that the most likely mechanism for long-range hole migration in brookite is thermal ionization into the valence band, followed by band-like conduction and subsequent retrapping. Hence, we predict brookite should be associated with the highest hole mobility of the three phases considered.
While the barrier to thermal ionization of hole polarons in anatase is still relatively small (0.24 eV), several pathways for hole hopping have smaller activation energies. Long-range diffusion of hole polarons in the [100] and [010] directions is possible by constructing diffusion paths consisting of the most favorable pathway (i), which has an activation energy of only 0.12 eV. An example of such a path is shown in Figure 5a. Longrange diffusion in the [001] direction requires combining elementary pathways of types (i) and (ii) (Figure 5b); hence, the activation barrier for the rate-limiting step in this case is 0.20 eV. The similarity of activation energies for ionization and hopping suggests that long-range hole migration in anatase may have a mixed character. At a very low temperature, the dominant process would be hole hopping, which would be highly anisotropic with diffusion in the [100] and [010] directions preferred. At higher temperatures, one would expect a mixture of hole hopping and thermal ionization/retrapping resulting in an intermediate hole mobility. We discuss this in more detail in the subsequent section.
The barrier to thermal ionization of hole polarons in TiO 2 -B is much larger (0.53 eV) and although activation energies for hole hopping are also large, they are all smaller than the barrier to ionization. As for anatase, one can again determine favorable pathways for long-range hole diffusion in each crystallographic direction and corresponding activation energies for the ratelimiting elementary steps. As shown in Figure 6a (Figure 6c). Since there is no barrier to hole self-trapping in TiO 2 -B and the barrier to thermal ionization is large, one expects hole polarons will readily form. As for anatase, the rate-determining activation barriers for long-range hole hopping are slightly smaller but comparable to the barrier for thermal ionization and so one expects migration to have a mixed character. However, the activation energies involved are much higher than for anatase and so TiO 2 -B is predicted to have the lowest hole mobility of the three phases.

DISCUSSION
Quantitative prediction of hole mobility would require a full kinetic Monte Carlo simulation; however, one can provide order of magnitude estimates using the most favorable diffusion pathways we have determined together with the Einstein relation (e.g., see ref 26). In anatase, at room temperature, the mobility associated with hole hopping between adjacent sites is of the order 10 −3 cm 2 /V s while the ratio of hopping events to ionization events is around 10 2 . Since the typical diffusion length of a band-like carrier could easily exceed that achieved by 10 2 hops between adjacent sites, this suggests that hole mobility at room temperature in anatase has a significantly mixed character (i.e., both ionization/retrapping and small polaron hopping). At nitrogen temperatures (77 K), hopping events are around 10 8 times more probable than ionization. Therefore, in this regime, mobility should be dominated by small polaron hopping and exhibit significant anisotropy. Since TiO 2 is usually intrinsically n-type, information on hole mobilities is scarce, particularly for the less common phases. The review by Bak et al. surveys some of the data available in the literature, including single-crystal and polycrystalline materials. 65 However, direct comparison to experimental data is challenging since it convolutes the intrinsic polaron mobility in the bulk material with effects of point and extended defects which are not easy to disentangle. Nevertheless, our estimates of mobility are comparable in magnitude to reported experimental mobilities and previous theoretical calculations. 26,65 The prediction that anatase should exhibit a mixed character for hole polaron mobility is consistent with time-resolved photoluminescence, photoconductance, and transient absorption spectroscopy studies of anatase single crystals, which attributed non-exponential decays (not seen in rutile) to the presence of multiple carrier trapping processes. 66 There is very little existing data for brookite or TiO 2 -B but our calculations suggest very distinct behavior that could be probed by electrical conductivity measurements. Hole mobility in brookite is predicted to be the highest of the three phases and involve thermal ionization from small polaron configurations followed by long-range band-like diffusion and retrapping. On the other hand, hole mobility in TiO 2 -B is associated with very The Journal of Physical Chemistry C pubs.acs.org/JPCC Article high activation energies (for both hopping and ionization). At room temperature, the small polaron hopping mobility is of the order 10 −7 cm 2 /V s with a ratio of hopping to ionization of around 10 2 . Therefore, as for anatase, we predict a mixed character of mobility at room temperature. At nitrogen temperatures, small polaron hopping again dominates, being around 10 10 times more probable than ionization, but unlike anatase, we predict fairly isotropic mobility. The results discussed above have implications for applications on TiO 2 in photocatalysis and photoelectrochemistry. Equilibrium anatase nanocrystals expose predominantly {101} surface facets. Since the most rapid diffusion pathways in anatase are along the [010] and [100] directions, we predict holes that are photogenerated in the bulk of the nanocrystal will readily migrate to the {101} surfaces where they may interact with adsorbates and facilitate reactions. 31 Diffusion to {001} facets on the other hand would be less facile, particularly at lower temperatures, where the mobility takes on a stronger hopping character. For TiO 2 -B, we expect much reduced hole mobility in general and migration of holes to surfaces of nanocrystals will be much less efficient than for the other phases. This would suggest that TiO 2 -B nanocrystals are less active toward photo-oxidation than the other phases.
We conclude this section with a discussion of the factors that may affect the accuracy of the predictions. One of the largest potential sources of inaccuracy is residual SIE that can lead to inaccurate results for the stability and localization of small polarons. Our approach is explicitly parameterized to minimize SIE, with linearity of the total energy with respect to fractional occupation satisfied to less than 0.05 eV. 14 We have also employed comparatively large supercells in our calculations, which should minimize undesirable finite size effects. The onedimensional model for electron transfer we employ in conjunction with MEHAM theory are more significant simplifications. More general approaches are available but come with increased computational costs. 67 While quantitative prediction of rates should not expected, this simpler approach should be sufficient for identifying the most favorable pathways for migration and order of magnitude estimates of electron hopping rates.

CONCLUSIONS
In this study, we have used hybrid DFT that has been explicitly parameterized to minimize self-interaction errors (ensuring accurate description of charge-carrier trapping) to model hole polaron migration, ionization, and trapping in the anatase, brookite, and TiO 2 -B phases of TiO 2 . We find that the hole polaron mobility increases in the order TiO 2 -B < anatase < brookite. There are also distinct differences in the character of hole polaron migration in each phase. In brookite, hole diffusion proceeds via thermal ionization of the small polaron, band-like diffusion of the hole, and subsequent retrapping. At room temperature, hole polaron diffusion in anatase and TiO 2 -B has a mixed character (combining short hops between neighboring oxygen sites with less frequent thermal ionization and retrapping). At nitrogen temperatures, the migration character in both these materials becomes polaron hopping-dominated, with anatase exhibiting significant anisotropy. The activation energies for polaron ionization and hopping in TiO 2 -B are very high (around 0.4−0.5 eV), suggesting very poor hole mobility. These predictions should serve as a guide for future experiments aiming to probe charge-carrier dynamics in TiO 2 , a very challenging problem due to the need to disentangle a number of effects due to point and extended defects. As well as having fundamental importance, we have discussed the implications for applications of TiO 2 in photocatalysis and photoelectrochemistry.