Introduction

There is an ever-increasing need for developing new materials with novel electronic functionalities in our technology-driven modern society. To this end, a recent focus has been to integrate computational and experimental techniques for the purpose of obtaining robust, predictive tools for gaining insights into structure–property relationships in functional materials in order to reduce the time and cost in the materials discovery process. In this connection, it is important to expand the search space of materials to include compounds containing elements from the last rows of the periodic table, which present an especially rich playground for hosting novel functionalities that arise through the interplay of effects of spin–orbit coupling, strong electron correlations and interactions of highly localized f-electrons with free-electron-like conduction electrons. Although widely known for nuclear energy applications, actinides and lanthanides have recently attracted attention also in connection with superconductivity,1 magnetism,2 and Kondo physics3 as well as energy-related applications.4,5

While a “1-to-1” approach (one system—one calculation—one experiment) has been the traditional strategy in lanthanide and actinide research, this approach is quite inefficient due to the large experimental cost of handling high-Z materials. Moreover, f-electron systems have proven notoriously hard to model within the first-principles density functional theory (DFT)6,7 framework, while dynamical mean field theory (DMFT) methods8,9,10 invoke semi-empirical correlation parameter U, which reduces their predictive power. Here, in order to address these challenges, we discuss a novel predictive strategy, which is based on adopting an integrated data and theory-driven approach. In this connection, we have built the f-electron structure database (fESD, http://correlatedmaterials-lanl.org) that combines experimental and DFT simulated data using state-of-the-art algorithms11,12,13 to train our machine learning models. We will show how the analysis of our database reveals new insights into the intricate electronic and structural properties of f-electron systems.

A unique feature of the fESD is that it incorporates high-quality simulated electronic structure data which capture the effects of electron-correlations and spin–orbit coupling in an all-electron environment. Many popular databases (e.g., AFLOW,14 Materials Project,15 Organic Materials Database16) use pseudo-potentials and plane-wave-based techniques for DFT simulations. In contrast, we use full potential, all-electron, relativistic linearized-augmented-plane-wave (LAPW) results,17,18 which is important for obtaining a predictive model for strongly correlated materials, computational cost notwithstanding.

Results and discussions

Database design and collection of data

Our database is primarily an electronic-structure database with query tools that allow searches of compounds with desired crystal and/or band structures. In addition, the database contains DFT simulated data with a variety of query capabilities. The database management system (DMS) is built on MySql and Java Development platforms. Additionally, a number of Python-based scripts are used for data parsing at various levels (e.g., crystallographic information files (CIF), DFT output). Work toward designing an application-programming interface (API), which will integrate various query tools and provide a Python-based library of functions for programmed access to the database, is currently in progress.

Inorganic Crystal Structure Database (ICSD)19 and Crystallographic Open Database (COD)20 are our initial sources, which currently contain 188,000 and 364,331 CIF files, respectively. There are 54,465 structure files in ICSD and 27,502 in COD with compounds containing lanthanide and actinide elements (e.g., outer shell f-orbital electrons). We started by downloading structure information such as crystal structures, lattice parameters, and atomic positions from these CIF files into our database. At present, the scope of our database is limited to ground state DFT-based electronic structures using available crystal structure data without sensitivity to pressure and temperature dependencies or the specifics of how the data was acquired (e.g., theoretically or experimentally). A richer labeling with more detailed attributes of the data, will be undertaken in future extensions of the database.

Data cleaning, verification, and crystal structure determination

While restructuring the data format, we developed a number of query tools with new search features which, to the best of our knowledge, are not available in other existing databases.14,15,16 In particular, we provide extended search capabilities such as space-group or lattice-system-based search, and search for a superstructure of desired symmetry. As a part of this development process, we performed the “data cleaning” step by employing supervised learning tools. This resolved ambiguities in CIF files resulting from missing information, and corrected the structure files needed for generating reliable electronic structure data, which are the core content of our f-electron database. The details of the generation of electronic and crystal structure data are described in the Methods section.

In our SQL database, some crystallographic features (e.g., lattice parameters, space groups) are missing or ambiguous in many compounds because of incomplete original CIF files that were parsed to create the database. These compounds cause convergence problems or inaccurate results in DFT simulations. We resolve this problem by adopting a machine learning approach for verifying and cleaning the incorrect information. Currently, our database hosts crystal data for approximately 82,000 f-electron compounds, out of which 8711 compounds contain missing or incorrect lattice systems. In this connection, we discuss three different supervised machine-learning algorithms to predict one of the seven possible “lattice systems” as follows. (1) Logistic Regression (LR),12 where one predicts categorical target variables using a logistic function involving linear combinations of feature values. (2) K-Nearest-Neighbor (KNN) algorithm,11 which predicts category of an unknown instance based on K “similar” instances. We use Euclidean distance as a measure of similarity and pick the majority of the three (K = 3) most similar materials to predict the lattice parameters of an incomplete system. (3) Multilayer Perceptron (MLP)13 in which a forward feeding network of perceptrons is used, as opposed to a single perceptron that is equivalent to a LR model. We use three layers of perceptrons with eight perceptrons in each layer.

The models were trained on a set of 275,926 labeled instances from the COD with accurate lattice constants (a, b, c, α, β, γ), volume, and the space group as training features. We ran a 10-fold cross-validation, and calculated the prediction accuracy, precision, and recall. Using the miss and hit between the predicted-class and true-class outcome, we constructed a representation widely known as the confusion (or error) matrix in the machine-learning community. As shown in Fig. 1, the confusion matrices for the three methods evaluated indicate clearly the superiority of the MLP algorithm, which essentially is a neural network approach. MLP yields an accuracy of 99.1% in predicting the crystal system, see Table 1. Although computationally more complex and expensive, we employed this algorithm for cleaning and obtaining high-fidelity structure information for our database.

Fig. 1
figure 1

Confusion matrices for different machine-learning algorithms. Confusion matrices for three methods: a Multilayer Perceptron (MLP), b K-Nearest-Neighbor (KNN), and c Logistic Regression (LR). MLP, which is a neural network approach, captures the correct crystal system with maximum accuracy

Table 1 Benchmarking the performance of three machine-learning techniques (LR, KNN, and MLP) discussed in the text

However, all three models exhibit some off-diagonal non-zero values in the confusion matrix indicating errors in prediction. Interestingly, we found such instances to be restricted to three unique classes, namely, trigonal, hexagonal, and rhombohedral, see Fig. 1a. In the crystallographic convention,21 there are two equivalent systems, i.e., the crystal system and the lattice system. Trigonal and hexagonal lattices belong to the crystal convention, while the rhombohedral and hexagonal are defined in the lattice system. We resolved such ambiguities of definition in the original CIF files by deploying our neural-network-based MLP tool (Fig. 1a), and thus obtained cleaner datasets for further ab initio simulations. In this way, we correctly predicted the missing crystal systems for 8711 compounds with the distribution given in Table 2. Note that our main purpose in determining crystal system information was to illustrate the viability of our machine learning tools. However, the correct coding of the crystal system and the space group for materials in the database can be of interest in searching the database for materials belonging to a given space group or a crystal system. For example, certain space groups and crystal systems have proven useful in successful searches of viable topological materials.22 Moreover, useful physical insights can often be obtained by comparing and contrasting the evolution of properties of materials within a space group or a given crystal system as well as by comparisons across different or related space groups or crystal systems.

Table 2 Distribution of the crystal systems, which were predicted correctly by using MLP tools

Perovskite heterostructure prediction

In the recent years, there has been considerable theoretical and experimental interest in creating layered compounds and superstructures by design due to their enormous potential for achieving various emergent functionalities (e.g., superconductivity,23,24 multiferroics25,26). A common strategy employed in synthesis efforts is to attempt the creation of a new superstructure with desired electronic properties by overlaying two or more compounds with different chemical compositions but similar crystal symmetries. In this connection, we investigated a set of perovskite compounds in our database with the goal of identifying compounds that would be most likely to form stable superstructures. We used lattice information and chemical intuition related to “ionic radii and Goldschmidt tolerance factors”27 as our primary screening tools for the initial narrowing-down of the search space for perovskite heterostructures. We then performed a data-mining query, and identified sets of paired superstructure combinations with the highest likelihood of forming superstructures. The next stage would then involve a more careful computational analysis of the small number of promising candidate systems so identified before it will be appropriate to encourage an experimental synthesis and validation cycle. Generally, we should expect that the prediction–synthesis–validation loop may well require more than one iteration to successfully identify a viable new material.

Figure 2 shows the results of a compatibility test for the ordering of superstructures of a number of double perovskites with chemical formula AA′BB′C3\({\mathrm{C}}_3^\prime\). Here, A/A′-site cations are either rare earths or actinides, B/B′ cations are transition metal elements, and C/C′ anions are mostly oxygen or halogens. Using the ionic radii28,29 of the two existing single perovskites ABC3 and A′B′\({\mathrm{C}}_3^\prime\) (see the inset in Fig. 2b), we assess the geometric stability of a possible double perovskite AA′BB′C3\({\mathrm{C}}_3^\prime\) by calculating the Goldschmidt tolerance factor,27,30 \(t\) = \({\textstyle{{r_A + r_C} \over {\sqrt 2 \left( {r_B + r_C} \right)}}}\), where rA, rB, and rC are the average ionic radii of (A, A′), (B, B′), and (C, C′), respectively; here, we have neglected other contributing factors such as the effects of octahedral tilting. We considered both ordered and disordered compounds on B and B′ sites, and following the method described by Vasala et al.,30 we treat all compounds with the same chemical formulae regardless of the ordering at the B-sites. Note that Goldschmidt tolerance factor is a good measure for testing the stability of AA′BB′C3\({\mathrm{C}}_3^\prime\) type perovskite compounds, where B or B′-site cations can be ordered or disordered.30 Considering each data point as a paired combination of a double-perovskite superstructure, we plot the associated tolerance factor t and the lattice-parameter-ratios \({\textstyle{{a1} \over {a2}}}\) and \({\textstyle{{c1} \over {c2}}}\) in Fig. 2a, where (a1, c1) and (a2, c2) are lattice parameters of ABC3 and A′B′\({\mathrm{C}}_3^\prime\) perovskites, respectively. Clustering of the data in Fig. 2a reflects the correlation between the type of anion/cation involved and the space group of the individual perovskites with the tolerance and lattice parameters. High tolerance factors can be achieved when \({\textstyle{{a1} \over {a2}}}\) and \({\textstyle{{c1} \over {c2}}}\) are close to unity, i.e., when the difference in ionic radii of A (rA − rA) or B (rB − rB) site cations is small. We also see the important role of space groups here such as two cubic perovskites usually result in a stable double perovskite in contrast to the cases where the space groups are mixed. This is to be expected because similar crystal structures of the two contributing perovskites lead to lesser cation size mismatch or smaller differences in ionic radii of the cations of individual sites and result in higher tolerance factors. This is why we would expect the most abundant double-perovskite systems to be formed from two cubic systems, as demonstrated in Fig. 2b. We see a linear trend in Fig. 2b where the data points are plotted in the \({\textstyle{{a1} \over {a2}}}\) and \({\textstyle{{c1} \over {c2}}}\) plane. The central line indicates that the most abundant double-perovskite systems would be formed from ABC3 and A′B′\({\mathrm{C}}_3^\prime\) cubic systems. Other trend lines in Fig. 2b show the cases in which mixed systems can be formed from cubic and other (non-cubic) crystal structures. For this reason, we allowed a mixture of anions to identify various possibilities in the stability trend lines. We generally see the trend of high stability (t) when both C and C′ anions are oxygen or hydrogen compared to the case when both anions are halogens (see Fig. 2a). This indicates that the most stable structures are formed when both anions are identical, which should be contrasted with the situation where anions are allowed to mix. In this way, we were able to cluster possible candidate perovskites with respect to anion mixing and distinguish them in terms of their stability.

Fig. 2
figure 2

Superstructure compatibility test for double-perovskites AA′BB′C3\({\mathrm{C}}_3^\prime\). Each data point is a paired combination of two single perovskite systems as shown in the inset to (b). Different colors imply different combinations of C and C′ as shown in the legend in (b). Space group for each C/C′ type is given in the parenthesis in the legend in (b)

We turn now to address the effects of octahedral distortion of the ideal cubic structure, which is very common in the perovskites.30,31 Octahedral distortion induces changes in the B-O-B' bond angle away from 180°, and affects the B-O and B-O' bond lengths to produce elongation or contraction along the c-axis. Figure 3a shows that there is a linear trend between the tolerance factor and the average value of the B-O-B' bond angle. Systems with higher 〈B-O-B'〉 bond angles possess smaller octahedral distortions, yielding higher tolerance factors that favor stability. An inspection of our data indicates that 〈B-O-B'〉 bond angle gets closer to being linear when the B/B' cation has a small number of d electrons. However, B and B' sites with Jahn–Teller active cations such as Mn3+ and Cu2+ lead to higher distortions with the bond angle 〈B-O-B'〉 ≈ 142°. Majority of perovskite structures we considered showed an elongation of the c-axis, which results in shorter average B-O and B'-O bond lengths, see Fig. 3b. We thus adduce that, for a given type of A-cation, the average bond length decreases linearly as the amplitude of the distortion decreases. This linear relationship, however, may not hold when two different A-site cations are involved, which accounts for the scatter in the data of Fig. 3b.

Fig. 3
figure 3

Descriptors for octahedral distortion in perovskite superstructures. a Tolerance factor (t) is plotted against the average bond angle 〈B-O-B'〉. t is seen to correlate with decreased octahedral distortion as the bond angle 〈 B-O-B'〉 gets closer to linear. b A plot of the average 〈B–O〉 bond length against the bond angle 〈B-O-B'〉. The bond length is seen to correlate with decreasing distortion (increasing bond angle). For various B' cations, average bond angles and bond lengths have been obtained by considering all six neighboring B cations. The inset in (b) shows neighboring BO6 and B'O6 octahedra where B, B', and O ions are marked with grey, yellow, and red spheres, respectively

It is also important to consider the effects of ordering of the B/B'-site cations. When the B and B' cations occupy the correct sites, we take the arrangement to be ordered. Following Vasala et al.,30 we consider two types of cation disorder, namely, anti-site (AS) defect and anti-phase boundary (APB). Here an AS defect is defined as the interchange of B and B' cations, while APB involves the boundary of two ordered domains where occupancies of B and B' sites are reversed. Specifically, we identified 105 ordered rock-salt and 50 disordered double perovskite structures of the form A2BB'O6 in the ICSD database. Figure 4a shows the difference in B and B' cation oxidation state (\({\mathrm{\Delta }}Z_B\) = \(\left| {Z_{B^{\prime}} - Z_{B} } \right|\)) as a function of the corresponding ionic radius difference (\({\mathrm{\Delta }}r_B\) = \(\left| {r_{B^{\prime}} - r_{B} } \right|\)). We observe that smaller ΔZB values favor the disordered arrangement. This makes physical sense because B and B' cations in this case are chemically similar, and therefore, these cations will tend to occupy various sites interchangeably. We can also understand this result in terms of the competing effects of electrostatic repulsion and entropy. When B and B' cations are similar, the system will have a tendency to become disordered to increase its entropy. On the other hand, in ordered systems which tend to have higher ΔZB values (top red dots in Fig. 4a), highly charged B' cations always have less charged B cations as nearest neighbors, so that the Coulomb energy in the ordered state is lowered compared to the disordered state.30,32 The difference in ionic radii, ΔrB, also provides a good descriptor for obtaining a handle on the viability of ordered vs disordered arrangements. It is known32,33 that a larger ΔrB enhances lattice strain, which favors ordered phases.

Fig. 4
figure 4

Descriptors of cation ordering in double perovskites. a Difference in B and B' cation oxidation states, \({\mathrm{\Delta }}Z_B\) = \(\left| {Z_{B^{\prime}} - Z_{B} } \right|\), is plotted against the corresponding difference in ionic radii, \({\mathrm{\Delta }}r_B\) = \(\left| {r_{B^{\prime}} - r_{B} } \right|\). A combination of smaller values of ΔrB and ΔZB is seen to favor disordered arrangements, while a higher value of these descriptors favors ordered arrangements. b Tolerance factor (t) vs ground state energy per unit cell volume. Ground state energy is seen to be higher for disordered cases even though the associated tolerance factors can be quite large

Keeping the aforementioned trends involving ΔZB and ΔrB in mind, we infer that perovskites such as La2CoFeO6ZB = 0 and ΔrB = 0.035 Å) are likely to be stable even in disordered arrangements.34 In contrast, a compound like La2NaIrO6ZB = 4 and ΔrB = 0.45 Å) will be expected to be stable only in ordered configurations35 (see Table 3). In this way, ΔZB and ΔrB can be good descriptors for identifying ordered vs disordered superstructures. Figure 4a shows that systems with ΔZB = 0 and ΔrB ≤ 0.1 Å are always disordered, while systems with ΔZB = 4 and ΔrB ≥ 0.1 Å are always ordered. Coexistence of ordered and disordered structures for ΔZB = 2 presumably reflects the competition between Coulombic and entropic effects. In order to gain a deeper understanding of the underlying energetics, we have computed ground state energies of various structures as shown in Fig. 4b. While the disordered configurations can achieve high tolerance factors (t > 0.85), they are seen to be energetically less favorable compared to the corresponding ordered systems in most cases, consistent with our analysis above based on the descriptors ΔZB and ΔrB. Notably, however, some ordered systems exhibit higher energy (>−75 Ry/a.u.3) and stability (t > 0.85), which is possible when Coulombic and entropic effects become comparable.

Table 3 Predicted double perovskites obtained via our superstructure compatibility test

Our analysis discussed above based on using ionic radii, bond angles, oxidation states, lattice parameter ratios, ground state energies, and Goldschmidt tolerance factors as descriptors correctly predicts several existing double perovskites, see Table 3 (top eight rows). These compounds have high tolerance factors and have indeed been synthesized at ambient pressure.30,34,35,36,37,38,39,40,41,42,43 We also predict four new double perovskites with high tolerance factors (bottom four rows in Table 3), which will be interesting to explore experimentally.

Electronic structure calculations and orbital-resolved electronic structure search

Interpretation of specific features in band structures can be complicated in the presence of various competing effects, particularly in the strongly-correlated f-electron systems of interest to this study. In this connection, contributions of various atomic sites and different orbital angular momentum channels to the Bloch wave functions associated with the band structure near the Fermi energy can provide key information for understanding the electronic behavior of the material. Keeping this in mind, we have implemented a unique search capability based on identifying different atomic orbitals contributing to various bands at k-points along the high-symmetry directions. As an example, Fig. 5 shows an advanced band structure query on Uranium Nitride (UN) that has dominant f-orbital contributions near the Fermi energy. Such information is helpful for gaining insight into the complex electronic, magnetic, and optical properties of f-electron materials, and for developing machine-learning models and predictive tools. Note that local dynamical correlation effects are missing in the simulated DFT results in our database. We address this aspect in the following section.

Fig. 5
figure 5

Web interface for band structures. As an example, we show the band structure of UN displaying f-orbital characters near the Fermi level. The different drop-down boxes at the top show that an energy window of −7 to 7 eV is chosen for uranium and that the f orbital characters are being considered

A search for strongly correlated actinides

We have computed frequency-dependent hybridization functions for all the compounds in our database to supplement our ground state DFT calculations. Our earlier work44 has shown that hybridization is a good descriptor for detecting localization trends, although an accurate description of localization phenomenon will require self-consistent DMFT calculations.9,10 However, our findings44 indicate that a high-throughput, local Green’s function based hybridization can at least qualitatively capture physically-relevant trends in the localization of the f-states. Our analysis with the present f-electron database using ≈350 Ce-based binary compounds shows that the maximum hybridization value near the Fermi energy decreases with increasing lattice spacing.

We comment briefly on technical aspects of our frequency-dependent hybridization function computations. In the Anderson impurity model,45 when an impurity electron is immersed in a sea of itinerant electrons, it hybridizes with the Bloch states of the surrounding electron sea.46 The interaction of the impurity cluster, which in our case are the f-electrons, with its surrounding bath can be monitored via the energy-dependent hybridization function Δ(E). Using Dyson equation, we can formulate an expression for Δ(E) as9,44,47

$$G_0^{ - 1} = \left( {E - E^{QI}} \right) - {\mathrm{\Delta }}(E) = \left( {E - H} \right) - {\mathrm{\Delta }}(E)$$
(1)

where G0 is the site projected Green’s function obtained from DFT and H is the hybridization-free impurity Hamiltonian corresponding to the quantum impurity energy EQI. A bigger overlap of f-states with the surrounding itinerant electrons will produce a larger hybridization function, which can thus serve as a descriptor for the localization/delocalization of f-electrons.

Figure 6 presents the calculated hybridization functions Δ(E) for 4f and 5f monopnictides LnX and AnX (Ln = Ce, An = Th, U, Pu, and X = N, P, As, Sb, Bi) with space group \(Fm\bar 3m\) and a structure similar to NaCl. We reproduce previously reported results on Δ(E) for CeX44 and the established trends in similar compounds,48,49 indicating the high fidelity of our calculations. Our analysis further identifies a number of similarities and trends in Δ(E) features in various 4f and 5f monopnictides as follows. (1) A distinct peak in Δ(E) is found between 2 and 3.5 eV below the Fermi level in both LnN and AnN. The size of this peak decreases monotonically (with small variations) as we scan from the top to the bottom of group 15 of the periodic table. The reduced width of Δ(E) peak in Sb and Bi compounds indicates a more localized nature of f-electrons, whereas a sharp corresponding peak in N compounds with two times larger area under the curve shows a greater degree of delocalization. (2) Figure 7 (top panel) shows that there is an overall increase in the value of the hybridization function at the Fermi energy in going from 4f to 5f monopnictides, reflecting the more extended nature of the 5f shells compared to their 4f counterparts. (3) Within the 4f and 5f shells, if we move from the left to the right of the periodic table, we see an increase in localization with an increase in f-electron count. (4) The bottom panel in Fig. 7 shows an anti-correlation between the degree of localization and lattice constant for both the LnX and AnX series. And, finally, (5) moving from the top to the bottom of group 15, we observe that a higher unit-cell volume produces a smaller Δ(EF) value. This can be understood from the fact that a larger ligand distance will reduce the interaction of the f-electron atom with the neighboring atoms, and thus yield a more localized band structure. Along this line, the 4f monopnictides exhibit higher volumes overall and increased localization compared to the corresponding 5f compounds.

Fig. 6
figure 6

Trends in f-orbital hybridization function. The calculated hybridization function Δ(E) for the series of 4f and 5f monopnictides, LnX and AnX (Ln = Ce, An = Th, U, Pu, and X = N, P, As, Sb, Bi). The ICSD identification numbers are shown in brackets

Fig. 7
figure 7

Volume vs hybridization function at the Fermi level. The value of hybridization function at the Fermi level, Δ(EF) (top), and the volume per unit cell (bottom) for the series of monopnictides considered in Fig. 6, shows an anti-correlation between the degree of localization and lattice constant

Although we have focused on identifying trends in terms of the hybridization functions, such an analysis could also be carried out based on an examination of features in DOS and/or band structures in conjunction with lattice symmetry information. We emphasize that any predictive search must first narrow down the search space by using domain-specific knowledge by using hybridization and/or other appropriate descriptors as a prelude to deploying more sophisticated techniques (e.g., charge self-consistent DMFT and/or structure relaxation). In this way, practical schemes can be obtained for robust theory-guided discovery of new correlated functional materials.

In conclusion, we have presented the design of our newly developed f-electron structure database and discussed the associated high-throughput tools for analyzing structural and electronic properties of materials, including query tools for identifying orbital characteristics of electronic states near the Fermi energy. fESD currently contains data on about 80,000 f-electron compounds. All structure data in fESD have been cleaned and corrected using machine-learning tools. We have included DFT-based high-quality electronic structure data using all-electron self-consistent computations on the large number of compounds in fESD as its core content. Potential for materials discovery via fESD is illustrated by considering a stability search of superstructures composed of two perovskite compounds. Using ionic radii, bond angles, oxidation states, ground state total energy, lattice parameter ratios, and Goldschmidt tolerance factors as descriptors, we show how stable double-perovskite superstructures can be predicted, and how insight into the roles of anion types and space groups involved for achieving higher structural stability can be obtained. We also demonstrate that our database can be used to gain a handle on the strongly correlated aspects of f-electron materials. For this purpose, we have carried out calculations of hybridization functions to supplement our DFT results that provide a descriptor for the strength and nature of interactions between the localized f-electrons and itinerant conduction states. In this way, a number of trends in the lanthanide and actinide based series of compounds are identified, including an anti-correlation between the hybridization strength at the Fermi energy and the volume of the unit cell. fESD is thus well equipped to spur the discovery of next generation f-electron-based materials with novel functionalities.

Methods

Generation of electronic and crystal structure data

Electronic ground states were obtained within the DFT framework using full-potential LAPW scheme as implemented in the WIEN2k package.18 Exchange-correlation effects were accounted for by using the generalized gradient approximation (GGA) with Perdew–Burke–Ernzerhof (PBE) functional.50 Calculations start with the transformation of CIF files into WIEN2k input files using the tool “cif2struct”. During this preparation stage, the space group symmetry of the system is correctly captured based on all atomic positions. Computations were carried out using a large plane wave cut-off (RKmax = 9.0), and a 10 × 10 × 10 k-mesh in the first Brillouin zone (BZ) in order to obtain well-converged energies and electronic structures.

Band structure calculations were carried out along the standard high-symmetry lines in the BZs. Density of states (DOS) was obtained over a uniform energy mesh and resolved into the basis of real orbitals. Hybridization functions Δ(E) for 4f and 5f electrons were calculated by solving the local Green’s function on the real axis.9,10 Since we are only interested in total Δ(E), we used a real harmonic basis with spin but without any symmetry to project the local Green’s function. In order to obtain an accurate estimate of Δ(E), the energy window was set to −20 to +10 eV around the Fermi level with a large k-mesh and a Lorentzian broadening of 0.25 eV. All calculations were maintained at the same level of numerical accuracy, so that trends over large datasets can be captured correctly.