Magnetic phase diagram of the infinite-U Hubbard model with nearest- and next nearest-neighbor hoppings

We study the infinite-U Hubbard model on ladders of 2, 4 and 6 legs with nearest (t) and next-nearest (t') neighbor hoppings by means of the density-matrix renormalization group algorithm. In particular, we analyze the stability of the Nagaoka state for several values of t' when we vary the electron density $(\rho)$ from half-filling to the low-density limit. We build the two-dimensional phase diagram, where the fully spin-polarized and paramagnetic states prevail. We find that the inclusion of a non-frustrating next nearest neighbor hopping stabilizes the fully spin-polarized phase up until |t'/t|=0.5. Surprisingly, for this value of t', the ground state is fully spin-polarized for almost any electron density 1 $\gtrsim \rho \gtrsim$ 0, connecting the Nagaoka state to itinerant ferromagnetism at low density. Also, we find that the previously found checkerboard insulator phase at t'=0 and $\rho$=0.75 is unstable against t'.


I. INTRODUCTION
The Hubbard model was first introduced in 1963 to explain itinerant ferromagnetism in transition metals 1−3 . Even though it was built as the simplest model capable of describing the behavior of correlated materials 4 , it has been proven to be useful in the interpretation of a wide variety of phenomena ranging from metal-insulator transitions 5,6 to high-temperature superconductivity 7 . In the past years this model has been brought back into focus due to its applicability in the description of ultracold atoms in optical lattices 8−10 .
Despite being a rather simple model when written down, most of the progress in the understanding of the Hubbard model has been made numerically, either by exact diagonalization of small clusters, mean field approaches or by using more sophisticated numerical techniques such as quantum Monte Carlo (QMC) or density-matrix renormalization group (DMRG). Only a few mathematically rigorous results regarding this model exist to date 11,12 . The Nagaoka's theorem 13,14 is one of these few exact results, making it an important and solid starting point to study the phase diagram of the Hubbard model. It states that, when the system has one electron less than half-filling, U → ∞, and the lattice satisfies certain connectivity conditions, the ground state of the system is a fully spin-polarized ferromagnetic state (FSP) and it is unique apart from the trivial spin degeneracy. These connectivity conditions require that the smallest loop must be no longer than 4 sites and the kinetic energy of the hole motion around this loop must not be frustrated. In the subsequent years to the Nagaoka's theorem, a lot of effort was put into trying to widen this isolated point of the phase diagram by relaxing some of the requirements of this theorem. For example, very recently, the condition regarding the loop size has been extended to cover larger loops 15 , proving that the two-dimensional honeycomb lattice (loop size = 6) also has a FSP ground state.
When the accumulated sign of the hoppings along the minimum loop is negative, the kinetic energy of the hole motion is frustrated and the Nagaoka's theorem is no longer valid. Nevertheless, it has been shown that the FSP state can survive in the presence of small enough values of frustrating hoppings 16 . On the other hand, high frustration can lead to an antiferromagnetic Néel order, as it happens in the isotropic triangular lattice with positive hoppings 17 . Surprisingly enough, in this case, the antiferromagnetic ground state is classical and has the maximum staggered magnetization possible 16,18 . This kind of classical antiferromagnetic state has also been found in the square lattice with a frustrating next-nearest neighbor hopping 18 .
Also, it is worth mentioning that the Nagaoka's theorem is only valid for finite lattices, where the one hole away from half-filling condition makes sense. As we get closer to the thermodynamic limit, clearly, this condition means that the electron density tends to half-filling. Because of this, there has been a long-standing question regarding the existence of the FSP phase in the thermodynamic limit at finite hole-doping. The Hubbard model on the square lattice with U → ∞ and varying electron density is the easiest model where this problem can be studied. Early calculations showed that the FSP phase was unstable against hole doping beyond the Nagaoka's theorem conditions 19,20 , although more recent ones obtain a critical hole density in the thermodynamic limit around ρ c ≈ 0.8 21−24 . Nonetheless, there is still much to uncover, as other recently published results suggest that two holes away from half filling the existence of the FSP state depends strongly on the boundary conditions and the sizes of the finite lattices 25 . This shows that the mechanisms responsible for stabilizing the FSP state away from the Nagaoka's theorem conditions are not yet fully understood.
DMRG calculations 23 show that, also in the square lattice, new interesting phases appear below the critical electron density, ρ c = 0.8. For example, a commensurate checkerboard insulator state emerges at ρ = 0.75, and leads to a phase separation region between ρ = 0.8 and 0.75; below this point the system behaves as a paramag-net. Another interesting study 26 performed with dynamical mean-field theory (DMFT) shows that the inclusion of a next-nearest neighbor hopping, whenever not frustrating the loop conditions within the Nagaoka's theorem, can stabilize ferromagnetic phases for smaller values of ρ (reducing the paramagnetic region).
The weak to intermediate coupling regime (i.e. the infinite-U condition is not longer fulfilled) has been widely studied 27−32 and it was found that for low densities and |t /t| = 0.5 there is a fully spin-polarized phase. The presence of this low-density ferromagnetism even at relatively small values of U is due to the strong particlehole asymmetry and the van Hove singularity near the bottom of the band.
In this paper, we study in more detail the problem of the stability of the Nagaoka ferromagnetism (FSP state) against hole doping, combined with the inclusion of a non-frustrating next-nearest neighbor hopping on the square lattice by means of numerical calculations using DMRG. We also follow the evolution and stability of the checkerboard insulator previously mentioned along with the phase separation.
This paper is organized as follows: in Sec. II we introduce the Hamiltonian and describe the details of the methods employed in the rest of the paper. In Sec. III we show and discuss the results for two-to six-leg ladders and build the two-dimensional (2D) phase diagram. In Sec. IV we present the conclusions.

II. MODEL AND METHODS
The object of our present study is the Hubbard model, which we can write down as where t ij is the hopping integral between sites i and j and U is the on-site repulsion. We take the repulsion between electrons U → ∞ to remain within the Nagaoka's theorem conditions; this means that we can never have two electrons in the same site. Taking into account only two different hopping terms, we can write our implementation Hamiltonian as where t connects the nearest neighbors on a square lattice and t the next-nearest neighbors (see figure 1). The presence of a hopping t = 0 makes the lattice no longer bipartite, breaking the particle-hole symmetry. Given the importance of the sign of the interaction in the Hamiltonian, it is worth mentioning that we have chosen it differently for t and t for later convenience. The new operatorsc † iσ =ĉ † i,σ (1 −n i,σ ) ensure the exclusion of the doubly-occupied states imposed by the infinite-U condi- tion; and its new commutation relations are responsible for the complications in diagonalizing the Hamiltonian.
Also, to ensure the validity of the Nagaoka's theorem connectivity condition we need to check that the accumulated sign around the smallest loop is positive, meaning that where we can see that the sign of t is irrelevant. So t has to be positive to fulfill the connectivity condition and we will take t = 1 as the energy unit from now on. Taking t > 0, we can be sure that the ground state of the system will always be the FSP state at one electron less than half-filling.
Having guaranteed this starting point, it is up to us to study in further detail the effect of the inclusion of nextnearest neighbor hopping terms t over the existence and stability of the FSP state upon further hole-doping. To do so, we use the density-matrix renormalization group algorithm 33 based on the matrix product state (MPS) representation 34 contained in the ALPS libraries 35,36 . We solve ladders with fully open boundary conditions for L y =2, 4, 6 (legs) and several L x (rungs) when possible. We vary t and the electron density ρ to map the phase diagram, using between m = 2000 and 9000 DMRG states to ensure the convergence of our results.
The first of our aims is to determine the critical value of ρ c up to which the FSP state survives. In order to accomplish this task, we need a reliable signature to help us determine whether this state is, or not, the ground state of the system for a given hopping t and electron density ρ. There are two complementary methods we can use that will also help us characterize the rest of the phase diagram. One is to determine the magnetization of the system through the spin structure factor where we can take k as a continuous variable given the open boundary condition. Obtaining the spin structure factor of the system gives out more information about the magnetic ordering, but it relies on the calculation of the spin correlations, which are considerably less accurate than the energy calculations. The other alternative depends on calculating the ferromagnetic magnetization of a system solely through energy calculations. This allows us to calculate the total spin of the ground state by exploiting the degeneracy of a ferromagnetic ground state. To do so, we calculate the lowest energy for all subspaces with different quantum number S z , and find S max . Then, the normalized ferromagnetic magnetization value for a given electron density can be calculated as M = 2S max z /N e . In the case of the FSP state, this magnetization gives out M = 1 and in the case of the paramagnetic state, M = 0. This method is expected to be more accurate but more time-consuming as we have to compute the minimum energy state for several subspaces. Nevertheless, there is an issue with the given M formula: when there is an odd number of electrons in the system, the minimum value of the spin projection S z is 0.5 (and not 0), meaning that the minimum value of M is not zero. To solve this, in these cases we subtract to the S max z one half, and one electron to the N e to be able to represent the FSP and paramagnetic limits, M = 1 and M = 0. The corresponding formula of the magnetization for the odd number of electrons then results in

A. Two-leg ladders
We start by studying the simpler 2-leg ladder systems as we expect that the main characteristics of the phase diagram do not change much upon adding more legs. To cover the phase diagram we use 2 × 20 and 2 × 30 ladders and take t from 0.0 to 0.5 at 0.1 intervals. Also, we briefly comment results for t > 0.5 up to t = 1.0. For each chosen value of t we calculate the ground state energy within all S z subspaces for all electron densities below half-filling. This allows us to compute the magnetization curves for all our t values as a function of the electron density ρ.
In the square ladder (t = 0.0), the FSP state is present at high densities and survives up to ρ c = 0.8, as can be seen by the classical magnetization M = 1 region in figure 2. This value was previously obtained in 2-leg ladder systems 22,23 and in the 2D limit 21,23 , proving that it does not scale with the number of legs. Below this critical value of the electron density, the magnetization M is lowered up until ρ = 0.75, where it reaches M = 0. At ρ = 0.75 lies the checkerboard insulator; this phase consists on plaquettes of four sites and three electrons (each plaquette has the same electronic density as the lattice) in a FSP state, but building an antiferromagnetic order between plaquettes (see figure 3). This exotic state can be seen through the spin structure factor as a set of two broad peaks centered at ± π 2 (shown in next section for larger ladders). Between this phase and the FSP there is a phase separation region that arises as a combination of the FSP state with ρ = 0.8 and the checkerboard insulator. This phase separation region can be seen in the spin structure factor as a combination of the peaks that belong to each of the phases, lowered and broadened by the mixture (also shown in the next section for larger ladders). This phase can also be characterized by the uneven charge distribution 23 , where a certain part of the system has the same density n = 0.75 and order as the checkerboard insulator whereas the rest of the system has ferromagnetic order with n = 0.8. Between ρ = 0.75 and ρ = 0.6 there are intermediate phases with non-zero magnetization which disappear in wider ladders. Below ρ = 0.6 the ground state is paramagnetic, and it is signaled by the null ferromagnetic magnetization. These results, obtained for 2 × 20 and 2 × 30 lattices, are in complete agreement with the previous DMRG study by Liu et al. 23 ; and they provide us a good benchmark to start analyzing the effect of introducing the next-nearest neighbor hopping t .
When t is turned on, we find that the FSP state region starts growing and ρ c is lowered, as can be seen in figure  2. For example, when t increases to 0.30, the critical value of the electron density for which the FSP can be found moves down to ρ c = 0.62. This enhancement of the stability of the FSP phase when including t has also been found using DMFT in 2D 26 , where they report that at t = 0.1 the critical value of the density is ρ c = 0.705. For the same value of t , we obtain ρ c = 0.775. The difference in our results may be a signal that ρ c depends on the number of legs when approaching to the 2D limit for t = 0 (unlike what happens for t = 0.0). To understand why the t stabilizes the FSP state, Park et al. 26 have solved the four-site plaquette with three electrons. They have shown that the existence of t lowers more the FSP state energy than the low-spin state energy because of the quantum interference of different hole paths. As a consequence, the gap between these two states increases with t .
It is noteworthy that, if the value of the next-nearest neighbor hopping is half of the hopping on the square lattice, t = 0.5, the FSP state survives for every value of electron density; that is ρ c → 0 (see green triangles in figure 2). This is a remarkable result, because it connects the Nagaoka ferromagnetic phase (an exact result near half-filling from the infinite-U limit) with the FSP state found at low density with relatively small U . The latter phase was studied among others by Taniguchi et al. 31 and they found this ferromagnetic phase for a rather small U in the low-density limit around t 0.5, where the Fermi energy is close to the van Hove singularity. They have also shown that with increasing U (until U = 5) the FSP phase is the ground-state of the system for a wider regime of values of t and ρ.
When increasing t above 0.5, the FSP state region starts going back, giving in to the paramagnetic phase. This is a consequence of the magnetic behavior of the system when t → ∞, in this case the ladder splits in two independent chains where the ground states is a paramagnet 37 . Also, for every value of t = 0 we find intermediate phases with interpolating ferromagnetic magnetization between the FSP and paramagnetic regions, but we expect them to shrink and disappear when adding more legs and going closer to the 2D limit. An important difference with the t = 0.0 case is that there is no checkerboard insulator from t = 0.1 on. This also means that there is no phase separation. Given that the checkerboard insulator phase is absent even at t = 0.1, we decided to follow its evolution more closely. Taking a fixed electron density ρ = 0.75, we made runs varying t at 0.01 intervals. In figure 3 we plot the results for t = 0.00, 0.05, and 0.06. Clearly, the four-plaquette structure can be seen up until t c = 0.05, but as soon as t is increased this structure disappears and the ferromagnetic magnetization rises.
Regarding the validity of these results for longer ladders, it can be seen in figure 5 that the results do not change upon adding more rungs. The 2 × 20 lattice already is a good representation of the 2-leg thermodynamic limit for every value of t . The previous discussion on the magnetization can be summarized into the phase diagram shown in figure 4. The FSP phase is painted in green, being the green squares the last points for which we find this kind of state, coming down from the one electron less than half-filling condition of the Nagaoka's theorem. In blue we show the paramagnetic phase, being the blue circles the last points for which we find a paramagnetic phase, with M = 0, when increasing the electron density from zero. In both of these cases it is usually easy to determine, by studying the dependence of the energy with the value of S z , if the ground state lies only on the S z = 0 subspace or if it can be found in all subspaces up to S max z . The blue triangles in line at ρ = 0.75 represent the points of the phase diagram for which we find the checkerboard insulator phase. In the middle of this phase and the FSP, we shadowed the area in which the phase separation exists. The yellow region corresponds to the intermediate phases with intermediate ferromagnetic magnetization. Closer to the FSP they present clear ferromagnetic peaks in the spin structure factor that evolve quickly into the paramagnetic state. Here, the tendencies of the energy as a function of S z make it more complicated to identify the exact value of S max z , and therefore to analyze the properties of the ground state in this region. As mentioned above, increasing t above 0.5 causes the FSP region to shrink again (and the paramagnetic one to grow).
With respect to the charge distribution, we find that both the paramagnetic and the fully spin-polarized ferromagnetic phases show almost homogeneous distribution. On the other hand, the intermediate phases show a minor variation of the charge distribution along the ladder which, however, does not resemble that of a separation of phases where only two well distinctive densities appear.

B. Four-and six-leg ladders
To shed some light on the 2D behavior of the system we extended our calculations to 4-and 6-leg ladders. Leaning on our results for the 2-leg ladders, we calculate the magnetization value in 4 × 10, 4 × 12 and 6 × 8 lattices around the transition points for all the same values of t as in the 2-leg case to see how they scale. Also, we calculate the spin structure factor and other correlations. Away from the transition points the convergence of the results behaves very well even for the bigger lattices. But, close to transition, the computational effort increases and it becomes difficult to determine precisely the S max z of the ground state of the system. This is the main reason why we had to limit the number of rungs and legs used in our lattices. Previous studies already show that, even at t = 0.0, the amount of DMRG states needed to obtain accurate results in big lattices is huge 23 .
We show on figure 5 the magnetization value M for two arbitrary values of next-nearest neighbor hoppings t = 0.2 (top panel) and t = 0.4 (lower panel), and several lattice sizes. From these results we can see the qualitative behavior remains the same for all the lattices. For example, in the top panel we can see that the critical value ρ c seems to move a little in comparison with the 2-leg lattice, but we have roughly the same critical density for the 4-and 6-leg lattices. Also, it seems that the region corresponding to intermediate phases shrinks as the FSP phase region grows, making the transition to the paramagnetic phase more abrupt, also seen in the lower panel. We have observed that this behavior, shown for the t = 0.2 and t = 0.4 cases, holds for almost every hopping value below t = 0.5; t = 0.0 being the special case. At t = 0.0, the checkerboard insulator prevails and prevents the FSP state to take over and, instead, the intermediate phases become paramagnetic. On the other hand, when t = 0.5 there is already no place to move the critical electron density down, so the ground state remains FSP for all values of the electron densities below half-filling. agree, as ρ c = 0.8 in DMRG and ρ c = 0.815 in DMFT. These similarities between our DMRG results and DMFT seem to indicate that, around the phase transition, the non-local physics are unimportant. Nonetheless, this could be a mere coincidence and more research is needed to elucidate the source of this agreement. We can then conclude that the 4-and 6-leg lattices provide us enough information about the scaling of the transition points to build the 2D phase diagram of the system until t = 0.5, at least up to a certain small error. Above t = 0.5 the picture is completely different. We found a more pronounced scaling of the critical densities and intermediate phases that made it not possible for us to extract valuable information. One possible reason for this to happen has to do with the large t limit. While in the 2-leg ladders the t → ∞ limit results in two isolated chains, in a 2D system it results in two sets of independent square lattices. The difference in these limits may be responsible for the behavior above t = 0.5. A more detailed study in wider ladders is needed and this is why we restrain our results to the region below t = 0.5.
In figure 6 we show, for the 4×10 lattice, several calculations of the spin structure factors obtained at t = 0.0. In the upper-left panel we can see the last density for which the ground state of the system is FSP, ρ = 0.8. This phase is signaled by a large peak centered at k = 0 with certain asymmetry given by the lattice. In the lower-left panel we can see the signature of the checkerboard insulator, ρ = 0.75. The structure factor is composed by four peaks situated at k = ± π 2 , ± π 2 that arise from the antiferromagnetic alternation of the plaquettes. Note that, in this case, the ferromagnetic magnetization is zero and S(k = 0) = 0, because the ground state exists only at S z = 0. In the upper-right panel we can see the phase separation at ρ = 0.775, signaled by a combination of the surrounding spin structure factors. It shows a low ferromagnetic peak and the checkerboard insulator weakened four-peak structure. Finally, in the lowerright panel we show the paramagnetic phase at ρ = 0.725 where no clear peak can be seen. It is important to be aware that, for this lattice size, all these states are one electron away from each other. The density step is With all this gathered information, now we can compute the full phase diagram (shown in figure 7) and compare it with the 2-leg case. We have decided to take the transition points as the average of the ones in 4-and 6-leg ladders, with error equal as the difference in these results. As we expected before, the 2D phase diagram resembles the 2-leg ladders one, but with less area left for intermediate phases, specially around t = 0.1, and a certain growth of the FSP phase. In contrast with t = 0, the intermediate phases shrink but they do not disappear in the four-and six-leg ladders. For these ladders, the intermediate phase is ordered ferromagnetically (peak in the structure factor at S(k = 0)) but its magnetization is not saturated. The checkerboard insulator, even though we clearly find it in the 4-leg lattices (as can be seen in figure 6), was much harder to find in the 6 × 8 lattice. We had to use m = 9000 DMRG states and a small Zeeman field applied to pin this order; and check that the arising checkerboard insulator is not an excited state. Given that in the 4-and 6-leg lattices we observed that this phase vanishes when t ∼ 0.02 (much smaller than in the 2-leg case), we conjecture that the phase boundaries for the checkerboard insulator phase in the 2D limit should be rather small (t c ∼ 0.02 at most). This is in consonance with the small spin gap of the checkerboard phase at t = 0, about 10 −3 t.

IV. CONCLUSIONS
We have used density-matrix renormalization group to study the phase diagram of the infinite-U Hubbard model on square ladders with nearest and next-nearest neighbor hopping amplitudes t and t , respectively; and always in the absence of kinetic frustration. We have found that, for all ladder sizes, the presence of a nonfrustrating next-nearest neighbor hopping amplitude stabilizes the fully spin-polarized phase. With increasing t , the fully spin-polarized phase region grows in the phase diagram until t reaches a value of one half the nearest neighbor hopping amplitude. For this particular value (t = 0.5) the ground state of the system is always a fully spin-polarized state, regardless of the electron density chosen below half-filling. We have connected in the infinite-U limit the fully spin-polarized state from the Nagaoka's theorem (valid for one hole over half-filling) with the low-density ferromagnet (also FSP) that arise due to the van Hove singularity in the bottom of the band. It would be interesting to further investigate the behavior of the single-particle spectral densities in the different phases, and in particular for t = 0.5 in the Nagaoka's theorem case and in the low-density regime, to unwrap the transition between different types of ferromagnetism. For t > 0.5 we need to explore wider ladders to uncover the 2D behavior, but our results indicates that beyond this point, the fully spin-polarized phase region starts to shrink as t increases.
With regards to the intermediate phases, we conjecture that the previously found checkerboard insulator phase only survives for small values of t in the thermodynamic limit. This may be a consequence of the proximity to the fully spin-polarized state. The checkerboard insulator phase can only exist at ρ = 0.75, where one hole lives in each four-site plaquette. But, ρ c quickly goes below 0.75 when t is included and the fully spin-polarized state prevails over the checkerboard insulator. Moving away from the t = 0 case, the intermediate phases are less interesting and have generally a ferromagnetic behavior which seems to connect continuously the FSP phase with the paramagnetic phase. In this region it is more difficult to analyze the spin and charge behavior so we cannot ensure the size of this region or characterize the nature of the transition in the thermodynamic limit.
Also, we expect these results to be of interest and to contribute to the on-going research in optical lattices, which is where these models with large repulsions within particles can be experimentally realized.
V. ACKNOWLEDGMENT