Oxygen Vacancies in Strontium Titanate: a DFT+DMFT study

We address the long-standing question of the nature of oxygen vacancies in strontium titanate, using a combination of density functional theory and dynamical mean-field theory (DFT+DMFT) to investigate in particular the effect of vacancy-site correlations on the electronic properties. Our approach uses a minimal low-energy electronic subspace including the Ti-$t_{2g}$ orbitals plus an additional vacancy-centered Wannier function, and provides an intuitive and physically transparent framework to study the effect of the local electron-electron interactions on the excess charge introduced by the oxygen vacancies. We estimate the strength of the screened interaction parameters using the constrained random phase approximation and find a sizeable Hubbard $U$ parameter for the vacancy orbital. Our main finding, which reconciles previous experimental and computational results, is that the ground state is either a state with double occupation of the localized defect state or a state with a singly-occupied vacancy and one electron transferred to the conduction band. The balance between these two competing states is determined by the strength of the interaction both on the vacancy and the Ti sites, and on the Ti-Ti distance across the vacancy. Finally, we contrast the case of vacancy doping in SrTiO$_3$ with doping via La substitution, and show that the latter is well described by a simple rigid-band picture.


I. INTRODUCTION
Strontium titanate, SrTiO 3 , is a perovskite-structure oxide with the ideal cubic P m3m structure at room temperature, and a band-insulating electronic structure due to the formal 3d 0 configuration of the Ti 4+ cations. In spite of its apparent simplicity, SrTiO 3 shows a wealth of interesting and sometimes technologically relevant properties, such as tunability of its high dielectric constant, 1,2 quantum paraelectricity, 3,4 and even superconductivity. 5 Although these properties have been known for many years, there remain many open questions 6 . For example, superconductivity occurs at unusually low doping levels 7 , indicating an exotic mechanism possibly related to its quantum paraelectricity. 8 At the same time, reports of a two-dimensional electron gas 9 and emergent magnetism 10 at surfaces and interfaces have rekindled interest in SrTiO 3 thin films for oxide electronics.
All of these phenomena require the existence of electronic charge carriers, which are usually introduced through oxygen vacancies (O V ), 11 substitution of Sr 2+ by a trivalent ion such as La 3+ , or of Ti 4+ by a pentavalent ion such as Nb 5+ . Interestingly, the resulting properties can be quite sensitive to the specific type of doping. [12][13][14][15][16][17][18][19][20][21][22] Regarding oxygen-vacancy doping in SrTiO 3 , perhaps the most pressing open question is the nature of the introduced charge, with different experimental measurements leading to apparently contradictory conclusions. On one hand, it is known that even very low concentrations of oxygen vacancies cause metallicity, with transport measurements indicating an increase in carrier density with increasing oxygen vacancy concentration. [23][24][25] On the other hand, there are multiple reports of optical absorption signals within the band gap, [26][27][28][29] such as for example a red luminescence at 2.0 eV, attributed to localized electrons forming Ti 3+ polarons trapped at isolated oxygen vacancies. 30 Computational studies aiming to clarify the physics of the O V state also show a range of conflicting scenarios, largely due to their different treatments of exchange and correlation effects. (For a detailed summary see Ref. 31.) Standard density functional calculations using the local density approximation (LDA) or generalized gradient approximation (GGA) predict a delocalized defect state at the bottom of the Ti-t 2g conduction band, consistent with the observed metallicity. [32][33][34] Use of B3PW hybrid functionals or the LDA+ U method, however, lead to a doubly occupied in-gap state, 0.77 eV (B3PW) or 0.11 eV (LDA+ U, with U = 5 eV) below the minimum of the conduction band, but do not capture the reported metallic conduction (although this can be reconciled by considering the formation of polarons). 35 Using spin-polarized GGA+ U calculations with U = 5.0 eV, the two electrons released by the missing oxygen are found to distribute between a localized magnetic in-gap state and a delocalized state in the conduction band, consistent with both sets of reported experimental behaviors. 31 Finally, we mention a recent DFT+ U + V study, including selfconsistent on-site and inter-site electronic interactions, which gives a good description of both stoichiometric and oxygen-deficient STO, with the details of the band structure and the vacancy formation energies agreeing well with experiments. 36 The authors also present a systematic study of the effects of various parameters such as structure and cell size, the used exchange-correlation functional, and the treatment of structural relaxations and spin polarization on the resulting properties.
These previous studies clearly indicate that electron interaction effects are important in describing the behavior of oxygen vacancies in SrTiO 3 , 37 in spite of the band-insulating nature of the stoichiometric host material. However, while such DFT+ U (+V) calculations are able to appropriately treat strong interactions between electrons in the transition-metal d states (and in prin-ciple the oxygen p states), this is usually achieved by introducing an artificial symmetry breaking resulting in long-range magnetic order. Furthermore, current implementations of the DFT+ U formalism are typically based on projections on site-centered atomic orbitals, and so leave electrons at the vacancy site uncorrelated. The importance of explicitly considering interactions on the vacancy site was shown in model calculations using a minimal three-orbital model representing the vacancy, coupled to a bath representing the t 2g bulk bands, 38 but a full first-principles description is still lacking.
Here, we study the electronic properties of oxygendeficient SrTiO 3 using a combination of density functional theory (DFT) and dynamical mean-field theory (DMFT). The latter accounts for all dynamic correlation effects between electrons on the same "site" and thus introduces genuine many-body effects in the electronic structure obtained from DFT. Similar to our previous work on oxygen-deficient LaTiO 3 , 39 we focus on the most important bands at the bottom of the conduction band, which can be expressed in a basis of maximally localized Wannier functions (MLWFs) 40,41 with predominant Ti-t 2g character plus one additional Wannier function located at the vacancy site. This results in a physically intuitive framework to study the distribution of charge between the vacancy site and the Ti-t 2g conduction bands.
We estimate the strength of the screened electronelectron interaction, both for the Ti-t 2g orbitals and for the vacancy level, using the constrained random phase approximation (cRPA). [42][43][44][45] We obtain a sizeable U parameter for the vacancy level, confirming that the corresponding correlation effects should not be ignored. We then find that the occupation of the vacancy site is indeed controlled by the strength of the electron-electron repulsion both on the O V as well as on the Ti sites, balancing the system between a doubly-occupied localized vacancy state and a singly-occupied vacancy with the other electron doping the conduction band. The latter case can also be viewed as a site-selective Mott insulator, in which the singly occupied in-gap state corresponds to the lower Hubbard band on the vacancy site.

II. COMPUTATIONAL METHOD
The key point in our analysis is the explicit inclusion of an electronic orbital corresponding to the oxygen vacancy point defect into the low-energy Hamiltonian constructed from the DFT results, together with the usual Ti-t 2g bands. This allows us to explore in detail the effect of the electronic correlations on the vacancy site, and is the main novelty of the methodology presented here.
To obtain accurate geometries and initial bandstructures, we perform standard DFT calculations using the projector-augmented wave (PAW) method, as implemented in the "Vienna ab-initio simulation package" (VASP) 46  rations of the PAW potentials used are Sr(3s 2 3p 6 4s 2 ), Ti(3s 2 3p 6 4s 1 3d 3 ) and O(2s 2 2p 4 ). To allow the system to accommodate the defect with a reasonable O V concentration, we use a 20-atom unit cell, 19-atom for the O V -defective system SrTiO 2.75 , corresponding to an oxygen vacancy concentration of 8.3%. 49 Our calculations for La x Sr 1−x TiO 3 (Section III F) are performed using larger 40-atom supercells. For calculations of geometries and band structures (Sections III A and III B), well-converged results are obtained by sampling the Brillouin zone with a 8×8×8 Γ-centered k-mesh and using a plane-wave energy cutoff of 800 eV. For the more demanding cRPA calculations (Section III D), the k-space sampling is done with a 4 × 4 × 4 mesh, and an energy cutoff of 500 eV is used. In Section III E (Ti-O V -Ti distance dependence) these parameters are 8 × 8 × 8 and 750 eV, and in Section III F (comparison with La doping), these are 4 × 4 × 4 and 700 eV. Where lattice relaxation is employed, structural degrees of freedom are relaxed until forces fall below a 10 −4 eV/Å threshold, with the symmetry constrained to disable the rotations of the oxygen octahedra that occur in SrTiO 3 below ∼100K. All calculations are performed with spin polarization excluded.
The low-energy correlated subspace for the DMFT calculations is then constructed using a basis of maximally localized Wannier functions (MLWF) 40,50 , employing the Wannier90 code 41 . We use the TRIQS/DFTTools package [51][52][53] to implement the DMFT calculations, av-eraging over both spin channels to enforce a paramagnetic solution. An effective impurity problem is solved for each inequivalent Ti site plus the vacancy site using the TRIQS/CTHYB solver 54 , while the different impurity problems are coupled through the DMFT selfconsistency. The O V site is treated at the same level as the Ti sites, as introduced in Ref. 39 (more details about this implementation are given in Section III C). The DFT+DMFT calculations are performed without full charge self-consistency. The local interaction is modeled using the Hubbard-Kanamori parametrization with spin-flip and pair-hopping terms included. 55 Within the Hubbard-Kanamori parametrization, the strength of the electron-electron interaction is described by the intraorbital Hubbard parameter U and the Hund coupling parameter J. Note that these parameters are different from the average interaction parameters U and J typically used in DFT+ U calculations. The double counting correction is computed within the fully localized limit according to Held 56 , and all calculations are performed at room temperature, β = (k B T ) −1 = 40 eV −1 . We use a fixed value of J = 0.64 eV on the Ti sites, whereas the values for the Hubbard U , both on the vacancy and the Ti sites, are varied to analyze the effect on the electronic properties. Full frequency spectral functions, A(ω), are obtained from the local Green's functions in imaginary time, G(τ ), using the Maximum Entropy algorithm 57 . The spectral weight around the Fermi energy, A(0), is calculated from the impurity Green's function as We calculate the screened Coulomb interaction within the low-energy correlated subspace using the constrained Random Phase Approximation (cRPA). [42][43][44][45] . In the cRPA method, a partial polarization function, P r , is calculated by excluding all possible electronic transitions taking place within the correlated subspace. The bare Coulomb interaction, v, is then renormalized through screening by the higher energy degrees of freedom through this P r , yielding the frequency-dependent The local interaction parameters of the Kanamori Hamiltonian, U and J, are then obtained from the static limit by calculating matrix elements of W r (ω = 0), with the MLWFs used as local basis orbitals within the DMFT calculation, and taking appropriate averages over orbitals.

A. Geometry optimization
Our calculated optimized LDA lattice constant for SrTiO 3 constrained to cubic P m3m symmetry is 3.862 A. (For a picture of the 20-atom ( √ 2, 2, √ 2) unit see Fig. 1(a)). This value agrees well with literature LDA calculations, 34 and it is not too far from the experimental value of 3.900Å . 58 In subsequent relaxations after removing one oxygen ( Fig. 1(b)), we keep the lattice constants fixed to our calculated values for stoichiometric SrTiO 3 . 59 The resulting structure and Ti-Ti distances are shown in Fig. 1(c). The removal of one of the oxygen atoms lowers the symmetry of the formerly cubic crystal to tetragonal P4/mmm, in which the unique axis corresponds to the Ti-O V -Ti chain (y axis in Fig. 1(b)), and divides the Ti sites in two inequivalent types: one next to the vacancy, with five Ti-O bonds, and one farther from the vacancy, with a complete oxygen coordination octahedron. (Note that we do not include the low-temperature antiferrodistortive rotations of the oxygen octahedra; our preliminary tests suggest that their effect on the vacancy state is rather weak.) The Ti-Ti distances in the plane perpendicular to the tetragonal axis are slightly modified, with the vacancy pushing the two closest Ti apart, increasing their distance by 2.5% with respect to the distance in stoichiometric SrTiO 3 . The next Ti-O-Ti distance on the same axis contracts by the same amount, due to the volume constraint. On the other inequivalent Ti-O-Ti chain along the y axis there is also a 1.6% contraction (expansion) of the Ti-O-Ti distance at the same (different) y position as the closest Ti-O V -Ti.

B. DFT Bandstructures
The calculated bandstructures for stoichiometric SrTiO 3 and defective SrTiO 3−δ are shown and compared in Fig. 2. The low-energy region around the gap is composed of a valence band of mostly O-p character, and a conduction band whose bottom part has mainly Ti-t 2g contributions, with some weight coming from the O-p orbitals, and minimal overlap with the Ti-e g and Sr-s bands at ∼3 eV above the gap. Valence and conduction bands are separated by a gap of 1.8 eV, strongly underestimating the experimental value of 3.25 eV. 60 The bands between approximately 0 and 3 eV can be expressed in terms of MLWFs centered on the Ti atoms and showing a strong t 2g orbital character with additional admixtures of O-p on the surrounding ligands, see Fig.2 The removal of one oxygen atom changes mainly the lower part of the conduction band region of the bandstructure, with the appearance of an additional band that crosses the Fermi energy, E F , making the system metallic, and accompanied by a lifting of degeneracies at the special points of the Brillouin zone. This band could accommodate, in principle, the two electrons released by the vacancy. However, its overlap with the bottom of the conduction band edge (Fig. 2) results in a partial transfer of charge into the Ti-t 2g bands. Therefore, a complete description of the low-energy behaviour of the system must include both the Ti-t 2g bands and this O V -induced band, and we proceed by incorporating it into the minimal ba- sis of Ti t 2g -like MLWFs in the DMFT treatment. We thus construct 12 + 1 MLWFs, using initial t 2g projections on the Ti sites plus one s-like projection centered around the vacancy site.
The relationship of this new band to the O V is clearly seen from the MLWF centered on the position of the missing oxygen: its real-space representation shows an approximately spherical orbital centered around the vacancy with tails reaching to the neighbouring ions ( Fig. 2(e)). Moreover, if one extracts a single MLWF centered on the vacancy site, the resulting Wannier band closely follows the new O V -induced Bloch band (orange thick line in Fig. 2(b)). As one can see from the MLWFprojected density of states (Fig. 2(c)), the O V -centered MLWF accounts for most of the weight of this additional band. The occupations of the Wannier states, as calculated in DFT, are 0.21, 0.28, and 1.02 for the two inequivalent Ti sites, next to and farther from the O V , and the O V site itself, respectively, in units of the electron charge.

C. DMFT results
Next, we perform DMFT calculations to investigate how an explicit local Hubbard-like interaction between the electrons affects the spectral properties and the charge distribution between the vacancy state and the Ti-t 2g bands. The input for the DMFT calculations is constructed from the tight-binding-like Hamiltonian expressed in the basis of the MLWFs, plus the Coulomb matrix that models the electron-electron interaction. As stated in Sec. II, we simplify the latter for the t 2g states by using the Kanamori form, for which only two sitedependent parameters have to be specified: U , the onsite intra-orbital Hubbard repulsion, and J, the Hund's coupling. For the O V site with only one orbital, there is only one parameter, U , describing the corresponding intra-orbital Coulomb repulsion. This procedure allows us to independently vary the strength of the local interaction on the different types of sites, including the vacancy, and hence to determine the evolution of the system when U (O V ) changes independently of U (Ti). We note that for other early transition metal perovskites, such as, e.g. LaTiO 3 , values of U (Ti) between 4-5 eV have often led to good agreement with experimental observations when using a minimal t 2g orbital subspace; 61,62 we therefore choose the range from 2-8 eV. On the other hand, given the lack of chemical intuition for choosing U (O V ), we have considered U (O V )= 0 and U (O V )=U (Ti) as reasonable limits. Fig. 3 shows our calculated DMFT site occupations, along with the site-resolved quasiparticle weight Z, and the corresponding spectral functions, for different choices of U on the Ti and the O V sites. The top row (Fig. 3(a)) corresponds to U (O V )= 0, that is treating the electrons on the vacancy site as "uncorrelated". In this case, the effect of increasing U (Ti) is to increase the O V occupation from < 1.5 for small U (Ti) to nearly 2, that is almost completely filled, for large U (Ti). For U (Ti) > 6 eV, a metal-insulator-transition (MIT) is observed, in which the Ti bands are totally depleted, and the system becomes a band insulator. This MIT is related to a shift of the O V (Ti) spectral weight to lower (higher) energies, until the corresponding overlap vanishes (see spectral function for U (Ti) = 8 eV). The quasiparticle weight Z is close to 1 for all sites, independent of U (Ti), indicating that the electrons remain uncorrelated, consistent with the transition from an uncorrelated metal to a band insulator. 63 Our DFT+DMFT result for higher U (Ti) and U (O V )= 0 is thus equivalent to the results obtained in previous studies within DFT+ U, 32 as discussed in Sec. I.
The remaining rows in Fig. 3 Increasing U (O V ) also strongly reduces the quasiparticle weight on the vacancy site, indicating strong local electronic correlations. Simultaneously, U (O V ) triggers a reduction of Z on the Ti sites, in particular on the site farther away from the vacancy. One can also observe a pronounced effect on the site-resolved spectral functions. In particular, increasing U (O V ) results in a clear gap opening in the spectral function on the vacancy site once it reaches half-filling. This indicates a site-selective Mott transition, that is a localization of one electron on the vacancy site while the doped Ti bands remain metallic, albeit with a quasiparticle renormalization on the farther Ti site of Z ≈ 0.7. Thus, this regime is characterized by a metallic conduction band doped with one electron per vacancy, and a split-off "in-gap" state con-taining one electron localized on the vacancy site. In the picture of the site-selective Mott insulator, this in-gap state corresponds to the lower Hubbard band of the vacancy site spectral function. Note that, in this regime of U (O V )=U (Ti) the filled region of the O V spectral function overlaps with the O-p bands, which start 2 eV below E F , but are excluded from the DMFT calculation.
Our results show that the specific treatment of electronic correlations on the O V site has a strong influence on whether the two electrons released by the missing oxygen are itinerant or localized. Several scenarios that have previously been discussed and reported for SrTiO 3−δ can be realized with particular choices of U (Ti) and U (O V ), from the band-insulating limit with localization of the vacancy charge in the band gap 64 (high U (Ti), zero U (O V )) to the case of a paramagnetic impurity, in which one electron remains trapped in the gap while the other is delocalized into the Ti d-bands (U (Ti)=U (O V )). 31 We note that the lack of charge self-consistency in our calculations may diminish the quantitative nature of these results; however, the qualitative picture should remain valid.

D. cRPA calculations
The results discussed in the previous section show that the localization/delocalization of the two electrons associated with the missing oxygen depends critically on the values of the interaction parameters, U (Ti) and U (O V ). In some cases, values for the interaction parameters can be estimated from previous experience or from comparison of certain calculated quantities with corresponding experimental measurements. In the present case, however, there is a lack of consensus among experimental studies (as described in Sec. I). Furthermore, due to the use of a non-standard orbital basis, in particular on the vacancy site, it is not a priori clear what values for U (Ti) and U (O V ) will provide the most realistic description of SrTiO 3−δ . Therefore, we now estimate the strength of the screened Coulomb repulsion corresponding to our orbital basis using the constrained random phase approximation (cRPA). [42][43][44][45] Within cRPA, the electronic degrees of freedom are divided into a "screening subspace" and a "target subspace", and excitations taking place exclusively within the target subspace are excluded from the screening (see Sec. II). Generally, the cRPA target subspace should be identical to the correlated subspace used for the DMFT calculation. However, there are also cases where a certain number of "uncorrelated" bands are included in the DMFT subspace, without considering a local Hubbardlike interaction for these bands. In such cases it might be appropriate to include these uncorrelated bands within the screening subspace in the cRPA calculation. 43,44 Thus, in order to gain additional insights, and to see how the calculated values depend on the specific subdivision of electronic degrees of freedom, we consider three different cases. In the first case, we construct only three t 2g -like Wannier functions for each Ti site. These orbitals then also define the cRPA target subspace, achieved through a projection of the Bloch states onto these Wannier orbitals. All other bands, including that related to the vacancy, act as the screening subspace. Following the notation established in Ref. 43 and used in several other studies, 45,65,66 we denote this case as (t 2g /t 2g ). Here, the first symbol denotes the orbitals/bands defining the cRPA target subspace, while the second symbol indicates the whole set of Wannier functions that has been constructed. The screened interaction parameters are then evaluated for the Wannier functions spanning the cRPA target subspace.
For the second scenario, (t 2g -O V /t 2g -O V ), a full Wannier representation of the low energy conduction bands in SrTiO 3−δ containing three t 2g -like Wannier functions per Ti plus an additional one located on the vacancy site is constructed, as described in the previous section. This whole set of Wannier functions is then used to define the cRPA target subspace. Thus, in this case the vacancy band is excluded from the screening and included in the target/correlated subspace.
The third, intermediate, case is denoted as (t 2g /t 2g -O V ). Here, the same Wannier functions as in the previous case are constructed, but only the subset of t 2g -like functions are used to define the cRPA target subspace. This essentially means that the effective interaction parameters calculated for the Ti-t 2g orbitals also include screening processes involving the vacancy-band, similar to the first case, but using exactly the same Wannier representation as in the second scenario. Fig. 4 shows the values of the screened intra-orbital interaction parameters U (Ti) and U (O V ), given by the corresponding averaged diagonal elements of the calculated U tensor for all three cases. The difference in U (Ti) between the (t 2g /t 2g ) and (t 2g /t 2g -O V ) cases is negligible. This similarity is expected, since the Ti-centered Wannier functions in the two sets have very similar spreads, due to the rather weak hybridization between the Tit 2g bands and the O V -band, and since the screening subspace is identical in the two cases. However, removing the vacancy band from the screening channel has a strong effect, as can be seen by comparing cases (t 2g /t 2g -O V ) and (t 2g -O V /t 2g -O V ), resulting in an increase of U (Ti) from 2.65 eV to 3.47 eV, and demonstrating the sensitivity of U on the specific screening channel. In addition, a sizeable interaction parameter of U (O V ) = 2.12 eV is obtained, which corresponds to approximately 60% of U (Ti). This clearly shows that interaction effects cannot be neglected for the vacancy orbital, and that, in addition, the presence of the vacancy band strongly influences the effective interaction parameters on the Ti sites. These findings thus further support our approach of including an explicit treatment of correlation effects on the vacancy orbital, reinforcing the crucial role played by both U (Ti) and U (O V ) already presented in Sec. III C.
The values of the averaged interaction parameters obtained for the case (t 2g -O V /t 2g -O V ), U (Ti) ≈ 3.5 eV and U (O V ) ≈ 2.1 eV, would correspond to DFT+DMFT results with a vacancy occupation slightly larger than 1 and a corresponding quasiparticle weight below 0.5, but still with a metallic spectral function on the vacancy site (Fig. 3). We note, however, that simply using the static (ω = 0) value of the screened interaction in a DFT+DMFT calculation with frequency-independent local interaction, might underestimate the corresponding interaction effects. Furthermore, recent work has shown that the random phase approximation can lead to an "overscreening" of the local interaction, 67 and thus an underestimation of U , in particular for strongly correlated systems. 68 The cRPA results should thus rather be viewed as providing a rough ballpark, or a lower bound, rather than definite values to be used in realistic DFT+DMFT calculations. 69 E. Influence of the Ti-OV-Ti distance As shown Fig. 1 and discussed in Sec. III A, structural relaxation using the standard LDA results in an outward relaxation of the two Ti atoms adjacent to the vacancy, and thus an elongation of the Ti-O V -Ti distance across the vacancy compared to the Ti-O-Ti distance in the ideal stoichiometric case. While this is consistent with previous work using LDA or GGA, 70,71 other calculations using hybrid functionals or LDA+U corrections have instead found a contraction of the Ti-O V -Ti distance, independent of other factors such as the size of the supercell. 72,73 As mentioned previously (see Sec. I), the latter calculations also find a tendency for charge localization, whereas LDA/GGA result in some charge spilling into the Ti bands. Therefore, we next manually tune the Ti-O V -Ti distance and calculate the resulting changes in electronic properties. Fig. 5 shows the total energy calculated using LDA and LDA+ U for 19-atom SrTiO 3−δ as a function of the Ti-O V -Ti distance (in units of the Ti-O-Ti distance in stoichiometric SrTiO 3 ). In these calculations, we only allow for a structural relaxation of those oxygen atoms that are situated next to the Ti in the same plane per-pendicular to the Ti-O V -Ti direction. Calculations are performed for three different choices of the static Hubbard U parameter. 74 The case with U = 0 corresponds to a standard LDA calculation and leads to an increase of the Ti-O V -Ti distance, consistent with the results presented in Sec. III A. Increasing U leads to a decrease of the relaxed Ti-O V -Ti distance, in accordance with the aforementioned trend in the literature, ultimately shortening this distance with respect to that in pure SrTiO 3 .
In order to further explore the relationship between the localization of the defect state and the Ti-O V -Ti distance, we perform LDA calculations for a SrTiO 3−δ structure in which we have systematically varied the Ti-O V -Ti distance. For each (fixed) Ti-O V -Ti distance, the oxygen ions are allowed to relax within the tetragonal symmetry constraint, while all other atoms and the lattice constant are kept fixed. Fig. 6 shows the evolution of the resulting DFT bandstructure as the Ti-O V -Ti distance is varied over a range that includes also the various equilibrium values obtained for different values of U within DFT+ U. We see that the Ti-O V -Ti distance affects primarily the position of the vacancy band, which is higher in energy, and therefore contains less electronic charge, at larger Ti-O V -Ti distances. We conclude that short Ti-O V -Ti distances tend to localize the charge around the defect, while longer distances tend to partially delocalize the charge.
Following the same procedure as in Sec. III C, we then perform a series of DMFT calculations for the different structures corresponding to different Ti-O V -Ti distances. We use U (Ti)= 4 eV and three different values for U (O V ) (0, U (Ti)/2, and U (Ti)). Our calculated site occupations are displayed in Fig. 7. We see that for all values of U (O V ), the occupation of the vacancy is decreased with increasing Ti-O V -Ti distance, while the Ti occupancies increase proportionally. The most pronounced change in occupation is observed for U (O V )=U (Ti) (Fig. 7, lower panel). In this case, the vacancy site occupation drops from close to 2 (for the shortest considered Ti-O V -Ti distance of 0.92) down to 1 (for distances equal or larger than 1.04). This case corresponds to the site-selective Mott-insulating state obtained previously in Sec. III C, in which the vacancy occupation is essentially locked to 1 due to the gap in the corresponding local spectral function. The strong dependence of the DMFT occupancies on U (O V ) for larger Ti-O V -Ti distances is of course consistent with the results already presented in Fig. 3, where the structure was relaxed within LDA, leading to a Ti-O V -Ti distance of 1.024 in these units.
These results indicate that the site-selective Mott insulating state, corresponding to a singly occupied vacancy site, is favorable for an elongated Ti-O V -Ti distance. We note that the lack of charge self-consistency does not allow us to extract reliable DFT+DMFT total energies as a function of the Ti-O V -Ti distance.

F. Comparison with LaSr substitution
As we mentioned in the Introduction, the question of how different doping sources affect the detailed electronic properties is an important one, with earlier studies suggesting that, while O V doping significantly alters the band structure as we have seen here, La Sr causes only Here, we test whether the DFT+DMFT treatment is consistent with these earlier findings by calculating the behavior when a small concentration of Sr ions is replaced by La. Fig. 8 shows the evolution of the DFT La x Sr 1−x TiO 3 bandstructure for three different values of x corresponding to La concentrations of less than or equal to half.
These results confirm that the main features of the SrTiO 3 bandstructure are robust against doping. The only, yet rather important, difference between the three cases is the increasing partial filling of the Ti-t 2g bands with increasing La concentration.
The DFT+DMFT spectral functions of the Ti-t 2g bands for these three La concentrations, shown in Fig. 9(a), confirm the metallicity of the system. Their respective integrals up to ω = 0 yield the electron occupation of the corresponding orbitals, and in our simulations all have values consistent with a homogeneous distribution of the excess charge of one electron per La ion into the Ti bands. Increasing the strength of the interaction parameter from U = 4 eV to U = 8 eV does not have a noticeable effect on A(ω). However, the corresponding quasiparticle weights, Z, (Fig. 9(b)) are reduced from the "uncorrelated" value of Z = 1, indicating some degree of electronic correlation in the corresponding bands. Z decreases, and becomes more sensitive to the value of U , as the concentration of La increases, consistent with the vicinity to the Mott insulating state at x = 1. Thus, our results confirm the validity of the rigid-band picture under La Sr substitution, while at the same time indicating moderate correlation effects, increasing with La concentration.

IV. SUMMARY AND CONCLUSIONS
We have presented a DFT+DMFT description of the long-discussed problem of oxygen vacancies in SrTiO 3 , focusing on a controlled and systematic treatment of electronic correlations on the vacancy site. This is achieved by using a minimal correlated subspace, which consists of the low energy Ti-t 2g orbitals plus an additional Wannier function located on the vacancy site. Our study reveals a strong influence of U (O V ) on the overall electronic structure of the system. Furthermore, our cRPA calculations support the importance of electronic correlations on the O V site, yielding a U (O V ) close to 60% of U (Ti).
In particular, we find a transition from a doublyoccupied (U (O V )=0) to a singly-occupied O V state (U (O V )=U (Ti)), accompanied by a charge transfer of one electron from the defect state into the conduction band. The latter state is equivalent to the scenario proposed by Lin and Demkov based on a minimal model of the vacancy levels, 7 and consistent with the spinpolarized GGA+U results of Hou and Terakura, 31 but without the need for an artifical spin-symmetry breaking. It can in principle reconcile apparently contradictory experimental observations of low temperature metallic conductivity in combination with localized ingap states. 24,25,28,29 Analysis of the effect of the distance between two Ti ions separated by an O V on the nature of the O V charge leads to two important conclusions. First, different treatments of electronic interactions at the DFT(+U) level lead to different Ti-O V -Ti distances. This in turn affects the relative O V and Ti-t 2g energy levels, with smaller Ti-O V -Ti distances favouring localization of the charge on the O V and larger distances promoting partial filling of the Ti orbitals with one of the electrons released by the missing oxygen. Second, this picture is largely maintained after a DFT+DMFT analysis, with the O V electron occupation also depending on U (O V ). Future DFT+DMFT calculations considering full charge selfconsistency and structural relaxations could provide further insight as to which Ti-O V -Ti distance and vacancy site occupation is indeed energetically preferred.
Lastly, we have compared O V -doped SrTiO 3 with the case of La Sr substitution, which represents another common way of electron doping the system. Our DFT+DMFT results, in line with previous DFT findings, show that La x Sr 1−x TiO 3 is much closer to the trivial doping case, with no significant change in the bandstructure other than an increasing uniform filling of the Ti-t 2g bands by the extra electrons donated by the La cations. This demonstrates that different routes to electron-doping in SrTiO 3 are certainly not equivalent, and can lead to rather different electronic structures.
Finally, we mention that the method that we have demonstrated for controlling and analyzing the strength of electronic correlation on a vacancy defect state within the DFT+DMFT formalism is applicable to other defects beyond the oxygen vacancies studied here. Our finding that the electronic nature of the vacancy state depends strongly on the strength of the local interaction, U (O V ), highlights the need for such a methodology. We hope that, in addition to contributing to the ongoing debate about the case of doped SrTiO 3 , our work will motivate similar studies on other materials in which correlated defects might play an important role.