Assessing GW Approaches for Predicting Core Level Binding Energies.

Here we present a systematic study on the performance of different GW approaches: G0W0, G0W0 with linearized quasiparticle equation (lin-G0W0), and quasiparticle self-consistent GW (qsGW), in predicting core level binding energies (CLBEs) on a series of representative molecules comparing to Kohn-Sham (KS) orbital energy-based results. KS orbital energies obtained using the PBE functional are 20-30 eV lower in energy than experimental values obtained from X-ray photoemission spectroscopy (XPS), showing that any Koopmans-like interpretation of KS core level orbitals fails dramatically. Results from qsGW lead to CLBEs that are closer to experimental values from XPS, yet too large. For the qsGW method, the mean absolute error is about 2 eV, an order of magnitude better than plain KS PBE orbital energies and quite close to predictions from ΔSCF calculations with the same functional, which are accurate within ∼1 eV. Smaller errors of ∼0.6 eV are found for qsGW CLBE shifts, again similar to those obtained using ΔSCF PBE. The computationally more affordable G0W0 approximation leads to results less accurate than qsGW, with an error of ∼9 eV for CLBEs and ∼0.9 eV for their shifts. Interestingly, starting G0W0 from PBE0 reduces this error to ∼4 eV with a slight improvement on the shifts as well (∼0.4 eV). The validity of the G0W0 results is however questionable since only linearized quasiparticle equation results can be obtained. The present results pave the way to estimate CLBEs in periodic systems where ΔSCF calculations are not straightforward although further improvement is clearly needed.


Graphic for TOC Introduction
A possible way to overcome the difficulties arising from the use of DFT methods to approach CLBEs in periodic systems is to make use of the GW formulation, earlier introduced by Hedin. 22 The GW method includes many body effects beyond the mean-field description of the electron-electron interaction in DFT via the so-called selfenergy, which is non-local and energy dependent and, in a sense, replaces the exchange correlation potential in DFT. 23,24 Reviewing the GW formalism here is beyond the scope of this letter and the interested reader is addressed to the excellent review of Aryasetiawan and Gunnarsson for detailed information 25 All available practical implementations of GW have a common feature, the need for large computational resources. This has led to different levels of approximations resulting in various flavors denoted as G 0 W 0 , GW 0 , eigenvalue self-consistent GW, or quasi-particle self-consistent GW. 26 All these approximations convert the KS particles into quasiparticles (qp) with a well-defined physical meaning. In the case of occupied states, the qp energies effectively represent ionization potentials whereas in the case of unoccupied states, they represent electron affinities. In the case of periodic insulators and semiconductors qp from GW calculations provide an accurate estimate of the fundamental band gap of these materials 27,28 and also of states associated to point defects such as F-centers in simple oxides. 29,30 Apart from the high computational cost of the GW calculations, one must also point out the dependence of the results with respect to the starting density, 31 when the results are not obtained from a fully self-consistent approach as is the case in G 0 W 0 , lowest level of this theoretical framework. In particular, it has been suggested that hybrid functionals should provide a better starting point. 32 The GW method has also been applied to molecular systems and small clusters. [33][34][35][36][37][38][39][40][41][42] In particular, the G 0 W 0 level has been applied to 100 molecules (the GW100 database) 43 with excellent performance for the vertical ionization potential; the same database has more recently been used to benchmark different implementations and levels of self-consistency of the theory. 44,45 The rather good success of GW in predicting the ionization potential of molecular systems strongly suggest that it may as well provide an estimate of the CLBEs in molecules and solids where, as above commented, ΔSCF calculations are cumbersome. The calculation of CLBEs obviously requires the presence of actual core levels in the calculation. Note, however, that many solid-state codes use pseudopotentials with a frozen core approximation. These make the direct calculation of CLBEs impossible although different approximations have been proposed to predict ΔCLBEs. 14 Moreover, one must point out that GW calculations for core levels may be very tricky since describing the self-energy at deep energies requires a full frequency method and solving the quasiparticle equation can be very complicated as discussed later on. The goal of the present work is precisely to investigate the performance of GW on predicting CLBEs using a series of simple molecules, where experimental data and ΔSCF results for several exchange-correlation functionals are available, as a convenient benchmark.

Computational details.
The ground-state calculations are based on density functional theory within the Perdew-Burke-Ernzerhof (PBE) 46 generalized-gradient approximation (GGA) for the exchange-correlation functional. The molecular data set has been taken from previous work. 17 The molecular structures have been fully relaxed using the PBE functional and a tight Tier 2 numerical atom-centered orbitals (NAO) basis set [47][48][49] using the FHI-aims code. 47 The quasi-particle calculations are performed with the Turbomole package using the def2-TZVP and def2-QZVP Gaussian basis sets 50 and extrapolated to the complete basis-set limit. 51 In addition to electron correlation, relativistic effects also play a role in determining the final value of calculated CLBEs.2 ,12,17 Within the GW method, these can be introduced through various formalisms. 52 , 53 Nevertheless, for sake of comparison, relativistic effects have been omitted. We note, however, that the relativistic contribution to the CLBEs is essentially atomic in nature and it increases with the atomic number. For the core levels studied in the present work contribution of relativistic effects to the CLBE vary from 0.13 eV for C to 0.75 eV for F, 17 increasing along the C−F series. Different levels of GW calculations have been carried out including the so-called one shot G 0 W 0 approach, where both the Green's function G and the screened Coulomb potential W are obtained from an initial electron density arising from a given density functional approach, and the quasi particle self-consistent GW (qsGW), where both G and W are iterated until self consistency is achieved. Clearly, qsGW results do not depend on the starting exchange-correlation functional and this has been numerically verified here. Consequently, qsGW results are taken here as benchmark. 54,55 Note also that two different implementations of the simpler G 0 W 0 level are used. Nevertheless, both implementations of G 0 W 0 and qsGW approaches use the full analytic expression from the reducible response function for the self-energy and hence the results are obtained from a full frequency method. The difference in the two G 0 W 0 approaches used lies only in how the quasi particle equation is solved. In one case, the quasi particle equation is iteratively solved, in the other, denoted as lin-G 0 W 0 , the quasi particle equation is linearized. For the HOMO level the linearized and solved version often give very similar results. For core levels this is not the case. Due to the much larger corrections for core levels, in all molecules studied here the G 0 W 0 @PBE self-energy has poles in the region where the quasi particle equation (

Results and discussion
The set of molecules used to assess the reliability of the different methods considered here consists of CH 4 , CF 4 , CO 2 , HCN, H 3 COCH 3 , H 3 COH, H 2 O, NH 3 , pyridine, and pyrrole. The calculations are carried out for the C, N, O, and F 1s core levels for which experimental data is available. For reproducibility, all molecular geometries are reported in Table S1 of the supporting information. In addition, Table S2 reports the PBE total energy at the optimized geometry corresponding to the tight Tier 2 NAO basis set. To provide insight into the quality of this basis set Table S2 also reports the PBE total energy, at the same geometry, obtained with Gaussian type orbital (GTO) basis sets of aug-cc-pCVTZ and aug-cc-pCV5Z level of quality, 58,59 computed using Gaussian09. 60 Inspection of Table S2 clearly shows that the tight Tier 2 NAO basis set quality is even higher than that of aug-cc-pCVTZ and only slightly below aug-cc-pCV5Z. Table 1 collects the CLBEs for the whole set of molecules as obtained from the different levels of theory together with experimental values taken from the literature. 61,62 Nevertheless, one must advert that with the current implementation, fully analytic frequency treatment and with only a smp parallelism in TURBOMOLE, the calculation for pyridine at the def2-QZVP is not feasible and the corresponding extrapolated result could not be obtained; here the def2-TZVP value is listed instead. A quick inspection of Table 1 shows that, on average, going from KS to lin-G 0 W 0 to qsGW, the calculated results are progressively approaching experiment with best results attained for the qsGW method. A deeper insight into the relative errors is provided in Table 2 showing that the error for each core level is different with a mean absolute error (MAE) larger than 25 eV for the KS-ε s estimates, about 3-9 eV for G 0 W 0 , and slightly less than 2 eV at the highest qsGW level. These results illustrate that qsGW significantly improves the agreement with experiment as compared to the KS eigenvalues, a result which is in line with previous findings concerning the first ionization potential (HOMO) in molecules. 43,63 Nevertheless, the accuracy reached by qsGW may not be sufficient to properly interpret XPS experiments. The present results illustrate the difficulties of the GW methods on describing CLBEs. We stress that describing the self-energy at deep energies requires a full frequency method, which makes the calculations quite costly.
Also, for CLBEs, G 0 W 0 starting from the PBE functional faces fundamental problems since the self-energy has poles in the region where the quasiparticle equation is solved giving rise to a multitude of possible solutions. Besides the obvious numerical issues this grossly complicates the physical interpretation. Using the linearized quasiparticle equation solves at least the numerical part of the problem, solutions become numerically stable and unique but accuracy remains poor. In practice, it is more suitable to predict the quasiparticle energies of valence states. 30 Using a hybrid functional as starting point may remediate this deficiency as shown in valence states comparing to experiment 31 and also to scGW calculations. 64 This is because the core levels are already deeper causing the poles of the self-energy to be deeper in energy as well. Alternatively one can start from the HF density but here the lin-G 0 W 0 leads to values, which, not surprisingly, are overestimating. The advantage is that in this case the QPE can be solved, the original poles are further away and come up in energy during the selfconsistency process. The average over overestimation is however roughly 6 eV (see Table S4); a significant improvement over HF-ε s again with the correct physical ingredients where final state effects are introduced by the many body terms of the GW approach but a deterioration with respect to G 0 W 0 @PBE0. In the case of the qsGW, the poles move deeper as well and, in addition, the contribution of the off-diagonal terms of the self-energy matrix elements make it less sharply peaked.
For CLBEs calculated using the PBE functional and employing a ΔSCF approach the MAE is of ~1 eV, this is larger than the 0.3-0.4 eV corresponding values for calculations at the Hartree-Fock level or using the meta-GGA TPSS functional. 17 In this sense, the MAE of ~2 eV for the qsGW results is really remarkable and puts this method as a good choice to estimate CLBEs, especially in periodic solids. Here one must point out that the qsGW result does not depend on the starting point, yet the fact that the accuracy reached is lower than ΔSCF with PBE0 or TPSS functionals 17 indicates that higher order terms are likely to be needed in the expansion of the screened potential W and the self-energy.
Finally, it is also important to note that, in most practical cases, one is not interested in the absolute CLBEs but in their shifts with respect to a given reference (ΔCLBE). Taking CH 4 , H 2 O, and NH 3 as references for the C(1s), O(1s), and N(1s), the values for the ΔCLBEs arising from KS-ε s , lin-G 0 W 0 @PBE, lin-G 0 W 0 @PBE0, and qsGW are reported in Table 3 and the statistic analysis of errors in Table 4. The later shows that MAE values of 0.7, 0.9, 0.4 and 0.6 eV are found for KS-ε s , lin-G 0 W 0 @PBE, lin-G 0 W 0 @PBE0, and qsGW respectively. Note, however, that the experimental ΔCLBE values seldom exceed 10 eV, which makes the relative errors still rather large. Another important feature emerges when comparing calculated and experimental ΔCLBE values. Figure 1 shows that while the trends are well reproduced, there is a significant dispersion with the best values corresponding to lin-G 0 W 0 @PBE0, and qsGW respectively.

Conclusions
A study has been carried out to assess the performance of different GW approaches in predicting core level binding energies (CLBEs). Results presented for a series of representative molecules allows us to establish some first conclusions.
First, the present calculations agree with previous studies 9, 10 showing that the Kohn-Sham orbital energies cannot be taken as a measure of the CLBEs and they do not represent an estimate of initial state effects. The Kohn-Sham orbital energies are always smaller than experimental CLBEs with differences with respect to experiment in the 20-30 eV range. Interpreting them as approximate initial states would hence imply unphysical positive relaxation energies. In the framework of DFT, initial state effects can be recovered by computing the molecular system with a core hole with the fixed electron density of neutral molecule as discussed in previous works 10 and will not be further commented here. The interpretation of the KS-orbitals as approximate Dyson orbitals hence seems more physical especially for core levels.
The most accurate results, as compared to experiment, are obtained at the quasiparticle self-consistent GW level (qsGW). Nevertheless, to achieve good agreement with experiment, extrapolation to the complete basis set limit is needed, especially for the absolute CLBE although it may be less problematic for the shifts. In any case, the rather good agreement with experiment evidences that GW is able to largely correct the failure of the KS and HF orbital energies covering a large part of relaxation energy, thus going in the direction of the right-answer for the right-reason.
We find that the qsGW calculated CLBEs are always larger than the experimental values, as expected due to underscreening. The mean absolute error with respect to experiment is reduced to about 2 eV; one order of magnitude smaller than for KS values and close to predictions from variational ΔSCF calculations which, for the PBE functional, are in the 1 eV range. Here it is important to point out that an absolute error of 2 eV on a quantity in the 300-500 eV range implies a percent error of 0.6% only, which is quite remarkable, and, actually, very similar to the percent error of 0.5% achieved when considering the qsGW prediction of the first ionization potential. 43 a somehow larger error (6.5%) is obtained in the case of lin-G 0 W 0 @HF. Note also that adding contribution from relativistic effects will further decrease the ΔSCF error to less than 1 eV. Yet, the error bar of XPS, especially when synchrotron radiation is used, can be as small as 0.1 eV, meaning that further developments are needed. Smaller errors of roughly 0.6 eV are found for the qsGW calculated core level binding energy shifts (ΔCLBE) although these are surprisingly larger than those arising from Kohn-Sham orbital energies. This is an indication that qsGW effects on same core level in different molecules are not equally taken into account.
The simpler G 0 W 0 and computationally more affordable level of the theory also leads to a significant overall improvement with respect to KS orbital energies. However, starting from the PBE density, the calculated CLBEs are still smaller than experiment evidencing limitations inherent to this level of approximation as commented in the previous section. The absolute mean errors (8.98 and 3.76 eV for the PBE and PBE0 starting points respectively) are significantly larger than the corresponding value for qsGW but smaller than the one arising from the direct use of KS orbital energies. A paired T-test indeed proves that the differences are significant.
To summarize, CLBEs derived from GW approaches represent a considerable improvement with respect to KS energies although the accuracy reached with the present implementations, even for the qsGW, is still lower than the one obtained from ΔSCF calculations. The simpler G 0 W 0 method applied on top of the PBE density also lead to results improved with respect to KS predictions but with too large errors. These errors are considerably reduced when starting the G 0 W 0 calculations from a density obtained from a hybrid functional. Similar considerations apply to the core level binding energy shifts with rather satisfactory results for the lin-G 0 W 0 @PBE0 and qsGW methods. For practical applications in computational materials science, G 0 W 0 on top of single point PBE0 density at the PBE optimized structure may provide a practical approach since this will be computationally less demanding than going to the qsGW level and lead to similar accuracy.

Supporting information
The Supporting Information is available free of charge on the ACS Publications website at DOI: xxxxxxxxxxxxxx Table S1. Cartesian coordinates of the molecules studied in the present work as optimized FHI-AIMS tight Tier 2 PBE. Table S2. Total PBE energy (eV) of the molecules studied in the present work with geometries as in Table S1 and corresponding to numerical atomic orbital (NAO) tight Tier 2 basis set using the aims code and using the FHI-AIMS code and to the aug-cc-pCVTZ and aug-cc-pCV5Z gaussian type orbitals (GTO) basis sets as obtained with the Gaussian09 code. Table S3. Total PBE energy (eV) of the molecules studied in the present work with geometries as in Table S1 as obtained with the def2-VPTZ GTO basis and different codes. Table S4. CLBEs as obtained from lin-G 0 W 0 @HF.    Table 4. Mean absolute error (MAE) for the ΔCLBEs reported in Table 3 and obtained from KS, lin-G 0 W 0 @PBE, lin-G 0 W 0 @PBE0, and quasi particle self-consistent (qsGW) methods.  1.-Calculated KS, lin-G 0 W 0 @PBE, lin-G 0 W 0 @PBE0, and quasi particle selfconsistent (qsGW) versus experimental core level binding energy shifts (ΔCLBE). The dashed line corresponds to a perfect agreement.