# Metal-insulator transitions, superconductivity, and magnetism in the two-band Hubbard model

###### Abstract

We explore the ground-state properties of the two-band Hubbard model with degenerate electronic bands, parametrized by nearest-neighbor hopping , intra- and inter-orbital on-site Coulomb repulsions and , and Hund coupling , focusing on the case with . Using Jastrow-Slater wave functions, we consider both states with and without magnetic/orbital order. Electron pairing can also be included in the wave function, in order to detect the occurrence of superconductivity for generic electron densities . When no magnetic/orbital order is considered, the Mott transition is continuous for (quarter filling); instead, at (half filling), it is first order for small values of , while it turns to be continuous when the ratio is increased. A significant triplet pairing is present in a broad region around . By contrast, singlet superconductivity (with -wave symmetry) is detected only for small values of the Hund coupling and very close to half filling. When including magnetic and orbital order, the Mott insulator acquires antiferromagnetic order for ; instead, for the insulator has ferromagnetic and antiferro-orbital orders. In the latter case, a metallic phase is present for small values of and the metal-insulator transition becomes first order. In the region with , we observe that ferromagnetism (with no orbital order) is particularly robust for large values of the Coulomb repulsion and that triplet superconductivity is strongly suppressed by the presence of antiferromagnetism. The case with , which has an enlarged SU(4) symmetry due to the interplay between spin and orbital degrees of freedom, is also analyzed.

## I Introduction

The single-band Hubbard model represents the simplest example to describe strongly-correlated systems, where the interplay between kinetic energy and Coulomb repulsion may give rise to a rich phase diagram, which includes insulating and conducting states, with possible superconductivity and/or spin/charge disproportionations leblanc2015 (); zheng2017 (). This model can be used to capture the low-energy properties of materials where spin and charge fluctuations involve predominantly one orbital (i.e., fluctuations among different orbitals are substantially quenched), as for example cuprate superconductors. In this regard, it is widely accepted that the single-band Hubbard model (or its strong-coupling limit, i.e., the so-called model) may grab the essential features of high-temperature superconductivity dagotto1994 (); imada1998 (). Still, there are many cases in which orbital fluctuations are relevant and give rise to important physical phenomena that cannot be captured within a single-band model. For example, the hybridization among different orbitals and the presence of the Hund coupling may produce appreciable effects at low temperatures, thus affecting both the normal and the superconducting phases. One of the most prominent examples is given by the iron-based superconductors, where all the orbitals of Iron atoms may play an important role in the conducting properties and the inclusion of multiband effects is necessary to correctly describe the relevant aspects of the electronic properties (e.g., the topology of the Fermi surface) lee2008 (); daghofer2010 (); yu2011 (); hu2016 (); fernandes2017 ().

Within multiband models, one key point that has been addressed in the past is to understand how the Mott metal-insulator transition (MIT) at integer fillings is affected by orbital degeneracy, inter-orbital Coulomb interaction, and Hund coupling. In this context, many works have been performed in the “symmetric sector”, namely disregarding any possible magnetic or orbital long-range order, in order to capture the correlation effects that are not spoiled by weak-coupling effects. This approach is justified by the choice of describing the physical picture that can be realized when magnetic and orbital order is suppressed by the presence of competing interactions, i.e., frustration (without including it explicitly in the model). For the single-band Hubbard model, this way of proceeding has been widely used within the Gutzwiller approximation brinkman1970 (); vollhardt1984 (), dynamical mean-field theory (DMFT) georges1996 (), slave-boson approaches kotliar1986 (), and variational Monte Carlo methods capello2006 (). For the -band Hubbard model, in the absence of the Hund coupling , it has been observed that the value of the Coulomb interaction , for which the MIT occurs at commensurate filling , reaches its maximum at half filling, i.e., for . This result has been obtained by using the Gutzwiller approximation lu1994 (), DMFT rozenberg1997 (); ono2003 (), and quantum Monte Carlo techniques gunnarsson1999 (). The presence of a finite term reduces the value of at half filling ono2003 (); han1998 (). Then, recent studies demedici2011a (); demedici2011b () highlighted the opposite trend for all the other (integer) fillings, where the presence of a finite increases (in this case, the existence of a correlated metal with tiny quasiparticle weight has been also emphasized demedici2011a (); demedici2011b (); fanfarillo2015 (); nomura2015 ()). One important issue that has been addressed in multiband Hubbard models is the nature of the MIT. Indeed, while in the single-band model different numerical methods brinkman1970 (); georges1996 (); capello2005 () established that the MIT is continuous at zero temperature, former studies of multiband models, based on the Gutzwiller approximation, suggested that the transition, at half filling, becomes first order whenever , while it remains continuous only at buenemann1997 (); buenemann1998 (). Similar results have been obtained more recently by means of the DMFT method ono2003 (); facio2017 ().

The analysis of the role of band degeneracy and Hund coupling in the development of superconductivity in multiband Hubbard models represents another topic of great interest, particularly relevant for iron-based superconductors. However, treating nonlocal pairing beyond perturbative approximations is particularly difficult. A recent DMFT study on a three-band Hubbard model highlighted the emergence of on-site (i.e., local) triplet superconductivity at finite doping for hoshino2015 (), in agreement with previous results obtained in the large limit, within an Hartree-Fock-Bogoliubov approach zegrodnik2012 () and the Gutzwiller approximation zegrodnik2013 (). Here, spin-triplet superconductivity is related to the emergence of local magnetic moments, which originate from the Hund coupling and are enhanced by an Ising anisotropy that suppresses fluctuations among different spin configurations. The presence of non-local pairing (i.e., with -wave symmetry) is much more difficult to assess within DMFT, since it would require a cluster extension, which is computationally heavy for multiorbital systems.

In addition to superconductivity, long-range magnetic order may be stabilized in a relatively large region of the phase diagram for . Within the two-band model, various calculations highlighted the existence of itinerant ferromagnetism for , which can be stabilized by the double-exchange mechanism for held1998 (); momoi1998 (); kubo2009 (); peters2010 (). In addition, recent DMFT calculations on the three-band model hoshino2015 () suggested the possibility to have antiferromagnetism close to half filling and ferromagnetism in a wide doping region at large values of the Coulomb repulsion. In the limit, the situation is delicate; in fact, the model with degenerate bands possess an enlarged SU symmetry, which is generated by spin and orbital degrees of freedom. Models with SU symmetry have been investigated within the strong-coupling limit, i.e., within the Heisenberg model hermele2011 (). In the square lattice for (corresponding to two electronic bands), a variety of numerical calculations suggested the presence of a spontaneous symmetry breaking in the ground state, for both one and two particles per site corboz2011 (); kim2017 (). Within the Hubbard model for , quantum fluctuations could be sufficiently strong to destroy magnetic/orbital order at small values of the Coulomb interaction even at half filling in the presence of a perfect nesting of the underlying Fermi surface (instead, for , the ground state has long-range magnetic order for any value of at half filling). In the weak- and intermediate-coupling limit, the cases with and have been considered in a generalized Hubbard-Heisenberg model at half filling, suggesting that for and small values of the antiferromagnetic coupling a -density wave state is stabilized assaad2005 ().

In this paper, we consider the two-band Hubbard model with degenerate bands on a square lattice, as the simplest case to investigate the role of the inter-orbital Coulomb repulsion and Hund coupling. We analyze the model by means of the variational Monte Carlo method. This approach, which works directly in two spatial dimensions, allows us to present a point of view that is complementary with respect to previous DMFT studies, which apply to infinite dimensions. First of all, we locate the MIT at commensurate fillings when no magnetic/orbital order is considered within the variational Ansatz. For the generic case with , the Mott transition appears to be continuous for ; instead, for , it is first order for small values of and turns to be continuous when the Hund coupling is increased. At half filling, the Mott transition is also accompanied by the stabilization of a sizable on-site triplet pairing, which survives in a wide region of doping around . A small singlet pairing with -wave symmetry is also observed in a narrow region close to for sufficiently small values of the Hund coupling. A finite singlet pairing can be stabilized also for , thus breaking the SU(4) symmetry in the variational wave function; in this case the MIT is first order; by contrast, when a fully-symmetric Ansatz is considered, the Mott transition becomes continuous. At quarter filling and close to it, neither triplet nor singlet pairing can be stabilized, indicating that superconductivity is not present around in the two-band Hubbard model with degenerate bands.

Symmetry-breaking states can be studied by including magnetic/orbital order within the variational Ansatz. At half filling, antiferromagnetic order is stabilized for small and intermediate values of the Coulomb interaction even for , suggesting that the SU(4) symmetry can be broken at small values of . At quarter filling, the metallic phase is stable for small values of , while the Mott insulator acquires both ferromagnetic and antiferro-orbital orders, in agreement with previous calculations held1998 (); momoi1998 (); kubo2009 (); peters2010 (). For , a wide region of ferromagnetism (without orbital order) is found for large values of the Coulomb repulsion. For small values of the Hund coupling, phase separation for may appear. By contrast, the region of stability of antiferromagnetism is limited to dopings close to . In the presence of magnetic order, triplet superconductivity is strongly suppressed close to half filling, coexisting with antiferromagnetic order.

The rest of the paper is organized as follows: In Sec. II, we introduce the two-band Hubbard model and the variational wave functions that are used within the Monte Carlo method. In Sec. III, we describe our results on the metal-insulator transitions, superconductivity, and magnetic/orbital order. Finally, in Sec. IV, we draw our conclusions.

## Ii Model and method

We consider the two-band Hubbard model defined by:

(1) |

where the kinetic term describes hopping processes of electrons within two degenerate orbitals:

(2) |

here () creates (destroys) an electron with spin on site and orbital and is the nearest-neighbor hopping amplitude on the square lattice. The interaction term includes four different contributions:

(3) | |||||

where is the electronic density per spin on site and orbital . These four terms represent the intra-orbital (inter-orbital) Coulomb interaction () and the spin-flip (pair-hopping) Hund term (). In the following, we set and georges2013 ().

Our numerical results are obtained by means of the variational Monte Carlo method, which is based on the definition of suitable wave functions to approximate the ground-state properties beyond perturbative approaches. In particular, we consider the so-called Jastrow-Slater wave functions that extend the original formulation proposed by Gutzwiller to include correlations effects on top of uncorrelated states gutzwiller1963 (); yokoyama1987 (). Our variational states are described by:

(4) |

here, is the density-density Jastrow factor, which is defined by:

(5) |

where is the electron density on site and orbital ; (that include also the local Gutzwiller term for and ) are pseudopotentials that are optimized for every independent distance . In the following, we will consider and . Moreover, the Fourier transform of the intra- and inter-orbital Jastrow terms will be denoted by and , respectively. The Jastrow factor has been shown to be crucial in describing a Mott insulating state within the single-band Hubbard model capello2005 (). As far as the two-band Hubbard model is concerned, the role of the Jastrow factor has been already highlighted in a variational Monte Carlo study of the orbital-selective Mott insulator tocchio2016a (). Then, is an uncorrelated state that is constructed from an auxiliary (quadratic) Hamiltonian:

(6) |

where is the kinetic term defined in Eq. (2), includes electron pairing and chemical potential:

(7) | |||||

and incorporate magnetic and orbital orders:

(8) | |||||

(9) | |||||

In Eq. (7), describes (on-site inter-orbital) triplet pairing, (nearest-neighbor intra-orbital) singlet pairing with symmetry, namely is its Fourier transform. In Eqs. (8) and (9), , , , and represent staggered and uniform parameters for magnetic and orbital orders. All these terms are further variational parameters that may be optimized in order to minimize the variational energy.

In the generic case with a finite Hund coupling , wave functions with no magnetic/orbital order can be constructed by considering only and . Notice that the latter one breaks the spin SU(2) symmetry whenever a triplet pairing is considered, without necessarily leading to a magnetic order. As far as the Jastrow factor is concerned, for the generic case with a finite Hund coupling, different intra- and inter-orbital pseudopotentials are allowed in Eq. (5), i.e., . For , a fully-symmetric wave function requires no pairing terms in Eq. (6), i.e., only can be retained in the auxiliary Hamiltonian, and a Jastrow factor that only involves the total electron density on each site, i.e., . Finally, states with magnetic and/or orbital order (with either staggered or uniform patterns) are easily obtained by also including and/or .

In order to assess the metallic or insulating nature of the ground state, we calculate the density-density structure factor , defined as:

(10) |

where indicates the expectation value over the variational wave function of Eq. (4). Indeed, a metallic phase has for , corresponding to the existence of gapless excitations, while an insulator is expected to have feynman1954 (); tocchio2011 (). Within our definition of the variational wave function, the metallic or insulating character can be also detected by looking at the small- behavior of the Jastrow factor, as shown within the single-band Hubbard model capello2005 (). In the two-band Hubbard model, the metallic phase is described by [and ], while the Mott insulating phase has instead (and ).

## Iii Results

In this section, we show our main results for the metal-insulator transitions at half filling () and quarter filling (), including the case where magnetic and orbital orders are prevented, and for superconductivity for densities between and . Most of the calculations are performed on the cluster with periodic (antiperiodic) boundary conditions along the () direction, in order to have a non-degenerate ground state for .

### iii.1 The metal-insulator transitions without magnetic/orbital orders

Let us start to study the MIT at commensurate electron densities, and , when no magnetic/orbital orders are allowed within the variational wave function. The results for and are reported in Fig. 1. For , , and , the Mott transition is first order, since two different wave functions, whose energies cross at , can be stabilized (in the vicinity of the MIT). While for small values of the Coulomb interaction, the best variational state is metallic with in the limit of , for large , the lowest-energy state is insulating with , see Fig. 2. This modification in the density-density correlations is triggered by the Jastrow factor, e.g., in the metal, while in the insulator, see Fig. 2.

We mention that the region where metastable solutions can be stabilized shrinks as increases, thus suggesting that the transition may become second order for a large enough value of the Hund coupling. Indeed, for , the MIT appears to be continuous, with no metastable solutions that can be obtained, see Fig. 1. Still, the small- behavior of the Jastrow factor is different for and , as in the single-band Hubbard model, where the Mott transition is continuous capello2006 (). Furthermore, our variational approach reproduces the well-known fact that decreases with increasing , since the Mott insulator with localized moments may take advantage of the Hund coupling ono2003 (); han1998 ().

Within the metallic regime, there is no appreciable gain when including superconducting pairing (either singlet or triplet) in the auxiliary Hamiltonian of Eq. (6); a similar result has been obtained in the paramagnetic solution of the single-band Hubbard model, where the metallic phase at half filling has vanishingly small pairing correlations tocchio2012 (); dayal2012 (). In addition, the intra- and inter-orbital Jastrow pseudopotentials are approximately equal for every distance, indicating that the correlation between two electrons on the same orbital is similar to the one between two electrons on different orbitals. By contrast, within the insulating phase, the intra-orbital Jastrow factor is larger than the inter-orbital one, implying that configurations with two electrons on the same orbital are penalized with respect to the ones with two electrons on different orbitals, as expected in the presence of a finite value of . Only for small values of , a (nearest-neighbor intra-orbital) singlet pairing with symmetry can be stabilized (see Fig 1), similarly to what occurs in the single-band Hubbard model at half filling capello2006 (); tocchio2012 (). Most importantly, a strong (on-site inter-orbital) triplet pairing is stabilized by the presence of a finite Hund coupling, giving a sizable gain in the variational energy with respect to the case with no pairing (see also Sec. III.3). Nonetheless, we must emphasize that the Jastrow factor with and is able to destroy the superconducting long-range order that is present in the uncorrelated wave function capello2006 (). Therefore, the presence of electronic pairing in leads to the existence of “preformed pairs” without phase coherence in the full correlated wave function of Eq. (4), as in the single-band Hubbard model. The relevant difference with respect to the latter case is that here “preformed pairs” do not form singlets with symmetry, but triplets with (on-site) symmetry.

We now briefly discuss the case with at . Here, whenever the variational wave function is taken to have a full SU(4) symmetry (i.e., by only considering the kinetic term in the auxiliary Hamiltonian and imposing ), the transition appears to be continuous (at ), with no metastable solutions in the energy optimization, see Fig. 3. However, by allowing different intra- and inter-orbital Jastrow factors in the variational optimization, another insulating solution exists, which is energetically favorable for , see Fig. 3. Then, this insulating state can be further improved by considering the electron (singlet) pairing (with symmetry) in the auxiliary Hamiltonian, further lowering the transition to . As before, the Jastrow factor prevents the existence of off-diagonal superconducting order.

Let us now investigate the case with , for which the results are shown in Fig. 4. In contrast to the half-filled case, here the Mott transition is always continuous and is marked by a progressive change in the small- behavior of the Jastrow factor, see Fig. 5. Remarkably, no gain in the variational energy is detected by allowing (on-site inter-orbital) triplet or (nearest-neighbor intra- or inter-orbital) singlet pairings, both in the metallic and the insulating phases. In addition, the intra- and inter-orbital Jastrow pseudopotentials are very similar, implying that the variational wave function remains fully symmetric not only for but also for . For , we find that , validating the fact that, in the absence of the Hund coupling, the maximum value of is obtained at half filling. Moreover, our calculations confirm the fact that, when restricting to the case with no magnetic or orbital order, the effect of the Hund coupling is to shift upward the MIT, as previously suggested by DMFT and slave-particle approaches demedici2011a (); demedici2011b (). In fact, the insulator with one electron per site does not have any substantial advantage from the presence of the Hund coupling, while the metallic phase may develop weak ferromagnetic correlations, thus gaining some energy. We would like to mention that, given the very gradual modification of the Jastrow factor (and correspondingly the density-density correlations), it is difficult to give a precise determination of when considering fully-symmetric wave functions (also for the case with and , see above); locating with high precision is however beyond the scope of this work.

### iii.2 The metal-insulator transitions with magnetic/orbital orders

The above picture for the metal-insulator transitions at and drastically changes when magnetic and/or orbital order is allowed within the non-interacting wave function, i.e., when also the last two terms of Eq. (6) are considered. At half filling, a finite (staggered) magnetic order can be clearly stabilized for (while no orbital order is detected). Notice that, in the case with , magnetic and orbital orders are related by SU(4) symmetry and, therefore, also an orbital order can be found. The optimized antiferromagnetic parameter of Eq. (8) is reported in Fig. 6, for and . In the former case, is significantly reduced with respect to the single-band model, which is also reported for comparison. In the presence of a finite , the triplet pairing is vanishing (or very small); however, a variational wave function with no magnetic order but a finite triplet pairing can be still stabilized as a local minimum at higher variational energies. Our results for are compatible with a finite magnetic order down to , with an exponentially small magnetization for . Given the smallness of the energy gain due to in the weak-coupling limit (i.e., ), we are not able to exclude the possibility that antiferromagnetism sets in at a (small) finite value of and not exactly at . Nevertheless, our variational calculations clearly support the existence of antiferromagnetism at half filling for intermediate values of . Moreover, since the SU(4) Heisenberg model with two (fermionic) particles per site is expected to be ordered kim2017 () and since a finite Hund coupling cooperates with the super-exchange mechanism to favor staggered magnetism, we foresee that magnetic order should survive for any value of up to .

For , no evidence for antiferromagnetic order is obtained, at least for . Instead, in the presence of a finite Hund coupling, a considerable energy gain is found in the strong-coupling regime by allowing both ferromagnetic and antiferro-orbital order, since virtual-hopping processes favor configurations in which two electrons on neighbor sites have parallel spins and reside on different orbitals kubo2009 (); peters2010 (). Indeed, for sufficiently large electron-electron repulsion, the best variational state is insulating with saturated magnetization (where ) and a finite in Eq. (9). By contrast, for small values of , a fully-symmetric metal with and no orbital order is found. No intermediate values of can be stabilized with orbital order. The results for are reported in Fig. 7, where a first-order phase transition between a metallic state with and no orbital order and an insulator with appears at .

### iii.3 Superconductivity and magnetism

In the single-band Hubbard model, several calculations have suggested that (singlet) -wave superconductivity emerges upon doping the Mott insulating state at half filling halboth2000 (); maier2005 (); eichenberger2007 (); gull2013 (); yokoyama2013 (); kaczmarczyk2013 (); deng2015 (); tocchio2016b (). Within the resonating valence-bond picture anderson1987a (); anderson1987b (); baskaran1988 (), this result can be explained by the existence of “preformed” electron pairs in the Mott insulator, where conduction is impeded by the strong electron-electron repulsion; then, phase coherence of mobile pairs emerge upon hole doping. In our variational picture, a necessary condition for having finite superconducting correlations is the presence of a finite pairing amplitude in the auxiliary Hamiltonian of Eq. 7. Indeed, in the single-band model, a finite BCS pairing with -wave symmetry can be stabilized for moderate and large values of tocchio2016b (). This picture becomes less robust in the multiband Hubbard model with degenerate electronic bands. For very small values of the Hund coupling (including ), a finite pairing amplitude with symmetry can be stabilized at half filling, see Fig. 3; however, drops to zero for very small dopings, i.e., around . In addition, for , no intra-orbital pairing can be stabilized in the wave function. We would like to mention that one way to recover a finite singlet pairing at reasonably large dopings is to break the symmetry between the inter- and the intra-orbital Coulomb repulsion, e.g., considering but still . In this case, orbital fluctuations are reduced (since configurations with two electrons on different orbitals are favored over the ones with a doubly-occupied orbital) and the resulting physical behavior can be assimilated to the one of two (weakly-coupled) single-band Hubbard models (one for each orbital). Therefore, a finite -wave pairing can be stabilized at finite dopings. We also mention that, in the opposite limit with , an on-site -wave pairing is present close to half filling, since doubly-occupied orbitals are favored over singly-occupied ones. Remarkably, these two kinds of pairings compete with each other and no singlet pairing can be stabilized away from half filling in the isotropic case with .

When no magnetic and orbital order are allowed in the variational wave function, a sizable interband triplet pairing is present in the vicinity of for and sufficiently large Coulomb repulsion , see Fig. 8. This is a consequence of the fact that, on each site, states are favored when ; a similar feature, with the developments of large local moments, has been also suggested by a recent DMFT study of the three-band model hoshino2015 (). However, in contrast to the latter work, which found that an Ising anisotropy in the Hund coupling is important to stabilize triplet superconductivity, we have evidence that a finite triplet pairing is present also in the isotropic case, which is modeled by the Hamiltonian of Eq. (3). It must be emphasized that, away from half filling, the presence of a finite electron pairing in the uncorrelated wave function implies a true long-range order, since the Jastrow pseudopotential has . As expected, the strength of triplet superconductivity is proportional to the Hund coupling, thus implying that the doping region in which enlarges with increasing , see Fig. 8.

When also magnetism is included in the variational wave function, superconductivity is largely suppressed. First of all, antiferromagnetic correlations are strong for electron densities close to half filling. Here, we can consider wave functions that contain both electron pairing and antiferromagnetism and optimize and together. The results are shown in Fig. 9 for . When is present, triplet pairing is strongly reduced close to half filling, leading to an antiferromagnetic metal with no pairing correlations. For , a tiny triplet superconductivity emerges around , where antiferromagnetism is still present, thus leading to a coexistence between these two order parameters. The pairing amplitude becomes much stronger when increasing the value of the Coulomb interaction, e.g., for , where displays a dome-like feature with a broad maximum at . However, in the presence of a finite Hund coupling also ferromagnetism becomes competitive in energy, especially when is large. A direct comparison between the superconducting state (with or without antiferromagnetic order) and the ferromagnetic one (with or without orbital order) allows us to draw the phase diagram of Fig. 10 for . In this case, the best variational state has antiferro-orbital order for , while a uniform ferromagnet exists in a wide region at finite electron densities and large . Instead, close to , a paramagnetic metal intrudes between these two ferromagnetic states. Our results are in qualitative agreement with previous variational kubo2009 () and DMFT peters2010 () calculations, which found the existence of uniform ferromagnetism at large values of the Coulomb repulsion for . Orbital order should survive in a tiny region close to quarter filling; however, even on the largest cluster that we considered (i.e., ) at (which is the closest available density to quarter filling that allows a direct comparison between ferromagnetic and paramagnetic states) the ferromagnetic wave function has a slightly higher energy than the paramagnetic one. For , phase separation is expected to take place close to the (first-order) transition between the paramagnetic-ferromagnetic transition. For larger values of , the ferromagnetic state can be stabilized also close to quarter filling, thus eliminating phase separation (not shown). The possibility to have triplet superconductivity inside the ferromagnetic region could be investigated by using Pfaffian wave functions spanu2008 (), in which pairing is considered for electrons with parallel spins. This kind of approach goes well beyond the scope of the present work.

## Iv Conclusions

We have considered the two-band Hubbard model with degenerate electronic bands by using variational wave functions and Monte Carlo techniques. At integer fillings with and , we have first investigated the metal-insulator transitions when both magnetic and orbital order are not considered. In this regime, our results for the location of the MIT, as a function of the Hund coupling , are qualitatively in agreement with previous DMFT and slave-particle approaches demedici2011a (); demedici2011b (). At half filling for , the transition is first (second) order for small (large) values of the Hund coupling, with a sizable triplet pairing within the Mott insulator (still, no superconducting long-range order is established at half filling, because of the strongly repulsive Jastrow factor). At quarter filling, the transition is second order with no finite pairing neither in the metallic nor in the insulating phase.

We have then included the possibility to stabilize magnetic and/or orbital order. At half filling, a clear evidence for antiferromagnetic order has been obtained for . In particular, the qualitative behavior of the magnetic parameter resembles the one of the single-band Hubbard model, where antiferromagnetic order sets in at ; therefore, our results suggest that the ground state for is antiferromagnetically ordered for any positive value of the Coulomb interaction . Triplet pairing is not present when a finite antiferromagnetic parameter is stabilized. At quarter filling, no sign of antiferromagnetic order is detected (up to ); instead for , the ground state shows a first-order phase transition from a metallic state for small values of the electron-electron interaction to an insulator with staggered orbital order and ferromagnetic correlations in the strong-coupling regime.

At intermediate electron dopings with , when both magnetic and orbital order are not included, a sizable triplet pairing is present for finite values of the Hund coupling and sufficiently large electron-electron interactions, i.e., when the Mott insulator at is doped. A similar trend has been recently found by DMFT calculations on the three-band Hubbard model hoshino2015 (). However, in our case, the Ising anisotropy in the Hund coupling is not necessary to obtain triplet pairing. We report that, at odds with the single-band Hubbard model, no sizable singlet pairing is instead present away from . When magnetic order is also considered within the variational wave function, triplet superconductivity is strongly suppressed by antiferromagnetic order close to ; furthermore, the region where superconductivity can be stabilized is also reduced by the presence of ferromagnetism, which is competitive in a wide range of densities for large Coulomb repulsions. The possibility to have a coexistence of triplet pairing and ferromagnetism could be considered by extending our variational approach to Pfaffian states, which is however quite expensive for multiband models and goes beyond the scope of this work.

###### Acknowledgements.

We thank M. Fabrizio and L. de’ Medici for useful discussions.## References

- (1) J.P.F. LeBlanc, A.E. Antipov, F. Becca, I.W. Bulik, G. Kin-Lic Chan, C.-M. Chung, Y. Deng, M. Ferrero, T.M. Henderson, C.A. Jiménez-Hoyos, E. Kozik, X.-W. Liu, A.J. Millis, N.V. Prokof’ev, M. Qin, G.E. Scuseria, H. Shi, B.V. Svistunov, L.F. Tocchio, I.S. Tupitsyn, S.R. White, S. Zhang, B.-X. Zheng, Z. Zhu, and E. Gull, Phys. Rev. X 5, 041041 (2015).
- (2) B.-X. Zheng, C.-M. Chung, P. Corboz, G. Ehlers, M.-P. Qin, R.M. Noack, H. Shi, S.R. White, S. Zhang, and G. Kin.-Lic. Chan, Science 358, 1155 (2017).
- (3) E. Dagotto, Rev. Mod. Phys. 66, 763 (1994).
- (4) M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys. 70, 1039 (1998).
- (5) P.A. Lee and X.-G. Wen, Phys. Rev. B78, 144517 (2008).
- (6) M. Daghofer, A. Nicholson, A. Moreo, and E. Dagotto, Phys. Rev. B81, 014511 (2010).
- (7) R. Yu and Q. Si, Phys. Rev. B84, 235115 (2011).
- (8) J. Hu, Sci. Bull. 61, 561 (2016).
- (9) R.M. Fernandes and A.V. Chubukov, Rep. Prog. Phys. 80, 014503 (2017).
- (10) W.F. Brinkman and T.M. Rice, Phys. Rev. B2, 4302 (1970).
- (11) D. Vollhardt, Rev. Mod. Phys. 56, 99 (1984).
- (12) A. Georges, G. Kotliar, W. Krauth, and M.J. Rozenberg, Rev. Mod. Phys. 68, 13 (1996).
- (13) G. Kotliar and A.E. Ruckenstein, Phys. Rev. Lett. 57, 1362 (1986).
- (14) M. Capello, F. Becca, S. Yunoki, and S. Sorella, Phys. Rev. B73, 245116 (2006).
- (15) J.P. Lu, Phys. Rev. B49, 5687 (1994).
- (16) M. J. Rozenberg, Phys. Rev. B55, R4855(R) (1997).
- (17) Y. Ono, M. Potthoff, and R. Bulla, Phys. Rev. B67, 035119 (2003).
- (18) E. Koch, O. Gunnarsson, and R.M. Martin, Phys. Rev. B60, 15714 (1999).
- (19) J.E. Han, M. Jarrell, and D.L. Cox, Phys. Rev. B58, 4199(R) (1998).
- (20) L. de’ Medici, J. Mravlje, and A. Georges, Phys. Rev. Lett. 107, 256401 (2011).
- (21) L. de’ Medici, Phys. Rev. B83, 205112 (2011).
- (22) L. Fanfarillo and E. Bascones, Phys. Rev. B92, 075136 (2015).
- (23) Y. Nomura, S. Sakai, and R. Arita, Phys. Rev. B91, 235107 (2015).
- (24) M. Capello, F. Becca, M. Fabrizio, S. Sorella, and E. Tosatti, Phys. Rev. Lett. 94, 026406 (2005).
- (25) J. Bünemann and W. Weber, Phys. Rev. B55, 4011 (1997).
- (26) J. Bünemann, W. Weber, and F. Gebhard, Phys. Rev. B57, 6896 (1998).
- (27) J.I. Facio, V. Vildosola, D.J. Garcia, and P.S. Cornaglia, Phys. Rev. B95, 085119 (2017).
- (28) S. Hoshino and P. Werner, Phys. Rev. Lett. 115, 247001 (2015).
- (29) M. Zegrodnik and J. Spalek, Phys. Rev. B86, 014505 (2012).
- (30) M. Zegrodnik, J. Spalek, and J. Bünemann, New J. Phys. 15, 073050 (2013).
- (31) K. Held and D. Vollhardt, Eur. Phys. J. B 5, 473 (1998).
- (32) T. Momoi and K. Kubo, Phys. Rev. B58, 567 (1998).
- (33) K. Kubo, Phys. Rev. B79, 020407 (2009).
- (34) R. Peters and T. Pruschke, Phys. Rev. B81, 035112 (2010).
- (35) M. Hermele and V. Gurarie, Phys. Rev. B84, 174441 (2011).
- (36) P. Corboz, A.M. Läuchli, K. Penc, M. Troyer, and F. Mila, Phys. Rev. Lett. 107, 215301 (2011).
- (37) F.H. Kim, K. Penc, P. Nataf, and F. Mila, Phys. Rev. B96, 205142 (2017).
- (38) F.F. Assaad, Phys. Rev. B71, 075103 (2005).
- (39) A. Georges, L. de’ Medici, and J. Mravlje, Annu. Rev. Condens. Matter Phys. 4, 137 (2013).
- (40) M.C. Gutzwiller, Phys. Rev. Lett. 10, 159 (1963).
- (41) H. Yokoyama and H. Shiba, J. Phys. Soc. Jpn. 56, 1490 (1987).
- (42) L.F. Tocchio, F. Arrigoni. S. Sorella, and F. Becca, J. Phys.: Condens. Matter 28, 105602 (2016).
- (43) R.P. Feynman, Phys. Rev. 94, 262 (1954).
- (44) L.F. Tocchio, F. Becca, and C. Gros, Phys. Rev. B83, 195138 (2011).
- (45) L.F. Tocchio, F. Becca, and C. Gros, Phys. Rev. B86, 035102 (2012).
- (46) S. Dayal, R.T. Clay, and S. Mazumdar, Phys. Rev. B85, 165141 (2012).
- (47) C.J. Halboth and W. Metzner, Phys. Rev. Lett. 85, 5162 (2000).
- (48) T.A. Maier, M. Jarrell, T.C. Schulthess, P.R.C. Kent, and J.B. White, Phys. Rev. Lett. 95, 237001 (2005).
- (49) D. Eichenberger and D. Baeriswyl, Phys. Rev. B76, 180504(R) (2007).
- (50) E. Gull, O. Parcollet, and A.J. Millis, Phys. Rev. Lett. 110, 216405 (2013).
- (51) H. Yokoyama, M. Ogata, Y. Tanaka, K. Kobayashi, and H. Tsuchiura, J. Phys. Soc. Jpn. 82, 014707 (2013).
- (52) J. Kaczmarczyk, J. Spalek, T. Schickling, and J. Bünemann, Phys. Rev. B88, 115127 (2013).
- (53) Y. Deng, E. Kozik, N. V. Prokof’ev, and B.V. Svistunov, Europhys. Lett. 110, 57001 (2015).
- (54) L.F. Tocchio, F. Becca, and S. Sorella, Phys. Rev. B94, 195126 (2016).
- (55) P.W. Anderson, Science 235, 1196 (1987).
- (56) P.W. Anderson, G. Baskaran, Z. Zou, and T. Hsu, Phys. Rev. Lett. 58, 2790 (1987).
- (57) G. Baskaran and P.W. Anderson, Phys. Rev. B37, 580(R) (1988).
- (58) L. Spanu, M. Lugas, F. Becca, and S. Sorella, Phys. Rev. B77, 024510 (2008).