Emergent Half Metal at Finite Temperatures in a Mott Insulator

Sustaining exotic quantum mechanical phases at high temperatures is a long-standing goal of condensed matter physics. Among them, half-metals are spin-polarized conductors that are essential for realizing room-temperature spin current sources. However, typical half-metals are low-temperature phases whose spin polarization rapidly deteriorates with temperature increase. Here, we first show that a low-temperature insulator with an unequal charge gap for the two spin channels can arise from competing Mott and band insulating tendencies. We establish that thermal fluctuations can drive this insulator to a half-metal through a first-order phase transition by closing the charge gap for one spin channel. This half-metal has 100% spin polarization at the onset temperature of metallization. Further, varying the strength of electron repulsion can enhance the onset temperature while preserving spin polarization. We outline experimental scenarios for realizing this tunable finite temperature half-metal.

S pin-polarized metals are of great technological importance, having applications in diverse fields ranging from spintronic devices to quantum computation 1,2 . Thus material realizations of half-metals continue to be a significant research direction on both theoretical [3][4][5][6][7][8][9][10][11][12] and experimental fronts [13][14][15] . Candidates for half-metals such as double perovskites [4][5][6][7] , Heusler alloys [8][9][10][13][14][15] , and manganites 16,17 , have therefore been extensively investigated, among others. Introducing vacancies in transition metal oxides such as NiO, MnO 18 , fluorinated BN, ZNO 19 , and ferrimagnetic inverse spinels such as Fe 3 S 4 20 have also been proposed as systems that can host half-metallicity. However, the degradation of spin polarization with temperature increase and experimental challenges in material growth 21 has long been barriers to producing spin current sources at room temperature in all half-metal candidates. Thus, predictions of half-metallic states in simple models are vital for uncovering novel mechanisms for half-metallicity and subsequent material realizations that can potentially circumvent the aforementioned difficulties.
In this context, we consider the ionic Hubbard model (IHM), with local repulsion U and staggered onsite "ionic" potential Δ, introduced initially for studying charge transfer effects in organic crystals 22 and for modeling ferroelectric perovskites 23 . In the half-filled IHM for U = 0, a finite Δ opens up a band-gap that promotes large occupation of the sub-lattice with −Δ potential. In the opposite limit, for Δ = 0 and large U, singly occupied sites are preferred on both sub-lattices. The simultaneous presence of the two agencies leads to a competition that causes a zero temperature band insulator to Mott insulator transition with increasing U, for a fixed Δ. This band to Mott insulator transition has been demonstrated in one dimension [24][25][26][27][28] and using dynamical meanfield theory (DMFT) on square 29,30 and Bethe lattice 29,31,32 . More importantly, the half-filled model in one dimension 33 and on square lattice using DMFT 34 have been shown to support zero temperature half-metallic state on introducing next nearest neighbor (nnn) hopping. However, transport response and survival of spin polarization of this half-metal, in presence of thermal fluctuations, has remained largely unaddressed.
In this article, by computing spin-dependent transport response at finite temperatures, we first establish a lowtemperature antiferromagnetic half-metal (HM 1 ), bracketed by a paramagnetic band insulator for U < U B and an antiferromagnetic Mott insulator for U > U M . The spin polarization of HM 1 is maximized as T → 0. Next, we reveal a ferrimagnetic halfmetal (HM 2 ) that emerges out of thermal excitations of the Mott insulating state for U values in the vicinity of U M . We show that, unlike HM 1 , the Mott insulating ground state ensures that halfmetallicity in HM 2 occurs only above an onset temperature via a thermally driven first-order insulator to metal transition. We demonstrate that at the onset temperature of metallization, HM 2 is fully spin-polarized. We find that the onset temperature is highly sensitive to correlation strength, acting as a control knob for tuning the temperature window of half-metallicity. Raising the onset temperature does not degrade the spin polarization in its vicinity. We explain the origin of HM 2 and map out the parameter space of its stability. We close with a discussion on the experimental realization of HM 2 in correlated oxides and cold atomic systems.
We define the IHM on a square lattice as follows, where, c y iσ (c iσ ) are electron creation (annihilation) operators at the site i with spin σ. t and t 0 are respectively the nearest and nnn hopping amplitudes. We choose t 0 =t<0 in our study. n iσ ¼ c y iσ c iσ is the number operator for spin σ at a site i and n i is the spin summed local number operator. Δ is the magnitude of "ionic" potential and takes positive (negative) values on the A (B) sublattice. U is the local Hubbard repulsion. μ denotes the chemical potential, and is adjusted to maintain half-filling. We employ a recently developed semiclassical Monte Carlo approximation (s-MC) [35][36][37] to investigate the finite temperature properties of IHM on large lattice sizes.
"Methods" sub-sections titled, Treatment of the interaction term, Solution strategy, and The nature of approximation cover the (s-MC) scheme in detail. In particular, in the "Methods"-The nature of approximation, it is discussed that s-MC is numerically inexpensive, can capture results beyond finite temperature mean-field theory, and has a reasonable quantitative agreement with Determinantal Quantum Monte Carlo (DQMC) over a wide temperature range 35,[37][38][39] .

Results
Transport response regimes. In Fig. 1, we show the spin-resolved resistivity, ρ σ (T), extracted from the low-frequency optical conductivity, as discussed in the "Methods"-Optical conductivity. We show the case of low U band insulator in Fig. 1a Fig. 1b and Fig. 1c that, as will be detailed below, support HM 1 and HM 2 respectively. Finally in Fig. 1d we provide an example of the large U Mott insulator. In all the plots, triangles and squares denote ρ ↑ (T) and ρ ↓ (T) respectively. We first compare the small U band insulator in Fig. 1a with the large U Mott insulator in Fig. 1d. In Fig. 1a, we find typical insulating nature (dρ σ /dT < 0) for both spin channels. We can see Fig. 1 Temperature evolution of spin-dependent transport. a-d Show spin-resolved resistivity, ρ σ (T) for correlation strength U values as indicated. ρ σ (T) is calculated in units of (πe 2 /ℏa 0 ), with e, ℏ and a 0 being the electronic charge, the Planck's constant divided by 2π and the lattice spacing respectively. The next nearest neighbor hopping t 0 =t and charge transfer energy Δ/t are chosen to be −0.2 and 1.0, respectively. The squares denote ρ ↓ and triangles represent ρ ↑ . The inset in b shows the low T(= 0.005t), spin-resolved density of states (DOS), N σ (ω) for U/t = 3.2, with a small arrow at the bottom indicating the Fermi energy. The up and down spin channels are also indicated by arrows. The inset in c shows the total resistivity ρ(T) for U/t = 4.8. d ρ σ (T) for the robust Mott insulator at U/t = 5.5. Results are shown for 32 2 system size with periodic boundary conditions. Standard deviation of the lowfrequency optical conductivity data is used to estimate error bars in the resistivity plots. Symbol sizes used are larger than the error bars.
in Fig. 1a that dρ/dT remains negative for all temperatures. This is because the correlation strength is insufficient to screen the one body potential and close the charge gap as temperature increases. Similar behavior was found in DQMC study of IHM for the t 0 ¼ 0 40 . For the chosen values of Δ and t 0 , this band insulating response occurs for U < U B (=2.8t), while the Mott insulator lies beyond U M (=4.75t). We see typical Mott insulating behavior representative of large U(=5.5t) in Fig. 1d with diverging resistivity as T → 0. Further, we see that the sign of dρ σ /dT for both spin channels changes from negative to positive around T/t~0.12, signaling a crossover from an insulator to a metal, mediated by thermal fluctuations. This temperature scale is the analog of T * in the standard half-filled Hubbard model, with zero Δ and t 0 , that governs the crossover from a paramagnetic insulator to a paramagnetic metal, previously seen in DQMC 41 and s-MC 35 .
In contrast to these limiting cases, in the regime U B < U < U M and U M < U < 5t, we find distinctly different thermal evolution of resistivity. Representative data for these two regimes are shown in Fig. 1b, c respectively. In Fig. 1b, spin resolved resistivity for U/t = 3.2 shows that dρ ↓ /dT is negative for T < 0.04t, while dρ ↑ /dT exhibits metallic behavior for all temperatures. At T = 0.04t, the of sign dρ ↓ /dT changes and ρ ↑ becoming equal to ρ ↓ . We define T/t = 0.04 as T P , the lowest temperature where ρ ↑ and ρ ↓ become identical. We will later show that below T P , we have a spinpolarized metal, whose polarization approaches 100% as T → 0. We refer to this half-metallic regime as HM 1 and T P define a deporalization temperature scale, above which spin polarization is completely lost.
In the resistivity data shown in Fig. 1c for U/t = 4.8, we find diverging ρ σ as T → 0, as is expected in a Mott state. However, at T/t = 0.02, dρ ↑ /dT switches sign while dρ ↓ /dT continues to remain negative up to T/t = 0.09. From the inset in Fig. 1c, showing the total resistivity as a function of temperature, we identify T/t = 0.02 as the onset of metallic response, where dρ ↑ /dT had changed its sign. We denote this onset temperature by T S . Also as for HM 1 , we denote T/t = 0.09 as T P , the temperature above which ρ ↑ (T) = ρ ↓ (T). Thus, for U ≳ U M , we find a window in temperature, T S < T < T P , where the Mott insulating ground state gives way to a finite-temperature half-metal, HM 2 . We will show that, within numerical accuracy, HM 2 is fully spin polarized near the onset temperature T S .
Density of states and spin polarization. The inset in Fig. 1b shows the spin-resolved density of states (DOS), N σ (ω), defined in the "Methods"-The spin-resolved DOS. The result is shown for U/t = 3.2 and at the lowest temperature of our calculation (T/t = 0.005). We find that only one "up" spin channel contributes to the weight at Fermi energy for very low temperature, as would be expected for usual half-metals. In contrast, lowtemperature DOS for U/t = 4.8 in Fig. 2a shows a finite charge gap that is unequal for the two spin channels. Consequently, we find that the temperature-induced filling of the charge gap is much more rapid for the "up" than for the "down" spin channel. We see this behavior in Fig. 2b, which shows the temperature evolution of N σ (0), the spin-resolved spectral weight extracted from the DOS, at ω = μ. The filling of the charge gap starts at T S (=0.02t) in a "spin asymmetric" manner, and the asymmetry persists up to T P (=0.09t). We discuss the origin of this spinasymmetric filling of the charge gap later in the paper. We add that the Mott state at low temperature is always characterized by an unequal charge gap for the spin channels for all U values investigated in our study. However, beyond U/t = 5, the temperature-induced filling of the charge gap occurs equally for both spin channels as seen in Supplementary Fig. 1a presented in Supplementary Note 1, and prohibits the formation of the HM 2 phase. The low-temperature s-MC DOS qualitatively agrees with DMFT calculations at T = 0 34 .
To quantify the spin polarization of the half-metals at finite temperature, we define the transport spin polarization as where J x σ is the spin resolved current operator along the x-direction. Here and later in the paper, the angular brackets indicate quantum as well as thermal averaging. Details of the calculation are provided in "Methods"-Optical conductivity. Fig. 2c, shows P(T) for U/t = 3.2, 4.8 and 4.9. For HM 1 at U/t = 3.2, we find 100% spin polarization or P(T)~1 at T = 0.005t. With temperature increase, the "down" spin channel that was gapped in the ground state, e.g. the inset of Fig. 1b, starts to fill up, causing the polarization to decay, eventually suppressing it to zero beyond T P (=0.04t). P(T) for HM 2 , at U/t = 4.8, is shown by circles. We find that P(T) = 0 in the insulating regime, for T < T S (=0.02t). It then abruptly jumps to 1 at T S , the onset of half-metallicity. With temperature increase beyond T S , P(T) drops until it reaches zero at T P = 0.09t, coinciding with the closing of the window of asymmetric charge gap filling, as seen in Fig. 2b. We also see that while T S increases from 0.02t for U = 4.8t to 0.06t for U = 4.9t, P(T) at T S , remains unchanged. We thus find an unexpected property that T S enhancement does not sacrifice spin polarization for HM 2 . We note that at T S , for all parameter points of HM 2 , the insulating to metallic spin channels resistivity ratio is~10 2 −10 3 . This translates into P(T)~1 or 100% spin polarization at T S , as can be seen from the expression of P(T) in the "Methods"-Optical conductivity.
Ferrimagnetic half-metal at finite temperature. To define magnetic order at finite temperature in two dimensions, we add a small SU(2) symmetry breaking magnetic field, as described in the "Methods"-The treatment of the interaction term. In Fig. 3a, we show the thermal evolution of the average sub-lattice magnetizations, S z α ðhn α" i À hn α# iÞ=2, α ∈ {A, B} for U/t = 4.8. We present the formulae for spin and sub-lattice resolved densities in "Methods"-The spin and sub-lattice resolved density. For T < T S , we find that S z A ¼ ÀS z B , indicating long-range antiferromagnetic order. However, for T S < T < T P , we see that jS z A j > jS z B j, while the staggered nature of magnetic order remains unchanged.
In Fig. 3b, we show the temperature evolution of the difference between the average sub-lattice local moments (δM = M A − M B ), where M α = 〈n α 〉 − 2〈n α↑ n α↓ 〉, and α ∈ [A, B]. We find that δM = 0 for T < T S , but has a finite value in the window T S < T < T P . Fig. 3c shows the net spin polarization of the system for U/t = 4.8 with triangles, δn = 〈n ↑ 〉 − 〈n ↓ 〉. The correlation between these δM and δn is apparent from Fig. 3b and c. The temperature dependence of the sub-lattice magnetizations in Fig. 3a follows from the relation S z A þ S z B ¼ δn=2. Due to the unequal magnitude of the staggered sub-lattice magnetizations, we refer to HM 2 as a ferrimagnetic half-metal. We also observe that, as a generic feature of HM 2 , δn, δM, and P(T) undergo an abrupt jump at T = T S , as seen here for U/t = 4.8. This behavior strongly indicates a thermal fluctuation-induced "first-order" transition from an antiferromagnetic Mott insulator to a ferrimagnetic half-metal. In Supplementary Fig. 1b reported in Supplementary Note 2, we provide typical hysteresis behavior for δM associated with the first-order transition, obtained through a cooling and heating protocol for HM 2 . Since our calculation uses fixed temperature steps that only bracket the actual critical temperature for the insulator to metal transition, we refer to the temperature scale of metallization as an "onset" scale. We find no discontinuity in the temperature evolution of sub-lattice magnetization for HM 1 as seen in Supplementary Fig. 1c presented in Supplementary Note 3.
To understand the above temperature-driven behavior, we first note that finite Δ causes a double occupation penalty of U + 2Δ (U − 2Δ) on the A (B) sub-lattice in the zero hopping limit. This differential penalty persists for finite hopping, albeit with a renormalized value. Supplementary Fig. 1d in Supplementary Note 4 shows typical temperature evolution of the sub-lattice double occupations 〈n α↑ n α↓ 〉 with α ∈ (A, B) for HM 2 for U = 4.8t. It clearly shows that, 〈n B↑ n B↓ 〉 > 〈n A↑ n A↓ 〉 at all temperatures. In Fig. 3a, we see that the larger cost of double occupation at A, initially triggers a finite magnetization S z A at T P ( = 0.09t), while S z B becomes non-zero only at a lower temperature. It also maintains jS z B j < jS z A j and consequently a finite δn following S z A þ S z B ¼ δn=2, below T P and down to T S . The sub-lattice local moments difference δM in Fig. 3b, also follows from this differential double occupation penalty. However, energetically, the large build-up of 〈n B↑ n B↓ 〉 becomes untenable at low temperatures. Below T S , the difference in the sub-lattice double occupations is reduced through a first-order transition, making the sub-lattice magnetizations equal, which quenches δn and limits the half-metallic state HM 2 to finite temperature. From Fig. 1c, we find that the conduction electron spin direction is identical to the net spin polarization δn, 'up' in this instance. It also translates into N ↑ (0) > N ↓ (0) for T S < T < T P , and can be interpreted as thermal fluctuation driven spinasymmetric filling of charge gap, as seen in Fig. 2b.
U-T phase diagram. We now discuss the U − T phase diagram, Fig. 3d, to investigate the stability of HM 2 and the tunability of T S . For T/t ≤ 0.02, U M (=4.75t) separates HM 1 and M-I phases. Figure 3d is a part of the full phase diagram shown in Supplementary Fig. 2 and discussed in detail in Supplementary Note 5. As seen in Supplementary Fig. 2, HM 1 starts at U = 2.8t and the T P (triangles) increases with U, reaching a maximum of (0.055t) close to U M . From static magnetic structure factor discussed in the "Methods"-Static magnetic structure factor, we find the HM 1 has an antiferromagnetic background and the magnetic transition temperature T N , coincides with T P . Supplementary Fig.  1c discussed in Supplementary Note 3 depicts typical temperature dependence of the sub-lattice magnetizations for HM 1 . The HM 2 phase emerges from the M-I in the range U M < U < 5t and above a U-dependent T S , indicated by diamonds. The ferrimagnetic order in HM 2 is destabilized at the corresponding T P values (squares), as seen for example, in Fig. 3a for U = 4.8t. Circles represent the antiferromagnetic transition temperature for the M-I phase beyond the HM 2 window or U > 5t.
We also see that T S increases rapidly with U for HM 2 . This behavior follows from the need of greater thermal energy to close the Mott gap at larger U values. However, as U grows, we also move deeper into the Mott state, progressively suppressing the effect of Δ. The increasingly dominant role of U manifests as a systematic reduction in magnitudes of δM and δn as is seen in Fig. 2 Thermal evolution of density of states (DOS) & polarization. a This shows the low T density of states, N σ (ω) for U/t = 4.8. The large arrows denote the spin channels. b It shows thermal evolution of N σ (0), the spectral weight at the chemical potential for U/t = 4.8. Here triangles denote N ↑ (0) and squares indicate N ↓ (0). The locations of the onset temperature of half metallicity T S and depolarization temperature T P are marked with black arrows. c Shows spin polarization P(T) as defined in Section 2 of the main text titled "Density of states and spin polarization", for indicated U values. Various T S and T P are demarcated with arrows color coded with the symbols. Lines are a guide to the eye. All results are shown for t 0 =t ¼ À0:2 and Δ/t = 1.0. Error bars in b and c are computed from the standard deviation in the spin-resolved DOS, and low-frequency optical conductivity data respectively. Symbol sizes used are larger than the error bars. Show the difference in sub-lattice local moments (δM) and the total spin polarization of the system (δn) respectively, for U/t = 4.8 4.9 and 4.95. In each case, the dashed lines between two temperature points indicate the location of the onset temperature for half-metallicity. d This shows the U − T phase diagram with a low-temperature antiferromagnetic half-metal HM 1 to antiferromagnetic Mott insulator (M-I) transition at U M = 4.75t. Triangles indicate T P , the depolarization temperature of half metallicity for HM 1 . The red shaded region refers to the finite-temperature ferrimagnetic half-metal HM 2 for U M < U < 5t and bounded by the onset temperature of polarization T S and the depolarization temperature T P . The T S (T P ) scales for HM 2 are demarcated by diamonds (squares). Circles refer to T N , the temperature above which the M-I state loses long-range antiferromagnetic order. All results are for t 0 =t ¼ À0:2, Δ/t = 1.0 and on 32 2 system. Standard deviations in the sub-lattice magnetization in a and the various ordering temperature scales in d are smaller than the symbol sizes; Error bars for U = 4.8t are generated from standard deviations the sub-lattice local moments in b and spin-resolved densities in c. Error bars for U = 4.9t and U = 4.95t are similar in magnitude to that for U = 4.8t in b and c. Fig. 3b, c. In these plots, we compare the temperature profiles of δM and δn for three U values. We find that while with U increase, the window of finite δM and δn move to higher temperatures, their magnitudes are systematically suppressed. As the sub-lattice moment magnitudes approach each other with increasing U, T P and T S merge into a single transition temperature T N , of the (M-I) phase beyond U = 5t. Finally, from Fig. 3d, we see that T S increases by 350% over the window in U for which HM 2 is stable. In the discussion section, we propose 4d transition metal element-based double perovskites as potential candidates for realizing HM 2 . We obtain a realistic estimate of the T S and T P using a hopping scale of 0.2-0.3 eV for double perovskite. It yields T S to be about 70-100 K and T P about 230-350 K in the middle of HM 2 in Fig. 3d at U/t = 4.875.
In concluding this section, we would like to briefly discuss the role of t 0 in stabilizing HM 2 . t 0 acts as a source of frustration that inhibits formation of long-range antiferromagnetic order and promotes metallic tendency. In fact, both half-metallic phases are absent for jt 0 =tj < 0:05. Increase in the magnitude of t 0 =t beyond 0.05, systematically shifts U M to larger U values to overcome frustration effects and consequently increases T S as seen in Supplementary Fig. 3a discussed in Supplementary Note 6. However, once U M gets too big for a large enough jt 0 =tj, the reduced impact of Δ is inadequate to support HM 2 , i.e., the temperature-driven gap-filling becomes spin symmetric, as previously described. We find that HM 2 is no longer realized for jt 0 =tj ≥ 0:3, at Δ/t = 1. In Supplementary Fig. 3b, presented in Supplementary Note 6, we show preliminary results from an exhaustive numerical search, from which we infer that HM 2 can be stabilized for Mott insulators with Δ/U~0.2-0.3, over a large window of U and Δ. Thus, HM 2 has a broad stability window in the U À Δ À t 0 parameter space, and Fig. 3d depicts only a fixed Δ À t 0 cross-section of this regime.

Discussion
Unlike all previous proposals where half metallicity is a ground state property, here we have demonstrated that a half-metal can emerge from a Mott insulating ground state. The insulating ground state ensures that the half-metal with 100% spin polarization occurs at a finite temperature and protects against temperature-induced depolarization effects. The correlation strength-dependent enhancement of the onset temperature naturally carries the fully spin-polarized metal to higher temperatures without depolarizing. We want to emphasize that the low-temperature half-metal HM 1 , has negligible δn and zero T S due to the reduced effect of U and, in this respect, is different from HM 2 . Lastly, HM 2 does not require interaction between large Hund's coupled core spins and itinerant electrons like in the double exchange mechanism; it is different from the Slater-Pauling rule-based, low-temperature half-metallic Heusler alloys and T = 0 half-metallicity in doped Mott insulators away from half-filling 42 . The experimental realization of this unique finite-temperature half-metal will be a significant step toward the goal of creating spin-polarized current sources at room temperature. With this in mind, we close with a discussion on the possible realization of HM 2 in solid-state and cold atomic experimental setups.
For correlated oxides representing the IHM, we consider epitaxial thin films of cubic double perovskites, X 2 ABO 6 . Here, A and B represent two species of 4d transition metal (TM) atoms, and X can be rare-earth, alkaline-earth, or alkali elements. To prevent high U and Hund's coupling effects, we do not consider 3d TM elements. In addition, we avoid 5d TM elements with strong spin-orbit coupling. In the 4d TM series, {Zr, Nb, Mo, Tc, Ru and Rh} starting with Zr, difference in the charge transfer energy (Δ) between successive elements range between 0.2 eV and 0.7 eV 43 . Also the local correlation strength U is moderate, ranging between 1eV and 3eV 44,45 . Hence, Δ/U~0.3 required for (HM 2 ), can be achieved. In an octahedral crystal field, the 4d orbitals are split by 3-4 eV 46 , into high energy e g and low energy t 2g levels. We expect that (A, B) combinations of 4d TM elements with a small Δ and relatively large crystal field splitting will facilitate the formation of partially filled t 2g bands. Suitably chosen X site element can achieve half-filled t 2g manifold and also can tune U by controlling the electronic bandwidth. Finally, nnn hopping is relevant for 4d TM elements 47,48 . Hence one can potentially identify a number of candidates that are insulators with small charge gap, such as Sr 2 RuMoO 6 49 . Based on the abovementioned literature we can crudely estimate U Re/Mo~3 eV, Δ~1−2 eV and a t 2g − e g splitting of 3 eV. These estimates yield Δ/U~0.3−0.6 eV, which is within the ballpark of the ratio needed for HM 2 for Sr 2 RuMoO 6 . We expect complete quenching of orbital angular momentum in a half-filled t 2g 4d orbital; however, one has to investigate the role of small Hund's coupling 50 , prior to identifying realistic material candidates.
While the model has not yet been realized on a square lattice in the cold atomic setup, the experimental techniques for such a realization exist. For example IHM has already been realized for fermionic cold atomic systems for hexagonal lattice 51 . This experiment investigated the band to Mott insulator transition by tracking suppression of double occupation over wide variations of Δ and U [0, 41t] and [0, 30t], respectively, in the units of nearestneighbor hopping strength t/h~174 Hz. This range of parameter variation should be easily possible for the square lattice as well. t 0 =t ratio can be tuned over a large window on shaken optical lattices 52 . Also, the highest T S /t~0.1 is within the ballpark of the experimental temperature scales of 0.2t 53,54 . We suggest that, δn would be the natural quantity to measure as a signature of HM 2 , analogous to spin polarization measurements for metallic (Stoner) ferromagnets in cold atomic systems 53 .

Methods
We briefly present the formal derivation of the s-MC approximation scheme, calculation method, nature of the approximation, and definitions of observables.
Treatment of the interaction term. To set up the semiclassical Monte-Carlo s-MC method, we first decouple the interaction term in Eq. (1) using standard Hubbard-Stratonovich (H-S) transformation, by introducing auxiliary fields (Aux. F.). For this we express the interaction term at each lattice site as the sum of square of the local number operator n i and the local spin operator S i as follows: here, the spin operator at the ith site (running over A and B sub-lattices) is S i ¼ 1 2 ∑ αβ c y iα σ αβ c iβ , with ℏ ≡ 1 andΩ is an arbitrary unit vector. The partition function for the IHM Hamiltonian, H ≡ H 0 + H 1 is Z = Tre −βH with H 0 and H 1 containing the one body and interaction terms respectively. The trace is taken over the occupation number basis. β = 1/T, is the inverse temperature where k B is set equal to 1. The window [0, β] is divided into M equally spaced slices separated by interval Δτ, with β = MΔτ. The slices are labeled by l. In the limit Δτ → 0 by using Suzuki-Trotter decomposition, we write e ÀβðH 0 þH 1 Þ % ðe ÀΔτH 0 e ÀΔτH 1 Þ M to the first order in Δτ. For a given imaginary time slice 'l', the partition function for the interaction term can then be written as, const: Here, we have introduced (Aux. F.)'s ϕ i (l) and ξ i (l) at each site and imaginary time slice. Following literature 55  slices run from l = M to the l = 1. We note that at this stage the partition function is exact, manifestly SU(2) symmetric and the {ϕ i (l), m i (l)} fields fluctuate in space and imaginary time. We now make the approximation of (i) retaining only the spatial dependence of the (Aux. F.) variables and (ii) use a saddle point value of the (Aux. F.) ϕ i (l) equal to iU〈n i 〉/2. This allows us to extract the following effective spin-fermion type Hamiltonian H eff where fermions couple to classical (Aux. F.) {m i }. For further details, we refer to our earlier work 35 .
In H 0 , we include a small SU(2) symmetry-breaking magnetic field term. The field avoids well-known Mermin-Wagner issues and allows for an unambiguous definition of long-range magnetic order in two dimensions at finite temperature. It also provides a global axis for the definition of up and down spin components. We have checked that the spin polarization in the band insulator and the HM 1 phases are zero within numerical accuracy, establishing that the small symmetry breaking term does not induce spurious spin polarizations. We solve the spin-fermion model by well-established cluster-based exact diagonalization coupled with classical Monte-Carlo (ED+MC) method 36,56 . Below we briefly outline the main steps of the ED+MC and refer to our earlier work detailing the solution methodology 35 .
Solution strategy. We start the calculation at high temperature T/t = 1 with random {m i } and {〈n i 〉} fields at each site. The {〈n i 〉} fields define the {ϕ i } variables, as mentioned above. We then sample the {m i } fields by sequentially visiting all lattice sites. For this, we first diagonalize H eff for the chosen (Aux. F.) configuration and then propose an update at a site. To decide the acceptance of a proposed update for m i at the ith site, we employ a standard Metropolis scheme. In this approach, the usual Boltzmann factor governs the acceptance probability, which depends on the energy from the exact diagonalization before and after the proposed update. This update process done sequentially at all lattice sites constitutes a Monte-Carlo system sweep. We conduct 6000 sweeps at each temperature. After every 20 sweeps, we perform a self-consistency loop to converge the {〈n i 〉} fields for a fixed {m i } configuration. We leave the first 2000 sweeps for equilibration and compute observables from the remaining 4000 sweeps, skipping every five sweeps to avoid self-correlation effects. We reduce the temperature in small steps and repeat the above protocol at each temperature. Finally, we average all observables over 20 separate runs with independent random initial choices of {m i } and {〈n i 〉} configurations at the starting (high) temperature.
Nature of the approximation. Earlier work 35 has established that at low temperatures, the Monte Carlo sampling leads to uniform m i configurations akin to unrestricted Hartree-Fock results. However, at higher temperatures where thermal fluctuations begin to dominate quantum fluctuations, the thermal sampling of the (Aux. F.)'s {m i } captures temperature effects considerably more accurately than a simple finite temperature mean-field theory, making s-MC a progressively superior approximation. s-MC based study of the half-filled Hubbard model in two and three dimensions 35,36 , and with nnn hopping 37 , have revealed that this approach can capture several results that are beyond simple finite temperature mean-field theory. These include non-monotonic dependence of antiferromagnetic ordering scale on U, finite T paramagnetic insulating phase with local moments, pseudo-gap to normal-metal crossover, and specific heat systematics with temperature. Moreover, these results have a reasonable quantitative agreement with DQMC 38,39 . s-MC is free of the usual fermion sign problem, analytical continuation issues and is numerically inexpensive. Recent s-MC based studies include finite-temperature properties of the Anderson-Hubbard model 57 , and temperature-driven Mott transition in frustrated triangular-lattice Hubbard model 58 .
Observable definitions. Here we define the various observables used in the paper.
Optical conductivity. The total d.c conductivity σ dc along the x-direction, is computed from the Kubo-Greenwood formula 59 for optical conductivity.
The spin resolved DOS. The spin-resolved DOS is defined as follows: N σ (ω) = ∑ γα |〈α, σ|ψ γ 〉| 2 δ(ω − ϵ γ ), where ϵ γ and jψ γ i are the eigenvalues and eigenvectors of H eff . Here α ∈ {A, B}, the two sub-lattices and σ refers to spin. Lorentzian representation of the above δ−function is used to compute DOS. The broadening of the Lorentzian is~BW/2N with BW being the noninteracting bandwidth and N is the total number of lattice sites.
The spin and sub-lattice resolved density. The average spin and sub-lattice resolved density is given by where, i is the site index. α ∈ A, B, and σ is the spin index. ϵ γ , jψ γ i and f(ϵ λ − μ) are as defined above. N is the total number of lattice sites. The sub-lattice resolved occupation is calculated by summing over densities for A and B, for each spin channel independently. The magnetization in Fig. 3a are constructed from these densities.
Static magnetic structure factor. S q as defined below, is computed for q = (π, π), from which T N is extracted for the various phases.
All symbols have the usual meaning as defined above and i, j run over all lattice sites.

Data availability
The data that support the findings of this study are available from the corresponding author upon request.

Code availability
Codes used to produce the findings of this study are available from the corresponding author upon reasonable request.