Ground state phase diagram of the doped Hubbard model on the four-leg cylinder

We study the ground state properties of the Hubbard model on a four-leg cylinder with doped hole concentration per site δ (cid:2) 12 . 5% using density-matrix renormalization group. By keeping a large number of states for long system sizes, we ﬁnd that the nature of the ground state is remarkably sensitive to the presence of next-nearest-neighbor hopping t (cid:2) . Without t (cid:2) the ground state of the system corresponds to the insulating ﬁlled stripe phase with long-range charge-density-wave (CDW) order, and short-range incommensurate spin correlations appears. However, for a small negative t (cid:2) a phase characterized by coexisting algebraic d -wave superconducting (SC) and algebraic CDW correlations. In addition, it shows short-range spin-and fermion correlations consistent with a canonical Luther-Emery (LE) liquid, except that the charge and spin periodicities are consistent with half-ﬁlled stripes instead of the 4 k F and 2 k F wave vectors that are generic for one-dimensional chains. For a small positive t (cid:2) yet another phase takes over, showing similar SC and CDW correlations. However, the fermions are now characterized by a (nearly) inﬁnite correlation length while the gapped spin system is characterized by simple staggered antiferromagnetic correlations. We will show that this is consistent with a LE liquid formed from a weakly coupled (BCS like) d -wave superconductor on the ladder where the interactions have only the effect of stabilizing a cuprate style magnetic resonance.


I. INTRODUCTION
The Hubbard model is the simplest model that captures essential features of strongly interacting electrons realized in solids containing transition-metal and/or rare-earth elements, characterized by a compromise between kinetic energy and strong local repulsion.It describes a tight-binding model for the band structure with hopping matrix elements t i j , supplemented by an on-site Coulomb repulsion U , where c † iσ is the electron creation operator on site i = (x i , y i ) with spin σ , and n i = σ c † iσ c iσ is the electron number operator.In this paper, we will consider a square lattice geometry compactified on a cylinder, specializing the hopping integral to t i j = t for nearest-neighbor (NN) and t i j = t for nextnearest-neighbor (NNN) sites, respectively.This is among the most studied models in the history of physics in the case of fermions [1][2][3][4][5][6].Its simple appearance is an illusion; after tens of thousands of papers it has been brought under mathematical control in one and infinite dimensions [7,8].However, due to spectacular progress in the computational methodology, recently some solid results were reported [9][10][11][12][13][14][15][16][17][18][19][20][21][22][23][24].In part this is due to the exponential growth of computational resources, having the effect that for instance quantum Monte Carlo (QMC) now yields results at temperatures low enough to be of physical relevance [22,23].On the other hand, new algorithms were developed inspired from quantum information where the difficulty with, e.g., the Hubbard model is to keep track of the many-body entanglement.The density-matrix renormalization group (DMRG) is the oldest of this family [25].While it was designed to study systems in one dimension (1D), it has been made to work well for quasi-1D systems having a finite spatial extend in the transverse direction, i.e., the ladder systems, which is now the standard setting to compute the physics in two dimensions (2D) [26].
We reported recently on the presence of a Luther-Emery (LE) liquid characterized by algebraic superconducting (SC) order coexisting with the (dual) algebraic charge-densitywave (CDW) order in a particular parameter regime of both the Hubbard and closely related t-J model, where the critical role of t was shown to control the pattern of stripes and, as a consequence, superconductivity [20,21].Here we will report on how the ground states of the Hubbard model look in a much larger regime of physically relevant parameters [27,28].Another aspect is that already the very first four-leg cylinder results showed a strong appetite to form "intertwined order," specifically of the spin-stripe variety [20][21][22][23]29].The stripes have a long history going back to the 1980s when it was discovered, that according to classic mean-field theory (Hartree-Fock), in the case of large U the doped systems develop rather complex textures: domain walls form in the antiferromagnetic spin background, where the carriers localize forming "rivers of charge" [30].According to the mean-field theory, these stripes are stable under the condition that there is one carrier per domain-wall unit cell, for the reason that a gap opens at E F rendering these stripe phases to be insulating.The outcomes of various methods characterized by qualitatively different systematic errors were compared and such insulating filled stripes ground states were established.This includes the DMRG results on the four-leg cylinders at t = 0 and U = 8t [17][18][19][20].
All along there were indications that various ground states of a quite different nature may be quite close in energy in Hubbard and t-J models [10,16,19].The t-J model has the advantage that it can be easily deformed, involving longer range exchange interactions and so forth.As a function of these extra parameters a wealth of different phases, "stripy" and otherwise, were found on the four-leg cylinder [10,20,21,31].The early ladder results showed also such stripes, with the difference that these turned out to be "half-filled": one carrier is associated with two domain wall unit cells [10,31,32].Meanwhile, such stripe instabilities were observed in high-T c superconductors of the so-called 214 family [4].Strikingly, these are also half-filled, giving further support for the approach in general.
Using the large bond dimension DMRG simulations we will report here a surprise.The nature of the ground states of the doped Hubbard system appears to be extremely sensitive to the next-nearest-neighbor hopping t .We can diagnose the various phases in a precise fashion, the reason being that these are 1D systems that will reveal the universal behaviors of such systems at long length scales.However, different from the canonical chain systems at short distances, the ladders are more like 2D systems.This will alter the "numbers" at long distances which in turn may reveal aspects of relevance to the 2D case, given that the convergence as function of ladder width is expected to be rather rapid.
The main result is a zero-temperature phase diagram which is surprisingly sensitive to a single parameter: t , Fig. 1.Not much happens upon varying the doping δ and U quite a bit: it looks similar for the very strong coupling t-t -J model in Fig. 2.However, varying t by a tenth of t has the effect of stabilizing completely different ground states.Using the 1D diagnostics, the "purple" regime for t 0 is identified as the commensurately pinned "Luttinger liquid" chargedensity/spin-density wave that coincides with the filled stripes [30] in the 2D language.The phases found at positive and negative t are in the first place characterized by coexisting algebraic SC and CDW order which are in dual relation: the signature of the LE phase, and we call these therefore the LE1 and LE2 phases, indications for which were already reported in [20] and [23,29].
In other regards the LE1 and LE2 phases are, however, entirely different.The LE1 phase as first reported in [20] is quite like the canonical LE phase realized on chains.The spin correlations are short ranged, associated with the spin gap characterizing the Cooper pairs, while they show a periodicity which is twice the charge periodicity.However, the latter is counted like half-filled stripes and it is no longer related to the 2k F and 4k F spin and charge periodicities hard-wired in the chain systems.The equal-time fermion propagator demonstrates that the gap in their spectrum is very small, if not even vanishing.In addition, the spin correlations are short ranged but characterized by a simple two-sublattice periodicity.
We will argue that this can be naturally explained by asserting that this LE2 phase is a smooth continuation from a weakly coupled (BCS) d-wave superconductor living on the ladder.It is easy to see that the vanishing fermion mass revealed by the equal-time fermion correlator is just anchored in the fact that the electronic modes of the ladder will touch the nodal direction where the d-wave gap is vanishing.Turning to the spins, it is well established that in weakly interacting 1D systems the random phase approximation (RPA) is quite accurate [33].According to RPA the massless spin fluctuations carry so little spectral weight that these are not seen in the equal-time correlations.Instead, it is a well established principle that spin fluctuations in a conventional d-wave superconductor are susceptible to strong enhancements at FIG. 2. Ground state phase diagram of the t-t -J model in Eq. ( 2) as a function of hole doping concentration δ and t .The black dots are the data points.Here t = 1 and J = 1/3.staggered wave vectors, the prevalent explanation for the magnetic resonance in near-optimal cuprate superconductors [34,35].In the conclusions we will further discuss potential relations between these findings and the situation in these experimental systems.

II. THE PHASE DIAGRAM: OVERVIEW
The methodology.The DMRG method [25] has acquired quite a merit for the computation of the ground states of strictly 1D (chain) systems.Although the computations become rapidly more demanding, it was early on realized that it can be applied to systems with a transversal extent, and the results converge relatively rapidly towards the 2D physics [26].Two-leg ladders are clearly too narrow, while it was very recently found that one has to employ very large bond dimensions to get fully convergent results on four-leg ladders.Here we employ such ladders of width (number of sites) L y = 4 with cylindrical boundary condition in this compact direction, and open boundary condition in the extended direction with a length up to L x = 96.We keep up to m = 20000 states in each DMRG block with a typical truncation error ∼ 10 −6 and perform around 60 sweeps, which leads to excellent convergence of our results.Further details are provided in the Appendix.
While we have access only to equal-time correlation functions it is a matter of principle that at long distances this system has to behave like a 1D system.These are exceptional in the regard that their IR fixed points are well understood [36].Information on their equal-time two-point correlation functions suffices to reconstruct the nature of these fixed points: the combination of the electron density, spin, pair, and single-electron two-point functions augmented by the bipartite entanglement entropy which we all compute suffices to find out the nature of the ground states.
As we will highlight, the precise properties of these fixed points may yet be different from the canonical truly 1D chain systems.The reason is that at short distance the finite extent of the ladder in the transversal direction does change the physics qualitatively.This short distance ("UV") physics sets the conditions for the long-wavelength properties which are much richer than in 1D systems.At short distance one encounters the vigorous, strongly coupled 2D quantum system.The fact that very large bond dimensions are needed in the DMRG simulations-much larger than in truly 1D systemsindicates that this is a densely many body entangled system where the (implicit) semiclassical language underlying the established 1D canon may well fall short.The open question is whether this will be eventually understood in terms of analytical theory and/or general emergence principle.
The models.Our main focus is on the Hubbard model Eq.(1).The main interest is in the regime of low doping relative to the Mott insulator realized at half filling.Given a total number of sites N = L x × L y there is n i = 1 electron per site when N e = N.The hole doping concentration is defined as δ = N h /N, where N h = N − N e is the number of doped holes.We focus on the doping regime 1/12 δ 1/8.In addition, the interest is in the intermediate-to-large interaction regime U > 8t where, throughout, we take t = 1 as unit of energy.At several instances it will be useful to compare the results for the Hubbard model with those for the t-t -J model, which captures the low energy physics of the Hubbard model in the strong coupling regime, U t.The c † iσ creates projected fermions (no double occupancy is created in the hopping process), while S i describes Heisenberg spins coupled by a superexchange interaction J ∼ t 2 /U .
A remarkable phase diagram.The main result of our simulation is the zero-temperature phase diagram summarized in Figs. 1 and 2. A first surprise is that it is rather insensitive to the interaction strength (Fig. 1, upper panel).In the intermediate coupling regime 8 U 12 it looks very similar.This is further amplified by the fact that it looks similar even for the t-t -J model (Fig. 2) which is associated with the large-U limit.A next surprise is that, at least in the doping range we consider, the phase diagram is rather insensitive to the doping (Fig. 1, lower panel), as we will discuss below, this is especially quite mysterious with regard to the "filled stripe" (FS) phase in the middle.
The profundity is in the extreme sensitivity to the nextnearest-neighbor hopping t , the reason to take it as the horizontal axis in the figures.As function of t we identify three phases with quite different thermodynamical identities, as we will show in the following details.The FS phase found for t 0 can be regarded as the 1D incarnation of an insulator.By varying t by the tiny amount of 0.1t we find instead LE phases that can be viewed as the 1D versions of SC phases.However, depending om the sign of t these are of a very different nature.Further emphasizing that these phases are very different is the fact that our simulations reveal that these phases are separated by phase separation (PS) regions.See the Appendix for details of PS and phase boundaries.
Let us now turn to the identification of these three phases, relying on the behavior of the equal-time correlations.In the sequence FS-LE1-LE2 these 1D phases deviate increasingly from the canonical 1D physics associated with chains.Let us discuss therefore these phases in this order.

III. THE FILLED STRIPE (FS) PHASE
A long time ago it was established that according to Hartree-Fock so-called stripe phases are formed in doped Mott insulators [30,37,38].These consist of antiferromagnetic Mott-insulating domains separated by magnetic domain walls on which the holes localize, forming in turn a periodic array.Mean field theory insists that these are most stable at a filling of one hole per domain wall unit cell: the "filled stripes."It was quite a surprise when quite recently state-of-the-art numerical calculations (including DMRG) showed these to be the ground state of the Hubbard model at a doping δ = 1/8 and t = 0 [17][18][19][20].To give an example, in Fig. 3 we show a result for the Hubbard and t-t -J model; these look quite similar everywhere in the purple regime of the phase diagrams in Figs. 1 and 2. To measure the charge order, we define the local rung density operator as n(x) = 1 L y L y y=1 n(x, y) with expectation value n(x) = n(x) .The charge density profile n(x) on am L x = 64 cylinder is shown in Fig. 3(a), revealing the ordering wave vector Q c = 2πδ consistent with filled stripes.
While it still remains challenging to precisely determine the long-distance behavior of the correlations in the Hubbard model, even though we have kept m = 20 000 states in the DMRG simulation (see the Appendix and [20] for details), it appears fairly easy to address this issue in the closely related t-t -J model, where the FS phase is also present [29].To find out the nature of the SC correlations we study the pair-field correlation, defined as ) is the spin-singlet pair-field creation operator and α = x, ŷ denotes the bond orientations.(x 0 , y) is the reference bond located at x 0 ∼ L x /4 and r is the distance between two bonds in x direction.We find for a characteristic set of parameters J = 1/3 and t = 1/12 that the pair correlator appears to be short ranged, (L x /2) ∼ e −L x /2ξ sc with correlation length ξ sc ∼ 18 (Fig. 3) (see the Appendix for more details).This evidence shows that the FS phase is an insulator.
We have also calculated the spin-spin correlation where (x 0 , y) is the reference site and r is the distance between two sites in the x direction.In first instance one would expect antiferromagnetism with the incommensurate modulation characterized by λ s = 2λ c , confirming the domain wall picture which is decaying algebraically like a Heisenberg spin chain.However, given the even width of the cylinder this eventually opens up a spin gap of the same kind as on a pure four-leg Heisenberg ladder.For the Hubbard model [18][19][20], the spin correlation indeed decays exponentially, F (r) ∼ e −r/ξ s with ξ s ∼ 6.5 (see the Appendix for details).The one of the t-t -J model decays in the same way, F (L x /2) ∼ e −L x /2ξ s , as shown in Fig. 3(c), where we take r = L x /2 for each L x cylinder for simplicity.
Although the literature emphasizes the 2D Hartree-Fock stripes [30,37,38], these "stripes on the ladder" may equally well be understood in terms of natural extensions of the strictly 1D theory.These are actually states that are just deformations of the ubiquitous Luttinger liquids.A striking property of strictly 1D systems is that the periodicities seen in charge and spin correlators are independent of the interaction strength.In the noninteracting Fermi gas one finds these to be algebraic, with ordering wave vectors for spin and charge set by the nesting vectors 2k F and 4k F , respectively, where k F is the Fermi wave vector.The only change as U is increased is in the exponents governing the algebraic decay.The reason is revealed by the U → ∞ limit, where the notion of the "squeezed space" was discovered [39]: consider any configuration of holes and spins and remove the sites where the holes are residing and pretend that the spin system living on this "squeezed" lattice is described by a Heisenberg model.This squeezing operation can be wired in by the so-called string correlators, and it was shown by DMRG that this works always at long distances, all the way to the free limit [40].Very recently this "squeezed space" universality was directly verified in cold atom experiments [41].
The meaning is simple: every charge is a holon, i.e., an electron bound to a kink in the antiferromagnetic spin system.One can view the Luttinger liquid and its descendants (like the LE liquid) as a 1D version of the filled stripes.The simplest way to embed this in 2D is by "putting holons on a row" [42].On a ladder geometry this leads automatically to a 2k F spin and 4k F charge periodicities.The next step is that at dopings like δ = 1/8 and δ = 1/12 the algebraic CDW order is actually at a higher commensurability with the lattice.In the absence of quantum fluctuations this holon lattice will be subjected to commensurate pinning to the lattice and this pinned crystal is just a "higher commensurate" Mott insulator (see [43] for a spectacular example revolving around black holes): the filled stripe phase.
The fate of an algebraic 1D solid in a commensurate background potential is enumerated by the sine-Gordon field theory [36], which reveals that when the pinning potential compared to the kinetic term enumerating the zero point quantum motions drops below a critical value the pinning ceases to exist, with the effect that the electronic density wave depins from the lattice.This is actually a main difference with the Luther-Emery (LE1, LE2) phases.Next to being such "floating charge density waves" these show automatically also superconducting correlations.
The last piece of information comes from the single particle equal-time correlator, G σ (r) = 1 L y y c † (x 0 ,y),σ c (x 0 +r,y),σ . (5) Consistent with the above interpretation, we find in the Hubbard model its correlation length ξ G ∼ 4 to be very short, of the order of the width of the ladder.The single-particle gap G ∼ 1/ξ G should be larger than the charge commensuration gap combined with the spin gap: given the charge-spin separation principle the electron fractionalizes into a spinon and a holon that can only be inserted at energies larger than the spin and charge gaps respectively.As seen in Fig. 3(b) it oscillates in a rather complicated, multiharmonic fashion where at least the periodicity of the CDW can be discerned.Given the very small correlation length this may well reflect the rather complex short distance physics, and we have not attempted to analyze in further detail.
In summary, the mystery associated with the filled stripes is, why is it so that a minute t suffices to destabilize it?As Hartree-Fock shows, commensuration is a muscular source of stability [30].In such a setting the effects of t would be secondary: it could make a difference when t becomes of order 1, but not for a t ∼ 0.1.Apparently the dense entanglement is changing these rules in a way that is beyond our present comprehension.

IV. THE FIRST LUTHER-EMERY PHASE (LE1) VERSUS FLUCTUATING STRIPES
The LE1 phase becomes stable already at a very small t , while it appears to be ubiquitous in a large range of U 's and dopings, only requiring that t is negative.We already identified it in an earlier study as a typical LE phase [20].Luther-Emery describes the universal SC-like state in 1D [36].It arises when attractive interactions are added to a Luttinger liquid.As a signature of the formation of singlet Cooper pairs, a gap opens up both in the single-fermion and spin excitation spectra.The equal-time pair correlation function in Eq. ( 3) should show an algebraic behavior, 1/x K sc .However, given that the IR fixed point in 1D is always strongly interacting, this goes hand in hand with CDW correlations characterized by a 4k F periodicity which also decays algebraically, 1/x K c .This CDW has to be now in the regime where the lattice commensuration is irrelevant.As a highlight, these exponents should be in a dual relation K sc = 1/K c .The LE1 phase exhibits nearly all these diagnostic features.
We observe the CDW modulation in the charge density correlator [see Fig. 4(a)].The charge density correlation decays algebraically, where the exponent K c can be extracted by fitting the Friedel oscillation, which is induced by the open boundaries of the cylinder, of the charge density distribution n(x) [20,21,44].Specifically, we use to fit the local density profile to extract the Luttinger exponent K c .Here, δn is the nonuniversal amplitude, φ is a phase shift, n 0 is the background density, and k F is the Fermi wave vector.An example is given in Fig. 4(a) for L x = 60 cylinder at doping concentration δ = 10% with U = 12 and t = −0.25.Due to the presence of CDW modulations, the SC correlations αβ (r) exhibit similar spatial oscillations with n(x).Following the procedure in Refs.[20,21], we find that the SC correlations in this phase always decay with a power law, whose exponent K sc , shown in Fig. 5 for a range of hole doping concentrations δ = 8.33%-12.5%, is obtained by fitting the results using Eq. ( 3).An example of the SC correlations at hole doping concentration δ = 10% is given in Fig. 4(c).The sure test for LE is whether the product of the exponents K c K sc = 1 is that we extract the exponents K c and K sc at various doping concentrations δ on cylinders of length up to L x = 72 and we find that the relation K c K sc ∼ 1 is satisfied within the numerical uncertainty (see Fig. 5).We also use another method to estimate the exponents K c and K sc from a series of cylinders (see the Appendix for details) and obtain reasonably consistent results, listed on the left side of each panel of Fig. 5.
As a further test, a key feature of the LE liquid is that it has a single gapless charge mode with central charge c = 1, which can be obtained by calculating the von Neumann entanglement entropy S = Tr ρ ln ρ, where ρ is the reduced density matrix of a subsystem with length x.For critical systems in 1 + 1 dimensions described by the conformal field theory (CFT), it has been established [45,46]  This establishes the unique signature of the LE phase insofar as the charge properties are concerned.However, it should also exhibit gaps in the single-fermion and spin response.We find that the equal-time spin-spin correlator in Eq. ( 4) decays exponentially with a finite correlation length ξ s ∼ 8.9 for U = 12, t = −0.25, and δ = 12.5%.In addition, we do observe that the modulation of the spin density is consistent with the Luttinger-liquid/stripe rule that the spin ordering wavelength is twice that of the charge, λ s = 2λ c .
Finally, the singlefermion correlator (5) exhibits also a finite correlation length.We find that G σ (r) ∼ e −r/ξ G with yet again a small correlation length ξ G as for the filled stripes.For U = 12, t = −0.25, and δ = 12.5% we find it to be ξ G ∼ 4.8, as compared to a spin correlation length ξ s ∼ 8.9.The charge is now massless and the single-fermion gap is bounded from below by the spin gap, ξ G /2 ξ s .However, we find that the single-particle correlation length is substantially smaller than 2ξ s .This may be interpreted as the single-particle correlation length reflecting a pairing gap, while the spin gap is substantially smaller reflecting the fact that the spinons may bind into a triplet excitation inside the SC state; see the next section.Anticipating the discussion of the LE2 phase, the most striking difference between the two Luther-Emery phases is in the behavior of the single-particle sector.Its large single-particle gap (i.e., G ∼ 1/ξ G ) leaves no doubt that the LE1 is a strongly coupled system, reminiscent of a local pair superconductor/CDW system.Quite differently from the LE2 phase, we have not managed to identify the oscillations seen in G(r) with any natural length scale, indicating that this is, as in the FS phase, associated with complicated short distance physics.
In fact, only in one regard is this LE1 different from the canonical 1D textbook version.The periodicity of the CDW in the latter case is again determined by continuation from the free fermion limit: it should correspond to 4k F , and we argue that this should count as filled stripes.However, the periodicitiy of the LE1 phase counts as half-filled stripes: the ordering wave vector is Q c = 4πδ of wavelength λ c = 1/2δ according to Fig. 4, indicating only half a hole per each CDW unit cell.
On the cylinder "2D-like" transversal space for motion is possible at short distances which is apparently quite sensitive to t .This short distance ladder physics should be viewed as a strongly coupled affair where many-body entanglement plays a crucial role given the fact that the LE behavior at long distances requires a very large bond dimension.These observations suggest that at short distances half-filled spin stripes are formed.We confirm the early finding that the short-range SC correlations are much stronger along the rungs than in the direction of the legs [20,21,29].As was pointed out very early [10,31], it may well be that some kind of resonating valence bond like state is realized on the half-filled stripes where spin singlets are exchanging with pairs, offering a potential explanation for the preferred on-stripe charge density of 1/2 charge per domain wall unit cell.At long distances this "fluctuating order" may then renormalize into the universal 1D LE physics, leaving behind as a fingerprint that the doubled charge-and spin periodicities are disconnected from the generic 2k F , 4k F wave vectors of truly 1D physics.Notice that in Ref. [29] it was pointed out that this is still consistent with the so-called generalized Luttinger theorem.

V. THE SECOND LUTHER-EMERY PHASE (LE2) VERSUS d-WAVE SUPERCONDUCTIVITY
For positive t we find yet another phase to be very robust.Like the LE1 phase it extends in a large regime of doping and U 's and persists up to large t as well.As we will now show, the charge properties show precisely the diagnostics of the LE phase.These are actually of the same kind as found in the LE1 phase: power-law SC pair-field correlations concomitant with the dual half-filled stripe algebraic CDW order; see Fig. 6.We already highlighted the similarity of the phase diagrams of the Hubbard and t-t -J models [29].The charge physics of the LE2 phases behave in both models in a quite similar fashion.We compare these in Fig. 6 for representative parameters: J = 1/3 and t = 0.18 for the t-t -J model, while U = 12, t = 0.25 for the Hubbard model.It clearly seen that both models have similar charge-and spindensity-wave oscillations along with exponentially decaying spin correlations.Moreover, both the SC pair-field (r) and charge correlations decay algebraically with corresponding exponents, e.g., K sc ≈ 0.92, K c ≈ 1.14, and K c K sc ≈ 1.05 ∼ 1 for the t-t -J model.The extracted central charge c ∼ 1 as shown in the Appendix, which further confirms the presence of the LE2 phase.A more detailed comparison of exponents and other physical quantities of the phases is presented in Table I.
This confirms that the LE2 phase is indeed of the LE kind.However, both the single-electron and spin correlations appear at first sight to violate this assignment.As we argued, the single-fermion correlations should be short ranged as is the case in the LE1 phase, revealing that Cooper pairs are formed.Taking representative parameters U = 12, t = 0.25, and δ = 12.5%, we find that within our resolution the single-particle gap ( G ∼ 1/ξ G ) is vanishing: the correlation length ξ G > 33 for a cylinder of length L x = 96, as shown in Fig. 6 (see the Appendix for more details).This is actually dependent on models.We find that in the t-t -J model the single-particle gap ( G ∼ 1/ξ G ) is finite but rather small, ξ G ∼ 10, compared to the LE1 case.Turning to the spin correlations, it becomes even more of a puzzle.A dramatic difference with the LE1 phase is that the spin periodicity observed through the equal-time spin correlator at short distances is no longer set by the generic λ s = 2λ c , but instead it reveals a simple two-sublattice staggered antiferromagnetism, with a wavelength λ s = 2, as shown in Fig. 6.This is independent of the details of the charge stripe periodicity, while we observe it throughout the LE2 phase independently of the parameters.
The next difficulty is that generically the spin correlation length is shorter than the single-electron correlation length.In the Hubbard case we find for instance that ξ s ∼ 5.9 for parameters U = 12, t = 0.25, and δ = 12.5%, where ξ G > 33 within the precision of DMRG.This same pattern repeats in the t-t -J, model where we find ξ s ∼ 3.8 while ξ G ∼ 10 for parameters J = 1/3 and t = 0.18.This appears to be at first sight completely unreasonable: the single particle sets a lower bound for the energy it costs to insert a spinon.The triplet gap should be at the least twice the spinon energy: it can be less because spinons may bind together in a triplet excitation, but it cannot be larger.For this simple reason it is expected to find a single-electron correlation length that is (much) larger than the spin correlation length, as mentioned above.
How do we understand these results?There is actually a natural explanation: in strictly 1D all singlet superconductors are effectively s wave since Cooper pairs cannot have angular momentum.This s-wave nature is therefore automatically hard-wired into the LE state.On the ladders there is, however, room for angular momentum, and, as we already alluded, the short-range SC correlations are clearly of the d-wave type.The sign of the order parameter flips upon going from the x to the y direction.How can this influence the long-wavelength physics?
Dealing with a relatively weakly interacting system, one can follow the canonical procedure for chain systems.One departs from the kinetic term of the free system, which is bosonized and then combined with the interaction terms that take a simple form in the bosonized representation, to solve subsequently in the interacting bosonic field theory.The free theory then leaves its imprint on the interacting system in the form of quantities like the 2k F and 4k F ordering wave vectors, as we already stressed above.Although it is not the standard procedure, one could view the LE state in a similar way as a descendent of the BCS superconductor taken as the free limit.It will fall short of explaining why the CDW and the superconductor are in a dual relation, but it does explain why LE is characterized by spin and single-particle gaps.
Let us consider how this works, departing from a d-wave BCS state on the ladder.The amplitude of the pair density appears to be quite different along the x and y directions according to the DMRG simulations, given the fact that the fourfold rotational symmetry of the square lattice is of course broken on the cylinder [29].Let us for the sake of the argument assume that an isotropic d-wave superconductor is formed.The fermion spectrum will then be characterized by nodal lines along the diagonals of the zone.Near half-filling one then expects a nodal point on the Fermi surface in the close proximity of (q x , q y ) = (±π/2, ±π/2) on the square lattice, as illustrated in Fig. 7(a).
On the four-leg cylinder there are four finite-size quantized momenta available in the y direction q y = 0, ±π/2, π: two of the four 1D modes in the x direction characterized by q y = ±π/2 will become massless at q x = ±π/2 because they intersect the nodal points [see Fig. 7(b)].The system can of course still be a superconductor since it is gapped at the other allowed momenta.The equal-time fermion propagator will now show an infinite correlation length at large distances since it is completely dominated by the presence of the massless points; this should be remembered from the bosonized interacting theory.This vanishing mass is not robust.Given that there is no C 4 symmetry axis on the cylinder, there has to be generically an s-wave admixture giving rise to a finite mass at the nodal point, while also this minimum gap may shift away from the nodal line.But this gap will be generically quite a bit smaller than the average d-wave gap.This explains the observation of the single-particle correlation length that may seem to diverge or stay quite finite pending the details of the model.This interpretation is directly confirmed by the oscillation observed in G σ (x).The massles "nodal point" should be in the close vicinity of q x = ±π/2, implying a wavelength λ G 4. This is precisely what is observed, and the DMRG result for G σ (r) is very similar to the result for the noninteracting d-wave model, as shown in Fig. 7(d).
How do we explain the spin correlations?It is an equally well established wisdom that the collective response of a not too strongly interacting 1D system is quite well approximated by the results of the RPA.One computes first the dynamical susceptibility of the free system χ 0 ( q, ω) (the Lindhardt function) and the effects of the interactions are included at the time-dependent mean field level χ ( q, ω) = χ 0 ( q, ω)/[1 − J q χ 0 ( q, ω)] with the interaction parameter J q depending on the nature of the response and interactions.The way this works in the d-wave SC state on the square lattice is well known since it has turned into the standard explanation for the so-called magnetic resonance observed in cuprate superconductors [34,35].A gap opens up in the Lindhardt spectrum associated with the SC gap in the single-particle momentum.This is at maximum at q = (π, π ), the momentum associated with a simple staggered order parameter.In the proximity of momenta like (π, 0) it will close, given the fact that node-tonode scattering is kinematically allowed.
A next ingredient is that upon switching on the enhancement factor J q a bound mode appears in the Lindhardt gap that is generically most strongly bound at the (π, π ) point for the reason that the density of electron-hole pair states is highest at these momenta.Similarly, the χ is strongly suppressed at the massless points for simple kinematical reasons.We computed this RPA susceptibility for the cylinder system using the electron dispersion relations we just discussed, finding that with regard to these features it is very similar to the square lattice.Practically, the simple staggered magnetic order becomes dominant once J (π,π ) exceeds 2.78t (see the Appendix for more details).
What does this mean for the equal-time spin-spin correlations?We computed the spin-spin correlation functions from the RPA dynamical susceptibility, finding for a reasonable choice of the interaction parameter J q an outcome that is virtually identical to the DMRG result in the LE2 phase; see Fig. 7(b).The way this works is simple.One anticipates that the massless spin excitations should dominate the spin correlator at large distances.However, for the reason that the associated density of states is strongly suppressed, this contribution is just not resolved.Instead, given that the dynamical susceptibility is dominated by the "resonance" at the staggered ordering wave vector, it dominates also the equaltime correlations.
The important message is that the DMRG results become comprehensible by asserting that at short distances the system renormalizes into a seemingly rather weakly interacting dwave superconductor living on the ladder.The main effect of the residual interactions appears to be in the spin channel, where enhancement factors ensure that a well developed "magnetic resonance" forms.However, this is not yet the whole story: this "nearly noninteracting" d-wave Luther-Emery view is in one regard quite misleading.The LE2 phase is still characterized by the "half-filled stripe" algebraic CDW order.As for the spin system one could look for a nesting type fermiology interpretation, but this turns out to fail completely [29].There is just nothing in the noninteracting d-wave band structure that relates to the ordering momentum of the charge density.In this regard the LE2 phase continues to be rooted in strong coupling physics.

VI. CONCLUSIONS
Resting on the fact that it has become possible recently to compute reliably the ground states of four-leg ladder Hubbard systems by DMRG, we have systematically investigated the ground state properties in the low doping regime.Given that these ladder systems renormalize into one-dimensional systems at large distances, one can use the well understood universal properties of 1D physics to diagnose the physics.At short distances the physics is like that of the strongly interacting 2D system that is apparently ruled by dense many body entanglement, as signaled by the need for very large bond dimensions in the DMRG.The outcome is in the form of LE liquids and the "commensurately pinned Luttinger liquid" (the filled stripes) that obey the 1D universality, although these rest on short distance data which are entirely different from the canonical chain systems.
Surprisingly, we find that the zero-temperature phase diagram is rather insensitive to the strength of the interactions when these are large enough, and to the doping when they are not too high.Instead, the ground states turn out to be exceedingly sensitive to the next-nearest-neighbour hopping t .We do not find a plethora of phases as is sometimes claimed: as function of t there are just three preferred phases.Why this is the case is not at all understood: standard product state mean field (Hartree-Fock) does not give any clue in this regard, and it has to be somehow rooted in the densely entangled nature of the short distance physics.
One part of the puzzle is why the simple commensurate pinning stabilizing the filled stripes is so extremely sensitive to t .When this occurs, the liquid phases that are formed are of a very different nature depending on the sign of t .The "LE1 liquid" for negative t shows the symptoms of a strongly interacting yet, in a way, conventional LE state.The behavior of correlations in the LE liquid suggests a highly collective fluctuating order physics: well developed half-filled spin stripes which are formed from Cooper pairs (the Luther-Emery charge density-pair correlation function duality), kept fluctuating by 1D zero-point motions.
Strikingly, the "LE2 liquid" setting at small positive t appears to be instead a continuation of a weakly interacting d-wave superconductor living on the ladder.The symptoms are in the first place the observation of a small (or even vanishing) gap in the single-particle correlations, G ∼ 1/ξ G , suggesting that the system remembers the nodal fermions of the 2D BCS d-wave superconductor.The spin correlations are even more informative: the two-sublattice spin correlations are seemingly uniquely explained by the generic fermiology motive of strong enhancement of staggered spin correlations in the d-wave SC state.This brings us to experiments.One may wonder whether the ladder is no more than a metaphor or whether there may be a more literal relation to the observations.In the mid 1990s the half-filled spin stripes were discovered in the 214 superconductors.It is well understood that the incommensurate static spin order in these systems follows precisely the stripe rule λ s = 2λ c , while the spin order parameter is quite large by any standard.Moreover, the doping dependence of the periodicity at least for δ 12.5% is linear, indicating a preferred density of half a hole per domain wall unit cell.This was also the origin of the idea of fluctuating order: the spin excitations measured by inelastic neutron scattering are suggestive that also in the absence of static order there are strong stripe correlations going on at higher energy.
However, in the 123 and 2212 families, especially near optimal doping, one may be closer to the truth using the fermiology language.The spin excitations are dominated by the magnetic resonance occurring at the staggered wave vector, which near optimal doping is only present in the SC state.This is naturally explained by the same RPA logic as we used to explain the magnetic correlations in LE2.
There has been much debate in the past trying to reconcile these seemingly very different forms of physics of these two families of cuprates that seem chemically so similar: all the action is supposed to take place in the same Cu-O layers.However, we now learn from the Hubbard ladders that tiny differences in a microscopic parameter-that, based on established wisdom, should not matter-can make a big difference, apparently causing the system to land in near opposite regimes of long distance physics.
Note added.Recently, we became aware of a work by Ponsioen et al. [47], who studied the same model on the 2D square lattice using infinite projected entangled-pair states (iPEPS) and proposed that the ground state of system at doping δ 0.14 is in a state with period-4 charge stripe but without SC for negative t , while it is in a state with coexisting uniform SC with long-range antiferromagnetic order but without charge stripe for positive t .This is different from this work in both negative and positive t sides, where we find coexisting quasi-long-range SC and CDW correlations but without (quasi-)long-range spin order.In contrast, full convergence in the "filled" stripe phase is known to be more challenging, partially due to the enlarged charge density oscillation periodicity and the enhanced quantum fluctuation near the phase boundaries as the filled stripe phase is relatively small in the phase diagram.To resolve this issue, we further implement higher SU (2) symmetry in the DMRG simulation and push the effective U (1) number of states up to m ∼ 36 000, which gives excellent convergence in our study on systems such as the L x = 96 cylinder.

Phase diagram of the Hubbard model at U = 8
The ground state phase diagram of the Hubbard model at U = 8 is shown in Fig. 8 for three different hole doping concentrations δ = 8.33%-12.5% on four-leg cylinders of length L x = 48, 40, and 32, respectively, where the phase boundaries are determined by the distinct pattern of charge density modulation of the adjacent phases.We keep m = 5000 U (1) states when scanning the phase diagram, which is large enough for determining the ordering vector Q of the stripy CDW and PS features.This is quite similar to the phase diagram in Fig. 1 of U = 12 shown in the main text, but with a slightly larger filled stripe phase.

The filled stripe phase
As mentioned in the main text, it is more challenging to obtain a fully converged ground state of the Hubbard model in the filled stripe phase than the LE phases and especially to directly establish the long-distance behavior of correlation functions, such as the superconducting correlation.To resolve this issue and directly determine the nature of the filled stripe phase, we focus on this phase with a characteristic set of parameters t = 0.05, U = 12, and δ = 12.5%.By keeping around m = 18 000 number of states, we are able to converge to the true ground state, which preserves all the symmetries of the Hamiltonian including the spin SU (2) rotational symmetry with S z i = 0. Examples are shown in Fig. 9 for the system on an L x = 64 cylinder.Consistent with the filled stripe phase [19,20], the charge density profile shows a modulation of wavelength λ c = 1/δ, and the spin-spin correlation decays exponentially, F (r) ∼ e r/ξ s with a finite correlation length ξ s ∼ 6.5.Both the single-particle and superconducting pairfield correlations are also short ranged: G(r) ∼ e −r/ξ G and (r) ∼ e −r/ξ sc with finite correlation lengths ξ G ∼ 4.0 and ξ sc ∼ 4.6, respectively.Again, they are consistent with the filled stripe phase in the t-t -J model discussed in the main text.
However, this filled stripe phase is only stable in a small region around or close to t = 0, which becomes increasingly fragile with the increase of U or the decease of the hole doping concentration.For example, the ground state of the system at t = 0 is in the filled stripe phase at U = 8, but this is no longer the case when U > 10.Even at U = 8, a tiny t ∼ −0.01 is large enough to drive the system out of the filled stripe phase.Examples of charge density profile n(x) of the Hubbard model on an L x = 32 cylinder are given in Fig. 10 with δ = 12.5% and t = 0 at different U .For each simulation, we start with a state with filled charge stripe pinned by an appropriate chemical potential, and then check the stability of the filled stripe in the following DMRG sweeps by turning off the pinning field.Our results show that the filled stripe disappears when U > U c ∼ 10.

Luttinger exponents
Here we present another method we used to extract Luttinger exponents K c and K sc and discuss the finite-size effect in the two methods.For a given L x cylinder, the CDW amplitude A cdw (L x ) deep inside the bulk can be extracted from n(x) = A cdw cos(Qx + θ ) + n 0 , where Q is the CDW ordering vector.Similarly, we can also select YY (L x /2) to represent the decay of the pair-pair correlation function.For a system with quasilong-range CDW and SC order, A cdw (L x ) and YY (L x /2) on a series of cylinders with different L x should follow The results of the finite-size scaling of A cdw (L x ) and YY (L x /2) in the U = 12 and t = −0.25 models with various doping concentrations are shown in Fig. 11.On the doublelogarithmic scale, all the results are nearly linear, consistent with the behavior discussed in the main text.However, as shown in Fig. 5 of the main text, we find some noticeable discrepancies between the exponents obtained by this method and those from the measurement on single long cylinder with given L x .In the limit of L x going to infinity, the K c and K sc obtained from these two methods should be equal.For finite systems, the A cdw of short cylinders might slightly deviate from the linear fit, and the Friedel oscillation n(x) clearly saturates in the bulk for the finite cylinders.From the small discrepancy between the exponents obtained by the two methods, we can estimate the error introduced by the finite-size effect.
As an additional support for the robustness of the LE physics, we investigate the t dependence of K c and K sc on the relatively short cylinder with L x = 32, U = 12, and δ = 12.5%, the same as the one used in plotting the phase diagrams.As illustrated in Fig. 12, the K c (K sc ) obtained keeping m = 5000 states is smaller (larger) than the extrapolated exponents.However, the trend of K c and K sc as t increases is still meaningful.In the LE1 phase, K c and K sc are nearly independent of t even close to the phase boundary.The properties of LE2 are also very stable: the increasing t term only slightly suppresses the CDW and enhances the SC ordering.

Single-particle correlation
One puzzling point is the slowly decaying single-particle correlation of the Hubbard model in the LE2 phase.Figure 13 U = 10 and δ = 10% on an L x = 60 cylinder.For comparison, the single-particle correlation G(r) in the LE1 phase at t = −0.25 is also given.Figure 13(b) shows the correlation length ξ G extracted by fitting the results using function G(r) ∼ e −r/ξ G for the LE2 phase.While a clear saturation tendency of ξ G vs L x is seen in the LE1 phase, it remains still unclear for the LE2 phase for cylinders of length up to L x = 64.

Entanglement entropy
As a support for our results in the main text, we have also calculated the von Neumann entanglement entropy of the Hubbard and t-t -J models in the LE2 phase, as shown in Fig. 14.We choose a characteristic set of parameters J = 1/3, t = 0.18, and δ = 12.5% for the t-t -J model on the L x = 96 cylinder.For the Hubbard model, to obtain fully converged entropy we choose a series of shorter cylinders with L x = 32, 48, 64 and parameter set t = 0.25, U = 12, and δ = 12.5%.The central charge can be extracted using [46]   contribution from a higher order oscillating term.For the Hubbard model, the extracted central charge c decreases from 2.3(2) to 1.8(1) as L x increases from 32 to 64.The estimated c(L x → ∞) = 1.2(1) using a second-order polynomial is shown in the inset of Fig. 14(a), which is roughly consistent with c ∼ 1. Consistently, the central charge c extracted from t-t -J model in Fig. 14(b) is c ∼ 1.05, which also agrees with c = 1.Therefore, we anticipate that the LE2 phases in both the Hubbard and t-J models are consistent with the LE liquid phases.

Phase separation
In the phase diagram, there are two additional regions for phase separation sandwiched by the filled stripe and LE phases, for which the evidence is provided in Fig. 15 for the Hubbard model.While both the half-filled charge stripe of wavelength λ c = 1/2δ and the filled charge stripe of wavelength λ c = 1/δ are clearly seen in the LE and filled stripe phases, respectively, the charge density profile n(x) in the phase separation region is clearly spatially inhomogeneous, which can be considered as a combination of both filled and half-filled charge stripes.This is true for both shaded regions in the phase diagram.To analyze the position of the magnetic resonation peak and the spin-spin correlation, it is instructive to start with the bare Lindhard susceptibility of the superconducting state: where f is the Fermi dispersion and we define k, q = ( k+ q k + k+ q k )/(E k+ q E k ).The three parts in χ 0 ( q, ω) are due to quasiparticle scattering, quasiparticle pair creation, and quasiparticle pair annihilation, respectively.To calculate the bare susceptibility χ 0 numerically, we replace 0 + by = 0.01t in the denominators.The calculated bare susceptibility χ 0 ( q, ω) with maximum response at q = (π, π ) is illustrated in Fig. 16(a).In the RPA approach, the interacting spin susceptibility is given by χ ( q, ω) = χ 0 ( q, ω) 1 − J ( q)χ 0 ( q, ω) , (A5) where J ( q) is the spin-spin response function assumed to be of the form J ( q) = −J[cos(q x ) + cos(q y )]/2.After J ( q) is turned on, a much sharper resonance peak at (π, π ) is expected at low energy, smaller than twice the maximum of the SC gap.In Fig. 16(b), we show the strong magnetic response χ ( q, ω) calculated at the J = 2.8t point.Finally, we calculate the equal-time spin-spin correlation F (r, t = 0) from q e iq x r χ ( q, ω)dω.At the noninteracting J ( q) = 0 point, the spin-spin correlation is generally incommensurate.When J ( q) is large enough, the period-2 oscillation originating from the sharp resonant peak at (π, π ) becomes dominate.In our calculation, the onset of commensurate period-2 oscillation happens around J ∼ 2.78t, which indicates a spin gap s ∼ 0.09t shown in Fig. 16(c).

FIG. 1 .
FIG. 1. Ground state phase diagram of the Hubbard model in Eq. (1) in (a) as a function of U and t at hole doping concentration δ = 12.5% where the dashed line labels t = 0, and in (b) as a function of δ and t at U = 12.Here t = 1 and the black dots are the data points.The phase boundaries are determined by the charge density profile from simulations by keeping m = 5000 states.See the Appendix for details.

FIG. 3 .
FIG. 3. (a) Charge density profile n(x) of the Hubbard and t-t -J models in the filled stripe phase.(b) Single-particle correlation G σ (r) of the Hubbard model at U = 12 and t = 0.05, whose correlation length ξ G ∼ 4 is fitted from the red squares.(c) Superconducting pair-field (L x /2) and spin-spin F (L x /2) correlations of the t-t -J model.Here J = 1/3 and t = 1/12.The hole doping concentration is δ = 12.5% for both models.

FIG. 4 .
FIG. 4. Ground state properties of the Hubbard model in the LE1 phase at U = 12, t = −0.25, and δ = 10%, measured on a L x = 60 cylinder.(a) Charge density profile n(x) fitted by the Friedel oscillation (solid line) using Eq.(6).(b) Spin-spin correlation F (r)(−1) r .The inset is the plot in semilogarithmic scale with exponential fitting |F (r)| ∼ e −r/ξs (red solid line).(c) SC pair-field yy (r) and singleparticle |G(r)| (inset) correlations in semilogarithmic scale.Note that only the red data points are used in the fitting.(d) Entanglement entropy S(L x /2) and the extracted central charge c (inset) at different doping concentrations.

FIG. 5 .
FIG. 5. Exponents K c and K sc in the LE1 phase of the Hubbard model at various doping concentrations δ as a function of cylinder length L x .Here t = 1, U = 12, and t = −0.25.The exponents on the left side of each panel are obtained from another method relying on A cdw and YY (L x /2), whose details are explained in the Appendix.The error bars represent the error from the fit of n(x) and YY .

FIG. 7 .
FIG. 7. (a) The Fermi surface of t = 1 and t = 0.3 model (details are in the Appendix).Dashed line: four quantized momenta q y = 0, ±π/2, π in the y direction.(b) The quasiparticle dispersion of the q y = π/2 mode.The inset is the zoomed-in illustration of the nearly vanishing gap of the dispersion.(c) Upper: The spin-spin correlation F (r) of the Hubbard model on an L x = 96 cylinder with δ = 12.5%, the same as Fig. 6(a2).Lower: The spin-spin correlation obtained from RPA calculation for the four-leg model shown in (a).(d) Upper: The single-particle correlation G(r) of the Hubbard model, the same as Fig. 6(a4).Lower: The single-particle correlation G(r) obtained from the noninteracting model shown in (a).

FIG. 8 .
FIG. 8. Ground state phase diagram of the Hubbard model at U = 8 as a function of hole doping concentration δ and t .The black dots are data points.

FIG. 9 .
FIG. 9. (a) Charge density profile n(x), (b) spin-spin correlation F (r), and (c) superconducting pair-field yy (r) (black) and single particle G(r) (red) correlations of the Hubbard model in the filled stripe phase on an L x = 64 cylinder at U = 12, t = 0.05, and δ = 12.5% by keeping m = 18 000 number of states.(d) Superconducting pair-field yy (r) and single-particle G(r) correlations of the t-t -J model at J = 1/3, t = 0.07, and δ = 12.5%.Dashed lines are guides for the eyes.

FIG. 10 .
FIG. 10.Charge density profiles n(x) of the Hubbard model at t = 0 and δ = 12.5% for U = 8-12.The filled stripe is no longer stable when U > 10.

FIG. 11 .
FIG. 11.(a) Finite-size scaling of A cdw as a function of L x in double-logarithmic scale for three doping concentrations.(b) Finitesize scaling of YY (L x /2) as a function of L x for the same systems.
FIG. 12.The Luttinger exponents K c and K sc as a function of t in two LE phases.All the exponents are calculated on L x = 32, U = 12, and δ = 12.5% systems keeping m = 5000 U (1) states.The dashed region indicates the PS and filled stripe phase.

FIG. 14 .
FIG. 14.(a) Von Neumann entanglement entropy of the Hubbard model on an L x = 64 cylinder, with t = 0.25, U = 12, and δ = 12.5%.The central charge extracted from the red data points is c ∼ 1.8.Inset: The central charge extracted from a series of cylinders with L x = 32-64.(b) Entropy of the t-t -J model at J = 1/3, t = 0.18, and δ = 12.5% in the LE2 phase.The extracted central charge is c ∼ 1.07.
To make direct connection with DMRG results, We employ the tight-banding dispersion( k) = −μ − 2t[cos(k x ) + cos(k y )] − 4t cos(k x ) cos(k y );(A3) here we use the parameter set inside the LE2 phase, t = 1 and t = 0.3.The chemical potential is μ = 1.0.We use a gap function with d-wave symmetry, ( k) = 0 [cos(k x ) − cos(k y )]/2 with 0 = 0.1.The lattice we used is L x × L y = 400 × 4; for simplicity, we consider periodic boundary conditions on both x and y directions.The quasiparticle dispersion is given byE k = √ 2 ( k) + 2 ( k).As discussed in the main text, k y can only take four discrete values, and the nodal point on the Fermi surface is in the close proximity of (π/2, π/2).The energy dispersion along the k y = π/2 line exhibits an almost invisible tiny gap ∼0.03t located at k x slightly larger than π/2.

FIG. 15 .
FIG. 15.Charge density profiles n(x) of the Hubbard model on an L x = 32 cylinder at U = 12 and δ = 12.5% at different t .

TABLE I .
Summary of the phases: "LE" denotes the Luther-Emery liquid, "FS" denotes filled-stripe phase.Corresponding Luttinger exponents K c , K sc , correlation lengths ξ s , ξ G , and central charge c for both Hubbard and t-J models are provided, whose values and errors are from the fitted lines in Figs.3, 5, and 6.