Electronic structure of (organic-)inorganic metal halide perovskites: the dilemma of choosing the right functional

Organic-inorganic metal halide perovskites (HaPs) are intensively studied for their light-harvesting properties. Owing to the interplay between strong electron-electron interaction and spin-orbit coupling (SOC), their quantitative theoretical description is still a challenge as evidenced by the wide variety of results available in literature. Here, various methodologies for computing their electronic structure are evaluated, also accounting for SOC. More specific, the GW approach as well as variants of the hybrid functionals PBE0 and HSE are at the center of our investigations. For both functionals, we explore methods to determine the mixing parameter $\alpha$, and for HSE, we investigate the impact of the screening-parameter $\omega$. An extensive investigation of PbI2, a precursor of many HaPs, leads to the conclusion that hybrid functionals with $\alpha$ tuned by the density-based mixing method are most suitable for obtaining band gaps comparable to $G_0W_0$ results. Moreover, this methodology is transferable to CsPbI3, and the same behaviour is expected for the entire family of lead-iodine perovskites.


Introduction
Organic-inorganic metal halide perovskites (HaPs) are in the focus of various optoelectronic applications, be it efficient solar-cells devices, lasers, or detectors. For instance, solar cells employing them as active layers, exhibit a power-conversion efficiency (PCE) of 25.5% [1], and perovskites/Si tandem cells have reached within less than five years a PCE of 29.5%, a higher value than for any other solar cell on the market. Perovskites have the chemical formula ABX 3 , where A is an organic or inorganic cation, like methylammonium (MA + ), formamidinium (FA + ) or Cs + ; B is a divalent metal cation, like Pb 2+ or Sn 2+ ; and X is an anion of the halogen group (Cl − , I − , Br − ). Among them, most suitable for solar-cell applications are lead-iodine based perovskites (APbI 3 ). Many current studies are focused on finding a replacement for the toxic element Pb, thereby preserving or even improving the physical properties. Double perovskites, combining a trivalent and a monovalent metal cation, have been proposed as an alternative, also thanks to good stability and an energy gap that is highly tunable in the visible range [2,3]. Other strategies to improve the material's stability consist in mixing organic and inorganic cations [4,5,6] or alternating layers of 3D and 2D HaPs [7]. Due to a growing complexity, resulting in large unit cells, one aims at accurate and, at the same time, computationally efficient approaches with high predictive power. The variety of methods to overcome the band-gap problem of semi-local density-functional theory (DFT), spans from employing hybrid functionals with varying amount of exact exchange to the GW approximation of many-body perturbation theory (MBPT), or a combination of both. All of them are computationally involved. Another particular challenge is to capture the strong spin-orbit coupling (SOC) effects [8] that significantly lower the band gap. The goal of this work is to provide insight into the performance of different methods and to determine a fully ab initio procedure to accurately compute the band gaps of lead-based perovskites and their precursor PbI 2 . To highlight the problem, we summarize in Figure 1 a collection of band-gap results for PbI 2 as computed with a variety of methods [9,10]. Obviously, the enormous spread is far from being satisfactory, suggesting a lack of predictive power. We will get back to these results step by step in Section 4. To shed light onto the role of different approaches, we are going to reexamine the performance of commonly  [22,23]) is indicated by the blue bar. Results of this work (computed with exciting) are depicted by stars, the others are collected from literature, i.e., references (a) [24], (b) [25], (c) [26], (d) [27], (e) [28], (f) [29], (g) [30]. These include computations with and without SOC, performed with semi-local and hybrids functionals, as well as G 0 W 0 calculations using different functionals as starting point.
used semi-local and hybrid functionals for the class of the perovskite materials. In this context, we mention other works to determine a reliable method to compute the electronic structure of HaPs. Among them, several focus on calculations based on MBPT [11,12,13,14], others employ HSE [15,16] or aim at improving hybrid functionals by tuning their parameters [17,18]. Here, we place a focus on the dependency of the results on the functionals' parameters. In particular, we apply the dielectric dependent hybrid (DDH) method [19] and the density-based mixing (DM) method [20]. We aim at finding optimal and transferable mixing parameters for two most popular hybrid functionals, i.e., PBE0 and HSE. Moreover, we investigate the dependence of one-shot GW (G 0 W 0 ) on different parametrizations as the starting point. All the calculations are performed using the full-potential all-electron computer package exciting [21]. The computational details are provided in the Appendix.

Hybrid functionals and their parameters
To put our work in a bigger context, we shortly summarize the status of hybrid functionals. First suggested by Becke in 1993 [31], they combine a (semi-)local functional with a fraction of Hartree-Fock (HF) exchange. This way, they partially compensate for the missing exchange-correlation (xc) discontinuity and the self-interaction problem. Part of their success in reasonably predicting band gaps is associated with the error cancellation between the effects of the HF method, which tends to highly overestimate them, and KS DFT that often underestimates them [32,33]. After the half-and-half hybrid functional [31], based on this idea, several such functionals were swiftly developed [34,35,36]. Among them, PBE0, proposed independently in two works [37,38], is often termed the first fully ab initio hybrid xc-functional. The PBE0 energy reads as follows: where the mixing parameter α is set to 0.25. This fraction of HF exchange was justified by Perdew and coworkers [36] by applying Görling-Levy perturbation theory [39]. Nevertheless, in Ref. [36], the authors suggest to optimize the mixing parameter for each system and property. PBE0 is frequently used for solids, though generally overestimating band gaps of typical semiconductors [40,41]. This overshoot is caused by the inclusion of the long-range (LR) tail of the Coulomb interaction that is expected to be effectively screened in periodic systems [42].
To overcome this problem, Heyd, Scuseria, and Ernzerhof presented a screened hybrid functional [43], known as HSE (first as HSE03), in which only the short-range (SR) part of the HF exchange is considered. The elimination of the LR part, besides improving the physical description of solid systems, is computationally convenient due to faster convergence with respect to the number of k-points [44]. To compute the SR HF exchange and substitute it by a fraction of SR PBE exchange, the Coulomb operator is split in a SR part and a LR part. To do so, one possibility is to make use of the error function and its complementary: where ω is the screening parameter. This parameter assumes different values in the different forms of HSE, such as HSE03 [43,45] and HSE06 [46]. As explained in the Ref. [45], HSE03 considers different values of the screening parameter for the HF-SR and the PBE-SR part (ω HF = 0.15/ √ 2 ≈ 0.106 a −1 0 and ω PBE = 0.15 · 2 1/3 ≈ 0.189 a −1 0 ), while HSE06 [46], uses the same value for both (ω = 0.11a −1 0 ) and α=0.25 (same as in PBE0). HSE06 shows great performance for small-and medium-gap semiconductors, while underestimating the gaps of wide-gap materials [40,41]. In this work, we focus on the most widely used hybrid functionals, PBE0 and HSE06 and evaluate two recently proposed first-principles techniques to optimize their parameters.

Dielectric dependent hybrid functionals
Marques and coworkers [20] proposed to adopt the inverse of the static dielectric constant as the mixing parameter for PBE0, α = ε −1 ∞ . The same idea is behind all the DDH methods developed later [47,48,49,50]. This relation is obtained by comparing the self-energy Σ of MBPT with the generalized KS equation that is solved in case of hybrid functionals [51,32]. For the large set of materials, investigated in Ref. [20], the DDH method shows an improvement in the calculation of energy gaps with respect to PBE0, reducing the average error from 29.42% to 16.53% as compared to experiment. This average error, however, does not reflect that PBE0 performs considerably well for intermediate-sized band gaps, while it tends to overestimate (underestimate) the gaps of narrow (wide-gap) band-gap materials that are characterized by strong (weak) electronic screening [33]. The application of this method can be somewhat ambiguous, since there are several methods to compute the dielectric constant. A common approach to access the dielectric function is based on linear-response theory in the framework of time-dependent DFT [52]. In this approach, one has to rely on the choice of an exchange-correlation kernel to evaluate the response function. Thereby, the random-phase approximation (RPA) [53] has been shown to perform well [54], due to error cancellation between the underestimation of the gap by, e.g. PBE, and the absence of electron-hole interactions [54]. Unfortunately, for complex materials, also RPA calculations can become expensive. The complexity in determining the dielectric constant motivated the authors of Ref. [20] to propose the DM method as an alternative. The idea of the DM method is to link the mixing parameter to a global estimator of the gap, obtained as the average of a local estimator over the unit cell. Local band-gap estimators, depending on |∇n|/n, have been proposed over the years in different contexts [55,56,57], and the idea of averaging it over the unit cell, to obtain a global estimator, was previously employed in the meta-GGA of Tran and Blaha (TBE09) [58]. The quantity proposed in Ref. [20] has the following form: In PBE0, the mixing parameter α andḡ are connected by the following linear relation, whereas in the case of HSE, the proposed relation is in the fourth power of the estimator: For the band gaps of of the materials set used in Ref. [20], average errors of 14.37% and 10.36% have been found for PBE0-DM and HSE06-DM, respectively, compared to 29.42% for PBE0 and 16.92% for HSE, thus showing a substantial improvement. Another evaluation of the method for a large set of materials can be found in Refs. [41,59]. A further advantage of this method, is that the estimator in Equation (4) is easy to compute. As found in Ref. [20], the functional used to compute the density has negligible influence onḡ, thus this quantity can be calculated at the end of a ground-state calculation performed with a (semi-)local functional. The same method has been investigated also by D. Koller et al. [49], who used a different definition ofḡ and did not find any satisfactory relation betweenḡ and the mixing parameter α. An advanced version of the method, suitable for interfaces, has been proposed more recently [60].
As derived by L. Hedin [61], the GW approximation is a powerful method of MBPT with high predictive power for many classes of materials [62]. This approximation yields the quasiparticle (QP) energies that can be obtained in, e.g., direct and inverse photoemission experiments. It is a common practice to formulate the GW approximation in a perturbative way, based on the mean-field solutions provided by DFT. This scheme, commonly known as single-shot GW or G 0 W 0 , leads to the following expression for the QP energies [63] in which Σ is the electron self-energy that plays the role of a generalized (non-local, energy dependent, and non-Hermitian) exchange-correlation potential, and Z nk is the QP renormalization factor: In these equations, ϕ 0 nk and ε 0 nk are the eigenfunctions and eigenvalues, respectively, of the independentparticle problem solved by employing the xc functional v xc . The solution of Equation (7) implies a dependence of the QP energies on the choice of the underlying KS eigenstates. The latter must provide a reasonable approximation for the QP states such that a perturbative treatment is valid. Otherwise, the single-shot QP corrections alone will not be sufficient to improve band gaps, electronic binding energies, or even the order of bands. Hybrid functionals are typically good starting points for the solution of the QP equation. For semiconductors and insulator, HSE06 is often the favorite choice, while G 0 W 0 on top of PBE0 tends to maintain the inherited overestimation of band gaps [64]. The starting-point dependence can be remedied by adopting self-consistent GW approaches [65,66]. We do not, however, apply them here.
We consider CsPbI 3 in its orthorhombic (space group Pnma) and cubic phase (space group Pm-3m), MAPbI 3 (space group Pm-3m), as well as bulk PbI 2 (space group P-3m1). Following the nomenclature of Ref. [67], we refer to the orthorhombic structure as γ-phase and to the cubic one as α-phase. The unit cells are shown Figure 2. For γ-CsPbI 3 , α-CsPbI 3 , and PbI 2 we adopt the experimental structural parameters [68,69,70]. For α-MAPbI 3 , the hydrogen positions cannot be resolved by X-ray-based crystallographic methods; therefore, we adopt the structure from Ref. [71]. The higher symmetry of the αphase (high-temperature phase) of HaPs is related to a dynamical disorder in the octahedral tilts. For a good theoretical description, without explicitly accounting for dynamical effects, it is convenient to consider lower-symmetry unit cells such as the γ-phase structure [72,73,68]. Therefore, we will compare the results obtained for γ-CsPbI 3 to experimental data, while we consider α-CsPbI 3 and α-MAPbI 3 to investigate the effects of the crystalline structure and the type of cation.  The four systems under investigation are made of lead and iodine, heavy atoms that are responsible for significant SOC effects. In Table 1, the fundamental band gaps computed with PBE and PBE+SOC are shown. The SOC-induced band-gap reduction is dramatic in all systems, ranging from 0.69 eV in PbI 2 to 1.13 eV in α-CsPbI 3 . Comparing the values for α-CsPbI 3 and γ-CsPbI 3 , the sensitivity of the band gap to the crystal symmetry and details of the atomic arrangement in the PbI 3 cage are apparent. Contrarily, the nature of the cation (organic or inorganic) plays a less decisive role in the electronic Ref. [73] band structure in the surrounding of the Fermi level. The differences between the band gaps of α-CsPbI 3 and α-MAPbI 3 , caused by the deviation of ∼0.01Å in the Pb-I bond length, are found to be within 100 meV, for both PBE and PBE+SOC.
For further analysis, we present in the top row of Figure 3 the band structures along selected highsymmetry paths. The bands computed with PBE (gray lines) and PBE+SOC (red lines) are aligned at the valence-band maximum (VBM) to facilitate comparison. In the perovskites structure, SOC effects in the valence bands are minimal, whereas in PbI 2 , especially at Γ, they are responsible for the splitting of several bands. Importantly, the reduction of the gap here is a consequence of the lifting of the degeneracy at the conduction-band minimum (CBm), while in PbI 2 the splitting occurs in the second unoccupied band. Also in this case, the difference between the band structures of the two α-phase perovskites is minimal, in contrast to that of the γ-phase. The energy gaps of the α-phase perovskites are located at R, in the γ-phase at Γ. In the α-phase perovskites, SOC lifts the triple degeneracy of the CBm, by splitting it into a single state (CBm) and a doubly degenerate state. In the γ-phase, SOC is not only responsible for splitting the doubly-degenerate CBm but also for changing the order of the states at Γ. Being found critically important, from this point on, we proceed with the discussion of the electronic properties, always taking SOC into account (unless noted otherwise). For a deeper analysis, we consider the partial density of states (bottom left of Figure 3). We observe (i) that in all four materials, the VBM (CBm) is dominated by iodine (lead) p-orbitals, and (ii) the organic/inorganic cation states do not contribute to the band-gap region. Concerning (i), in all considered materials, the VBM also exhibits a contribution from Pb-s orbitals as evident from Figure S1 of the Supporting Information where the orbitalresolved DOS are presented. The orbital character is also clear from the KS wave functions presented in the bottom right of Figure 3. Notably, the CBm orbitals do not reveal the typical p-orbital shape but rather the relativistic p 1/2 and p 3/2 forms [74], reflecting the strong SOC. The latter and the composition of the band-edge states, which PbI 2 shares with APbI 3 , provide already justification for taking PbI 2 as a representative compound for exploring the electronic structure with different methodologies. Finally, we will demonstrate that all conclusions drawn for precursor PbI 2 can be transferred to the lead-iodide perovskites.

Hybrid functionals and one-shot GW
In this section, we discuss the performance of standard methods (either DFT or MBPT) in predicting the electronic band gaps for the studied materials. In Table 1  no results for hybrid functionals and GW are provided since the organic/inorganic cation has no direct influence on the band gap in 3D lead-iodide perovskites. For the high-temperature phase of CsPbI 3 (αphase), no experimental reference is available, however, as for the other APbI 3 materials, we can expect the gap to be close to the one of γ-CsPbI 3 [73]. All the theoretical results summarized in Table 1, turn out to be significantly lower than the experimental value of γ-CsPbI 3 . As shown for MAPbI 3 [73], the reason is that the local octahedral environment is not symmetric and shows structural similarity to the γ-phase. This is why we will not use the α-phase to evaluate the different methods. By comparing the calculated gap (including SOC) of γ-CsPbI 3 with the experimental references (Table 1), we observe that G 0 W 0 @HSE, giving a value of 1.72 eV, reproduces experiment best. Moreover, PBE0 overestimates the gap by only 130 meV (∼ 7%). HSE and G 0 W 0 @PBE lead to an underestimation by 420 meV (∼ 25%) and 350 meV (∼ 21%), respectively. For PbI 2 , the best agreement with experiment is achieved by G 0 W 0 @PBE, showing a slight overestimation by 135 meV (∼ 5%). Also HSE and G 0 W 0 @HSE perform well with values being by 225 meV (∼ 9%) too low and by 210 meV (∼ 8%) too high, respectively. PBE0 overestimates the gap by 370 meV (∼ 14%). For both materials, G 0 W 0 @PBE0 overshoots slightly more than PBE0 (by 70 meV in PbI 2 and 130 meV in γ-CsPbI 3 ). In Figure 4, the energy gaps from Table 1 (considering SOC) are plotted. The dashed lines indicate how the results obtained with the hybrid functionals and G 0 W 0 on top of them change with respect to the mixing parameter α. It should be noted that comparison between experimental and computed values bear some uncertainties due to, e.g., crystal structure, presence of defects or thermal effects. Moreover, in our calculations, electron-phonon coupling is not taken into account [78]. Therefore, we also focus on comparing different computational methods. From Figure 4, some trends for the two materials can be observed: (i) The ascending order of the gap computed with the different methods is the same. (ii) The

PbI2
-CsPbI3 γ Figure 4: Energy gaps (considering SOC) of PbI 2 (left) and γ-CsPbI 3 (right) as a function of mixing parameter α (see also Table 1). The dashed lines are guides to the eye that reflect the known linear behavior.
differences between the methods are comparable. These observations justify transferability of the methods within materials of this class. Summarizing, G 0 W 0 @HSE appears overall to perform best for predicting the band gap of the materials under investigation. However, it is computationally very demanding, which leads us to explore methods to tune the parameters of the hybrid functionals PBE0 and HSE.

Dielectric dependent hybrid method
The key quantity of the DDH method is the electronic dielectric constant ε ∞ . As discussed in Section 2.1.1, RPA@PBE gives a good estimate. As SOC effects are decisive for describing the electronic properties in HaPs and their precursors, it is important to include SOC also when computing the dielectric response. Such calculations can be computationally expensive for complex materials such as perovskites. Therefore, we mimic the effects of SOC through a scissors operator, taken as the difference between the energy gaps obtained with PBE+SOC and PBE. From the band structures in the top panel of Figure 3, we can immediately observe that the scissors approximation is fully valid for PbI 2 , since the lowest unoccupied band is rigidly shifted, while the same is not true in the case of HaPs. To nevertheless justify and validate the use of the scissors shift within the DDH method, we compare the energy gaps and the mixing parameters obtained when employing PBE+SOC and PBE+scissors, respectively, for the calculation of the dielectric constant within the independent particle (IP) approximation. The results are shown in Table 2. As expected, for PbI 2 , they differ by a small amount only, leading to the same values for the optimized mixing parameter and for the energy gap. In the HaPs, instead, the usage of a scissors shift leads to a bigger dielectric constant. The differences in the mixing parameter α and thus in the energygaps are, however, small for all the HaPs, i.e., α being 0.01 bigger when SOC is accounted for explicitly. The energy gaps of γ-CsPbI 3 and α-CsPbI 3 obtained by the two methods are within 0.05 eV. From this analysis we conclude that, for determining the mixing parameter and further the energy gap, the scissors operator is legitimate. Moreover, through this analysis we confirm that SOC is fundamental for the DDH method, since by adopting PBE only, the final gap is up to 200 meV (for γ-CsPbI 3 ) bigger than the one obtained from PBE+SOC. For the reasons just pointed out above, for tuning the PBE0 mixing parameter with the DDH method  Figure S2 in the Supporting Information. (PBE0-DDH), we compute the dielectric constant employing RPA@PBE, and we make use of a scissors correction to account for SOC effects. For PbI 2 , the energy gap obtained with PBE0-DDH is equal to the one computed with G 0 W 0 @PBE, which underestimates the experimental gap by ∼5%. Moreover, PBE0-DDH improves over PBE0 that shows an overestimation by ∼14% compared to experiment. Also for γ-CsPbI 3 , the result of PBE0-DDH is comparable to that of G 0 W 0 @PBE (difference within 50 meV); however PBE0 can already reproduce the experimental gap well with a difference of only ∼7%, considerably better than PBE0-DDH (difference ∼18%). From the results in Table 2, we observe that the organic/inorganic cation has a smaller effect on the mixing parameter than the crystal structure. This is in accordance with the electronic structures of the different compounds, analyzed in Section 4.1.1. The effect of the structure is, however, moderate with α ranging from 0.12 for α-CsPbI 3 to 0.15 for PbI 2 and γ-CsPbI 3 . By transferring the mixing parameter obtained for PbI 2 to compute the gap of α-CsPbI 3 the error is 130 meV only (Table 2.1.1), and for α-MAPbI 3 we expect it to be even smaller since the α values are even closer.

Density-based mixing method
To the best of our knowledge, the DM method was not applied to HaPs before; neither were SOC effects taken into account. As we see from Table 3, the effect of SOC onḡ is negligible. We attribute this finding to the fact thatḡ is obtained by an average over the unit cell (Equation (4)), as well as to the minor differences between the total electron densities of PBE and PBE+SOC, as SOC mainly affects the conduction bands.
The DM method applied to PBE0 (PBE0-DM) leads to a larger overestimation of the gaps with respect to PBE0. In contrast, the DM method in combination with HSE (HSE-DM) improves over HSE for the gaps of both materials. For PbI 2 , the overestimation is only 40 meV (∼1.6%), making this method the best choice overall. For γ-CsPbI 3 , HSE-DM underestimates the gap by 130 meV (∼8%) which is comparable with PBE0 (difference of 130 meV). From Table 3, we observe that the parameters obtained for the four materials are in a small range (from 0.28 to 0.30 for PBE0-DM and from 0.36 to 0.38 for HSE-DM), a finding that, again would justify the 4.3 Screening parameter ω and G 0 W 0 calculations  Figure S2 in the Supporting Information. use of the parameters found for PbI 2 for the other compounds. Moreover, the organic/inorganic cation has more influence than the crystalline structure, since the difference in the parameters between α-CsPbI 3 and α-MAPbI 3 is bigger than between α-CsPbI 3 and γ-CsPbI 3 . As this difference is, however small, we again conclude that the organic cation does not have a major impact on the total electron density and its gradient.

Screening parameter ω and G 0 W 0 calculations
In this section, we investigate the impact of ω, the second parameter of HSE, on the band gap of PbI 2 . Moreover, we discuss the starting-point dependence of the G 0 W 0 gaps and how they are influenced by the choice of the parameters α and ω. To recall, in the limit α = 0, both PBE0 and HSE are identical to PBE. Moreover, by construction, HSE satisfies two limits: (i) For ω=0, it coincides with PBE0, and (ii) for ω → ∞ it becomes equivalent to PBE as follows from Equation (2). The left-top panel of Figure 5 visualizes the linear dependence of the PBE0 and HSE band gaps as a function of α (the symbols indicate the values from Tables 1, 2, and 3). The color scheme shows how the gaps change with respect to ω for different values of α. The latter can be seen more explicitly in the middle-top panel, in which for specific values of α (0.25 and the ones obtained from the tuning methods), the results for different values of ω are plotted. The α-dependence gets less steep for bigger values of ω, which reflects the definition of HSE, since for ω → ∞ all curves go asymptotically to the PBE value. These results are combined in the right-top panel of Figure 5, in which the color map encodes the energy gaps as a function of α and ω. The isolines represent those combinations that reproduce the experimental, DDH, and DM band gaps, respectively (see also Tables 1, 2, and 3). The DDH method provides a mixing parameter only for PBE0. Through the PBE0-DDH isoline we extract the α value for HSE (ω = 0.11 a −1 0 ) which is α = 0.28. We now perform the same analysis for G 0 W 0 . The results are included in Tables 2 and 3 and, together  with the values from Table 1, are plotted in the bottom of Figure 5. Also in this case, the gaps change linearly with respect to α, but in a much narrower range. Likewise, the middle panel shows the dependence of the results on ω for selected values of α; while the color map on the right highlights the combined dependence on both parameters. From this analysis, one can deduce that the dependence of the energy gaps on α and ω in G 0 W 0 is similar to that at the DFT level, but restricted to a smaller region. While on the DFT side, the gap changes in a range between 1.2 to 4.4 eV, on the G 0 W 0 side, only between 2.2 to 3.8 eV. Overall, the starting-point dependence with respect to the hybrid parametrization is significant, even if the range is smaller.

Conclusions
We have systematically investigated the performances of several methods of DFT and MBPT to compute the band gaps of PbI 2 and APbI 3 . To avoid mere benchmarking against experimental results, which comes with uncertainties, we also consider comparison between different theoretical approaches. We have verified that SOC is fundamental, irrespective of the method. This makes the calculations numerically challenging. Overall, G 0 W 0 @HSE performs overall best but is computationally remarkably expensive. Among the approaches that involve the tuning of parameters, PBE0-DDH gives results comparable to G 0 W 0 @PBE, but is also cpu-intensive, partly owing to SOC. Instead, in the DM methods, SOC has no major effect on the estimatorḡ. Based on the PBE electron density to computeḡ, makes this method accessible for complex materials such as APbI 3 . At the DFT level, HSE-DM gives overall the best estimation of the energy gap, with a performance comparable with G 0 W 0 @HSE. Both tuning methods lead to similar mixing parameters for the four investigated materials, with the biggest difference found for PBE0-DDH, where α=0.12 for α-CsPbI 3 and α=0.15 for γ-CsPbI 3 and PbI 2 . This means that the atomic species Pb and I have larger effects on the screening parameter and on the density gradient (|∇n|/n) than the crystalline structure and the (in)organic cation. This reflects the fact that the band gap region is dominated by p-Pb and p-I states, which also explains the transferability of the method from PbI 2 to APbI 3 . Investigating the dependency of the PbI 2 energy gap on different combinations of α and ω, we find that both parameters have significant impact in both hybrid functionals. This also applies when G 0 W 0 is computed on top, however, less pronounced. Our results suggest that one can transfer our findings from the studied precursors to the respective HaPs compounds.

Appendix A: Computational details
All the DFT and MBPT calculations are performed using the full-potential all-electron computer package exciting [21]. The code employs the linearized augmented planewave plus local orbitals ((L)APW+lo) basis [79,80,81,82] to expand the Kohn-Sham wavefunctions. In the (L)APW+lo method, the unit cell is partitioned into two regions: muffin-tin (MT) spheres around the atomic nuclei of radius R MT , in which the basis functions are atomic-like functions, and the interstitial region between the MT spheres, in which the basis functions are planewaves. The MT radii are chosen to be R Pb = 2.9 a 0 , R I = 2.9 a 0 , R Cs = 2.9 a 0 , R C = 1.1 a 0 , R N = 1.0 a 0 , and R H = 0.9 a 0 . The large sphere sizes of Pb and I avoid any core leakage, which is highly important when SOC is taken into account. The basis-set size is determined by the number of planewaves used. As the planewave cutoff, G max , depends on the muffin-tin radii, it is common practice to express the cutoff parameter as the dimensionless product R MT · G max , where R MT is the radius of the smallest sphere. In the DFT calculations (with PBE and hybrid functionals), k mesh and R MT G max are chosen such to guarantee a numerical precision of the band gap within 20 meV for PbI 2 , α-CsPbI 3 , and α-MAPbI 3 and within 50 meV for γ-CsPbI 3 . All G 0 W 0 band gaps are converged up to 100 meV. In Table 4, the values of G max and R MT G max used for the investigated systems are summarized. The electronic structure is computed employing the xc-functionals PBE, PBE0, and HSE as well as G 0 W 0 on top of them. The k-grids are shown in Table 5. The convergence of the energy gaps with PBE0 turned out slower than with PBE, HSE, and G 0 W 0 . Nevertheless, for PbI 2 , the same grid is used in HSE and PBE0 calculations such to investigate the asymptotic behavior of HSE → PBE0 for ω → 0. SOC is included if specified. In case of PBE, it is treated via the second variational scheme [83], where in the corresponding term of the Hamiltonian, H SO ∝ ∂v 0 (r)/r∂r (σ · L), v 0 is the spherical component of the effective potential. For the considered set of materials, the SOC treatment within the second-variation scheme requires almost all KS functions available to achieve precise results. The corresponding number of the empty states is significantly higher than that required to converge the hybrid and G 0 W 0 calculations. Therefore, to maintain the same level of precision with and without SOC, we use in all these calculations the higher number of empty states, i.e., those required for SOC, which are 300 for PbI 2 , 1000 for γ-CsPbI 3 , 460 α-CsPbI 3 , and 1000 for α-MAPbI 3 . Calculations with hybrid functionals, consist of a nested loop [84,44]. In the outer loop, the non-local exchange is computed by employing a mixed-product basis [85,86]; in the inner loop, the generalized KS matrix equation is self-consistently solved by updating in each step only the local part of the effective potential. SOC is included self-consistently in the inner loop, again through second variation. Since the gradient of a non local potential is not trivial to compute, the effective potential in H SO employs PBE [87]. This is justified since this contribution is small.  Table 5: k-mesh used for the calculation of the electronic structure and the dielectric constant (with and without SOC). In the G 0 W 0 calculations, always the same grid as in the underlying DFT calculation is employed.
Regarding HSE, the SR part of the HF exchange is computed by the difference of the total HF exchange and the HF-LR contribution as done by Schlipf and coworkers [44]. The authors show that by doing so, in the limit of q + G → 0, the singularities that would occur in the full HF exchange and HF-LR cancel out by Taylor expanding the exponent coming from the LR contribution, thus leading to the constant π/ω 2 . This approximation was also used by us in Ref. [88]. However, this expression implies numerical difficulties to study the asymptotic behavior of HSE for ω → 0 caused by the necessity to employ extremely dense k-point grids. To overcome this problem, we have derived an alternative expression to compute the SR part of the bare Coulomb potential in the limit of q + G → 0. More details are given in Appendix B. In our G 0 W 0 calculations [89], SOC effects are incorporated by applying the second-variation procedure after the QP corrections to the Kohn-Sham energies have been computed. This is an approximate but computationally efficient way for treating SOC in GW . As discussed in Ref. [90], this scheme may lead to wrong predictions in materials where SOC induces band-inversion such as topological insulators. In other systems like the ones investigated here, it is expected to produce results consistent with those from a more rigorous treatment.
Calculations of the dielectric functions are performed using both the independent-particle approximation and the RPA kernel [91]. To evaluate the dielectric constants, the parameters are chosen such that the mixing parameter α of PBE0 is determined with a precision of 5 · 10 −3 . The k meshes used to compute ε ∞ are shown in Table 5. In the RPA calculations, the parameter gqmax that governs local-field effects is chosen 2 for all systems. The number of empty states is 10 for PbI 2 , 20 for α-CsPbI 3 , and α-MAPbI 3 and 100 for γ-CsPbI 3 . These values are reduced with respect to those used for the electronic structure since SOC is not taken into account. The parameterḡ of the DM method (Equation (4)) is computed from the PBE density. We verified that, using the parameters of Tables 4 and 5, α is determined with a precision of 5 · 10 −3 . The crystal structures displayed in Figure 2 and the KS wavefunctions in Figure 3 have been produced with the software package VESTA [92].
Appendix B: Treatment of the singularity in the Coulomb potential The Fourier transform of the short-range part of the Coulomb potential employed in the HSE functional is given by It is easy to show [44] that lim q+G→0 v C (q + G) = π/ω 2 .
However, one quickly faces numerical instabilities when studying the parametric dependence of HSE results for small values of ω. Therefore, the case of q + G → 0 and ω → 0 requires a special treatment.
In this work, we estimate the limit by isotropic averaging in a small region around the Γ point of the Brillouin zone (BZ). Our goal is to compute the integral where V k = Ω BZ /N k is a small volume surrounding Γ, Ω BZ is the volume of the BZ, and N k is the total number of the k-points used to sample the BZ. To be able to compute the integral analytically, we replace V k with a sphere of radius Using the notation q ≡ q + G, q ≡ |q + G|, and β = 1/4ω 2 , Equation (11) turns into It can be shown, that for big values of ω (small values of β), the integral correctly recovers the value π/ω 2 from Ref. [44]. Our approach allows us to study the parametric dependence of the HSE results on both α and ω in the entire parameter space.

Supporting Information
Supporting Information is available from the Wiley Online Library or from the author.