Correlation driven metallic and half-metallic phases in a band insulator

We demonstrate, using dynamical mean-field theory with the hybridization expansion continuous time quantum montecarlo impurity solver, a rich phase diagram with {\em correlation driven metallic and half-metallic phases} in a simple model of a correlated band insulator, namely, the half-filled ionic Hubbard model (IHM) with first {\em and} second neighbor hopping ($t$ and $t'$), an on-site repulsion $U$, and a staggered potential $\Delta$. Without $t'$ the IHM has a direct transition from a paramagnetic band insulator (BI) to an antiferromagnetic Mott insulator (AFI) phase as $U$ increases. For weak to intermediate correlations, $t'$ frustrates the AF order, leading to a paramagnetic metal (PM) phase, a ferrimagnetic metal (FM) phase and an anti-ferromagnetic half-metal (AFHM) phase in which electrons with one spin orientation, say up-spin, have gapless excitations while the down-spin electrons are gapped. For $t'$ less than a threshold $ t_1$, there is a direct, first-order, BI to AFI transition as $U$ increases, as for $t'=0$; for $t_4<t'<\Delta/2$, the BI to AFI transition occurs via an intervening PM phase. For $t'>\Delta/2$, there is no BI phase, and the system has a PM to AFI transition as $U$ increases. In an intermediate-range $t_2<t'<t_3$, as $U$ increases the system undergoes four transitions, in the sequence BI $\rightarrow$ PM $\rightarrow$ FM $\rightarrow$ AFHM $\rightarrow$ AFI; the FM phase is absent in the ranges of $t'$ on either side, implying three transitions. The BI-PM, FM-AFHM and AFHM-AFI transitions, and a part of the PM-FM transition are continuous, while the rest of the transitions are first order in nature. The PM, FM and the AFHM phases have, respectively, spin symmetric, partially polarized and fully polarized electron [hole] pockets around the ($\pm\pi/2$, $\pm\pi/2$) [($\pm \pi, 0$), ($0. \pm \pi$)] points in the Brillouin zone.


I. INTRODUCTION
The ionic Hubbard model (IHM) is an interesting extension of the Hubbard model and is perhaps the simplest model for a correlated band insulator which has been of a fair amount of interest to the condensed matter community in the last two decades or so. The model was first proposed for describing the neutral-ionic transition in charge-transfer organic chains [1][2][3][4] . The IHM has also been invoked for understanding polarization phenomena in perovskite materials [5][6][7][8][9] and covalent insulators such as FeSi and F eSb 2 10 . Recently this model has also been simulated in ultra-cold atom experiments 11 .
Specifically, the IHM is the Hubbard model on a bipartite lattice with a staggered potential (∆) added. The staggered potential ∆ doubles the unit-cell and, in the absence of U , gives rise to a gap in the single particle excitation spectrum such that at half-filling the system is a band insulator (BI). But when the Hubbard interaction U is introduced, it opposes the larger occupancy of particles on one sublattice and holes on the other sublattice (favoured by ∆) and prefers single occupancy on each site. Thus there is a competition between the staggered potential and the Hubbard interaction and one would expect interesting phase transitions as a function of U/∆. In particular, at half-filling, the ground state of the IHM is a band insulator (BI) for U << ∆ and a Mott insulator (MI) for U >> ∆. In a Phys. Rev. Lett. 12 co-authored by two of us, we showed using single site dynamical mean field theory (DMFT) that in the simplest IHM, with only nearest neighbour hopping t, a sufficiently large Hubbard U can close the gap in the single particle spectrum all the way to zero, resulting in an intervening correlation induced metallic phase provided the system is constrained to be paramagnetic. However, if anti-ferromagnetic (AF) order is allowed for, a direct transition from the BI phase to an anti-ferromagnetic insulator (AFI) phase occurs, preempting the transition to the correlation induced paramagnetic metallic (PM) phase [13][14][15][16][17] .
There have been fairly extensive studies of the physics of intermediate values of U (i.e., of order ∆) in the context of the IHM with only nearest neighbour hopping in one dimension [18][19][20][21][22][23][24] which suggest that a bond-ordered (BO) insulating phase (and not a PM phase) separates the BI and MI phases. In two and higher dimensions the existence and nature of the intermediate phase is still controversial despite several attempts based on various approximate methods of solving the (nearest neighbour hopping) IHM 15,17,[25][26][27] . Determinantal Quantum Monte Carlo 27 studies suggest that the intermediate phase at half-filling is a paramagnetic metal, whereas an insulating bond ordered phase is observed in cluster DMFT (+lanczos) 17 studies. The variational cluster approach to the IHM at halffilling 15 leads to the same result as obtained by the single site DMFT 13,14 approach allowing for AF order mentioned earlier, that the system undergoes a direct transition from a BI to an AF Mott insulator [except for a sliver of AF halfmetal phase seen inside the AFI phase close to the transition point 13 ]; however, when magnetic ordering is suppressed, single site DMFT 12,14 shows an intermediate metallic state. Note that there is in principle a metallic point at the critical U corresponding to the transition from the BI to the BO phase even in the 1D IHM. An interesting ferrimagnetic half-metallic phase has been predicted in the IHM doped away from half-filling 16 . A half metallic phase has been suggested in the IHM on a bilayer honeycomb lattice 28 .
It is an obvious possibility that the introduction of frustration might inhibit the formation of AF order and stabilize the correlation induced PM phase. However, an explicit calcu-lation which allows for AF order and nevertheless shows a stable PM phase in such a context has been missing. In this work, we present and discuss the results of such a calculation for the IHM on a square lattice with the inclusion of a second neighbour hopping which can frustrate the AF order. We demonstrate that for t /t sufficiently large, the system does indeed support a stable correlation induced PM phase for a range of U . Furthermore, as a bonus, other interesting phases now appear for ranges of U that depend on the value of t /t, such as a ferrimagnetic metal (FM) phase [which has nonzero values of the uniform as well as staggered magnetization], and an anti-ferromagnetic half-metal (AFHM) phase [in which, apart from the AF order, the single particle density of states (DOS) of one of the spins is gapped where as that of the other spin is gap-less at the Fermi level], before the system eventually turns into an AF insulator as the strength of U is increased. The calculations on the IHM in the presence of second neighbour hopping t which we discuss in this paper have been carried out mainly using single site DMFT 29 combined with the hybridization expansion continuous time quantum monte carlo (HYB-CTQMC) method 30 as the impurity solver. As mentioned earlier, the goal is to get the strongly correlated stable paramagnetic metallic phase to become stable by frustrating the AFM order. The DMFT calculations have been supplemented by the much simpler (unrestricted) Hartree-Fock (HF) calculations, whose results we discuss as well. Effects arising from the inclusion of a second neighbour hopping on AF order have been studied extensively earlier in the context of the Hubbard model on various lattices 31,32 , but to the best of our knowledge, its effects on the IHM have not been explored much except for one recent work on the onedimensional IHM 33 .
The main findings of our study are summarized in the phase diagram in Fig. 1, drawn for ∆ = 1.0t . For 0.0 < t < 0.05t, the system does indeed undergo a direct, first order transition from the paramagnetic BI to the Anti-Ferro (AF) insulator (AFI) as U is increased. However, for 0.05t < t < 0.5t, the system undergoes a transition from the paramagnetic BI to the paramagnetic metal (PM) as U is increased from zero. On increasing U further, magnetic order turns on, but the nature of the phases and phase transitions encountered differ for different ranges of t /t . For 0.05t < t < 0.1t, AF order together with a single particle excitation gap in one spin channel turn on abruptly at a threshold U , resulting in a first order transition from the PM phase to an AF half-metal (AFHM). On further increasing U the gapless spin channels develops a gap as well, and the system undergoes a continuous transition from AFHM to an AF (Mott) insulator (AFI). For 0.1t < t < 0.36t, the magnetic order that turns on is initially a ferrimagnetic order , with both staggered and uniform magnetization turning on simultaneously; hence the PM phase undergoes a transition into a ferrimagnetic metal (FM) phase. As depicted in 1, we find that this transition is continuous for 0.1t < t < 0.18t, and first order thereafter. As U increases further, the uniform magnetization decreases, though the staggered magnetization keeps increasing. Coincident with the decay of the uniform magnetization to zero, a single particle excitation gap opens up in one spin channel, resulting in a continuous transition into an AF half-metal (AFHM). On further increasing U the gapless spin channel develops a gap as well, and the system undergoes a continuous transition from the AFHM into the AFI phase as before. For larger values of t , 0.36t < t < 0.46t, the ferrimagnetic metallic phase disappears and the paramagnetic metal phase directly undergoes a first order transition into AFHM phase, followed by a continuous transition to the AFI phase. For even larger values of t , there is no AFHM phase, and the system undergoes only two transitions; first from BI to the paramagnetic metal and then from the metal to the AFI. Even the BI phase ceases to exist for t > 0.5t . Amazingly, as we show in the following section, all the same phases show up in a simple unrestricted Hartree-Fock (HF) study of the t − t IHM, although the parameter values where the transitions occur , as well as some aspects of the topology of the phase diagram differ between the HF and the DMFT+CTQMC calculations. Needless to say, we expect the latter to be more accurate.
The rest of this paper is organized as follows. In Section II we present the model and its physics within a simple unrestricted HF approximation. In Section III we present and discuss the phase diagram and different observables calculated using DMFT+CTQMC. These are the central new results of our paper. In Section IV, we discuss (the mainly quantitative) differences in the phase diagrams obtained by the two methods. We end this paper with a concluding discussion in section V. Some details of the calculation methods used are presented in Appendix A for the HF approximation, and in Appendix B for DMFT+CTQMC.

II. IHM HAMILTONIAN, SIMPLE LIMITS AND HF MEAN FIELD THEORY
In this section we present the Hamiltonian for the IHM, its physics in simple limits, and possible phases within the simplest unrestricted HF mean field theory. The Hamiltonian for the IHM can be written as Hereĉ † iσ creates electron with spin σ on site i of a 2-d square lattice which we regard as made of two square sublattices A and B whose sites are nearest neighbour of each other; t is the nearest neighbour hopping, t is the second neighbour hopping 36 , U denotes the on-site coulomb repulsion and +∆(-∆) denotes the "ionic" potential for the A(B) sublattice. The chemical potential (µ) is chosen so that the average occupancy ( n A + n B ) /2 ≡ (n A + n B )/2, is 1 corresponding to the "half-filling" constraint.
A. Simple limits of the IHM Some simple limits of the above model are relatively easy to understand. In the atomic limit (t = t = 0), when U/2 < ∆ the ground state of the system has two electrons on every B site and no electrons on the A sites, resulting in the maximal charge density wave order, and with a single particle excitation gap of ∆ − U/2. This is the atomic limit of the (correlated) BI phase. In the other case, when U/2 > ∆, the ground state of the system is the (atomic limit of the ) Mott Insulator (MI) phase, with 1 electron at each site, and a gap for single particle excitations equal to U/2 − ∆. Clearly at U = 2∆ the system has gapless single particle (electron or hole) excitations, and can of thought of as "metallic".
For U = 0 the model can be diagonalized exactly, in terms of two bands of electron creation and destruction operators with quasi-momentum or wave-vector labels k ≡ (k x , k y ). The (spin independent) band dispersion relations in the full Brillouin Zone (BZ) are given by When t = 0 the two bands are separated in energy by a band gap (E gap ) of 2∆ along the square contour corresponding to k y = ±π ± k x . So at half filling the lower band is completely filled and the upper band is empty, resulting in a band insulator (BI). But when t = 0, the eigenstates for K = (0, ±π), (±π, 0) and (±π/2, ±π/2) are no longer degenerate. The bottom of the upper, conduction band is + K = ∆ at K = (±π/2, ±π/2) whereas the top of the lower, valence band is − K = 4t − ∆ at K = (0, ±π), (±π, 0). Hence, for t = 0 and ∆ > 2t the two bands have an indirect band gap E gap = −4t + 2∆. As long as E gap > 0 the system (at half-filling) continues to be a Band Insulator (BI). As one increases t from 0, the band gap monotonically decreases from 2∆ and become zero at t = ∆/2, whence we get a BI to (band) metal transition.

B. HF Mean Field Theory of the IHM
For U = 0 the IHM is no longer exactly solvable. A simple approximate solution can be obtained using the (unrestricted) Hartree-Fock mean field theory (HF-MFT), where, as described in detail in Appendix A, at the simplest level one approximates the Hamiltonian as an effective quadratic Hamiltonian allowing for mean field order parameters corresponding only to the sublattice and spin resolved occupancies n ασ ≡ n ασ , or, equivalently, to the staggered occupancy (δn ≡ (n B − n A )/2), the staggered magnetization (m s ≡ (m A − m B )/2) and the uniform magnetization (m f ≡ (m A +m B )/2) where m α is the magnetization on any site belonging to sublattice α (see Appendix A for details). The resulting band dispersion relations, now spin dependent, are given bỹ clearly one can interpret these results in terms of an effective spin dependent staggered potential∆ σ = ∆−U (δn+σm s )/2 and an effective spin dependent uniform potential −σU m f /2. The effective band gap (Ẽ gap,σ ) for electrons with spin σ, determined by the difference between the bottom of the conduction band (˜ + kσ ) and the top of the valence band (˜ − kσ ), similarly to the U=0 case, isẼ gap,σ = −4t + 2∆ σ , which is interaction and spin dependent. The order parameters are then determined self consistently by populating these effective non-interacting bands as per Fermi-Dirac statistics.
Phase diagram of the t − t IHM in the t -U plane at half filling using Hartree-Fock Mean Field Theory at zero temperature for ∆ = 1.0t. For t < 0.06t, the system undergoes a direct transition from BI to AFI phase with increase in U . For 0.08t < t < 0.24t, the system goes from a BI phase to AFI phase via an intermediate AFHM phase upon increase of U , whereas for 0.24t < t < 0.26t, the system shows three transitions, from BI to FM, from FM to AFHM and eventually from AFHM to AFI as U is increased. For 0.26t < t < 0.5t, the PM phase intervenes between the BI and the other phases, whereas for t > 0.5t there is no BI phase, and the transition is from the PM into other phases. Full lines in the phase diagram indicate first order transition lines, whereas dashed lines indicate continuous transitions.
The details are provided in Appendix A. The phases are characterised both by the order parameters that are nonvanishing, and the nature of the band dispersions. Amazingly, the resulting phase diagram obtained within the HF approximation, depicted in Fig. 2, shows all the same phases as obtained using the much more sophisticated and accurate DMFT+CTQMC calculations, albeit for quantitatively different ranges of the parameter values. Since qualitative features of this phase diagram can be easily understood given the simplicity of the HF approximation, we discuss key aspects of these results next (while relegating the details to Appendix A), before presenting a discussion of the key aspects of our DMFT+CTQMC calculations.

C. HF Phase Diagram of the IHM
First consider the HF phase diagram (Fig. 2) in the paramagnetic regime. In this case the system does not have any magnetic order and hence both the staggered magnetisation, m s , and the uniform magnetisation, m f , are zero. This implies that∆ ↑ =∆ ↓ . This spin symmetry is also reflected in the single particle excitation band gapẼ gap,σ . The resulting phase is adiabatically connected to the noninteracting BI phase, hence we use the same label for it. As U increases, the density difference between the two sublattices, δn, decreases because U does not prefer holes or double occupancy. As an effect the gap in the single particle excitation spectrum,Ẽ gap , decreases. One might expect that as one increases U further, E gap might go to zero at a certain U , whence one would get a band insulator to metal transition. On the contrary, when t = 0, one finds that the system never goes to a paramagnetic metallic state within the HF theory, i.e.,Ẽ gap is nonzero for all U . However, more accurate calculations using single site DMFT do lead to a metallic state, as pointed out first by Garg et al 12 . Later work [13][14][15][16][17] showed that this metallic solution is not stable against AF ordering 13 . As is clear from Fig. 2, this result continues to hold even when t = 0, as long as t < 0.26t for ∆ = 1.0t. One of the main new results of the work presented in this article is that for sufficiently large t (eg., t > 0.26t for ∆ = 1.0t as in Fig. 2),Ẽ gap does go to zero within the paramagnetic sector, in both the HF and DMFT methods, leading to a robust paramagnetic BI to Paramagnetic metal (PM) transition, and a stable correlation-induced PM phase. Now we discuss results from the unrestricted HF theory where the spin symmetry breaking is allowed. Panels (a) and (b) of Fig. 3 show the magnetic order parameters m s and m f at zero temperature within the HF approximation as functions of U for different values of t . As U increases, the staggered magnetisation m s turns on at a threshold value U c , where U c is a function of t . At smaller t values, the magnetic transition is a first order AF transition, with m f = 0 for all values of U . As t increases, the discontinuity in m s across the transition decreases. We find that the transition is continuous for 0.215t < t < 0.26t, and then again becomes first order as t increases beyond 0.26t, indicating the presence of multicritical points at (U = 2.82t, t = 0.215t) and (U = 2.16t, t = 0.26t). The net magnetization m f seems to turn on for the first time for t ∼ 0.24t when the magnetic transition is continuous. For larger values of t , m f also shows a first order jump indicating that the magnetic transition has become first order in nature. The first order nature of the transition in m s as a function of U in Fig. 3 is consistent with earlier results for t = 0 on the square lattice 13,25 , and differ from the results in 14 . Similar changes in the order of magnetic transition with increase in t have been observed in the t − t Hubbard Model 32 . Note however, that the FM-AFHM and the AFHM-AFI transitions are continuous, and can be identified by kinks (changes of slope) in the m s versus U/t curves. Clearly, these lead to other multi-critical points in the phase diagram.
We note that U c has an interesting non-monotonic dependence on t . U c initially decreases with increasing t , reaching a minimum around t = 0.26t, and then increases again. As we know from earlier studies on the IHM with t = 0 and bilayer systems, the value of U c , at which the BI-AFI transition occurs, decreases with decreasing ∆, the non interacting band gap 13,34 . Since the effective band gap,Ẽ gap , decreases as t increases, it is natural to expect that U c should decrease with increasing t . Another aspect of the next neighbour hopping is that t introduces frustration against AFM ordering 31,32 , which effectively should increase the U c required for ordering. We believe that these two competing effects of t are responsible  for the observed non-monotonic trend of U c with increasing t . This physics is seen in the phase diagrams obtained both within the HF and the DMFT methods. Additional and interesting aspects of the magnetically ordered phases emerge from an examination of the self consistent HF band dispersions, and their gaps and overlaps. Fig. 4 shows the minimum of the upper (conduction) band, (˜ + min,σ ), and the maximum of the lower (valence) band, (˜ − max,σ ), measured relative to the chemical potential µ and for each spin channel, as a function of U for t = 0.3t. When the AF order parameter m s is nonzero, the effective staggered potential is spin dependent (∆ ↑ =∆ ↓ ). In fact, for a positive m s , corresponding to the A site being preferentially occupied by up spins, the effective staggered potential for ↑ spins, (∆ ↑ ), is smaller than that for ↓ spins, (∆ ↓ ). So the effective band gap for ↑ spins, (Ẽ gap,↑ ) is less than that for ↓ spins,Ẽ gap,↓ . Because of this spin dependence of the effective band gap, if one dopes the system one  gets a region of half metallic (HM) phase over a range of U and doping 16 .
As illustrated in Fig. 4, for sufficiently large values of (t > 0.26t for ∆ = 1.0t), and for values of U beyond a threshold, when the system is still spin symmetric, the upper and lower bands overlap for both the spin channels with the Fermi level being inside the bands, indicating a paramagnetic metal (PM) phase, with electron pockets near the conduction band minimum, and hole pockets near the valence band maximum. As U increases further, the spin symmetry is broken. There is a small range in U past the magnetic transition for which the down spin hole pocket disappears as its valence band is below the Fermi level while the down spin electron pocket and the up spin electron and hole pockets are still present. This corresponds to the ferrimagnetic metalic phase where m f = 0. On further increasing U , the bottom of the down spin conduction band goes above the Fermi level, and there is a broad range of U where the ↑ spin gap is zero while the ↓ spin gap is finite. This parameter regime corresponds to an AF ordered half-metal (AFHM) phase. As one increases U , δn decreases but m s increases and at a certain U , δn = m s . If one increases U further,Ẽ gap,↑ starts to increase from its negative value, and eventually it again becomes positive, resulting in a transition from the AFHM phase to an anti-ferromagnetic insulator (AFI) even within HF theory. Fig. 5 shows the band dispersions along specific paths (mentioned in the caption) in the Brillouin Zone for five representative values of U chosen to correspond to all the five phases that arise for t = 0.3t as U is increased, which again clearly brings out the above mentioned features of these phases. Effective band dispersions, calculated using HF-MFT, plotted along the path in the Brillouin zone corresponding to (0, 0) → (π/2, π/2) → (π, π) → (0, π) → (0, 0) (vertical doted lines indicate the break points) for different values of U at ∆ = 1.0t and t = 0.3t. Zoomed in plots of the dispersions near zero frequency are presented in the lower panels. (a) For U = 0.4t, the upper band has a minimum at (π/2, π/2) and the lower band has a maximum at (0, π). The chemical potential lies within the gap. Hence the system is in the (paramagnetic) BI phase. (b) For U = 2.2t, both the upper and the lower bands cross the chemical potential. The upper band has an electron pocket at (π/2, π/2) while the lower band has a hole pocket at (0, π). This is the paramagnetic metallic phase. (c) For U = 2.6t, the band structure is different for ↑ and ↓ spins. The upper bands for both the spins cross the chemical potential, whereas only the lower ↓ spin band crosses the chemical potential. Hence we get a ferrimagnetic metallic state. (d) For U = 3.4t, the ↑ spin bands cross the chemical potential whereas the ↓ bands are gapped out, corresponding to the AFHM state. (e) For U = 5.0t, bands for both the spins are gapped out, and the system is in an AFI phase.

D. Single particle DOS within HF theory
All of the features of the phases and phase transitions discussed above can be seen more directly in the single particle density of states (DOS) for real frequencies, which can be calculated exactly within the HF approximation. Specifically, Fig. 6 shows the results of our calculations for the spin dependent single particle DOS averaged over the two sublattices, is the single particle DOS for sublattice α and spin σ. Here G ασ (ω) is the Green's function on sublattice α, ω + = ω + iη with η → 0. The top panels in Fig. 6 show the DOS for t = 0.3t and ∆ = 1.0t for the same values of U as in Fig. 5. For U = 2.2t, in the PM phase, the system has spin symmetry (with m s = m f = 0), the DOS for the two spin channels is the same and finite at the Fermi level (corresponding to ω = 0). For U = 2.6t in the FM phase, where both the magnetic order parameters are non-zero, the DOS is spin asymmetric, but finite at the Fermi level for both the spin channels. For U = 3.4t, A ↓ (ω) has a gap around the Fermi level, while A ↑ (ω = 0) is finite. Because of the half-filling constraint, this necessarily implies that m f = 0; so this is the AFHM phase. For U = 5.0t gaps in the DOS are present in both the spin channels, corresponding to the AF (Mott) Insulating phase. The lower panels of Fig. 6 show the DOS at the Fermi level, A σ (ω = 0), as a function of U for four illustrative values of t , bringing out the sequence of phases and phase transitions that occur as U increases consistent with the complete phase diagram in Fig. 2. The existence of the BI (for smaller U values) as well as AFI (for larger U values) phases where the DOS in both spin channels is zero, and the AF half-metal where A ↓ (ω = 0) = 0 but A ↑ (ω = 0) = 0, for different ranges of U is clearly evident. For t = 0.3t the PM phase, where A ↑ (ω = 0) = A ↓ (ω = 0) = 0, and the FM phase, where A ↑ (ω = 0) = A ↓ (ω = 0) = 0, appear as well. For t = 0.52t the PM phase is present even at U = 0, and there is no BI phase. Although HF-MFT is very useful for the insights it provides into the nature of the various phases as discussed above, it is a mean field theory that neglects quantum fluctuations entirely. Hence one expects that it is likely to overestimate the regions of stability of the phases with broken symmetry. Needless to say, the inclusion of the quantum fluctuations, especially in the strong correlation regimes, pose considerable challenges, and in this paper we have restricted ourselves to doing this within the DMFT approximation, which correctly includes at least all the local quantum fluctuations, although the challenge of solving the resulting self-consistent impurity problem necessitates the further approximation of CTQMC. We discuss these DMFT calculations next.

III. DMFT PHASE DIAGRAM
As mentioned earlier, the central result of our paper is the low temperature phase diagram of the half-filled IHM obtained within DMFT+CTQMC for the 2d square lattice in the presence of a second neighbour hopping term. Fig. [1] shows the phase diagram for ∆ = 1.0t. Recalling the salient features of this phase diagram, we note that for t < 0.05t, increasing U leads to a first order transition from the paramagnetic BI to an AFI, similar to the t = 0 case. However, as anticipated in our motivation for studying the t − t IHM , for t > 0.05t the first transition induced by increasing U is a continuous transition from the BI phase into the correlation induced paramagnetic metal (PM) phase. As U is increased further, two new phases not present in the t = 0 case (except at one value of U for each ∆ or with doping 13,16 ), namely the AFHM and the FM phases, intervene before the eventual transition into the AFI phase for large enough U . The AFHM phase is seen for a broad range of values of the second neighbour hopping, 0.05t < t < 0.46t, and abuts the AFI phase, whereas the FM phase is seen only for 0.12t < t < 0.36t , and abuts the PM phase. A direct transition from the PM to the AFHM phase occurs for 0.05t < t < 0.12t and 0.36t < t < 0.46t, which is first order in nature. So is the PM to FM transition, except for the range 0.12t < t < 0.18t, where it is continuous. Both the FM-AFHM and the AFHM-AFI transitions are continuous in nature. For 0.46t < t < 0.5t, the system undergoes only two transitions with increasing U , a continuous transition from the BI to the PM phase, followed by a first order transition into the AFI. For t > 0.5t there is no BI phase. In this range of t , the system goes through only one transition as U increases, namely, a first order tran- sition from the PM to the AFI phase. We note that indeed, as mentioned earlier, the values of various transition lines within the DMFT+CTQMC method are very different from those obtained from HF theory. Due to quantum fluctuations captured within the DMFT+CTQMC method, the magnetic order sets in at larger values of U as compared to the HF theory as ex- pected. Furthermore, the width of the paramagnetic metal phase and the ferrimagnetic metal in U − t plane is enhanced while the AF half-metal phase is diminished in width in the phase diagram from the DMFT+CTQMC method compared to that from the HF method.
In the rest of this section, we discuss the details of the DMFT results based on which this phase diagram has been constructed. While the identification of magnetic order is relatively easy in the DMFT method, it is much harder to distinguish insulating from metallic or half metallic behaviour, as the CTQMC produces imaginary time data, and precise spectral information in real frequency is therefore hard to generate. Although the frequency dependent spectral functions can not be obtained in a simple reliable manner from DMFT+CTQMC method, as we show later, one can look for features of metallic or insulating behaviour in momentum distribution function n kσ , which is obtained by summing over Green's functions at imaginary frequencies. The insights gained from the HF method are also helpful in this regard. Furthermore, note that the use of the CTQMC methods constrains the calculations to be done at a non-zero temperature. The calculations presented in our paper are at βt=50.0, which we believe is at a low enough temperature that we expect our results to be close to the zero temperature results.

A. Staggered density and Magnetization
The staggered magnetization m s and the uniform magnetization m f calculated within the DMFT+CTQMC are shown as functions of U/t in Fig. 7 for ∆ = 1.0t and inverse temperature β t = 50.0 for various values of t . For t = 0 the AF order turns on with a strong first order transition at U c ∼ 4.2t. When t is turned on, first the threshold value U c decreases, reaches a minimum around t = 0.18t ≡ t , and then increases again with further increase in t . Furthermore, the discontinuity in m s at the first order transition decreases with increasing t , and the system undergoes a continuous transition for the range 0.12t < t < 0.18t, where both m s and m f turn on smoothly at U c . For t ≥ 0.2t, again the staggered magnetization and the uniform magnetization turn on with a jump. For t = 0, at half-filling the system has particle-hole symmetry leading to m A = −m B , whence m f = (m A + m B )/2 = 0, i.e., the magnetic order is purely AF. When t is turned on m f becomes non-zero due to the absence of particle-hole symmetry in the system, although it remains rather small for a range of intermediate values of t , as shown in panel (b) of Fig. 7 for ∆ = 1.0t. For t < 0.12t, m f is zero for all values of U . For 0.12t < t < 0.36t, there is a small range U 1 < U < U 2 in which m f is non zero and the system has ferrimagnetic order, whereas for U > U 2 the system has pure AF order, and for U < U 1 the system is paramagnetic. For t > 0.36t, the uniform magnetization is again vanishingly small for all values of U . In this regime the system undergoes a transition from a para-magnetic phase to the AF phase and there is no ferrimagnetic phase. Fig. 8 shows the spin resolved as well as the mean staggered densities given by δn σ = (n Bσ − n Aσ )/2 and δn = (δn ↑ + δn ↓ )/2 respectively, for inverse temperature βt = 50 and ∆ = 1.0t. δn is a decreasing function of U for any value of t . An abrupt, first order drop in δn σ is seen across those values of U c where the magnetic transition is strongly first order, i.e., for small and large values of t . For intermediate values of t where the magnetic transition (PM to FM transition) is continuous or weakly first order, δn seems to change smoothly across the transition, but shows a slope change across the continuous AFHM to AFI transition.

B. Momentum Distribution function
At zero temperature, whether a system is metallic or insulating is ideally characterized by the single particle density of states (DOS) or spectral function, calculated from the imaginary part of the single particle real frequency Greens function (see equation 4 and the related discussion). For a metal, the DOS at the Fermi energy is non zero, while for an insulator the DOS at the Fermi energy is zero. But within the CTQMC method, the Green's function is calculated only at imaginary times and at a finite temperature, and hence its Fourier transform is available only for imaginary Matsubara frequencies, iω n ≡ i(2n+1)πk B T where n is an integer. In order to obtain the real frequency DOS at low frequencies close to the Fermi energy (ω = 0), an analytic continuation of the imaginary fre- quency Green's function is required. Perhaps the best known way of doing the analytic continuation is to use the maximum entropy method (MEM) but it has an in-built bias, involving an initial guess for the DOS, which makes it hard to get reliable real (low) frequency data to sharply distinguish between a metal and an insulator. As an alternative procedure that avoids analytic continuation altogether, in this paper we use the momentum distribution functions n α,k,σ = 0 −∞ dωA ασ (k, ω) where A ασ (k, ω) is the momentum resolved spectral function for sublattice α and spin σ, to distinguish between the metallic and insulating nature of the system. Within the CTQMC n α,k,σ is calculated by doing a summation over the Green's function at Matsubara frequencies as In a paramagnetic insulating phase at half filling, adiabatically connected with the noninteracting BI phase, where the lower band is fully occupied while the upper band is vacant, with the chemical potential lying in the gapped region between the two bands, the spin-resolved momentum distribution function n k,σ ≡ 1 2 α n α,k,σ will have a constant value of 1/2 −π 0 π k x π 0 −π k y n k, ↑ −π 0 π π 0 −π n k, ↓ (0, 0) ( π 2 , π 2 ) (π, π) (0, π) (0, 0) throughout the Brillouin zone (BZ). In a paramagnetic metallic state, adiabatically connected to the state in which the two (non-interacting) bands cross the Fermi level, we will have (interaction-renormalized) electron pockets in the conduction band leading to a larger (close to but less than 1.0) n k,σ and hole pockets in the valence band, leading to a smaller (close to but greater than 0.0) n k,σ , in such a way as to keep the total particle number fixed at half filling ((n A + n B )/2 = 1). One can similarly use the spin dependence of the momentum distribution functions combined with the information on the magnetic order parameters to distinguish the FM and AFHM phases, as discussed below. Fig.s 9, 10 and 11 show false color plots of n k,σ over the full BZ for t =0.2t for three values of U at ∆ = 1.0t. For U = 3.4t (Fig. 9), the momentum distribution function is spin symmetric and has clearly visible renormalized electron and hole pockets as mentioned above, reflected respectively in the peaks at k = (±π/2, ±π/2) and the dips at k = (0, ±π), (±π, 0) in n k,σ , corresponding to the PM phase. Note also the particle-hole asymmetry because of the presence of t . As U is increased and the magnetic order sets in, there develops an asymmetry in (n k,σ − 0.5) for the ↑ and ↓ spin channels, as shown in Fig. 10 for U = 3.8t, with the ↑ spin channel having bigger electron and hole pockets, and also a larger net occupancy, compared to the ↓ spin channel. This is clearly a FM phase. On further increasing U , as shown in Fig. 11 for U = 4.6t, only n k,↑ has Fermi pockets while n k,↓ has the uniform value of 0.5 ev- erywhere in the BZ. Thus the down spin channel is insulating with a finite band gap while the up spin channel is a metal with both conduction and valence bands crossing the Fermi level, while the net occupancy in each channel is the same. This clearly corresponds to the AFHM phase of the system. The complete picture about the multiple metal-insulator transitions and crossovers in the IHM, shown in the phase diagram of Fig. 1 requires the momentum distribution functions calculated within DMFT+CTQMC as illustrated above. A quick way of obtaining this picture is to study the variation in n k,σ at just the centers of the electron and hole pockets as functions of the parameters of the IHM. Fig. 12 shows (n kσ − 1 2 ) at k = (π/2, π/2) (red and green curves), and (π, 0) (magenta and blue curves) as a function of U for several different representative values of t , chosen to show the different sequences of phases and phase transitions that are possible. For example, consider the panel for t = 0.2t, for which all five phases that we have discussed exist for different ranges of U , separated by four phase transitions (see Fig. 1). In the insulating phases (n k,σ − 1/2) is zero and the onset of metallicity is indicated by (n k,σ − 1/2) deviating from zero. One can see that the first transition is continuous and from the BI to the PM phase, with spin symmetric occupancies, across the blue curve in the phase diagram of Fig. 1. This is followed by the second transition, into the FM phase where all the (n k,σ − 1/2) are still nonzero but become different for ↑ and ↓ spins. The third transition from the FM to the AF half-metal, where (n k,σ − 1/2) is different from zero for only one spin channel, whereas the other spin channel is insulating. The final transition at large U , is a first order transition to the AFI phase, with n k,σ = 1/2 again for both the spin channels. One can similarly see that the occupancy deviations from 0.5 shown in the other panels are consistent with Fig. 1. Note that a sharp distinction between the metal and the insulator can be made strictly only at T = 0 based on whether the spectral functions have a gap around the Fermi level or not, which is reflected in the distribution function n k,σ as well. However, CTQMC calculations are possible only at finite temperatures; therefore, even at the lowest possible temperatures accessible numerically, one can clearly distinguish between metallic and insulating phases only if the gap in the insulating phase is much larger than T . Close to the metal-insulator transition point, where the gap becomes very small, numerically it becomes very hard to locate the transition, as it is really a crossover. But within our numerical accuracy, we can still clearly see the presence of a correlation induced paramagnetic metallic phase, ferrimagnetic metallic phase and AF half-metal in the IHM, which is one of the main achievements of this work.

IV. CONCLUDING COMMENTS
In conclusion, we have presented a detailed phase diagram of the half-filled IHM in the presence of a second neighbour hopping t within DMFT+CTQMC. We have demonstrated that in this simple model of a band insulators (BI), due to the (nk,σ − 1/2) at k = (π/2, π/2) (red and green curves) and k = (0, π) (magenta and blue curves) vs U for eight different values of t = 0.08t, 0.1t, 0.2t, 0.34t (top set) and t = 0.0t, 0.05t, 0.25t, 0.6t (bottom set). In the metallic phases (nk,σ − 1/2) is greater than zero at k = (π/2, π/2), suggesting electron pockets, and is less than zero at k = (0, π), suggesting hole pockets. For t = 0.2 one can clearly see the transition from the BI to the PM phase, followed by a transition into the FM phase, where nk,σ becomes different for ↑ and ↓ spins. This is followed by a transition into an AF half-metal where (nk,σ − 1/2) is different from zero for only one spin channel. Finally at large U , nk,σ = 1/2 for both the spin channels, corresponding to the AFI phase. Other panels similarly show the phases and transitions for other values of t , and also how the width in U of the different phases changes with t .
frustration of the anti-ferromagnetic (AF) order induced by t , one can stabilize a correlation induced paramegnetic metallic (PM) phase intervening between two insulating phases, namely the BI and the AF-Mott Insulator (AFI) for relatively large values of t /t, whereas the PM phase occurs as a precursor to an AF half-metallic phase, or a ferimagnetic metallic phase for smaller values of t /t. We believe that it is interest-ing how e-e interactions can dynamically close the gap in the band insulator and result in the formation of stable metallic phases which can be paramagnetic or ferrimagnetic, and even more interestingly, anti-ferromagnetic half-metallic. The question of a whether there can be such a stable correlation induced metallic phase intervening between the BI and the AF Mott Insulator phases in the half-filled IHM has been floating around in the community for almost a decade and had not been answered to the full extent earlier. In this work, we have used the momentum distribution function calculated within the DMFT+CTQMC for the analysis of metal-insulator transition, which, we think, is more reliable then using analytically continued single particle DOS. Our work presents a clean and complete answer to the above question by analysing the IHM at half-filling in the presence of a second neighbour hopping, and furthermore, has uncovered a rich set of additional phases and phase transitions not known hitherto. We showed that, while for t below a threshold t 1 an increase in U induces a BI to AFI transition as in the t = 0 case, there is a range t 1 < t < t 2 where the correlation induced PM phase and then a new AFHM phase intervene between the BI and the AFI phases as U increases. In the next range of values of t /t, t 2 < t < t 3 , the system undergoes four transitions involving five phases as the strength of U is increased. First the paramagnetic BI undergoes a continuous transition into a PM, which at U c becomes a ferrimagnetic metal (FM) with the onset of uniform as well as anti-ferromagnetic order. As U is increased further, a gap opens up in the single particle excitation spectrum of one of the spin channels while the other spin channel remains gapless, resulting in the formation of a half-metal which has pure AF order (the AFHM phase). On further increasing U , both the spin channels attain gaps in their single particle excitation spectra and the system becomes an AFI. For larger values of t , t 3 < t < t 4 , the FM phase does not appear, and the transition from the BI to the paramagnetic metal is followed by a transition into the AFHM phase which eventually transforms to an AFI as U increases. For even larger values of t , t 4 < t < t 5 , the system undergoes only two transitions as U increases with a broad PM phase intervening between the BI and the AFI phases. For t > t 5 there is no BI phase, and there is only one transition, between the PM and the AFI phases. The sizes and the end points of the different ranges of t mentioned above vary with ∆, and we hope to study their systematics in future work.
Furthermore, in our work we also showed that a similar phase diagram can be obtained in a simple unrestricted HF theory, although within the simple HF theory the region of stability of the PM phase is suppressed compared to the DMFT results due to two reasons, both arising presumably from the missing quantum fluctuation in the HF calculations. First, the transition values U c for the onset of magnetic order are lower in the HF theory in comparison to DMFT; second, the AFM half-metal phase is too broad in the HF theory. Within the HF theory, unlike in the DMFT results, once t is large enough to stabilize the correlation induced PM phase, the AFHM phase seems to always intervene between the PM and AFI phases, and a direct transition from the PM to the AFI is never seen.
It will definitely be very interesting to investigate the ex-istence of these correlation driven metallic and half-metallic phases as well as the nature of various phase transitions involved, experimentally in real materials as well as in cold atom experiments where IHM has already been simulated on a 2d honeycomb lattice 11 , and theoretically, using the more sophisticated and computationally demanding cluster DMFT and other methods 37 . We hope very much that our work presented here stimulates further research in these directions.
Appendix A: Hartee-Fock theory In the HF approximation we break up the interaction term within a mean-field approximation as Un i↑ni↓ U [ n i↑ n i↓ +n i↑ n i↓ − n i↑ n i↓ ], (A1) Furthermore, in this work we make the simplifying ansatz that n iσ is the same, n Aσ or n Bσ , for all A or B sites respectively. Transforming to creation (and destruction) operators labelled with wave vectors k according tô where α = A, B is the sub-lattice index, we can write the (operator part of the) HF hamiltonian aŝ where A kσ ≡ ∆ − µ + k + U n Aσ and B kσ ≡ −∆ − µ + k + U n Bσ , with k ≡ −4t cos k x cos k y and k ≡ −2t(cos k x + cos k y ). Suppose that the one particle eigenstates of the Hamiltonian are |γ ± kσ >= (P ± Akσĉ † Akσ + P ± Bkσĉ † Bkσ )|0 > with eigenvalues ± kσ , which are clearly the effective one-particle band energies. Then the eigenvalue equation is The band energies at wave-vector k and for spin σ are therefore given bỹ where L ± kσ ≡ − k /[A kσ −˜ ± kσ ]. The particle occupancy at site α for the state (|γ ± kσ >) is given by n ± αkσ =< γ ± kσ |ĉ † αkσĉ αkσ |γ ± kσ >= P ±2 αkσ . In the ground state, only the single particle states with negative energies are occupied, Hence we get the following self consistency relations for the four occupancies at site α and for spin σ: For each set of parameter values (U/t, ∆/t, and t /t), we numerically solve the above four equations and the half filling constraint equation to determine the four (spin and site resolved) occupancies and the chemical potential µ. The HF order parameter results we presented and discussed in section II are obtained from these occupancies using the relations m A = n A↑ − n A↓ , m B = n B↑ − n B↓ n A = n A↑ + n A↓ , n B = n B↑ + n B↓ m s = (m A − m B )/2, m f = (m A + m B )/2 δn = (n A − n B )/2 (A7) We note here that to ensure the correctness of our HF calculations we have benchmarked them against 1) the U c Vs t plot for t − t Hubbard model in Vollhardt et al 32 . Furthermore, the analytical form of band gap in the non interacting t − t IHM on a square lattice is 2∆ − 4t . The values of the gap calculated from the numerically calculated DOS from our code match well with this analytical form of the gap.