How good are recent density functionals for ground and excited states of one-electron systems? Open Access

Density-functional theory^1,2 (DFT) has become one of the workhorses of computational chemistry, material science, and related fields, as modern density-functional approximations (DFAs) only require a reasonable amount of computational effort while providing a level of accuracy sufficient for semi-quantitative predictions on a broad range of systems.^3–5 Although hundreds of DFAs have been proposed, thus forming the infamous zoo of DFAs,⁶ new DFAs continue to be developed with the aim to find more universally applicable DFAs that combine suitable levels of accuracy and numerical effort.

New DFAs can be constructed along various strategies.^5,7–9 The traditional route to construct DFAs is to start from first principles and to impose known limits and constraints; this is the way along which many well-known functionals such as PBE,¹⁰ TPSS,^11,12 and SCAN¹³ have been constructed.

Semi-empirical fitting is another route for constructing DFAs. Here, the general idea is to introduce flexibility in the functional form by introducing several independent DFA components that are weighted by parameters, which are optimized against some training dataset. Classical examples of semi-empirically fitted functionals include B3LYP,¹⁴ Becke’s 1997 functional¹⁵ (B97), and several refinements thereof such as the HCTH [Hamprecht–Cohen–Tozer–Handy] functionals by Handy and co-workers,^16,17 as well as the Minnesota family of DFAs by Truhlar and co-workers^18,19 that has been reviewed by Mardirossian and Head-Gordon.²⁰ 

In reality, the classification of functionals into ones built solely from first principles vs ones formed by semi-empirical fitting is not always clear: for instance, the TPSS exchange functional is parameterized to yield the exact energy for the hydrogen atom’s exact ground state density,^11,12 while the SCAN functional¹³ includes parameters that are fit to data on noble gases. Modern semi-empirical functionals,^21–26 in turn, typically employ a combination of the two approaches by restricting the fits to known constraints.

DFAs from either route are widely used, given their suitable numerical accuracy and reasonable computational effort. However, the functionals obtained from the two routes tend to exhibit different behavior. For instance, while semi-empirical DFAs often deliver excellent descriptions of the total energy, they may fail to reproduce electronic densities of the same quality: a famous article by Medvedev et al.²⁷ initiated an intense debate about this in the literature;^27–32 it was even pointed out that any general mathematical measure of density error is too arbitrary to be universally useful.³³ DFAs built on physical first principles, in contrast, often yield steady performance in a variety of applications, but may not achieve the same level of accuracy as tailored functionals for specific types of systems.

One of the most important limitations of present-day DFAs, regardless of their design, is the self-interaction error (SIE): an artificial interaction of the electrons with themselves. This error is related to density delocalization error and the fractional electron problem,^34,35 and leads to incorrect dissociation limits³⁶ and barrier heights,³⁷ for instance. Recent avenues for circumventing SIE in DFAs involve determining the electron density with another method, such as Hartree–Fock^38,39 or multiconfigurational wave function theory.^40,41 Other types of approaches have also been proposed in the literature. To solve the self-interaction problem, Perdew and Zunger⁴² (PZ) proposed an orbital-by-orbital self-interaction correction (SIC),

where E^KS is the Kohn–Sham (KS) energy functional,² and the self-interaction error (SIE) is defined by

Here, n_iσ is the electron density of the i-th occupied orbital with spin σ, and E_J and E_xc denote the Coulomb and exchange-correlation energy functionals, respectively. The idea behind PZ-SIC is that the self-interaction error defined by Eq. (2) vanishes for the exact functional,⁴² and thereby, the Perdew–Zunger functional of Eq. (1) is a better estimate for the total energy than the uncorrected Kohn–Sham DFA E^KS; indeed, the PZ functional is exact for one-electron systems such as the $H2+$ molecule with approximate DFAs.

Despite the simple logic used to construct the PZ-SIC functional, the PZ-SIC method turns out to be quite complicated. The introduction of the explicit orbital dependence in Eqs. (1) and (2) breaks the unitary invariance of the energy functional,⁴³ requiring costly unitary optimization of the orbitals (see Ref. 44 for discussion). Even though the resulting method is known to correct charge transfer errors and barrier heights, it does not lead to improved atomization energies with GGA and meta-GGA functionals in general.⁴⁵ 

Continued research has illuminated other important theoretical aspects of PZ-SIC. First, the orbital dependence in Eqs. (1) and (2) has been recently shown to require the use of complex-valued orbitals for proper minimization, as real-valued orbitals can be shown to correspond to high-order saddle points.⁴⁶ When complex-valued orbitals are employed, the total energy is lowered, and PZ-SIC does lead to improved atomization energies for some GGA functionals; however, more accurate atomization energies can be obtained at significantly smaller cost with several standard DFAs.⁴⁷ 

Second, the orbital dependence in Eqs. (1) and (2) has also been shown to lead to the existence of several local minima in the orbital space.⁴⁶ This problem has been recently shown to persist also in a related SIC method⁴⁸ based on the use of Fermi–Löwdin orbitals (PZFLO-SIC), where various choices for the orbital descriptors lead to distinct local electronic minima.⁴⁹ The existence of such local minima is a significant and underappreciated aspect of PZ-SIC and PZFLO-SIC calculations, as finding the true ground state may require extensive sampling of the space of the various possible localized electronic configurations or bonding situations.

Despite their theoretical shortcomings, PZ-SIC and PZFLO-SIC have been found useful in many applications,^50–52 and we are positive that several of the aforementioned issues in PZ-SIC and PZFLO-SIC can be addressed by developments in the related theories by changing the way the self-interaction correction is applied. One possible way to achieve improved results would be to revisit DFAs based on the requirements of SIC calculations.⁵³ It is known that present-day DFAs yield poor estimates for the noded electron densities that are involved in SIC calculations.^54,55 Sun et al.⁵⁴ demonstrated that the ground and excited state densities of the hydrogen atom (as well as of $H2+$⁠, see below) lead to large relative errors in the exchange-correlation energy compared to the exact values, but we are not aware of any self-consistent calculations on this issue.

Following the recent discussion in the literature on the accuracy of DFAs on the electron densities of small atoms and ions^27–33 and motivated by the obvious connection of one-electron errors (OEEs) to the PZ-SIC and PZFLO-SIC methods, in this work, we will analyze the OEE of various functionals for the 1s ground state as well as the 2p and 3d excited states of hydrogenic ions Z^(Z−1)+, whose exact energies are well known to be given in atomic units by

where Z is the atomic number, n ≥ l + 1 is the principal quantum number, and l is the angular momentum.

As mentioned above, calculations of ground and excited states of the hydrogen atom and the 1σ_g ground state and 1σ_u excited state of $H2+$ have been discussed by Sun et al.⁵⁴ with non-self-consistent electron densities, while the 1s ground states of hydrogenic mononuclear cations as well as the 1σ ground states of hydrogenic diatomic cations have been discussed recently by Lonsdale and Goerigk⁵⁶ using self-consistent calculations. The novel contribution of this work is to address (highly) excited states with noded electron densities of hydrogenic cations self-consistently. Importantly, similar to the 1s ground state, the 2p and 3d excited states (as well as the analogous higher excited states like 4f) are the lowest states of the corresponding symmetry, and the ground-state Kohn–Sham scheme is applicable to such excited states as well as shown by Gunnarsson and Lundqvist.⁵⁷ 

We pursue thorough density-functional investigations of the 1s, 2p, and 3d states of hydrogenic ions in benchmark-quality Gaussian basis sets, especially suited for this purpose with a selection of 56 popular DFAs, including the recently developed, highly sophisticated machine-learned DeepMind 21 (DM21) local hybrid functional.⁵⁸ 

The layout of this work is as follows. The computational details are presented in Sec. II, and the results are given in Sec. III. A summary of our findings and an outlook for further investigations are given in Sec. IV. Atomic units are used throughout, unless specified otherwise.

We only use free and open-source software (FOSS) in this work, following the philosophy discussed in Ref. 59. PySCF⁶⁰ is an electronic structure code for all-electron calculations by using Gaussian-type orbitals (GTOs). As we are targeting one-electron states of specific symmetry (s, p, or d states), following Gunnarsson and Lundqvist,⁵⁷ we truncate the basis set in all calculations to contain functions only of the pursued symmetry: calculations on the 1s/2p/3d state only include the basis functions of the corresponding symmetry (s, p, or d functions, respectively) from the chosen parent basis set. This procedure has two important features: the 2p and 3d excited states become the ground state in the reduced-basis calculation, and the computational requirements are smaller since fewer integrals need to be calculated in the reduced basis than in the original basis set.

The one-electron guess—which is exact for one-electron systems and thereby is also expected to be accurate for calculations employing DFAs as well—is used in all calculations.⁶¹ To ensure that the SCF procedure converges to the global minimum instead of a saddle point, the following procedure was used. First, a regular SCF calculation was performed with PySCF with default settings; direct inversion in the iterative subspace (DIIS) is used to accelerate these calculations.^62,63 Next, convergence to saddle point solutions was checked: cases where the SCF converged to a final energy higher than that of the initial guess were restarted, with new calculations employing iterative diagonalization with level shifting⁶⁴ instead of DIIS to converge to the ground state. All calculations reported in this work are fully converged to a threshold of 1 × 10⁻⁷ E_h.

For the GTO basis sets, we use the family of hydrogenic Gaussian basis sets⁶⁵ (HGBS-n) that have been designed for high-accuracy calculations on atoms and small molecules. A special feature of the HGBS basis sets is that the basis for atomic number Z is determined by a universal even-tempered basis set for the ions Y^(Y−1)+ for Y ∈ [1, Z], whereas augmented hydrogenic Gaussian basis sets (AHGBS-n) add further functions for describing the Z = 0.5 one-electron ion.⁶⁵ The parameter n controls the relative precision of the hydrogenic Gaussian basis; (A)HGBS-n reproduces the exact total energies of the one-electron ions to an approximate relative accuracy of 10⁻ⁿ.⁶⁵ The motivation of this approach in Ref. 65 was that a many-electron atom experiences a screened nuclear charge that can be rewritten in terms of a radially dependent effective charge Z^eff = Z^eff(r) with the asymptotic limits Z^eff(0) = Z and either Z^eff(∞) = Z_∞ with the asymptotic limit Z_∞ = 0 for Hartree–Fock and DFT or Z_∞ = 1 for the exact effective potential.⁶¹ 

Another feature of the HGBS basis sets is that the functions of various angular symmetries are determined independently of each other, which facilitates the formation of polarized counterparts of the basis sets that are essential for studying molecules and excited states, as additional shells are added to the basis like lego blocks.⁶⁵ Following Ref. 65, the basis set with p ≥ 1 polarization shells and accuracy n is denoted as (A)HGBSPp-n. The definition of polarization shells varies by atom (see Ref. 65 for discussion); however, as we only include the functions of the pursued symmetry in each calculation, the choice of the polarization level of the (A)HGBSPp-n basis set does not matter as long as the original basis contains functions of the highest targeted angular momentum for the targeted atom, that is, d functions in this work.

For the reasons listed above, the hydrogenic Gaussian basis sets of Ref. 65 are ideally suited for the present study—as will be demonstrated in Sec. III by benchmarks with functions from the polarization consistent (pc-n) basis sets⁶⁶ and their augmented versions⁶⁷ (aug-pc-n)—and, as will be discussed in Sec. III A, we will take the exponents from the AHGBSP3-n basis sets in this work. All basis sets were taken from the Basis Set Exchange.⁶⁸ 

The Libxc library⁹—which implements over 600 DFAs—is used in PySCF to evaluate the DFAs. The library provides access to a vast variety of DFAs, of which 49 were chosen for this work; see Table I for the complete list of investigated functionals. Our selection includes functionals from the first to the fourth rung of Jacob’s ladder,¹²³ that is, local density approximations (LDAs), generalized gradient approximations (GGAs), meta-GGAs, as well as global and range-separated hybrid functionals. In addition, we consider six hybrids of rSCAN and r²SCAN with varying fractions of Hartree–Fock exchange discussed in Ref. 120; these functionals were defined in the PySCF input files as weighted combinations of r(²)SCAN exchange and Hartree–Fock exchange + 100% r(²)SCAN correlation. The DM21 functional was also chosen for this study; we use the original implementation in PySCF of Kirkpatrick et al.⁵⁸ This brings up the total to 56 functionals for this study. An unpruned (300,590) quadrature grid is used in all calculations, including the non-local correlation component in B97M-V, ωB97M-V, and LC-VV10.

As the total energies scale as E_n ∝ Z² according to Eq. (3), the results will be analyzed in terms of absolute relative errors (AREs). Hydrogenic estimates show that the approximated exchange-correlation energy scales like Z in the large Z limit,¹²⁴ meaning that the relative errors should tend to zero like 1/Z. The ARE for a given state of a given ion is given by

The information in the AREs is analyzed with two more error metrics. The mean state error (MSE) measures the overall functional error over all ions by averaging the ARE over all ions,

The overall error (OE) for a functional is obtained by further averaging the MSE over all considered states (1s, 2p, and 3d),

Before pursuing density-functional calculations, we analyze the basis set truncation errors (BSTEs) for the one-electron cations in the polarization consistent and hydrogenic Gaussian basis sets. We aim for a mean BSTE smaller than 5 × 10⁻⁵ E_h for the whole benchmark set ranging from H⁰ to Kr³⁵⁺ to ensure that our results are converged close to the complete basis set limit.

Unrestricted Hartree–Fock (UHF) is exact for one-electron systems and thereby gives the exact energy $EnUHF$ in the studied basis; therefore, the difference between $EnUHF/basis$ and the exact analytical energy [Eq. (3)],

is a variational measure of the BSTE for the state with given n of the studied hydrogenic ions.

The calculated mean BSTEs for a variety of polarization consistent and hydrogenic Gaussian basis sets are shown in Fig. 1; additional results can be found in the supplementary material. Unsurprisingly, uncontracting the (aug-)pc-n basis sets—yielding the unc-(aug-)pc-n basis sets—results in a noticeable decrease in the BSTE, because the contractions were determined in Ref. 66 with the BLYP functional^71–72 that suffers from SIE for the 1s state, while the p and d functions in the basis set describe either polarization effects or the occupied p or d orbitals in the screened neutral atom. Although the large uncontracted polarization consistent basis sets exhibit satisfactory performance for the 1s state, they result in much larger errors for the 2p and 3d states; this error is again caused by the p and d orbitals in the neutral atom being screened by the core electrons, which results in the lack of tight p and d basis functions that are necessary for the 2p and 3d states of the one-electron ions.

In contrast, the primitive (not contracted) hydrogenic Gaussian basis sets of Ref. 65 show uniform accuracy for the 1s, 2p, and 3d states, and as can be observed in Fig. 1, the targeted mean BSTE threshold is roughly achieved already with the AHGBSP3-7 basis set. The AHGBSP3-9 basis sets yield errors below the desired threshold for all states and is, therefore, chosen for all the remaining calculations of this study.

Although this analysis was limited to Hartree–Fock calculations, we note that the basis set requirements of Hartree–Fock and DFT are known to be similar.¹²⁵ Furthermore, reliable reference energies for DFAs can be obtained with fully numerical methods,^137–128 and exploratory calculations presented in the supplementary material confirm that the BSTEs in the AHGBSP3-9 basis are also small for DFAs.

We begin the analysis with a graphical study of the results of the SPW92, PBEsol, revTPSS, MN15-L, BHandH, and DM21 functionals in Fig. 2. This limited set of functionals contains LDA, GGA, and meta-GGA functionals from first principles (SPW92, PBEsol, and revTPSS, respectively), semiempirical functionals (MN15-L and DM21), as well as hybrid functionals (BHandH and DM21).

As will be discussed in Sec. III B 2, revTPSS is the most accurate meta-GGA functional for the 1s state. In Fig. 2, revTPSS is outperformed by DM21 only for He⁺; otherwise, revTPSS affords much lower errors than the five other functionals in the figure. In contrast, the performance of DM21 is inconsistent. DM21 has lower errors for light ions than for heavy ions, but the curve is kinked for the light ions. DM21’s curve becomes smooth for heavy ions, but DM21 is also less accurate for heavy ions. MN15-L also shows a kinky behavior with lower errors for light ions; these non-systematic features of DM21 and MN15-L can be tentatively explained by their semiempirical character; the curves for the first principles functionals are smoother.

The functional errors for the 2p state are shown in Fig. 2(b). The performance for the 2p state is strikingly different compared to the 1s state shown in Fig. 2(a). The plots for the 2p state in Fig. 2(b) show more structure and curve crossings. The behavior of DM21 is qualitatively different from that of the other functionals: DM21 shows large relative errors for light atoms and lower relative errors for heavy atoms, while most of the other functionals shown behave similarly to each other. The only other exceptions to this are the SPW92 and MN15-L functionals that show dips at Z ≃ 3 and Z ≃ 4, respectively; the two functionals are thus oddly more accurate for some values of Z than others.

The errors for the 3d state are shown in Fig. 2(c). BHandH has small errors for all ions for the 3d state. The behavior of DM21 and MN15-L again differs qualitatively from the other functionals. Although DM21 shows less variation for the 3d state than for the 2p state, MN15-L does the opposite: MN15-L has large errors for light ions and becomes nearly as accurate as BHandH for Z ≃ 22; however, the relative error increases again for heavier ions.

The MSEs and OE for all studied functionals are shown in Table I. Although Table I contains all of the data used in the present analysis, additional tables showing the rankings of the functionals in terms of the errors for the 1s, 2p, and 3d states as well as in terms of the overall error can be found in the supplementary material.

Clearly, the performance of all LDAs is practically identical. This suggests that the functional error for LDAs is limited by the simple functional form. Although LDAs show larger errors than GGAs and meta-GGAs for the 1s state, they perform better than GGAs and meta-GGAs for the 2p and 3d states.

PBEsol¹⁰⁴ is the GGA that yields the smallest errors for the 2p and 3d states, as well as the smallest overall error. The XLYP GGA has a lower error than PBEsol for the 1s state. Although XLYP and even PBEsol are better for 1s states than any LDA, they have higher OEs than any LDA because of their considerably poorer performance for the 2p and 3d states. Analogous findings apply also to all other studied GGAs.

The best meta-GGA for the 1s state is revTPSS,^99,100 closely followed by rSCAN¹²¹ and r²SCAN¹¹⁸ (see Table I or the supplementary material). The best meta-GGA in terms of overall error is MN15-L.⁹⁷ 

Hybrid functionals have better accuracy as they contain some Hartree–Fock exchange that is free of self-interaction. The best hybrid GGA functionals in terms of overall error are BHandH⁷³ and QTP17.¹⁰⁹ BHandH has low MSEs for all states and has the best overall performance, which can be understood by its composition of 50% of Hartree–Fock exchange and 50% LDA exchange + 100% Lee–Yang–Parr correlation. QTP17 has the second-best performance for all states; it, too, contains a mixture of Hartree–Fock (62%) and LDA exchange (38%).

The best functionals of each rung in terms of overall error are SPW92, PBEsol, MN15-L, and BHandH. The corresponding error distributions are summarized in comparison to DM21 in Fig. 3. Interestingly, the performance of the DM21 functional appears similar to that of MN15-L.

Following Medvedev et al.,²⁷ the calculated OE for all functionals and ions is plotted against the publication year in Fig. 4. As is clear from this plot, the improvement in one-electron error is not fully systematic and features a significant amount of spread and some notable outliers such as M11-L, GAM, and TASK. However, in the recent decade, hybrid functionals have become better overall. As an example, the various QTP functionals dominate the bottom right of the figure; these functionals are closely related in functional form and contain high amounts of Hartree–Fock exchange, which decreases the one-electron error. Unsurprisingly, hybrid functionals based on rSCAN and r²SCAN perform well, and the functionals with large fractions of Hartree–Fock exchange share the bottom right of the figure with the QTP functionals.

The study of Lonsdale and Goerigk⁵⁶ employed an uncontracted aug-cc-pVQZ basis set^129,130 with the following nucleus-dependent quadrature grids for the study of 1s states of hydrogenic cations: (45,770) for H and He, (50,770) for Li–Ne, (55,770) for Na–Ar, and (60,770) for K–Kr. However, it appears that K and Ca were excluded from Ref. 56 (see Figs. 3 and 10 in Ref. 56).

Our trends and absolute values for the MSE for the 1s state are in satisfactory agreement for the subset of functionals studied in both works, although we did identify basis set incompleteness issues in some results of Ref. 56 as discussed in the supplementary material. The largest basis set incompleteness effects are observed for the M11-L and M06-L Minnesota functionals, which are known to converge remarkably slowly to the basis set limit.¹³¹ 

Lonsdale and Goerigk⁵⁶ only studied one LDA functional (SVWN); we considered more LDAs and found them to have similar performance. Lonsdale and Goerigk⁵⁶ found OLYP to be the best GGA functional for the 1s state; we also considered XLYP and found it to yield a considerably lower MSE for the 1s state than OLYP. Lonsdale and Goerigk⁵⁶ included a broader set of hybrid functionals separating global, range-separated hybrids and double hybrids; however, our main motivation is the connection to self-interaction corrected methods where hybrid functionals are typically not used. We found rSCAN50/r²SCAN50 to be the best hybrid functional for the 1s state, while Lonsdale and Goerigk⁵⁶ determined TPSSh and SCAN0 to be the best hybrids. All rSCAN- and r²SCAN-based hybrid functionals, i.e., rSCANh/r²SCANh, rSCAN0/r²SCAN0, rSCAN50/r²SCAN50, and TPSSh, have a good performance for the 1s state. The revTPSS functional is found in our work as well by Lonsdale and Goerigk⁵⁶ to be the best non-hybrid meta-GGA functional for the 1s state.

We used exactly solvable hydrogenic cations in their 1s ground state and 2p and 3d excited states to determine the self-consistent one-electron error for 56 density functionals including the novel DM21 of Kirkpatrick et al.,⁵⁸ employing the methodology of Gunnarsson and Lundqvist⁵⁷ for the excited state calculations. In accordance with an earlier finding by Sun et al.⁵⁴ apparently based on non-self-consistent calculations for the hydrogen atom and molecule and one LDA functional, we find for 36 hydrogenic cations that all LDAs perform better for the excited 2p and 3d states than any of the tested GGAs and meta-GGAs. The performance of various LDAs appears to be almost identical since the calculated errors are nearly indistinguishable, suggesting that the errors are limited by the simple functional form used in LDAs. Sun et al.⁵⁴ pointed out that larger errors for excited states are a necessary consequence of orbital nodality.

The revTPSS functional is the best-performing meta-GGA for the 1s state, tightly followed by the rSCAN and r²SCAN functionals. MN15-L shows a better overall performance than LDAs for all states; however, the performance of MN15-L is non-systematic like that of DM21.

Hybrid functionals like BHandH and QTP17 have the best overall performance as they explicitly include some fraction of Hartree–Fock exchange. Moreover, both BHandH and QTP17 are mixtures of Hartree–Fock and LDA exchange, leading to good observed accuracy.

DM21 turns out to be only close to exact for the 1s state OEE from H⁰ to B⁴⁺ (see Fig. 2). For He⁺ to B⁴⁺, DM21 also shows good performance for 2p and 3d states. However, over the whole range of investigated species H⁰ to Kr³⁵⁺, DM21 exhibits various trend changes and overall inconsistent performance. Thus, one might improve the next generation of the DM21 functional by including more one-electron cations in the training sets for various elements in the periodic table. This might increase the consistency of promising properties of such kind of machine-learned functionals.

We found PBEsol to be the most accurate GGA functional for the 2p and 3d states. PBEsol is also the most accurate GGA functional overall. These findings are interesting to contrast with that of Lehtola, Jónsson, and Jónsson,⁴⁷ who showed that PBEsol is one of the few functionals whose accuracy improves when PZ-SIC is applied with complex orbitals. Therefore, the development of novel DFAs with reduced one-electron error could be useful for PZ-SIC calculations, as the reduced one-electron errors [Eq. (2)] would affect the numerics of the PZ correction [Eq. (1)] and might alleviate well-known issues with PZ-SIC and PZFLO-SIC discussed in Sec. I.

Note added in proof After the acceptance of this paper, we became aware of a preprint by Lonsdale and Goerigk [132] that includes discussion on excited states of hydrogenic cations.

See supplementary material for exploratory finite element studies of basis set truncation errors in density-functional calculations. Sorted rankings of the functionals by errors for the 1s, 2p, and 3d states as well as the overall error. Bar plots of the errors for all studied functionals. Comparison of the 1s data to the study of Lonsdale and Goerigk⁵⁶ with additional basis set incompleteness studies. Tables of the calculated total energies for the 1s, 2p, and 3d states for all studied functionals.

We thank Jens Kortus for their valuable comments on the manuscript. S.S. has been funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) - Project No. 421663657 - KO 1924/9-2. This work was part of the OpenSIC project. We thank ZIH Dresden for computational time and support.

The authors have no conflicts to disclose.

Sebastian Schwalbe: Conceptualization (equal); Data curation (lead); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Software (equal); Supervision (equal); Visualization (lead); Writing – original draft (equal); Writing – review & editing (equal). Kai Trepte: Conceptualization (equal); Data curation (supporting); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Software (equal); Validation (lead); Visualization (supporting); Writing – original draft (equal); Writing – review & editing (equal). Susi Lehtola: Conceptualization (lead); Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Software (equal); Supervision (lead); Validation (equal); Visualization (supporting); Writing – original draft (equal); Writing – review & editing (lead).

The data that support the findings of this study are available within the article and its supplementary material.