Quantum chemical insights into hexaboride electronic structures: correlations within the boron p-orbital subsystem

The notion of strong electronic correlations arose in the context of d-metal oxides such as NiO but can be exemplified on systems as simple as the H2 molecule. Here we shed light on correlation effects on B62− clusters as found in MB6 hexaborides and show that the B 2p valence electrons are fairly correlated. B6-octahedron excitation energies computed for CaB6 and YbB6 agree with peak positions found by resonant inelastic x-ray scattering, providing a compelling picture for the latter. Our findings characterize these materials as very peculiar p-electron correlated systems and call for more involved many-body investigations within the whole hexaboride family, both alkaline- and rare-earth compounds, not only for N- but also (N ± 1)-states defining e. g. band gaps. The hexaborides are a family of materials that exhibit an exotic range of low-temperature properties the physics of which is still under debate. Here, the authors adopt ab initio multi-configurational wavefunction theory to investigate the strong correlations of p-electrons in CaB6 and YbB6 and making a comparison with resonant inelastic x-ray scattering data.

T he physical properties of MB 6 hexaborides are remarkably diverse and have received persistent interest 1,2 , as a thorough understanding of this peculiar family of materials could not be achieved yet. Some of its members possess properties relevant for applications, generally determined by specificities of their electronic structures: CaB 6 is a semiconductor 3 , YB 6 turns superconducting for temperatures lower than 8 K 4 , LaB 6 is widely used as thermionic electron emitter 5 , EuB 6 displays colossal magneto-resistance 6 , the heavy-fermion system Ce 1-x La x B 6 has been intensively studied in the context of multipolar phases and magnetically hidden orders 7,8 , while SmB 6 is believed to host correlated topological states 9 . However, even for the seemingly simplest member of the family, CaB 6 , there are important electronic-structure features that are not at all clear, for example, the size and nature of its band gap: density-functional computations yield a vanishing band gap [10][11][12][13][14] while experiment indicates a gap of~1 eV 3 . The band structures of SmB 6 and YbB 6 are also a topic of active debate, especially in relation to their surface states 9,15-17 , for SmB 6 , even predictions made for the symmetry of the Sm 3+ -ion ground-state term are recently contradicted by experiment 18,19 .
In this context, we argue that one aspect constantly neglected in the study of these compounds is correlations within the B 2p electronic subsystem. To address such physics, we here employ ab initio multi-configurational wavefunction theory. For CaB 6 and YbB 6 , representative divalent hexaborides, the ground-state wavefunctions turn out strongly multi-configurational, with large admixture of higher-lying configurations. Trying to understand the nature of B 2p N-particle states, a rich N-particle excitation spectrum is evidenced in the ab initio calculations. These computed p-p excitation energies compare well with peak positions in B K-edge resonant inelastic x-ray spectra (RIXS) reported by Denlinger et al. 3,20 . Analysis of the many-body wavefunctions in terms of site-centered valence orbitals suggests that the ab initio data can be mapped onto a six-site, half-filled effective Hubbard model. Unfolded to the extended lattice, this might provide a convenient frame for further insights into the electronic properties of the B 2p subsystem. The computationally more expensive option is using quantum chemical methods to determine the nature of (N ± 1) quasi-particle states 21-23 defining e. g. band gaps.

Results and discussion
Correlated electronic structure of CaB 6 . The crystal structure of cubic MB 6 hexaborides is CsCl-like, with Pm-3m symmetry 1,2 . A given B 6 octahedron is surrounded by eight M nearest-neighbor sites defining a cube (see Fig. 1). M can be an alkaline-earth atom, Y, or a rare-earth ion (e.g., La, Ce, Sm, Yb). Interestingly, even for CaB 6 , a conventional closed-shell p-electron compound at first sight, standard band-structure calculations and experiment provide conflicting results on the fundamental gap: computations within the local density approximation (LDA) of density functional theory, for instance, indicate a metallic ground state 10,14 while a gap of~1 eV is found experimentally 3 . One obvious question is then what kind of correlations are missed in local density approximation (LDA) based density functional theory in this class of materials, but the answer is not evident. From post-LDA GW calculations, both finite 24 and vanishing 13,14 band gaps were reported. These contradictory sets of data suggest that the results are sensitive to specific details in the numerical implementation. Moreover, the GW scheme is not designed for dealing with short-range correlations; the latter turn out to be massive, as discussed below.
For an isolated B atom, the Aufbau principle implies an [1s 2 ] 2s 2 2p 1 valence configuration; on a B 6 2− octahedral unit, there are 6·3 + 2 = 20 electrons associated with the B 2s and 2p shells. But those 20 electrons can now be regarded as occupying symmetryadapted orbitals distributed over several B sites of a given octahedron. In cubic geometry (i. e., regular octahedron), there are a 1g , e g , and t 1u symmetry-adapted functions associated with the six s atomic orbitals and a 1g , e g , t 1u , and t 2g linear combinations related to the set of 18 p atomic orbitals.
To gain insight into the B 2s/2p electronic structure, we here rely on ab initio quantum chemical methods 25 as implemented in the ORCA program package v5.0 26 . The material model is a pointcharge embedded B 6 2− cluster, for both CaB 6 and YbB 6 (see Fig. 1). Details on the technicalities of this embedding process are given in the Methods section. The quantum chemical investigation was initiated as a preliminary Hartree-Fock (HF) calculation that provides starting orbitals for a much more sophisticated post-HF treatment: complete-active-space self-consistent-field (CASSCF) computations and subsequent multi-reference configuration interaction (MRCI) with single and double excitations 25 . To make the highly demanding MRCI calculations feasible, we enabled the use of symmetry, according to the D 2h subgroup; B 1s orbitals were not correlated.
Starting from the closed-shell HF solution, basic information on the nature of the low-lying excited states and on the nature of most important correlation effects was gathered from a series of CASSCF calculations in which different numbers of orbitals were considered as active. In quantum chemical terminology, "active" is the set of orbitals within which all possible determinants are generated on the basis of a given number of electrons. These are chemically or physically relevant orbitals, to e. g. bonding or spectroscopy, respectively; less relevant orbitals are assigned either double (core/semi-core levels) or zero (virtual levels) occupancy. By analyzing different active orbital spaces, the six molecular orbitals shown in Fig. 2(a) were found to be strongly correlated; in O h point group symmetry, relevant to a B 6 2− unit in cubic hexaborides, these orbitals transform according to the A 1g , E g , and T 1u irreducible representations. We found that for the aspects discussed here a six-orbital active space (a 1g + e g + t 1u ) accommodating six electrons (for lower-energy orbitals on the B 6 2− octahedron, double occupancy is imposed in CASSCF) already provides a balanced description; we refer to this active space as CAS (6,6).
Using this active space, a total of one septet, five quintets, 15 triplets, and eight singlet states were considered in a stateaveraged 25 CASSCF calculation. Results computed this way for CaB 6 are listed in Table 1. Electronic configurations are here expressed in terms of symmetry-adapted molecular-like orbitals (a 1g , e g , t 1u ); this is common practice when discussing the electronic structure of octahedral atomic units and also the usual way output numerical data are printed by quantum chemical programs. A first observation is the pronounced multiconfigurational structure of the computed wavefunctions. The leading ground-state configuration, for example, has a weight of only 61% in the many-body ground-state wavefunction ( 1 A 1g state) and a similar situation is observed for many of the excited states; the remaining weight has to do in each case with other electron configurations. The multi-configurational character of most of the states in Table 1 is also reflected in the natural-orbital occupation numbers, obtained through diagonalization of the density matrix 25 :for the 1 A 1g ground state, for instance, the a 1g , e g , and t 1u occupation numbers are 1.85, 1.70 (x2), and 0.25 (x3), respectively; at the HF level the occupation of the anti-bonding (see Fig. 2(a)) t 1u levels is 0 since the HF ground-state wavefunction is 100% a 1g 2 e g 4 . Also illustrative for the effect of electronic correlations: the (restricted open-shell) HF 1 A 1g -5 T 2g splitting, for instance, is less than 1 eV, much smaller than the CASSCF/MRCI values provided in Table 1.
An instructive exercise is re-expressing the CAS(6,6) groundstate wavefunction in terms of site-centered orbitals. To achieve this, the symmetry-adapted a 1g , e g , and t 1u composites were transformed to a set of six symmetry-equivalent functions using the orbital "localization" module in ORCA 26 . Two of these symmetry-equivalent site-centered orbitals are displayed in Fig. 2(b). Their asymmetric shape clearly indicates p-s mixing: dominant p but also s character, with substantial weight at the opposite B site. In this basis, the ground-state wavefunction amounts to 54% determinants with single orbital occupation, i.e., j "#"#"#i, 36% j "#"# 02i + …, 7% j "# 0202i + …, and 3% j020202i + … configurations, where ", #, 0, and 2 stand for spinup, spin-down electron, zero and double orbital occupation, respectively. As such, the wavefunction can in principle be mapped onto an effective six-site Hubbard-type model (t-U or t-U-V). Estimating hopping (t) and Coulomb repulsion (U, V) integrals [27][28][29] defining such effective models is beyond the scope of this work. But the quantum chemical wavefunctions and excitation energies reported here provide useful ab initio benchmarks for tuning the relevant parameters. Remarkably, p-electron Hubbard-U values as large as 9 eV were estimated for e. g. carbon and phosphorus allotropes 30,31 .
Good insight into correlation effects on the B 6 2− unit is attained by inspecting the wavefunction describing two opposite B sites (the hopping integral is larger for opposite-corner sp hybrids). To this end the following numerical test can be made: allow all possible occupations for two sp hybrids on opposite corners of the octahedron (c.f. those pictured in Fig. 2(b)) but restrict the occupation of each of the other four sp hybrids to 1 (without re-optimizing the orbitals). The weight of configurations with 0, 2 and 2, 0 orbital occupation for the opposite-corner orbitals is only 19%, indicative of sizable correlations. For comparison, in the H 2 molecule at equilibrium H-H interatomic distance the weight of those ionic configurations in the HF wavefunction (uncorrelated limit) is 50%, being reduced to about 44% in the exact wavefunction 32 . Support to our computational results is provided by the good agreement between calculated N-particle excitation energies and excitation energies extracted from B K-edge RIXS experiments 3 . This is apparent by comparing numerical (MRCI level of theory) and experimental values in the last two columns of Table 1 and Fig. 3(a). In other words, the numerically obtained excitation energies fit the measured excited-state relative energies; this good agreement allows for a one-to-one assignment of the main features for less than~4 eV energy loss in the RIXS spectra. But the interpretation of the electronic ground state itself and of the lower N-particle charge excitations resulting from this assignment is peculiar, with no correspondence in the extensive literature available on MB 6 hexaborides.
That computations beyond the CAS(6,6) level are necessary is visible from the large corrections brought through the MRCI treatment to some of the excitation energies, in particular, to the relative energies of singlet and triplet states at 2.4-3.5 eV and of the 7 A 1u septet, see Table 1. These corrections amount to up to 1.5 eV and point to configuration-interaction (CI) effects that are not yet captured with a six-orbital active-space CASSCF. This is evident, e. g., from the substantial reduction of the weight of the leading configuration in the CI wavefunction of each of those states. Moreover, we found that more approximate post-CASSCF CI schemes such as spectroscopy-oriented CI (SORCI) or secondorder perturbation theory (PT2) treatments such as N-electron valence PT2 correction (NEVPT2) yield inaccurate excitation energies (see Supplementary Table 1). Relative shifts of this magnitude when improving the description of electronic correlations (by appropriately enlarging the active orbital space in CASSCF and/or post-CASSCF correlation treatment) are typically found in d-electron oxides, e.g., for d 6 states in LaCoO 3 33 .
Correlated electronic structure of YbB 6 . For comparison, similar quantum chemical computations were carried out for the N-particle B 6 2− electronic structure in the rare-earth hexaboride  Weights of leading configurations are provided within brackets. Computational results based on complete active space self-consistent field (CASSCF) and more sophisticated multi-reference configuration interaction (MRCI) are given in comparison to experimental resonant inelastic x-ray scattering (RIXS) peaks.
YbB 6 . The relative energies of the first twelve excited states and leading electronic configurations are listed in Table 2. In agreement with the computational results obtained for CaB 6 , the ground-state wavefunction is of strong multi-configurational character, with its leading configuration amounting to only 60%. Also in this case, the calculated MRCI energies can be compared with B K-edge X-ray spectra; MRCI N-particle excitation energies are overlaid for this purpose onto the experimental plots of Denlinger et al. 20 in Fig. 3. It is seen that the relative energy of the well-defined first peak is rather well reproduced in both compounds. At higher energies, the level structure is denser; qualitatively, there is reasonably good correspondence between the peak found experimentally in each of the two materials at 2-3 eV and the excitations computed in that energy range. The B p-p excitation energies computed for YbB 6 are slightly lower as compared to CaB 6 . The likely explanation is that the larger effective charge of Yb ions in YbB 6 (2 or slightly larger than 2 according to experimental estimates 17,34-36 , a representation that is also adopted for the surroundings of the B 6 octahedron in our calculations) stabilizes the t 1u B-cluster orbitals as compared to the e g (also a 1g ) components and reduces this way the e g → t 1u excitation energies-the t 1u composites have stronger antibonding character, are more extended, and therefore sense more effectively the nearest-neighbor cations.

Conclusions
While strong correlation effects within a whole manifold of electron configurations were already pointed out for the case of d-ion clusters 37-41 we here describe similar physics occurring on pelectron octahedral clusters as found in MB 6 hexaborides. In particular, we highlight the role of correlation effects on B 6 octahedral units. For divalent CaB 6 and YbB 6 our quantum-chemically computed valence-shell excitation energies agree well with peak positions in earlier B K-edge RIXS measurements 3,20 , providing an interpretation of the RIXS spectrum. Shifts to lower energies of the N-particle excited states in YbB 6 are ascribed to a relative stabilization of the anti-bonding t 1u cluster orbitals as a result of larger cation charges in the Yb compound. Given (i) the limitations of density functional theory in describing the electronic structure of hexaboride compounds, documented in the case of CaB 6 in literature extending over a few decades 3,[10][11][12][13]24 , and (ii) more recent controversies on the electronic structures of e. g. SmB 6 9,18,19 and YbB 6 15-17 , an extension of the present investigation of N-particle excitations to quasi-particle band structures is highly desirable. The foundations of quasi-particle band-structure calculations with ab initio wavefunction-based methods were laid 22,23,42 and such an approach provides better control over the underlying approximations as compared to mainstream computational schemes based on density functional theory. The failures of the latter in hexaborides appear to be related to the correlated nature of the B 6 2− electronic wavefunctions.

Methods
As material model, a B 6 2− cluster was embedded in a field of M (with M = Ca, Yb) and B point charges (see Fig. 1 and Supplementary Table 2 for the cluster coordinates). The point charge field was created using the EWALD program 43,44 . For the lattice parameters, experimental data were used: lattice vectors a = 4.1537 Å and an x-axis B atom position x B = 0.201 for CaB 6 [45][46][47] and parameters of a = 4.1792 Å and x B = 0.202 for YbB 6 17,48 . To account for covalency effects, an initial charge value of +1.8 e was chosen for Ca ions in the embedding, which corresponds to an ionic charge of −0.3 e for B. In contrast, in YbB 6 , initial effective charges of +2.0 e for Yb 17,34-36 and −0.33 e for B were assumed. For both embeddings, these charges were slightly adjusted by EWALD to ensure convergence of the Madelung sum. The final point charge values do not deviate more than 0.01 e from the initial values. Charge neutrality was ensured in the overall system by slightly adjusting a set of  Weights of leading configurations are provided within brackets. Computational results based on complete active space self-consistent field (CASSCF) and more sophisticated multi-reference configuration interaction (MRCI) are given in comparison to experimental resonant inelastic x-ray scattering (RIXS) peaks.
about 150 B point charges at the boundary of the whole array of approximately 40 000 point charges. To eight M nearest neighbors and to six adjacent B 6 octahedra around the reference B 6 2− quantum cluster (see Fig. 1) we assigned capped effective core potentials (cECPs); the pseudopotentials of Kaupp et al. 49 , Dolg et al. 50 . and Bergner et al. 51 . were employed for Ca, Yb and B, respectively. Their main purpose is to prevent electron density from the quantum cluster being artificially polarized by nearby point charges. Consequently, for the quantum cluster, an all-electron correlation-consistent polarized valence quadruple zeta quality (cc-pVQZ) basis set 52 was used in conjunction with the CASSCF method described above.

Data availability
Supplementary Data (excitation energies of CaB 6 using further methods, cluster and point charge field geometries) is attached and Orca calculations outputs are deposited at