Fingerprints of an Exotic Orbital-Selective Mott Phase in the Block Magnetic State of BaFe$_2$Se$_3$ Ladders

Resonant Inelastic X-Ray Scattering (RIXS) experiments on the iron-based ladder BaFe$_2$Se$_3$ unveiled an unexpected two-peak structure associated with local orbital ($dd$) excitations in a block-type antiferromagnetic phase. A mixed character between correlated band-like and localized excitations was also reported. Here, we use the density matrix renormalization group method to calculate the momentum-resolved charge- and orbital-dynamical response functions of a multi-orbital Hubbard chain. Remarkably, our results qualitatively resemble the BaFe$_2$Se$_3$ RIXS data, while also capturing the presence of long-range magnetic order as found in neutron scattering, $only$ when the model is in an exotic orbital-selective Mott phase (OSMP). In the calculations, the experimentally observed zero-momentum transfer RIXS peaks correspond to excitations between itinerant and Mott insulating orbitals. We provide experimentally testable predictions for the momentum-resolved charge and orbital dynamical structures, which can provide further insight into the OSMP regime of BaFe$_2$Se$_3$.

The intuitive origin of the block magnetic state admits multiple descriptions. Within Hartree-Fock treatments [18], the "block" structure arises from magnetic frustration because the system is located in parameter space between a ferromagnetic state, induced by doubleexchange at large Hund coupling, and a staggered antiferromagnetic state, induced by superexchange at small Hund coupling. Alternatively, recent efforts based on an itinerant perspective predicted that the metallic orbitals are the "drivers" while the localized spins are the "passengers", with the block structure arising from Fermi surface nesting [25]. From the experimental perspective, Resonant Inelastic X-ray Scattering (RIXS) and X-ray Photoemission (XPS) studies report the presence of both localized and itinerant carriers in BaFe 2 Se 3 , also suggesting OSMP physics [20,21]; however, up to now, no sufficient theoretical study has been carried out to support these indirect experimental claims of an OSMP ground state in BaFe 2 Se 3 . Our primary goal is to fill this gap and provide evidence from a theoretical perspective that indeed the "123" ladder is in an OSMP state.
One of the first steps towards unveiling the characteristic excitations of an OSMP state is to calculate its single-particle spectral function A(q, ω) and its intraand inter-orbital dynamical spin S(q, ω), charge N (q, ω), and orbital L(q, ω) structure factors. Recently, theoretical predictions for S(q, ω) were presented [26] in a block antiferromagnetic (AFM) OSMP. Using determinantal quantum Monte Carlo and the maximum entropy method, the low temperature A(q, ω) was also reported employing an approximation to a multi-orbital Hubbard model [27]. However, N (q, ω) and L(q, ω), crucial for RIXS experiments, have not been addressed yet. Therefore, here for the first time, we use the density matrix renormalization group (DMRG) technique to compute the momentum-resolved charge and orbital response functions (see Figs. 2a and Figs. 2b), as well as the singleparticle spectral function (at zero temperature), of the block magnetic OSMP of a multi-orbital Hubbard chain.
There are multiple reasons for focusing on a chain geometry as opposed to a ladder. First, a three-orbital Hubbard chain is already equivalent to a three-leg oneorbital ladder; thus a three-orbital Hubbard ladder maps onto a six-leg one-orbital ladder, which is challenging even with DMRG. Second, in the usual "snake" geometry of DMRG, interorbital hoppings mutate into longdistance hoppings in the one-orbital chain analog, here involving sites effectively separated by eight lattice spacings. Such long-range hoppings compromise the accuracy of DMRG. Finally, RIXS in the eV scale usually addresses local excitations, and there should not be much difference between chains and ladders, as recent efforts in neutron scattering using both chains and ladders showed [26].
We compare our results against previously gathered RIXS data for BaFe 2 Se 3 [20], which measures magnetic, charge, and orbital excitations simultaneously [28][29][30][31][32]. By calculating the orbital response functions of the competing paramagnetic metal (PM) and ferromagnetic (FM) insulator states, we show that block OSMP has a characteristic two-peak structure that is distinctive and in striking agreement with RIXS results on BaFe 2 Se 3 [20]. Moreover, we identify the observed dd peaks as orbital excitations between localized and itinerant orbitals. Our study strongly suggests that the ground state of BaFe 2 Se 3 is indeed an OSMP with block magnetic order.
Results - Figure 1 plots the local charge occupation (panel b) and fluctuations (panel c) of each orbital vs the total electronic density, which demonstrates that the OSMP is stable and robust against variations in the hole and electron doping [33]. First, we focus on the results for 4 electrons per site (n = 4), shown by the vertical dashed line in Fig. 1 (b,c) (note that n = 4 in a three-orbital model is the analog of the realistic n = 6 in a five-orbital model). At this density, orbital c has one electron per site, and thus becomes Mott localized, while the two other orbitals remain fractionally occupied, and thus metallic. The reason for this exotic behavior relies on the crystal-field splitting and different bandwidths of the orbitals: as the Hubbard U interaction opens a gap and shifts energy states "up and down" relative to the gap, it becomes energetically favorable for orbital c to have a half-filled lower Hubbard band, which is formed by moving electrons from the other bands into c. The corresponding intra-orbital charge fluctuations δN 2 γ = n 2 γ − n γ 2 in orbital c are significantly suppressed, compatible with a localized state, while charge fluctuations in orbitals a and b remain finite, compatible with a metallic state. The inset of Fig. 1b establishes the magnetic order at this filling, and plots the spin structure factor at n = 4, displaying the characteristic peak at q = π/2 that is in agreement with the experimentally observed non-trivial block-type AFM order. Note that orbital c (representing the d xy orbital, see methods) has the standard characteristics of a Mott phase at n = 4, where charge degrees of freedom are "frozen" (localized) and accompanied by well-formed local magnetic moments, m c ∼ 0.98µ B . In contrast, the large local charge fluctuations of the other orbitals suggest a metallic behavior typical of itinerant electrons. The total on-site local moment at n = 4 is m tot 1.97µ B , a robust value slightly larger than in experiments [13,15,19].
The coexistence of localized and itinerant carriers is also evident when examining the electronic density of each orbital vs the global filling near the commensurate value n = 4 (Fig. 1b). The linear behavior of n a and n b suggests a band-like picture, while the robust plateau in n c extending over a wide range of doping (n = [3 : 4.5]) indicates a Mott gap in the single-particle spectral func-  [26]. (b) Average local occupation of each orbital vs overall electron density n, where n = 4 represents the filling of 4 electrons per site (using DMRG, 32 sites), which is a realistic situation when a three-orbital model is used. The inset shows the total static spin-structure factor with a peak at q = π/2, representing the block-type AFM order at n = 4, with two spins "up", followed by two "down", and then a repeated periodicity, as illustrated in (a) upper orbital. (c) Average local charge fluctuations on each orbital vs n. tion of orbital c. Since all these results provide ample evidence that we correctly capture the block-magnetic properties of BaFe 2 Se 3 , henceforth we fix the filling n = 4 to address the features of the block OSMP in the RIXS dynamical excitations, the focus of our publication.
Dynamical spectra of the block OSMP - Figures 2d and 2e show the orbital-resolved single-particle spectral functions. Orbitals a and b have a finite quasiparticle weight at the Fermi level E F , as also shown in the density-of-states in Figure 3a, confirming once again the itinerant nature of electrons in these orbitals. Moreover, FIG. 2. Sketches of (a) charge and (b) orbital excitations in a block OSMP state. (c) is a sketh of "orbital fractionalization" discussed in the FM insulator context below. Actual DMRG data shown correspond to the orbital-resolved (d-e) single-particle spectra, (f-g) charge excitations, and (h-i) orbital excitations spectra within the block magnetic OSMP state. A a/b (k, ω) has reduced but nonzero quasi-particle weight at the Fermi level (EF ) while A c (k, ω) has a characteristic Mott gap ∆c ∼ 1.5 eV. The charge excitations N a/b (k, ω) have a gapless continuum at finite ω while N c (k, ω) has an excitation gap of ∆ N c ∼ 1.5 eV. L x/y (k, ω) represents inter-orbital excitations between orbitals a/b and c (d xz/yz ↔ dxy) that are gapped with an excitation gap ∆ L ∼ 0.7-0.9 eV. L z (k, ω) are inter-orbital excitations between orbitals a and b (dxz ↔ dyz). (j-k) Orbital excitations spectra within the FM insulator competing state. L x/y (k, ω) has a gap ∼ 0.9 eV and L z (k, ω) is gapless with a robust q = π peak indicating quasi-long range staggered orbital ordering. Spectral functions are calculated using ∆ω = 0.05 eV, broadening η = 0.1 eV, up to 1200 DMRG states, 24 sites (d-g), 16 sites (h-k), and 8 DMRG sweeps.
we find that the orbital c has a gap of value ∆ c ∼ 1.5 eV [34] in A c (k, ω) (Figs. 2e and 3a), compatible with our analysis in Figs. 1b and 1c. The spectral function demonstrates again the coexistence of gapped Mott (localized) and gapless itinerant carriers and agrees well with an earlier study of a similar model [27]. Note also the presence of a pseudogap in the itinerant d xz/yz orbitals (not captured by earlier work [27]). It is likely that this pseudogap originates from correlations between the block ordered spins of Mott orbital c and itinerant electrons in orbitals a and b (Fig. 3a), since they are not decoupled from one another. We also highlight that the Fermi momentum of orbitals a and b is approximately q F π/4, leading to scattering at 2q F π/2 that is comparable to the spin-structure factor peak at q = π/2 corresponding to the block AFM order [25]. It is conceivable that the block phase of localized spins of the Mott orbital c is driven by the Fermi nesting of itinerant orbitals, as in manganites and heavy-fermion systems [25]. In fact, re-cent experimental [19] and theoretical [26] work have confirmed that the magnetic excitations in BaFe 2 Se 3 cannot be fully described using an effective Heisenberg model, highlighting the role of the other degrees of freedom in this material that we are focusing on here.
To arrive to our main conclusions, we now calculate the charge (Figs. 2f and 2g) and orbital (Figs. 2h and 2i) dynamics of the block-antiferromagnetic OSMP. Since the low-energy spin-dynamics have already been reported in the literature, we only focus on the high-energy charge and orbital dynamics that are more relevant to RIXS measurements at energy losses above 1 eV. As with A a/b (k, ω), the charge excitations of the itinerant orbitals display a gapless continuum (Fig. 2f) because charge fluctuations can propagate freely in a metallic system (at fixed-energy transfer ω many states allow for charge fluctuations without an energy cost). In contrast to this, the charge excitations of the Mott orbital c display a gap ∆ N c ∼ 1.5 eV (Fig. 2g), and N c (k, ω) shows incoher-ent charge excitations, compatible with doubly occupied orbitals d xy (i.e. c) configurations with an energy proportional to the Mott gap (Fig. 2a). We can contrast our results with the much studied one-orbital two-leg ladder Hubbard model, where charge excitations are gapped in the half-filled Mott phase but display a gapless continuum in the metallic phase away from half-filling [35]. Our results show the simultaneous presence of localized and itinerant fermions in our multi-orbital model. Consider now the q-resolved orbital excitations in the block OSMP. Figure 2h (L x/y (k, ω)) shows the excitations from (to) the itinerant orbitals a/b to (from) the localized orbital c. Figure 2i displays a gapless continuum in orbital excitations, L z (k, ω), between the itinerant orbitals a and b (see illustration Fig. 2b). These dynamical spectra can be understood using the density-of-states (Fig. 3a), and poles of the single-particle spectral function of the OSMP (Fig. 3b). Figure 3a shows that electron scattering from P 2 → P 4 requires an energy transfer ∼ 0.9 eV (vertical arrow in Fig. 3b). This scattering creates a gapped response peak in the orbital excitations at q 0 and ω 0.9 eV (Fig. 2h). The charge gap of orbital c in N c (k, ω) can also be understood via the scattering between points P 1 and P 4 of the density-of-states, shown by another vertical arrow in Fig. 3b. Additionally, the gapless continuum in L z (k, ω) is very similar to the itinerant charge excitations, N a/b (k, ω) because only the (almost degenerate) itinerant orbitals contribute in the calculations of L z (k, ω).
Dynamical spectra of the FM insulator -To identify features of orbital excitations in the block OSMP that are unique, we also show the q-resolved orbital excitations in the competing FM insulator known to exist at U/W = 4, [24] see Figs. 2j and 2k. In this FM phase, which is not yet known to be stable experimentally, theoretical calculations show that the single-particle spectral functions of all the orbitals have a Mott gap [27]. Figures 2j and 2k show that the corresponding L x/y (k, ω) and L z (k, ω) display clearly sharper excitations in the FM phase in comparison to the block OSMP. Remarkably, L z (k, ω) is gapless with a ω 0 peak at q = π, implying quasi-long-range antiferro-orbital order in the ground state, as discussed before [27]. In fact, the L z (k, ω) of this FM phase resembles the S (k, ω) of the spin-1/2 Heisenberg chain with a two-spinon continuum. By analogy, we conjecture that the L z (k, ω) of the FM phase may denote the existence of fractionalized orbital excitations [36] (see illustration Fig. 2c). Orbital and spin excitations were also studied recently in the context of Kugel-Khomskii models [37][38][39]. Here, we show novel results of orbital fractionalization in a FM insulator using a general multi-orbital fermionic Hamiltonian.
Discussion -Our main result is in Fig. 4, where we compare the calculated orbital dynamical response vs RIXS Fe-L 3 edge experimental data for BaFe 2 Se 3 at zero-momentum transfer [20]. Relating the full RIXS in- tensity to dynamical structure factors is not straightforward [31,40]. For example, for the single-band Hubbard model, the RIXS intensity can be directly related to dynamical structure factors only in limiting cases. To further complicate matters, the experimental effort [20] did not report the exact momentum transfer of the experiment. While the scattering geometry employed nominally probes q = 0 excitations along the chain, there were difficulties in aligning the chain orientation during beam time [41] and some finite q excitations may have been mixed into the spectra. Additionally, RIXS measures excitations of many different channels (e.g. spin, charge, and orbital), making it difficult to differentiate between those distinct channels. However, comparisons with the dynamical response functions can still provide at least qualitative insights into the dominant features and energy scales of RIXS experiments. Moreover, a merit of our effort is that we calculate excitations for each channel separately, and can, therefore, provide detailed predictions for future experiments. Until additional data become available our focus on RIXS experiments at q = 0 momentum transfer [20] is the natural avenue to pursue.
At zero energy transfer, there is an "elastic" peak commonly observed in RIXS experiments. We note that magnon contributions to the RIXS data generally occur at energy losses ω < 0.2 in the Fe-based superconduc- tors, and, therefore, are hidden in the elastic contribution in the RIXS experiment [19,20,26]. For this reason, this low-energy region will not be our focus. To identify the unique features of orbital excitations in the block OSMP (U/W = 0.8), we also present results for the PM (U/W = 0.04) and FM (U/W = 4) competing phases. Overall, we find a close agreement between our dynamical spectra in the block OSMP at q 0 and the two-peak structure observed in experiments, involv-ing localized dd excitations of the iron 3d orbitals. The OSMP two-peak structure (Fig. 4b) is distinct from that in the FM phase, where the two peaks are almost identical in position and width. With regards to the PM phase results (Figs. 4a), there are some similarities with the spectra shape of the OSMP phase. The presence of these similarities are natural since according to Ref. [24], at U/W = 0.04 the population of orbital c is already near 1, although the block magnetic phase is still not fully developed. One might argue that the PM regime data could also fit the experiments with appropriate broadening, particularly considering that high-frequency results involving interorbital excitations do not typically change abruptly across phase transitions; however, neutron scattering experiments tell us that the material is magnetically ordered in a block state. Thus the PM state cannot be chosen as the correct state even though the RIXS data might be consistent with our predicted dynamical spectra.
Our model suggests that the experimentally observed peaks represent excitations between localized d xy and itinerant d xz/yz orbitals. They are gapped with an activation energy 0.7 − 0.9 eV in the calculations, and ∼ 0.35 − 0.45 eV in the experiments (transitions involving e g orbitals, outside our model, should appear at higher energies). The difference in absolute gap numbers could be fixed by tuning our model parameters still within the OSMP phase (or by using ladders instead of chains, which is technically much harder, as discussed before). Since "fine tuning" is not our goal, but a qualitative understanding of results, we consider our (already costly) simulation results sufficient for our main qualitative conclusion: the block OSMP is the one that most closely resembles the experimental data. In addition, in the real samples, there is broadening due to phonons and non-crystallinity effects that could make it challenging to define the magnitude of the gap accurately. Regardless, it is clear from experiments that there is a gap in the RIXS q = 0 response. However, the paramagnetic phase has a gapless q = 0 response in L x/y (q 0, ω) (Fig. 4a), even-though there is a precursor to a two-peak structure, establishing another difference with the OSMP results. Moreover, while the ferromagnetic insulator shows a gapped L x/y (q 0, ω) response, this is with the presence of a single peak at ω 1 eV. In summary, the gapped two-peak response is only evident in the block OSMP phase.
Summary and Conclusions -Our results suggest that the experimentally observed gap is not a conventional semi-conducting gap, but instead originates from the inter-orbital excitations of a magnetically ordered OSMP. Figure 4b shows that peak A occurs at approximately 0.9 eV, for the zero momentum transfer, as a result of vertical (∆q = 0) scattering across E F from the itinerant d xz/yz orbitals to the d xy Mott orbital (Fig. 3b, P 2 → P 4 ). The shoulder/peak labeled as B represents scattering from the localized d xy band below E F to the itinerant hole pocket bands d xz/yz above E F , with zero momentum transfer q/π = 0 (Fig. 3b, P 1 → P 3 ). These peaks cannot be described using a simple weak-coupling framework.
The orbital d xy charge excitations in our calculations generate a response at ω = 1.4 eV of intensity ∼ 100 times smaller than the orbital excitations (Figs. 2f-2i). This is because local charge fluctuations are suppressed significantly in the Mott orbital, namely the probabilities associated with a doubly occupied or empty site configuration are small compared to configurations where sites are half-filled. Additionally, the features in N a/b (k, ω) and L z (q, ω) originate from itinerant carriers that are sensitive to the incident x-ray energy ( ω in ) of the RIXS experiments. It is known that localized excitations do not shift with ω in , while itinerant carriers produce a response that shifts linearly with ω in , becoming part of the fluorescence at large ω in [42]. Therefore, the L z (k, ω) peak at the low-energy transfer ω (Fig. 4b, dashed) will shift with the incident energy. The same is true for itinerant charge excitations N a/b (k, ω). In fact, RIXS experiment on BaFe 2 Se 3 also find (fluorescence) peaks that shift with incident energy and merge with the localized excitations (A and B of Fig. 4b), suggesting the existence of both localized and itinerant degrees of freedom at the experimental low temperature, as also found for the OSMP regime in our calculations.
To understand better the characteristics of the states at peak positions A and B, we also show real-space L x/y orbital correlations vs energy transfer at fixed distances (1 st to 4 th nearest neighbor) in Fig. 5. In the block OSMP the high energy states at peak A have positive L x/y correlations up to the 4 th nearest neighbor in realspace, while states at peak B have negative L x/y correlations only up to the nearest neighbor. In Fig. 5, we explicitly show that the peak A of the block OSMP represents states with ferro-orbital ordering, while B represents states with short-range anti-ferro-orbital fluctuations. However, the ground state of block OSMP does not have an orbital ordering unlike the ground-state of the competing ferromagnetic insulator at U/W = 4, which has uniform magnetic [24] and staggered L z orbital order (Fig. 2h). For the U/W = 4 FM case, the corresponding L x/y is gapped and has excitations only at ω ∼ 1.0 (Fig. 2h). In real-space, these high-energy states show short-range positive correlations, persisting only up to the 2 nd nearest neighbor (not shown).
Conclusions -We calculated the momentumresolved orbital and charge dynamics of an orbitalselective Mott phase using the three t 2g orbitals, and compared our results to the available RIXS data for BaFe 2 Se 3 at zero momentum transfer [20]. We find localized dd excitations that produce particular peaks A and B that are very similar to those observed experimentally, providing theoretical support that BaFe 2 Se 3 is indeed an orbital-selective Mott insulator. Moreover, we predict the q-resolved charge and orbital dynamical spectra that can be measured by RIXS in the future. Recent RIXS experiments display a very similar response for BaFe 2 Se 3 and BaFe 2 S 3 , implying that an orbital-selective phase is present in the BaFe 2 S 3 high-pressure superconductor as well [1,43]. We encourage the measurement of the band structure and density-of-states of these materials using angle-resolved photoemission and scanning tunneling microscopy, to confirm our predicted band(s) crossing of the Fermi surface and Mott like band(s) below the Fermi surface. We also encourage future q-resolved RIXS measurements for BaFe 2 Se 3 and BaFe 2 S 3 compounds. Finally, the novel orbital fractionalization proposed in the FM insulator phase defines a new avenue of research that will be pursued soon.

Methods
Model -We use a multi-orbital Hubbard model composed of kinetic energy and interaction terms as H = H K + H I . The kinetic energy part of the Hamiltonian is where c † iγσ (c iγσ ) creates (destroys) an electron at site i, orbital γ, and spin σ. The first term represents nearest-neighbor electron hopping from orbital γ to γ with a hopping amplitude t γγ . We denote the relevant orbitals d xz , d yz , and d xy as a, b, and c, respectively. The second term contains the orbital-dependent crystal field splitting. The parameters are (eV units) a = −0.1, b = 0.0, c = 0.8, t aa = t bb = 0.5, t cc = 0.15, and t ac = t bc = t ca = t cb = −0.10. The non-interacting bandwidth is W = 4.9t aa . This set of parameters is known [24] to produce bands that emulate iron-based superconductors, with hole pockets at q = 0 and an electron pocket at q = ±π [5]. Our reported results are all at zero temperature.
The interaction term, in standard notation, is The first term is the intraorbital Hubbard repulsion U . The second is the interorbital repulsion between electrons at different orbitals, with U = U − 2J H . The third term is the Hund's coupling J H , and the last term represents the on-site interorbital hopping of electron pairs (P iγ = c iγ ↑ c iγ ↓ ). Explicit definitions for all the operators in the model are in the supplementary material sections I and II. The rich physical properties realized by this model were discussed extensively in Ref. [24]. Most of the data presented here was gathered at J H /U = 0.25 (as used extensively before [24,44,45]) and U/W = 0.8 where the block-type AFM OSMP is known to be stable [24]. However, we also show results for the paramagnetic metal (PM, small U ) and ferromagnetic insulator (FM, large U ) phases to highlight unique features of the block OSMP by contrast.
Observables -To characterize the OSMP, we calculate the dynamical response functions using DMRG within the correction-vector formulation in Krylov space [46,47]. The single-particle photoemission spectral function is obtained using O i = c iσγ . The intraorbital particle-hole (charge) excitations arise using O i = σ n iσγ − n iσγ , where we explicitly subtract the ground-state contribution to measure only the fluctuations. Finally, the orbital excitations are obtained using These operators are derived from the L = 2 angular momentum operators in the t 2g orbital basis, see supplemental material section II. Code Availability -Computer codes used in this study are available at https://g1257.github.io/ dmrgPlusPlus/.
Data Availability -The data that support the findings of this study are available from the corresponding author upon request.
Supplemental: Fingerprints of an Exotic Orbital-Selective Mott Phase in the Block Magnetic State of BaFe 2 Se 3 Ladders

OPERATORS AND OBSERVABLES
In this section, we define all the operators that are used in the main text: where γ is the orbital index, σ and σ are the spin index, and σ κ are the Pauli matrices with κ = {x, y, z} being Cartesian components. The average occupation and charge fluctuations ( Fig.1 of main text) are defined as where L is the number of sites in the lattice. The static spin-spin correlations (inset of Fig. 1) is calculated using where S i = γ S iγ . In this supplemental, we also show results for n γ k defined as c jγ = c jγ↑ + c jγ↓ , (S4)

Spectral functions and sum rules
To characterize the OSMP, the dynamical response functions (shown below) are calculated using DMRG where the local operator O i can represent any degree of freedom of the model. In general, these functions are Fourier transformed into the crystal momentum domain to calculate the momentum-energy resolved spectra that is relevant to experiments: Note that within DMRG, the site j is fixed to the center of the lattice (d = L/2 − 1) to reduce the edge effects and computational cost, and therefore the modified Fourier transform becomes Additionally, we use open boundary conditions in the DMRG simulation and therefore the quasi crystalmomenta are defined as The dominant orbitals of an iron atom in iron-based superconductors are the five 3d orbitals, corresponding to well-known L = 2 orbital angular momenta. The corresponding operators {L x , L y , L z } are written in the basis of L z = {−2, −1, 0, 1, 2}, forming 5 × 5 matrices: DMRG is a real-space algorithm, and therefore we must write these matrices in the orbital basis using the transformation where {|x 2 − y 2 , |z 2 , |xz , |yz , |xy } are the five iron 3d orbitals. The operators {L x , L y , L z } represented in the orbital basis are: In the main text, we employ the iron t 2g orbitals, i.e., we drop the contribution from |x 2 − y 2 , |z 2 orbitals. After this approximation, we obtain the local orbital angular momentum operators used in the main text (site index understood) (S10)

ADDITIONAL RESULTS
In this section, we show tests performed to ensure the quality of the presented results. We first show n γ (k) that is in agreement with previous studies with a similar model [27]. A significant change in n a/b (k) indicates metallic behavior of these orbitals. On the contrary, n c (k) shows little variation in k that is associated with a gapped orbital c.
Furthermore, all spectral quantities defined by equations S5 and S6 must satisfy a sum rule. In general, it can be shown that integrating the spectral function over momentum k and energy ω gives a quantity related to a unique static observable. As an example, we use the single-particle spectral function A c (k, ω) of orbital c below (above) the Fermi energy (E F ) representing the filled n c (k) of the localized orbital c has only a slight variation in k, suggesting a gap in the single-particle states of orbital c. These results are in agreement with previous studies using Quantum Monte Carlo [27].
Im ψ 0 |c ic 1 ω + H + E g + iη c † jc |ψ 0 e −ik(i−j) , where c is the orbital index and c jc = c jc↑ + c jc↓ . To obtain the sum rule of this quantity, we first sum over the momentum k. This sum simply results in Lδ i,j , giving us the local response that is the single-particle density of states, shown in Figure S2. Further integration over ω of the imaginary part (Lorentzian poles) of the electron (hole) part simply gives the total electron (hole) density of orbital c: This is explicitly done within our calculations by integrating over the density of states (Fig. S2). The integration of the electron (blue curve) and hole (red curve) portions lead to approximately 1.0 that is consistent, within the accuracy of our results, with calculations of the local density from the ground-state. Note that a singly occupied orbital is one of the characteristics of a Mott insulator.
We also performed finite-size scaling on the density of states at the Fermi energy E F 4.3 eV. Figure  S3 shows that the L → ∞ extrapolated quasi-particle weight, A(ω = E F ), of orbital c is an order of magnitude smaller than the almost degenerate orbitals a and b. The near zero weight of orbital c at the Fermi energy is consistent with a Mott phase, providing further evidence of the presence of an OSMP as the ground-state. Finite-size scaling of the orbitalresolved electron part of the density of states at the Fermi energy (EF ). The quasi-particle weight of the Mott orbital approaches 0 (more accurately, 0.003) with increasing system size while the weight remains finite for the itinerant orbitals (0.04 and 0.07). Note that Ae(EF ) of the itinerant orbitals a and b is an order magnitude larger than the Mott orbital c, further emphasizing that orbital c is a Mott insulator.

REPRODUCING DATA USING DMRG++
The full open source code, sample inputs, and corresponding computational details can be found at https: //g1257.github.io/papers/86/.