|Home | About | Journals | Submit | Contact Us | Français|
Constrained-occupancy delta-self-consistent-field (ΔSCF) methods and many-body perturbation theories (MBPT) are two strategies for obtaining electronic excitations from first-principles. Using the two distinct approaches, we study the O 1s core excitations that have become increasingly important for characterizing transition-metal oxides and understanding strong electronic correlation. The ΔSCF approach, in its current single-particle form, systematically underestimates the pre-edge intensity for chosen oxides, despite its success in weakly correlated systems. By contrast, the Bethe-Salpeter equation within MBPT predicts much better lineshapes. This motivates one to reexamine the many-electron dynamics of X-ray excitations. We find that the single-particle ΔSCF approach can be rectified by explicitly calculating many-electron transition amplitudes, producing X-ray spectra in excellent agreement with experiments. This study paves the way to accurately predict X-ray near-edge spectral fingerprints for physics and materials science beyond the Bethe-Salpether equation.
X-ray absorption spectroscopy (XAS) is a powerful characterization technique in physics, chemistry, and materials science, owing to its element specificity and orbital selectivity. With the help of density-functional theory (DFT), the interpretation of XAS is greatly facilitated by simulating spectral fingerprints for hypothetical structures from first-principles. Satisfactory X-ray absorption spectra have been simulated across a wide range of systems from small molecules [1, 2] to condensed matter [3–9] and even complex interfaces .
In this letter, we introduce a novel first-principles theory for simulating X-ray absorption spectra extending the constrained-occupancy delta-self-consistent-field (ΔSCF) method [6, 7, 9]. The ΔSCF method assumes a fixed core-hole potential and employs only single-particle orbitals for finding transition amplitudes. It is generally thought this theory cannot capture certain many-electron effects as incorporated in many-body perturbation theories, such as the Bethe-Salpeter Equation (BSE) [11–14]. However, if the one-body transition amplitude in the ΔSCF method is recast into a many-body formalism, we show that the method can in fact capture many-electron effects in X-ray excitations. Furthermore, we outline how the excitation spectrum can be enriched by considering multiple electron-hole (e-h) excitations, permitting extensions beyond the BSE, on the level of the Mahan-Noziéres-De Dominicis (MND) theory[15, 16].
For our examples, we choose the O K edges (1s → np transitions) of transition-metal oxides (TMOs) to illustrate the utility of the many-body ΔSCF method. The study of the O K edges for TMOs is fueled by the quest for next-generation energy materials, for rechargeable battery cathodes [17–20], fuel cells [21, 22], water-splitting catalysts [23–25], and transparent conductive layers . Many of these materials are TMOs with complex chemical properties due to their d orbitals. Additionally, XAS has been employed to understand electron correlations inherent in TMOs, such as metal-insulator transitions [27–29], high-Tc superconductivity [30, 31], and emergent phenomena at perovskite interfaces . XAS can also serve as powerful guides to advance theories for correlated electron systems, including the DFT+U method , dynamical mean-field theory , and exact diagonalization approaches [35, 36]. In almost all of the aforementioned examples [17, 19–21, 23–28, 30–33, 35, 36], there are measurements at the O K edge. Despite its utility, very few studies have simulated this absorption edge for TMOs from first-principles.
Five TMOs are selected for benchmarking: the rutile TiO2, VO2, and CrO2 as well as the corundum α-Fe2O3 and the perovskite SrTiO3. They vary greatly in structure, band gap, or magnetism. The rutile VO2 (> 340K) and CrO2 are metallic, whereas TiO2 and SrTiO3 are insulating. CrO2 is ferromagnetic (FM) while Fe2O3 is antiferromagnetic (AFM). The O K edges from previous experiments [23, 27, 37–39] are shown in Fig. 1(a). These spectra are angularly averaged except for CrO2, where the polarization is perpendicular to the magnetization axis . The TM-3d-O-2p hybridization manifests as sharp double peaks around 530 eV, which result from the t2g-eg splitting in the octahedral field. The intensity ratio of the two peaks often serves as a diagnostic tool. Understanding gained through these simple TM compounds will find great utility in interpreting XAS for more complex contexts, e.g., liquids, interfaces, heterogenous materials, where other characterization techniques such as X-ray diffraction are less effective. For instance, the O K-edge peak intensities have been used to deduce oxygen redox chemistry in rechargeable battery cathodes , and yet the correlation between experimental spectra and first-principles modeling of structures remains to be established, in the absence of accurate predictions of the near-edge lineshape.
The X-ray absorbance σ(ω) is obtained from Fermi’s golden rule
where ε and R are the photon polarization and many-body position operator respectively. |Ψi is the initial (ground) state and |Ψf is any excited eigenstate of the system, while Ei and Ef are their energies. Predicting X-ray absorption requires: (a) an accurate prediction of the excitation energies Ef and (b) a reliable approximation to the amplitude Ψf|ε·R|Ψi. We first contrast the single-particle ΔSCF core-hole approach and the BSE using the five TMOs, and identify a systematic failure in the current ΔSCF method. We argue that the failure is due to the omission of many-electron response to the core hole, which motivates us to develop an advanced ΔSCF approach by using many-electron wavefunctions instead of single-particle orbitals.
In the ΔSCF core-hole approach the core-excited atom is treated as a single impurity with one electron removed from the excited core level. Depending on whether or not the X-ray photo-electron is included, it is termed as an excited-electron and core-hole (XCH) or full core-hole (FCH) calculation [6, 7, 9]. One places the core-excited atom in a sufficiently large supercell and then finds the equilibrated electron density using constrained-occupancy DFT. This equilibrated state, with the presence of a core hole, is referred to as the final state, and that of the pristine system as the initial state. The working approximation is to use single Kohn-Sham orbitals for finding the transition amplitude
where f’s are the unoccupied orbitals in the final state (with tilde) and ψh is the core orbital in the initial state. S is a many-body overlap that reflects the excited-state response of the remaining electrons other than the photoexcited one, and is normally treated as a constant in the sudden approximation . The excitation spectrum Ef − Ei is approximated based on differences in DFT orbital energies, f − εh, with some further realignment reflecting the local atomic context and the experimental energy scale .
To account for strong electron correlations in TMOs, the DFT+U theory  is employed where the DFT energy is captured by the Perdew-Burke-Ernzerhof (PBE) functional, and the on-site Coulomb potential U are from Ref. . We interpret the DFT+U orbital energies as quasiparticle (QP) energies and perform FCH rather than XCH calculations so as not to favor any particular occupation. More numerical details can be found in the Supplemental Material .
Strikingly, the FCH approach systematically underestimates the intensity of the near-edge peak associated with the t2g manifold for all selected TMOs (Fig. 1(b)). The t2g peak at ~ 529.2eV of CrO2  suffers from the most severe underestimation. It becomes a weak, broad feature in the FCH simulation as compared to the strong, sharp peak in experiment. The t2g-eg peak intensity ratios of VO2 and Fe2O3 are also too low - both are predicted as 0.7, compared with 1.7 and 1.0 as measured, respectively.
where ψc’s are the initial-state orbitals and c iterates over empty orbitals only. The exciton amplitude and energy Ef can be solved from the BSE as in Ref. [11–14]. The core-level BSE calculations in this work employ the OCEAN code [14, 45].
As is shown in Fig. 1(c), the BSE substantially improves on the O K-edge line shapes. The “t2g” peak intensity is retrieved for each investigated TMO, particularly for CrO2. The simulated t2g-eg intensity ratios are almost as measured for TiO2, α-Fe2O3, and SrTiO3 (below 537 eV). These BSE spectra are also in better agreement with experiments than the previous results using the initial-state rule [46, 47].
In fact, the core-hole effects predicted by the ΔSCF method in these oxides are counterintuitive. Typically, excitonic effects tend to sharpen the absorption edge due to e-h attraction [5, 11, 12, 43]. In Fig. 2, we show the core-hole effects from ΔSCF by comparing the initial- and final-state spectra. The core-hole attraction does redshift the spectra by over 1 eV in both the TM 3d pre-edge and the 4sp  region (near 544 eV). However, the FCH core-hole effect substantially reduces the t2g peak intensity. Similarly underestimated pre-peak intensity was encountered before [49–51] but no satisfactory explanation has been provided to date.
With this comparison, we can now explain why the one-body ΔSCF formalism tends to underestimate the peak intensity. The single-particle approximation in Eq. (2), Sf|ε·r|ψh, does not contain any information about the initial state of the valence electrons, when S is a constant. The spectrum is found by summing up single-particle transitions in the final state in which the valence electrons have already equilibrated with the core hole. The final-state picture alone does not reflect the full process of X-ray excitations: first a core hole is generated by an X-ray photon, then the valence electrons begin to relax in presence of the core-hole potential. This is why the final-state spectra predicted by the one-body formalism can be subject to spurious charge transfer effects, if the electron density happens to transfer onto the optically active orbitals due to core-hole attraction. In the above examples of O 1s excitations, a substantial amount of electron density is transferred into the 2p character local to the excited O (0.55 e− for TiO2 and 0.72 e− for CrO2), which blocks some bright transitions from 1s → 2p and leads to the reduced pre-edge intensity. Conversely, a many-body perturbation theory, such as the BSE, considers the evolution of many-electron systems using diagrammatic approaches, which is capable of capturing the correct time-ordering and dynamics in photoexcitations. This can be seen in Eq. (3), the BSE begins with empty initial-state orbitals and considers their hybridization under the core-hole effect, hence retrieving the edge sharpness.
To alleviate the deficiencies of the one-body ΔSCF method, we seek a better approximation to Ψf|ε·R|Ψi so as to incorporate the many-electron response as the system evolves from the initial to final state. We abandon the single-particle approximation used in Eq. (2) and develop an approximate many-electron wavefunction. To avoid solving an intractable full configuration-interaction problem and leverage the current ΔSCF approach, we assume that the valence Kohn-Sham (KS) orbitals are non-interacting. Consequently, a single Slater determinant comprising the KS orbitals is a consistent many-electron eigenstate. This approximation is analogous to the independent-electron model with a single core-hole considered in the MND theory[15, 16]. It may fail to capture certain many-body correlations and we will defer this discussion until later.
Within the independent-KS-orbital approximation, the ground state |Ψi can be constructed as a Slater determinant built from the N lowest occupied valence orbitals of the pristine supercell and the core orbital of the impurity atom. Here, N is the number of valence electrons in the supercell. Similarly, the final state |Ψf can be constructed as a Slater determinant of arbitrary N+1 KS orbitals in the final-state picture. To enumerate the many-electron final states, we generalize the final state index f to represent any excited configuration of these N+1 electrons: f now represents a unique (N + 1)-tuple, f = (f1, f2, ···, fN+1), where fi’s label distinct final-state orbitals that the N+1 electrons can occupy.
The single-particle orbitals necessary for constructing |Ψi and |Ψf can be readily obtained from respective ΔSCF calculations performed over the supercell. However, |Ψi and |Ψf are represented with the initial- and final-state orbital basis respectively. To eventually access Ψf|ε·R|Ψi, we perform the orbital transformation from the initial to final state: |i = Σj |ψjψj|i = Σj ξij|ψj, in which i and ψj are the ith final-state and jth initial-state orbitals in the supercell, respectively. ξij is the transformation coefficient. The resulting matrix element has the same form as in Eq. (3) but the amplitude of a given final state f becomes a single determinant of the transformation coefficients 
Within the independent-KS-orbital approximation, the energy of |Ψf can be found by summing up the energies of the occupied KS orbitals in the final state
The energy difference Ef − Ei for producing spectra may employ a similar alignment scheme as in the single-body approach, since the initial(ground)-state energy Ei is fixed. In the following discussion, it is more meaningful to define the relative energy Ωf = Ef − Eth, where the threshold energy, , corresponds to the lowest-energy final-state configuration (1, 2, …, N + 1).
Determinant expressions similar to Eq. (4) were also obtained in previous work [52–54] but they are rarely applied in a solid-state context from first-principles. Thus it is of great interest to examine whether the many-electron ΔSCF formalism in Eq. (4) and (5) can reproduce the correct lineshapes for the investigated TMOs.
Evaluating the many-electron seems formidable at first glance because a solid contains many electrons, which leads to a combinatorially huge number of final states. To facilitate the calculation, we regroup the final-state configurations according to the number of e-h pairs excited. For example, we denote a single configuration with one core-excited e-h pair as f(1) = (1, 2, ···, N, fN+1), which is used in producing the one-body ΔSCF XAS and fN+1 is the single unoccupied orbital in the final state. Based on this configuration, we can define a double configuration with one more valence e-h pairs as f(2) = (1, ···, i − 1, i + 1, ···, N, fN, fN+1), where fN+1 > fN > N, and so forth (Fig. 3).
Despite a large number of excitation configurations, not all of them contribute equally to the near-edge XAS, particularly near the absorption onset. First, for an insulator, Ωf will increase proportionally with the number n of e-h pairs excited across the band gap (Eg): Ωf ≥ (n − 1)Eg. This significantly reduces the number of necessary final-states for producing the first few eVs of XAS for an insulator. Second, the many-electron transition amplitude may decrease rapidly with increasing number of e-h pairs excited, as we will show next.
We reexamine two extreme cases: the insulating TiO2 and the metallic CrO2. Fig. 4 shows the XAS calculated from the many-body formalism as in Eq. (4) and (5). The orbitals used for calculating the transformation coefficients ξij in are obtained from the same set of ground-state and FCH calculations  as previously described. For simplicity, only the states at the Γ-point of the supercell Brillouin zone are used for producing the XAS. Remarkably, the simulated XAS lineshapes (Fig 4) with the many-electron formalism are in excellent agreement with experiments. In particular for CrO2, the edge sharpness is completely retrieved (Fig 4(b)) despite the disappearance of the first peak in a single-electron theory.
We find calculated XAS for TiO2 converges at the order of f(1). As expected, the contributions from f(2) (Fig 4(c)) appear at higher energies due to the sizable band gap, and are much weaker due to reduced wavefunction overlap. However, it is more challenging to achieve numerical convergence in the metallic CrO2. While the first peak is retrieved mainly from f(1) excitations, the correct peak-intensity ratio can only be reproduced when the next-order, f(2), is taken into account, which substantially increases the intensity of the absorption feature ~ 4 eV above onset (Fig 4(b)). Meanwhile, we find f(3) can be neglected up to the first 8 eV (Fig 4(d)).
Both the BSE (Fig. 1) and the many-electron formalism (Fig 4) yield very similar spectra for TiO2, but the many-electron formalism produces a spectrum for CrO2 in closer agreement with experiment than the BSE. The experimental peak-intensity ratio is 1.7; 1.6 from the many-electron formalism; but only 1.3 from the BSE. We attribute this improvement the explicit inclusion of many-electron response to the core-hole potential via the transition amplitude [Eq. (4)]. Since the Slater determinants are exact solutions to the independent-electron model, the many-electron formalism has taken into account all possible many-electron processes in the MND theory. In addition to the e-h ladder diagrams in the BSE, the MND theory also considers non-ladder diagrams, such as diagrams with crossing Coulomb lines and bubble diagrams that account for the orthogonality catastrophe [15, 16, 52]. Instead of adopting a diagrammatic approach, the many-electron formalism considers theses processes via incorporating response of each electron to the core hole in the determinant amplitude . Therefore, the many-electron formalism goes beyond the two-particle correlation in the BSE. For insulators like TiO2, the many-electron formalism produces results very similar to the BSE because the non-ladder diagrams are not important near the absorption edge. However, this is not the case for CrO2 since multiple e-h-pair production is more likely in a metal and the non-ladder diagrams become significant for determining the near-edge lineshape. The significance of the multiple electron response can also be appreciated from the fact that the CrO2 spectrum converges at a higher-order f(n).
At last, we briefly discuss the validity of the independent-KS-orbital approximation. On the f(1) level, there is a one-to-one correspondence between the many-electron configurations and the empty single-particle final-state orbitals, and the latter are typically good approximations to QP band topology and wave-functions [43, 55]. This approximation may break down if there are one or more valence e-h pairs. For example, valence e-h pairs in the f(2) configurations of an insulator may hybridize to form bound excitons; the two excited electrons in the f(2) configurations may interact through strong on-site Coulomb repulsion. In the case of CrO2, however, the strong metallic screening can largely reduce e-h binding energy, maintaining the independent-KS-orbital approximation. As for strongly correlated effects, it will be worthwhile to further consider embedded models based on the current determinant formalism.
We have shown that the BSE and a newly developed many-electron ΔSCF approach are highly predictive for O K-edge fingerprints of TMOs. The systematic underestimation of the peak-intensity ratio within the original one-body ΔSCF approach is attributed to the absence of many-electron response in this formalism. We have demonstrated how to rectify these shortcomings by (1) extending the wavefunctions in the original ΔSCF to a many-electron form (2) calculating the determinant form of the transition amplitude. This many-electron formalism is transferable and not at all peculiar to TMOs. We leave the discussion of shakeup effects, improvement of efficiency, and other examples to near-future work.
Theoretical and computational work was performed by Y. L. and D. P. at The Molecular Foundry, which is supported by the Office of Science, Office of Basic Energy Sciences, of the United States Department of Energy under Contact No. DE-AC02-05CH11231. W. D. was supported by the Joint Center for Artificial Photosynthesis, a DOE Energy Innovation Hub, supported through the Office of Science of the U.S. Department of Energy (Award No. DE-SC0004993). We acknowledge fruitful discussion with Chunjing Jia (Y. L.) and Bill Gadzuk (J. V.). We also would like to acknowledge the referees for giving very patient, professional, and thoughtful comments in the review process. Computations were performed with the computing resources at the National Energy Research Scientific Computing Center (NERSC).