DFT with corrections for an efficient and accurate description of strong electron correlations in NiO

An efficient and accurate description of the electronic structure of a strongly correlated metal-oxide semiconductor like NiO has been notoriously difficult. Here, we study the capabilities and limitations of two frequently employed correction schemes, a DFT+U on-site correction and a DFT+1/2 self-energy correction. While both methods individually are unable to provide satisfactory results, in combination they provide a very good description of all relevant physical quantities. Since both methods cope with different shortcomings of common density-functional theory (DFT) methods (using local-density or generalized-gradient approximations), their combination is not mutually dependent and remains broadly applicable. The combined approach retains the computational efficiency of DFT calculations while providing significantly improved predictive power.


Introduction
Optically transparent and electrically conducting transitionmetal oxides (transparent-conducting oxides, TCOs) are important functional materials for a variety of optoelectronic devices [1][2][3][4]. Understanding their electronic structure and being able to model it accurately is an important * Author to whom any correspondence should be addressed.
Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. requirement to advance these research fields. Unlike other TCOs like TiO 2 or ZnO, which can be treated rather well by electronic-structure calculations, [5,6] NiO is notoriously difficult to study with state-of-the-art approximations within the density-functional theory (DFT) [7,8]. For many materials, DFT has become a reliable tool for calculating various properties, with a good balance of computational efficiency and predictive power. To conquer strongly correlated oxides, one typically must go beyond DFT due to its meanfield approximation, which leads to an underestimation of localization of electrons in states of valence d or f shells. However, many-body approaches that are able to describe strongly-correlated materials accurately are often still computationally too demanding for their broad application [9,10].
A solution that does not add significantly to the computational cost of a common DFT calculation in the local-density approximation (LDA) or generalized-gradient approximation (GGA) has been proposed: combining parameterized Hubbard-model Hamiltonians with DFT (so-called DFT+U(+J)), one can improve the description of both electron-electron Coulomb and exchange interactions by adjusting the Hubbard-model parameters (U and J, respectively). However, with such a correction, semi-empiricism in the form of the two Hubbard parameters [11] (or one effective parameter U eff = U − J) [12] is introduced, which require a system-dependent choice. Although rigorous ways of determining the required parameters have been proposed, [13][14][15] the resulting parameters enable to reproduce often not all the desired target quantities [16]. For example, electronic band gaps often turn out to be still severely underestimated with such rigorously obtained Hubbard parameters. Thus, in many studies this approach is not followed. Instead, Hubbard parameters are adjusted by comparing results to a variety of reference material properties (such as band gap, lattice or spin structure, etc). Doing so, one usually finds that different parameters are needed for reproducing different observables. Such an approach has very limited predictive power. Any good agreement with respect to some reference observable is obtained by introducing deviations in other aspects of the electronic-structure description, which may lead to qualitatively nonphysical findings. Furthermore, in the case of NiO, quantitative agreement with the experimental band gap cannot be achieved with such a Hubbard approach acting solely on Ni 3d states.
Here we demonstrate how the electronic band structure of the difficult case of NiO can be modeled in quantitative agreement with experimentally observed characteristic features of the electronic structure of NiO. Two shortcomings of LDA and GGA DFT need to be overcome that can be addressed separately, namely the underestimated localization of O 2− ligand states at the valence band edge (VBE) [17][18][19] and the underestimated on-site correlation effects of Ni 2+ ions [20,21]. To do so we investigate a combined approach of DFT+1/2 [22,23] and DFT+U. It retains the computational efficiency of DFT methods and is in general applicable to any strongly correlated oxide.
In this work, all ground-state properties, i.e. total energies and atomic forces, are calculated using PBE. For calculating electronic band structures we apply and combine two correction methods by including orbital relaxations due to electron excitations approximately: DFT+U [11,12] and DFT+1/2 [22,23].
Hubbard U corrections are carried out using the simplified, rotationally invariant Hubbard Hamiltonian of Dudarev et al, with effective U values discussed in detail in the text [12]. This DFT+U method corrects the on-site correlations at the Ni 2+ ions. Note that this approach is qualitatively equivalent to the one by Liechtenstein et al [11]. The latter gives slightly larger band gaps for the same U and J values, within 0.1 eV, and the (unphysical) crossing of oxygen dominated states at the Γ point (discussed later in the text) appears already at smaller values (U = 5 eV, J = 1 eV).
To determine U values by means of linear-response calculations [28], the change of the on-site number of electrons N I with respect to a small change in the on-site potential V J , χ IJ = ∂NI ∂VJ is monitored. U is obtained from fitting the difference of linear-response functions calculated self-consistently, χ, and non-self-consistently, χ 0 , as U = χ −1 − χ −1 0 , following [28]. In this way we obtained U = 5.15 eV, see figure S2 in the supplemental information.
To correct the localization of the ligand states at O 2− ions in transition-metal oxides, we choose the DFT+1/2 self-energy correction [22,23]. It is an efficient adaption of Slater's halfoccupation rule for solids, and it is implemented for VASP pseudopotentials [29]. For further details, e.g. how the cutoff radius for each species is obtained by maximizing the influence on the electronic band gap, see [22,30] and figure S3 in the supplemental information.

Reference single particle electronic-structure
Below the Néel temperature of 523 K [31,32], an antiferromagnetic (AFM) spin structure stacked along the [111] direction (AFM type II) is observed in NiO. We describe this by a non-cubic unit cell with two formula units (figure S1(a)). The AFM state is more stable than a ferromagnetic or a nonmagnetic (non-polarized) state in our spin-polarized PBE calculations, in line with previous LDA [33] and PBE [34,35] calculations.
The calculated local magnetic spin moment of a Ni ion is µ(Ni) = ±1.25 µ B . This is smaller than experimentally observed, µ = 1.66 µ B to 2.2 µ B [31,[36][37][38][39][40]. AFM NiO correctly turns out to be a semiconductor with a band gap. However, with PBE the smallest indirect band transition is between the L point in the valence band (VB) and a point between P1 and Z in the conduction band (CB) and has a value of 0.79 eV (see figure S1(b)) for the definition of highsymmetry points within the Brillouin zone). The smallest direct band gap, located between P1 and Z, has a value of 1.01 eV. Both PBE values significantly underestimate the commonly agreed upon experimental value of 4.3 eV [41].
A clear classification of NiO as Mott-Hubbard or chargetransfer insulator remains difficult. In addition to the size of the band gap, broad consensus has been reached on the orbital contributions to the states near the band edges, e.g. from comparing (inverse) photoemission spectroscopy data at different energies [41,42]. Ni 3d states split into three contributions: the upper Hubbard band, which dominates the conduction-band edge (CBE), the lower Hubbard band, located about −8 eV below the Fermi level, [18,41], and a third band, which hybridizes with O 2p states to form the VBE. Thus, an excitation into Ni 3d e g states is evident, whereas the correct model to describe this excitation is a subject of ongoing discussion [42][43][44][45][46][47][48][49][50][51][52][53]. In the following we focus on comparing our results to known data for these quantities: energies and orbital contribution of electronic bands.
The DFT DOS is plotted in figure 1(a) (see figures S4 and S5 for a detailed view of the band edges and individual orbital contributions, respectively). The CB is dominated by Ni 3d − e g states, with small admixture of O 2p states. The VB is formed by Ni 3d − t 2g and e g orbitals, as well as O 2p states. O contributions dominate lower lying levels (around −5 eV and −6.6 eV) but are present at the band edge, too. Ni contributions dominate at the VBE and form another lower peak around −1.6 eV. In the AFM spin-structure, these features belong to Ni atoms of opposite magnetic moment. These findings are in line with earlier DFT results [18,34]. Compared to the above recalled experimental findings, however, deviations occur in terms of band gap size, orbital character, and location of the lower Hubbard band.

DFT+U
We now come to the effect of an on-site Hubbard correction, applying U eff on the Ni 3d states [12]. Its motivation is to split lower and upper Hubbard bands, i.e. shifting empty and filled 3d states towards higher and lower energies, respectively, and thus effectively increasing the band gap. Applying U = 5 eV (PBE+U(5 eV), see figures 1(b) and 2(b)) apparently fulfills this purpose [54][55][56][57] and is physically motivated [12,13]. µ(Ni) is increased to 1.68 µ B , which is within the range of experimental data. The lower Hubbard band is shifted considerably downwards, leading to a Ni 3d − e g dominated peak around −6 eV. The 3d − e g contributions now overlap with the smaller fractions of Ni 4s and 3p and O 2s and 2p states that are found in this energy region using PBE. This lower Hubbard band is in accordance with experiments [41]. The CB is shifted upwards. However, the magnitude of the shift is not sufficient to reproduce experimental findings. The DOS is centered ≈3.7 eV above the VBE, i.e. lower than the experimental peak observed around 5 eV above the VBE [41]. The indirect (direct) band gap is only increased to 3.09 eV (3.42 eV), i.e. lower than in experiments (4.3 eV). Due to the downward shift of Finding that the band gap increases with U motivates a choice of larger U values around 10 eV [18,42,58,59]. The changes induced by PBE+U(10 eV) are displayed in figures 1(c) and 2(c). All trends discussed comparing PBE with PBE+U(5 eV) are intensified. µ(Ni) is further increased to 1.91 µ B . The VBE is further depleted of Ni states, leading to almost pure O bands (2p with small amounts of 2s). This disagrees with the experimentally observed VB character and the interpretation of NiO as not being a pure charge-transfer insulator [44][45][46]. Finally, although the band gap is further increased, the magnitude of this increase is significantly less by increasing U from 5 to 10 eV, as compared to the increase from 0 to 5 eV. The (direct) band gap for U = 10 eV is 3.7 eV (4.03 eV), still underestimating the experimental value. This can be explained by analyzing the band structure: Although the unoccupied Ni 3d − e g states (with little dispersion) are shifted further up, the tail of a single O dominated band (mixed character of O s and p and Ni s and p states), which has a much larger band width, is not shifted upwards. Thus, the CBE is this band tail now, with important consequences: the location of the lowest energy of the CB is moved to the Γ point, leading also to the smallest direct band gap being located at Γ ( figure 2(c)). This state, however, has O character ( figure 1(c)). This is problematic for two reasons: First, this state prevents the sufficient increase of the band gap via an on-site U correction applied to Ni 3d orbitals, since it is not directly affected by such a correction. Second, the orbital character of this band, namely O like, would indicate a charge transfer towards O. Thus, this feature is nonphysical, as pointed out in [34].
From these results, it is evident that the electronic structure of NiO cannot be well described with a scheme where a potential parameter U, acting only on Ni 3d orbitals, is introduced: For a moderate value around 5 eV the band gap remains underestimated by 28%. Larger values lead to a depletion of Ni 3d density at the VBE and are unable to increase the band gap sufficiently. This also means that methods which build on such results inherit these qualitative problems of the DFT+U approach (see D. Discussion).

DFT+1/2
Underestimation of band gaps is a notorious problem of DFT methods. Even without the added complexity of magnetic transition-metal elements, the calculated Kohn-Sham energy gap differs from the true band gap by the discontinuity term [60][61][62]. Recently, a scheme has been extended from atoms and molecules to solids, which implements the halfoccupation approximation to overcome the band-gap problem [22]. While the DFT+1/2 method has shown promise for a variety of materials, including the TCOs TiO 2 and ZnO, its application to compounds of TMs with fractionally filled valence d shell, such as NiO, does not give the same good agreement with experimental observations [63].
We analyze the effects of PBE+1/2, starting with the Ni sublattice. In PBE calculation, Ni 3d orbitals dominate both band edges. In such a case, the self-energy correction acts equally on the VBE and the CBE and is expected to leave the band gap unchanged. We confirm that the 1/2 method has no effect on the Ni 2+ ions ( figure S3(a)). Consequently, this method is not fit to improve the position of Ni 3d states.
In contrast, the contribution of O (mostly 2p) is more pronounced at the VBE than the CBE, as expected for an ionic compound with oxygen ions of a nearly formal O 2oxidation state. Adding the 1/2 correction to the O potential affects the band gap (with an optimal cut-off radius r cut (O) = 2.1 Bohr, see figure S3(b)). The direct (indirect) band gap obtained with this PBE+1/2 method is 1.00 eV (1.13 eV) ( figure 2(d)). Thus, PBE+1/2 severely underestimates the experimental band gap value, even though it slightly increases in value compared to the PBE result.
We analyze the effects of the 1/2 approach by comparing the PBE+1/2 DOS to the PBE result ( figure 1(d)). The main effect is a shift of O 2p (and 2s) states to lower energies. This is visible in the two O dominated peaks between −4 and −8 eV, but also at the VBE. This is in line with the effect that a self-interaction correction (SIC) to DFT has for NiO [18]. This is expected, as both correction schemes can be shown to have the same effect on Kohn-Sham states [22]. However, the DFT+1/2 approach has a practical advantage of providing a clear recipe for optimizing the empirical parameter (r cut ), whereas in SIC [5,6] an empirical parameter needs to be chosen with respect to agreement with reference data. This is difficult for a material like NiO, where multiple corrections are necessary. In any case, the location of the lower O peaks agrees very well with the outcome of the LDA+SIC70 approach in [18], validating the choice made in the latter.
Another important detail is visible when comparing the band structures of PBE and PBE+1/2: The unoccupied O dominated band with large band width is shifted upwards by the 1/2 correction, away from the empty Ni 3d − e g states. This indicates that a larger band gap increase is possible by the application of a Hubbard U, compared to the case where the correction on O is omitted.

Combination DFT+1/2+U
We have demonstrated so far that i) the on-site correlation correction of the DFT+U method shifts Ni 3d states qualitatively into the correct direction and that ii) the self-energy correction of the DFT+1/2 method increases the gap between empty and filled O 2p states. Both correction schemes effectively increase the band gap, but neither is able to provide a convincing electronic-structure description. Since both corrections act on different aspects, this cannot be expected.
Therefore, we study the combination of both corrections. The DOS of PBE+1/2+U(5 eV) is shown in figure 1(e). First, the Ni dominated lower Hubbard band that forms due to the application of U = 5 eV is retained, but it is shifted to slightly higher energies (≈−5 eV). This is in excellent agreement with the computationally expensive calculation presented in [34], using hybrid-DFT data as input for a G 0 W 0 calculation or the GW+DMFT result in [10]. However, compared to experimental results, this feature would be expected at an even lower energy (≈−7 eV).
In the VB region, the too strong depletion of Ni 3d states in PBE+U results is corrected: Already for PBE+U(5 eV), the Ni 3d states at the VBE have less density than a second peak ≈−2 eV below, which is not in line with the shape of the experimental spectrum [41]. In the combined PBE+1/2+U(5 eV) approach, states with strong Ni 3d − t 2g contributions, admixed with O 2p states, are observed as peak with the highest density in the VB region.
The CB is centered around 5 eV above the Fermi level, in agreement with the experimental spectrum [41]. The artificial overlay of the empty O dominated band with the Ni 3d − e g states at the CBE is removed ( figure 2(e)). Our computed indirect (direct) band gap for DFT+1/2+U(5 eV) is 4.6 eV (4.77 eV), in quite good agreement with the experimental value. Possible reasons for the remaining small difference of ≈0.3 eV are i) an uncertainty that may be attributed to experimentally obtained band gaps, in particular for oxidic systems where point defects like O vacancies can easily occur and affect the measured band gap or ii) further effects ignored in our present study, e.g. a temperature broadening of levels from structural dynamics at finite temperatures [30,64]. In any case, recalling that all investigated features of the electronic structure calculated with DFT+1/2+U(5 eV) are computed at DFT numerical cost, their agreement with those from established experiments is striking.
From here, two avenues are open to further sophisticate our method. On one hand, for a well studied material, such as NiO, it is straight forward to keep on treating U as empirical parameter. In contrast to works where the localization of the O ligand states was not corrected, with the combined DFT+1/2+U approach, we are able to fine tune this parameter to reproduce all relevant quantities with just one U value. Doing so and accepting 4.3 ± 0.1 eV as target band gap that should be reached by a corrected DFT calculation, we slightly lower the U parameter to 4.5 eV. This small adjustment of U does not change any of the features discussed for DFT+1/2+U(5 eV) qualitatively. Quantitatively, µ(Ni) = 1.80 µ B remains within the range of experimental values and the direct (indirect) band gap reduces to 4.35 eV (4.50 eV). Thus, a practical electronicstructure description is achieved, which is able to reproduce all the addressed experimental observations quantitatively.
On the other hand, we can revisit the physical origin of the U parameter and optimize its value by evaluating the linear response of the electronic charge in the corrected state (Ni 3d) to the added potential parameter U [28]. We obtain a value of U = 5.15 eV in a (2 × 2 × 2) cell, containing 32 atoms ( figure S2). U values from this procedure are size dependent but moderate systems provide already good estimates. We obtain a final extrapolated value of U 5. To analyze whether the two corrections mutually affect each other, we repeat the optimization of the U parameter based on the optimized DFT+1/2 calculation in the 32-atoms cell. We obtain U = 5.13 eV, i.e. the change compared to optimizing U based on pure DFT is negligible. Analogously, optimizing the cutoff radius for the 1/2 correction on oxygen leads to the same result when performed based on DFT+U(5 eV).
Both correction approaches, empirical tuning of the U parameter to reach the experimental band gap and optimizing the U parameter by a linear-response calculation are equally successful, provided that the DFT+1/2 correction is applied on O sites (with the properly determined correction radius). The difference in U is only 0.7 eV, which leads to rather small changes in most of the addressed electronic-structure features. µ(Ni) changes only by ∆µ = 0.02 µ B and the band gap by ∆E g = 0.36 eV. The unphysical shift of O states below the Ni 3d states, letting the O states form the CBE and impeding the sufficient increase of the band gap in pure DFT+U methods, does not occur in either method.
Since U may also change the equilibrium lattice constant of the cubic NiO crystal, we optimize the latter with both U = 4.5 and 5.2 eV. Compared to the DFT-PBE value, deviations are negligible (∆a 0 remains within 0.1%) and, consequently, so are changes in the electronic structure with respect to these structural effects.

Discussion
We conclude by comparing the results obtained with the DFT+1/2+U approach with those of various other theoretical studies of the electronic structure of NiO, see table 1.
Only recently, a many-body approach without dependence on a starting point has been employed for NiO [10]. Compared to experiments, the electronic band gap is slightly underestimated by 0.3 eV. The magnetic moment is in line with our results and experiments. The Ni 3d contributions are in line with our DOS. Spectral features and contributions to band edges agree with experimental findings and our results. The contributions to the CB at Γ are in line with our results (mixed character of O 2s and 2p and Ni 4s, however, in their projected contributions s-type contributions are dominant over p-type contributions). The CBE is, however, located in the region between P1 and Z (figure 2(e)), where the CB is dominated by Ni 3d − e g states in our results as well as at Z in [10]. Unfortunately, no band diagram is reported in [10]. When the dependence on a starting point is not avoided, results of many-body calculations vary strongly (see G 0 W 0 results in table 1): One one hand, LDA+U or PBE lead to an underestimation of the band gap. On the other hand, results that are in very good agreement with experiments are observed with a HSE03 starting point. This is in line with the performance of the underlying DFT methods (the underestimation of the band gap is reduced comparing PBE > LDA+U > HSE03). All features of the DOS computed by DFT+1/2+U agree well with the G 0 W 0 corrected many-body results with a starting point that leads to good agreement with experiments [34].
A dependence on the underlying DFT method is also observed when coupling DFT with a DMFT description of the correlated states. It has been shown that the term to account for double counting of electron correlation for the respective states drastically changes the qualitative features of the electronic structure of NiO [65]. Since determining the correct double counting is not straight forward, in practice this means that qualitative features of DFT+DMFT calculations are inherited by the underlying DFT description. For example, in [35,66] the lowering of the band tail at Γ is reported in the embedded DMFT result. Based on our analysis above, this is an error of the underlying DFT+U description. Thus, an improved single-particle description as presented here by means of the DFT+1/2+U approach would be beneficial for more accurate future many-body studies.
Hybrid DFT itself can be viewed as intermediate between many-body approaches and DFT calculations in terms of computational cost. In terms of accuracy, their performance depends on the construction, in particular the amount of exact exchange. For a given functional, the proper amount is a system dependent parameter [67,68]. In practice, this is usually neglected and fixed parametrizations are frequently used, that attempt to balance the error for a wide range of materials. Nevertheless, some popular hybrid functionals perform well for NiO. The band gap values of HSE and B3LYP are in good agreement with experiments. For the latter, the parabolic state at Γ does not cross the manifold of Ni 3d states, in line with our results as well as the many-body results in [34]. Note, that the dependence on the exact-exchange fraction is evident by a variation of 0.3 eV comparing the band gap values obtained with HSE03 and HSE06, respectively. Last but not least, a semi-local functional, mBJLDA, has been demonstrated to perform well for NiO [69,70]. The band gap of 4.1 eV and the main features of the DOS are similar to our results and, thus, agree with experimental findings and many-body results [69].

Conclusion
Strongly correlated materials, such as NiO, have long been an elusive problem for an accurate electronic-structure description by DFT methods. Here we have shown that a combination of two approaches, the DFT+U on-site correction for Ni 3d states and the DFT+1/2 self-energy correction for O states, which are unable to provide satisfactory results when applied individually, provide an accurate and efficient electronic-structure description. In combination with the DFT+1/2 method, a rigorously obtained U parameter agrees well with an empirically chosen optimal value. Thus, the inherent ambiguity of choosing U that exists in the DFT+U method is overcome. Thus, the DFT+1/2+U approach is in principle widely applicable, also for compounds where an empirical choice of U cannot be easily made due to a lack of available high quality reference data. The DFT+1/2+U approach is a computationally efficient single-particle method whose electronic properties agree well with more expensive theoretical approaches and experimental results. It can be applied to more complex atomistic systems, including atomic or extended defects, surfaces, etc. Thus, it enables practical computational studies of complicated problems in material science. In addition, approaches to study many-body effects in a material like NiO [18] can benefit from this accurate singleparticle description as a potentially improved starting point over DFT+U.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.