Performance of the G 0 W 0 Method in Predicting the Electronic Gap of TiO 2 Nanoparticles

Using a relativistic all-electron description and numerical atomic centered orbital basis set the performance of the G 0 W 0 method on the electronic band gap of (TiO 2 ) n nanoparticles (n=1-20) is investigated. Results are presented for G 0 W 0 on top of hybrid (PBE0 and a modified version with 12.5% of Fock exchange) functionals. The underestimation of the electronic band gap from Kohn-Sham orbital energies is corrected by the quasiparticle energies from G 0 W 0 method which are consistent with the variational ΔSCF approach. A clear correlation between both methods exists regardless of the hybrid functional employed. In addition, the vertical ionization potential and electron affinity from quasiparticle energies show a systematic correlation with the ΔSCF calculated values. On the other hand, the shape of the nanoparticles promotes some deviations on the electronic band gap. In conclusion, this study shows ( i ) a systematic correlation exists between band gaps, ionization potentials and electron affinities of TiO 2 nanoparticles as predicted from variational ΔSCF and G 0 W 0 methods, ( ii ) that the G 0 W 0 approach can be successfully used to study the electronic band gap of realistic size nanoparticles at an affordable computational cost with a ΔSCF accuracy giving results that are directly related with those from photoemission spectroscopy, ( iii ) the quasiparticle energies are explicitly required to shed light on the photocatalytic properties of TiO 2 and ( iv ) that G 0 W 0 approach emerges as an accurate method to investigate the photocatalytic properties of both nanoparticles and extended semiconductor.


Introduction
Nanoparticles and nanoclusters have been widely investigated due to the broad variety of applications in biomedical, optical and electronic fields. 1,2In general, nanoparticles show remarkably high catalytic performance compared with the bulk phase due to their large surface area and quantum confinement effect. 3Particularly, semiconductor nanoparticles display specific electronic properties such as the lowest excitation ¾often referred to as band gap as in the corresponding bulk material¾ that make them especially attractive for applications in photocatalysis. 4For instance, to use sunlight for photocatalytic water splitting, leading to an almost inexhaustible sustainable energy source, requires materials with band gaps in the UV-VIS range.
Not surprisingly, the size and morphology of the nanoparticles are important factors in determining their electronic structure and, hence, their potential photocatalytic activity. 5perimentally it is difficult to discern between the effect of size and shape on the electronic properties of nanoparticles.On the other hand, computational modeling provides an unbiased approach to analyze the influence of these factors on the electronic properties.Using appropriate models, it is possible to represent different morphologies for a given composition or to vary the composition for a given morphology.The effect of shape and size of nanoparticle in defining the corresponding band gap and concomitant photoactivity has been illustrated in recent work on bottom-up 6,7 and top-down 8 models of TiO2 nanoparticles including explicitly over one thousand atoms.In these works, the band gap of the nanoparticles has been studied using density functional theory (DFT) methods.In particular, the electronic (or fundamental, Egap) and the optical (Ogap) band gap 9 have been considered.The former can be measured from photoemission techniques whereas the latter is accessible through optical spectroscopy.The charged states either of cationic (free extra hole, h + ) or anionic (free extra electron, e -) nanoparticles are associated to Egap whereas the generation of an excited electron-hole (e --h + ) pair corresponds to Ogap (more details can be found in Ref. [6] and references therein).Note that, in absence of excitons, Ogap and Egap of a bulk solid coincide.[8] In the framework of DFT, the use of the Kohn-Sham (KS) one-electron energies to approach either Egap or Ogap is a common practice, even although values predicted from Generalized Gradient Approximation (GGA) exchange correlation potentials are too small, 10 even incorrectly describing antiferromagnetic insulators such as NiO as metals. 11,12In the case of nanoparticles, the band gap is calculated as the difference between the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO) KS energies (ΔEH-L) and it is also strongly underestimated (up to 50%) at the standard GGA level.The use of hybrid approaches leads to values closer to experiment 13 although the amount of Fock exchange or the choice of the parameters defining more sophisticated range-separated hybrid functionals appear to constitute serious issues. 14,15In fact, a 12.5% of Fock exchange has been found to reproduce the band gap of anatase and rutile polymorphs of TiO2 whereas a 25%, as in PBE0, nicely reproduces the experimental gap in ZnO bulk structure. 16Recently, Perdew et al. 17 provided a theoretical basis to the well-known empiric fact that hybrid functionals as PBE0 yield realistic generalized KS gaps for semiconductor materials.In the case of the Egap, a more reliable method which is directly applicable to nanoparticles, is to use a ΔSCF approach which implies taking differences of total energy of a given N electron system and that of the system with N+1 or N-1 electronic system.This difference is the total energy upon adding (or removing) an electron, and is therefore the electron affinity (or ionization energy).However, this approach is not directly applicable to extended solids since it results in charged unit cells requiring some to introduce some charge compensation technique.
An alternative to ΔSCF methods to calculate the Egap is provided by the many-body perturbation theory GW approach suggested long ago by Hedin 18 and only applicable in practice thanks to the availability of more powerful computational and methodological developments. 19The GW acronym stems from the way the self-energy is calculated, as it is given by the product of the Green function G and the screened Coulomb interaction W. GW allows one to calculate the quasiparticle (qp) energies and their values correspond directly to electron removal or addition energies, being ideally matching those calculated using ΔSCF method and requiring the calculation of the N electron system only.4][25][26][27] In general, the evaluation of the qp energies often leads to a good correlation with experimental values for the energies of the valence excitation and gaps in several semiconductors 21,[28][29][30][31] and also in molecules. 324][35][36][37] The present study is in the line of previous work using the G0W0 approach to interpret photoelectron spectroscopy experiments involving (TiO2)n clusters with up to 10 units, 38 to explore two crystalline phases of dye-sensitized TiO2 clusters, 39 and to study the properties of rutile TiO2 nanoclusters. 40In the present work, TiO2 nanoparticles containing up to 20 units are chosen as models to investigate the trends in Egap as predicted using ΔEH-L from KS orbital energies, ΔSCF and G0W0 methods.Using a relativistic, all-electron, description with numerical atomic centered orbital basis set and modified hybrid functional, we show that a systematic correlation exists between band gaps, ionization potentials and electron affinities of TiO2 nanoparticles as predicted from ΔSCF and G0W0 methods, which could open the way for an accurate description of the electronic properties of oxide nanoparticles of realistic size at an affordable computational cost.

Models and Computational Details
All calculations are performed using the first-principles electronic structure theory based on DFT explicitly including all electrons.A numerical atom-centered (NAO) orbital basis set is used, as implemented in the Fritz Haber Institute ab initio molecular simulations (FHI-aims) program package. 41A light grid and tier-1 basis set is used, which has a quality similar to a valence triple-ζ plus polarization Gaussian Type Orbitals (GTO) basis. 8Additionally, light grid and tier-2 basis set is also used to better assess the accuracy of tier-1 basis set (see Supporting Information).The convergence threshold for geometrical relaxation of all nanoparticles is set to 10 -4 eV Å -1 .The presence of a transition element like Ti requires the inclusion of relativistic effects to ensure a correct convergence during the relaxation steps, zero atomic order regular approximation (ZORA) 42,43 is hence used in the calculations.Different hybrid functionals including a fraction of non-local Fock exchange are employed.These are PBE0 (25% Fock) and a modification of PBE0 containing 12.5% Fock and hereafter referred to as PBEx. 44The reason for using PBEx is that it has been proven to properly reproduce the main features of the electronic structure of stoichiometric and reduced anatase and rutile. 39Moreover, in view of the strong dependence of the G0W0 method with respect to the initial density, 45 hybrid functionals should provide a better starting point. 46set of (TiO2)n nanoparticles with n running from 1 to 20 units (Figure 1) is selected from the recent works 6,47,48 with the structures determined by global optimization using interatomic potentials 6 and, in each case, refined at the PBE, PBEx and PBE0 level using the NAO basis set commented above.The electronic band gap of the (TiO2)n particles has been estimated from different approaches but using always the minimum energy structure (Figure 1) consistent with the exchange-correlation potential used.The first one consists simply in taking the difference between the HOMO and LUMO orbital energies and will be denoted as ΔEH-L.The second method consist in making use of total energy differences and involves the vertical ionization potential (IPv) and vertical electron affinity (EAv): Egap = IPv -EAv.This approach is usually referred to as ΔSCF, as it requires the variational self-consistent energy of the neutral, cation and anion.For a given DFT method, the ΔSCF approach provides the best possible results, as it implies variational energies.
The third approach explored and benchmarked with respect to the ΔSCF results is the G0W0 approach.This makes use of a perturbative expansion to include the many-body effects through the many-body self-energy Σ.In the GW approximation 49 the self-energy is calculated as: ∑ (,  % , ) = ) *+ ∫ ′(,  % , ′)(,  % ,  % + ) 23 (1)   where (,  % , ′) is the one-particle Green´s function and (,  % , ) is the screened Coulomb interaction.The GW self-energy can be used to perturbatively correct the KS eigenvalues by means of the linearized quasiparticle equation: In the G0W0 or one-shot GW, the self-energy is calculated only once, whereas a more rigorous approach would require a fully self-consistent evaluation of Σ.Since the  ) 67 quasiparticle (qp) energies in Eq. ( 2) are evaluated perturbatively on top of a preliminary single-particle calculation, the G0W0 approach strongly depends on the starting point; 19,40,41 an issue that has been studied in detail in recent articles 50 with the general conclusions that hybrid functional provide better results.
In the following we refer to G0W0 based on PBEx and PBE0 as G0W0@PBEx and G0W0@PBE0, respectively.For a detailed account of the all-electron implementation of GW in FHI-aims we refer the interested reader to the original article of Ren et al. 51 It must be mentioned that performing GW calculations in an all-electron code such as FHI-aims has the advantage that possible pseudopotential errors are avoided, and eventually it can be used in periodic systems where the use of the variational ΔSCF approach requires additional approximations.As discussed extensively in the literature, the pseudopotential derived errors do not affect ground-state DFT energy but may become significant in GW calculations if there is significant spatial overlap between core and valence wave-functions. 52This issue was recently solved by Maggio et al. 53 reporting a careful comparison of GW qp energies obtained using Gaussian type orbitals and plane waves.In addition, the compact and inherently local nature of the NAO basis functions leads to a more rapid convergence with the number of basis functions. 54

Results and Discussion
We first focus on the trend of the Egap as a function of the number of TiO2 units as predicted by the ΔEH-L, ΔSCF and G0W0 approaches.The vertical ionization potentials and electron affinities obtained from ΔSCF and G0W0 methods are also discussed and compared with other theoretical and experimental works.Finally, the effect of the nanoparticle shape will be also investigated for the case of larger TiO2 nanoparticles, in this work containing 18 to 20 units.It must be noted that light grid and tier-2 basis set is used to evaluate the basis set convergence at G0W0@PBE0 level in (TiO2)n nanoparticles with n= 1-11 units, obtaining similar results than those obtained using light grid and tier-1 basis set (further details are given in Supporting Information).

Electronic band gap.
Figure 2 shows the Egap calculated following G0W0, ΔSCF and ΔEH-L procedures and using the PBEx and PBE0 hybrid functionals.Overall, the trends are systematically consistent for each method, with values from PBE0 being larger than those obtained from PBEx, although the effect is quite large for ΔEH-L (Figure 2c) and much less for G0W0 and ΔSCF, already indicating the higher accuracy of the latter approaches.More in detail, relative to G0W0 and ΔSCF, Egap associated to ΔEH-L is systematically underestimated and particularly sensitive to the exact exchange energy contribution associated to the hybrid functional (see Figure 2c).The percent of exact exchange energy contribution induces a significant increase of the band gap according to Egap(PBE0)/Egap(PBEx) ≈ 0.6 ratio.Note that, for a given functional, much larger particles are required to reach the bulk value depending on the method of calculation of the DFT energy. 8The ΔSCF method is more reliable and, actually, for a given DFT method provides the best possible estimate since all energy values used to estimate ionization potential and the electron affinity are variationally obtained ¾ i.e. from separate well-converged self-consistent calculations.Therefore, for either PBE0 or PBEx, the ΔSCF values are taken as the appropriate approach for the DFT level.
Let us now focus on the performance of the many body perturbation theory (MBPT) based on G0W0 method.Here, the Egap values are fully consistent with those arising from ΔSCF calculation to the point that a systematic correlation between both methods exist (see Figure 3a).
This result opens the possibility to use G0W0 for larger systems like extended solids where the variational ΔSCF is not suitable.G0W0 approach emerges as the appropriate computational technique to describe the electronic properties of semiconductors.This is consistent with previous works using GW techniques to approach the properties of bulk TiO2. 55In the case of bulk TiO2, however, in Ref. 49 the accuracy of G0W0 was established by comparing to the experimental value only whereas here the comparison between G0W0 and ΔSCF for a significant number of cases provides a more solid basis.This is further evidenced by the linear fit shown in Figure 3.It is noted that the oscillations observed in the smallest nanoparticles (n= 1-3) are related to the nanoparticle structure; similar oscillations have been previously reported. 42Although the Egap goes down, it is far from bulk values (see Figure 2).This is not surprising, since TiO2 nanoparticles of ~20 nm diameter composed by more than 10000 units are required to reach the bulk behavior. 8omparing the electronic band gaps obtained by G0W0 to those obtained from ΔSCF and ΔEH-L leads again to systematic correlation, as displayed in Figure 3.The straight line is particularly noted for the case of Egap(G0W0) versus Egap(ΔSCF) in Figure 3a, where the Egap values calculated at PBEx and PBE0 levels are grouped in the same fitting, showing that the functional (exact exchange energy contribution) affects equally in both methods.As it has been discussed above, the effect of the functional is much more pronounced when Egap is calculated via ΔEH-L.This is the reason why two different linear fittings are observed in Figure 3b.Although the functional promotes the shifts in the straight line, the deviation between both approaches is systematically consistent for PBEx and PBE0, as evidenced by fitting parameters in Figure 3. Similar linear correlations have been previously reported for a wide range of materials. 56,57 to here, the present analysis provides compelling evidence that the electronic band gap calculated using G0W0 is consistent with the ΔSCF approach, thus overcoming the problems arising from the use of Kohn-Sham orbital energies which are largely affected by self-interaction errors.
Nevertheless, one can argue that the success is due to an error cancellation consequence of a systematic error in the calculation of vertical ionization potentials and electron affinities.Results in the next subsection show that this is not the case.

Vertical Ionization potential and electron affinity.
In MBPT it is shown that the energy of the GW eigenvalues, usually referred to as quasiparticles (qp) energies, corresponds directly to electron removal or addition energies. 19,58erefore, vertical ionization potential (IPv) and electron affinity (EAv) can be obtained directly from the qp energies as IPv = -qpHOMO and EAv = -qpLUMO, respectively.Table 1 compiles the IPv and EAv values for the set of explored TiO2 nanoparticles as obtained from G0W0 and ΔSCF methods.For comparison, the corresponding potentials for periodic TiO2 determined through combined experimental-theoretical approaches 59,60 are included.Regarding the experimental values for (TiO2)n nanoparticles, photoemission data are only available for negatively charged clusters 61 and, therefore, we are not able to directly compare with the results for neutral clusters.However, our results are compared to computational analyses reported previously by Chiodo et al. 62 for small (TiO2)n nanoparticles with n = 3-10 using all-electron calculations with relativistic effects in consistency with our results.These authors reported ranges of 7.8-9.5 eV and 2.5-3.2eV for IPv and EAv, respectively.Our results reported in Table 1 are fully consistent with those reported by Chiodo et al. 56 Additionally, EAv for clusters from 3 to 10 units were experimentally determined to be in the range 2.6-3.5 eV. 55o further analyze the trends in IPv and EAv values as predicted from ΔSCF and G0W0 calculations, Figure 4 presents a linear correlation between the two sets of calculated IPv and EAv values.The amount of exact exchange energy contribution in the PBE0 and PBEx hybrid functional affects the results in the same way which is consistent with results in Figure 3a.In general, G0W0 predicts IPv and EAv values systematically slightly below those calculated using ΔSCF (see Table 1).To end this subsection, note that the values of IPv/EAv are significantly higher/lower than the corresponding potentials for periodic TiO2 anatase as it was already reported. 54These differences are due to the quantum confinement effect in the TiO2 nanoparticles.
Influence of the shape of TiO2 Nanoparticles on band gap.
For a given (TiO2)n nanoparticle size, distinct isomers exists within a narrow range of energy per unit. 6Figure S1 shows the set of isomers investigated and the electronic band gap (ΔEH-L, ΔSCF and G0W0) are reported in Tables S3 and S4 at PBEx and PBE0 level, respectively.It is therefore interesting to analyze the performance of the G0W0 approach on determining the Egap of the different isomers.With this idea in mind, the electronic properties of two different isomers, obtained during the global optimization search, 6 are investigated for the three largest particles considered in the present work which are those with n = 18, 19 and 20 (Figure 5).For a given size (n), the energy differences for the different isomers are almost the same for the PBEx and PBE0 methods; the n=19 nanoparticle shows the largest energy difference followed by nanoparticles with n=20 and n=18.
The energy differences may seem large but correspond to ~0.2 eV/unit only.
The results reported in Table 2 clearly show the influence of the nanoparticle shape, in agreement with previous findings reported by Cho et al. 6 The present results show that the minimum energy structures with a different atomic structure show smaller band gaps than the most stable ones, with differences in Egap as large as 1 eV.It is worth pointing out that these observations are consistent with experiments reported by Ba-Abbad et al. 63 Depending on the shape and size of the synthesized nanoparticles, these authors observed Egap variations of 0.15 eV.Nevertheless, even if the shape has a marked influence in small nanoparticles, this is notably reduced when increasing the size of the particles.Therefore, small band gap deviations in large nanoparticles composed by hundreds of TiO2 units are expected independently of the shape.Finally, note that the difference between Egap(DSCF) and Egap(G0W0) for the higher energy isomers is similar to that corresponding to the global minima, as expected, indicating the performance of the G0W0 method in predicting Egap of TiO2 nanoparticles barely depends on the particle size.

Figure 1 .
Figure 1.Structure of global minima of the (TiO2)n nanoparticles, n = 1-20 used in this work.Blue and red spheres correspond to titanium and oxygen atoms, respectively.

Figure 2 .
Figure 2. Calculated electronic band gap of (TiO2)n nanoparticles using (a) G0W0, (b) ΔSCF and (c) ΔEH-L methods at PBEx (orange) and PBE0 (purple) levels.The dashed lines in (c) correspond to the experimental band gap of bulk anatase (3.2 eV) and rutile (3.03 eV) systems which are provided for comparison.All these values are compiled in Figures S1 and S2 at PBEx and PBE0, respectively.

Figure 3 .
Figure 3. Linear relationship between Egap(G0W0) and (a) Egap(ΔSCF) and (b) Egap(ΔEH-L).The smallest nanoparticles, n= 1, 2 and 3 units are not considered in the linear fitting due to their oscillation effects.Fitting equations are shown in the inset.For simplicity, the Egap term is removed from the axis legends and linear regressions.All these values are compiled in Figures S1 and S2 at PBEx and PBE0, respectively.

Figure 4 .-
Figure 4.-Linear relationship between (a) vertical ionization potential (IPv) and (b) vertical electron affinity (EAv) obtained from G0W0 and ΔSCF calculations.Results for the smallest nanoparticles, n= 1, 2 and 3 units are not considered in the linear fitting due to their oscillation effects.The linear regression is also displayed.All these values are compiled in Figures S1 and S2 at PBEx and PBE0, respectively.

Figure 5 .
Figure 5.-(TiO2)n nanoparticles with n = 18, 19 and 20 units.The nanoparticles shown at the top are the most stable energetically and they are compared with the nanoparticles below, which have different shape and higher energies.Positive values of ΔE are associated with energies above the global minima, which are taken as reference.

Table 1 .
53,54e of IPv and EAv using ΔSCF and G0W0 methods at PBEx and PBE0 level.Potentials for periodic TiO2 anatase are also included.53,54Allvalues are in eV.

Table 2 .
-Electronic band gap of the global minima and one higher energy isomer of the (TiO2)n nanoparticles with n=18-20 as predicted from ΔEH-L, ΔSCF and G0W0 calculations carried out with the PBEx and PBE0 hybrid functionals.Numbers correspond to the most stable structures located at the top in Figure5.Those located between parenthesis correspond to the bottom structures in Figure5.ΔE has the same meaning as in the caption of Figure5.All values are in eV.