A continuous Mott transition between a metal and a quantum spin liquid
Ryan V. Mishmash, Ivan Gonzalez, Roger G. Melko, Olexei I. Motrunich, Matthew P. A. Fisher
AA continuous Mott transition between a metal and a quantum spin liquid
Ryan V. Mishmash,
1, 2
Iv´an Gonz´alez, Roger G. Melko,
4, 5
Olexei I. Motrunich, and Matthew P. A. Fisher Department of Physics, University of California, Santa Barbara, California 93106, USA Walter Burke Institute for Theoretical Physics, California Institute of Technology, Pasadena, California 91125, USA Centro de Supercomputaci´on de Galicia, Avda. de Vigo s/n, E-15705 Santiago de Compostela, Spain Department of Physics and Astronomy, University of Waterloo, Ontario, N2L 3G1, Canada Perimeter Institute for Theoretical Physics, Waterloo, Ontario, N2L 2Y5, Canada Department of Physics, California Institute of Technology, Pasadena, California 91125, USA (Dated: September 24, 2018)More than half a century after first being proposed by Sir Nevill Mott, the deceptively simplequestion of whether the interaction-driven electronic metal-insulator transition may be continuousremains enigmatic. Recent experiments on two-dimensional materials suggest that when the insu-lator is a quantum spin liquid, lack of magnetic long-range order on the insulating side may causethe transition to be continuous, or only very weakly first order. Motivated by this, we study ahalf-filled extended Hubbard model on a triangular lattice strip geometry. We argue, through useof large-scale numerical simulations and analytical bosonization, that this model harbors a contin-uous (Kosterlitz-Thouless-like) quantum phase transition between a metal and a gapless spin liquidcharacterized by a spinon Fermi surface, i.e., a “spinon metal”. These results may provide a rareinsight into the development of Mott criticality in strongly interacting two-dimensional materialsand represent one of the first numerical demonstrations of a Mott insulating quantum spin liquidphase in a genuinely electronic microscopic model.
I. INTRODUCTION
Strongly correlated electronic systems may have in-sulating phases that originate entirely from electron-electron interactions. These insulators, and their phasetransitions to metallic phases have a long history reach-ing back into the pioneering work of Mott [1, 2]. How-ever, despite decades of study, metal-insulator transitionsdriven by strong correlations—Mott’s namesake—remainrather poorly understood. Central to this difficulty is thefact that Mott transitions exhibit strong quantum fluc-tuations, which can inherit correlations from both theadjacent metallic and insulating phases. Thus, the na-ture of the Mott transition may depend crucially on theproperties of each of these phases.Conventional insulating phases, such as those withmagnetic long-range order, appear to predominantly giverise to first-order Mott transitions, as has been observedin a number of experimental systems in the past [3–7].The reason for first-order behavior is simple: The prop-erties of both the spin and charge sectors change qualita-tively at the transition, the former developing magneticlong-range order and the latter localizing to form an in-sulating state. In contrast, systems that harbor uncon-ventional, exotic insulating phases showing no symme-try breaking down to zero temperature—so-called quan-tum spin liquids [8–11]—offer a promising playgroundfor finding the long-sought-after continuous Mott tran-sition. For example, one beautiful possiblity is that thespin sector on the insulating side may be described by aspinon Fermi surface coupled to a U(1) gauge field [12](the so-called “spin Bose metal” [13], hereafter referredto as simply the “spinon metal”). In this case, the be-havior of the spin correlations would be qualitatively un-changed [14] upon crossing the transition, making the na- ture of the transition determined entirely by the chargesector. Thus, as proposed in Refs. [15–17], perhaps theelectronic Mott transition in d spatial dimensions can bein the ( d + 1)D XY universality class, the same as ob-tained for bosons [18]!Fortunately, this sort of physics is more than justa theorist’s dream, as recently several experimentalgroups have found strong evidence for spin-liquid be-havior proximate to a Mott transition in two separatequasi-two-dimensional triangular lattice organic materi-als. In 2003, a putative spin-liquid phase in κ -(BEDT-TTF) Cu (CN) was discovered [19], which is insulatingat ambient pressure with no apparent long-range orderbut can indeed be driven metallic by application of mod- t t’ V V U V V x = 1 x = 2 U/t
Mott transition:
U/t = 1 . C2S2 metal C1S2 spin liquid C0S0 VBS . FIG. 1.
Schematic of the half-filled extended Hubbardmodel on the two-leg triangular strip and its phasediagram.
Top: Our electronic model contains electron hop-pings t and t (cid:48) in addition to repulsive Hubbard interactions upto fourth neighbor [see Eqs. (1)-(2)]. As shown, we view thetwo-leg triangular strip as a 1D chain and attack the problemwith DMRG and bosonization. Bottom: The phase diagramof our model as a function U/t for the chosen characteristicparameters (see text). a r X i v : . [ c ond - m a t . s t r- e l ] J un erate pressure [20]. More recently, Itou et al. [21] founda spin-liquid candidate in EtMe Sb[Pd(dmit) ] . Fur-ther experiments indicated the existence of highly mo-bile gapless spin excitations in both compounds [22, 23],although the precise nature of the spin excitations in κ -(BEDT-TTF) Cu (CN) at the lowest temperatures isstill highly controversial [24]. These findings suggest thatthe spinon metal is likely a good starting point for un-derstanding the spin-liquid behavior observed in thesetwo materials [12, 14]. In addition, the pressure-inducedMott transition from the metal to the spin liquid is ob-served to be either only very weakly first order [20], orperhaps even continuous [25–27].Motivated by these experiments, we consider a modelof interacting electrons on a half-filled triangular lattice“strip” geometry (see Fig. 1), which we solve using large-scale density matrix renormalization group (DMRG) cal-culations. By increasing the strength of the repulsiveelectron-electron interactions, we drive the ground stateof the system from a metallic Fermi liquid-like phase toan insulating phase identified as the electronic spinonmetal [12, 28] via an intervening continuous Kosterlitz-Thouless-like quantum phase transition. Our realizationof this spin liquid phase constitutes perhaps the firstnumerical demonstration of a Mott insulating quantumspin liquid in an interacting microscopic model involv-ing itinerant electrons that is beyond the strictly one-dimensional (one-band) limit [29]. Furthermore, we areable to characterize this exotic phase in a very thor-ough fashion. Further increasing the electron interactionseventually drives the system into a spin-gapped valencebond solid (VBS) insulator—the phase realized by theeffective Heisenberg spin model that our half-filled elec-tronic model approaches at strong repulsion. Our cal-culations thus represent a direct quasi-one-dimensional(quasi-1D) analog of tuning a two-dimensional (2D) half-filled Hubbard-type model from a metal to a quantumspin liquid to a conventional ordered phase via increas-ing overall electron repulsion [30–33], a result with clearpotential relevance to the Mott physics observed in theorganic spin liquid materials [20, 27, 34]. II. EXTENDED HUBBARD MODEL ON THETWO-LEG TRIANGULAR STRIP
The most appropriate microscopic model for thetriangular-lattice organic materials is a Hamiltonian con-sisting of electron hopping plus moderately strong, pos-sibly extended [35, 36], Coulomb repulsion. As is well-known from some 30 years of research on the high-temperature cuprate superconductors [10], such a modeldoes not succumb easily to either exact analytical fieldtheory nor direct numerical simulations in two dimen-sions due to the fermionic “sign problem”.Recently, some of us have proposed a novel approachto the 2D limit of such models through a sequence ofstudies on quasi-1D ladder geometries, which have the significant advantage that they can be solved exactlywith DMRG [37–39]. Sheng et al. used this line of at-tack to extensively study an effective spin model appro-priate for the “weak” Mott insulating regime of the or-ganic materials [14, 31] and indeed found exceptionallystrong evidence that quasi-1D descendants of the spinonmetal exist as the ground state over a large region of thephase diagram [13, 40]. The low-energy degrees of free-dom of this exotic spin liquid are modeled as mobile andcharge-neutral spin-1/2 fermionic spinons coupled to aU(1) gauge field. In 2D, these gapless spinons give riseto a spin structure factor with power-law singularities re-siding on an entire “Bose surface” in momentum space.However, in quasi-1D the Bose surface is reduced to a setof points, so that quasi-1D descendants of the 2D spin liq-uid are dramatically recognizable on ladders, making thequasi-1D approach very fruitful [13, 40–42].Inspired by these recent developments and restrictingourselves to the two-leg triangular strip for numericaltractability, we consider the following extended Hubbardmodel (see Fig. 1): H = − (cid:88) x,α (cid:2) t c † α ( x ) c α ( x + 1) + t (cid:48) c † α ( x ) c α ( x + 2) + H . c . (cid:3) + 12 (cid:88) x,x (cid:48) V ( x − x (cid:48) ) n ( x ) n ( x (cid:48) ) , (1)where c α ( x ) destroys an electron at site x with spin α = ↑ , ↓ , n ( x ) ≡ (cid:80) α c † α ( x ) c α ( x ) is the electron numberoperator, and we take the system to be half-filled withone electron per site.In the usual on-site Hubbard model, we would have V ( x − x (cid:48) ) = U δ x,x (cid:48) . However, inspired by the resultsof Ref. [28], we allow for longer-ranged repulsion in our (cid:239) (cid:239) (cid:239) (cid:239) q/ π E l ec t r on / s p i non d i s p e r s i on " ( q ) = − t cos( q ) − t ! cos(2 q ) − µk F − k F k F − k F k F = − k F − π / FIG. 2.
Electron/spinon bands on the two-leg trian-gular strip.
In the noninteracting
U/t = 0 limit, the groundstate of our model for t (cid:48) /t > . Hamiltonian. For concreteness, we take the followingmodel potential: V ( x − x (cid:48) ) = U , | x − x (cid:48) | = 0 κU e − γ | x − x (cid:48) | , ≤ | x − x (cid:48) | ≤ , | x − x (cid:48) | > . (2)The reasoning for considering such longer-ranged repul-sion in the model Hamiltonian is twofold. First, suchterms are well-motivated by recent ab initio calcula-tions [35, 36], which indicate a substantial long-rangedtail in the effective screened Coulomb repulsion appro-priate for κ -(BEDT-TTF) Cu (CN) . Second, on thetwo-leg ladder, such terms fight the spin-gap tenden-cies present in the metallic phase of the t - t (cid:48) - U Hubbardmodel (i.e., our model with κ = 0; see, for example,Refs. [43–45]), thus at least allowing for the possibilityof a direct, continuous transition between a spin gaplesstwo-band metal and two-band spinon metal spin liquid.Guided by the weak and intermediate coupling analysis ofRef. [28], in what follows we choose characteristic param-eters t (cid:48) /t = 0 . κ = 0 .
5, and γ = 0 .
2, leaving the singledimensionless ratio
U/t to control the overall strength ofelectron repulsion.
III. MOTT METAL-INSULATOR TRANSITIONAND REALIZATION OF THE ELECTRONICSPINON METAL
We first sketch the low-energy effective theory describ-ing the putative metal to spinon metal transition andthen present strong numerical evidence that this exoticscenario is indeed realized. In the absence of interactions(
U/t = 0), our model for t (cid:48) /t > . α S β denotes aLuttinger liquid with α gapless charge modes and β gap-less spin modes [43]—is stable in our extended Hubbardmodel, Eqs. (1)-(2), in the presence of infinitesimal U/t .At half-filling, there is an allowed eight-fermion umklappterm in our two-band system (see Fig. 2). Bosonizing(see, e.g., Refs. [29, 46–48]) this interaction gives H = 2 u cos(4 θ ρ + ) , (3)where θ ρ + is the density field for the overall charge mode,i.e., δn ( x ) = 2 ∂ x θ ρ + /π is the coarse-grained electron den-sity. Assuming the C2S2 metal is stable against openingof a spin gap [28], then the fixed-point Lagrangian L C2S2 involves four gapless bosonic modes, one being θ ρ + (seeAppendix B and Ref. [28] for details). For free electrons,the scaling dimension of the eight-fermion umklapp termis ∆[ H ] = 4 >
2, so that H is strongly irrelevant atweak coupling. However, increasing U/t in our micro-scopic model will feed into “stiffening” θ ρ + in L C2S2 , thusdecreasing ∆[ H ]. Eventually ∆[ H ] = 2, beyond which the umklapp is relevant so that u grows at long scales pin-ning θ ρ + into one of the minima of the cosine potential in H . The resulting phase is a remarkable C1S2 Luttingerliquid, which is precisely the electronic spinon metal [13],The remaining “charge mode” does not transport chargealong the ladder but rather represents local current loopfluctuations; it encodes long-wavelength fluctuations ofthe spin chirality as discussed in Ref. [13].The critical theory describing the C2S2 → C1S2 metal-insulator transition is a sine-Gordon-like theory [49], witha technical complication arising because θ ρ + is coupled tothe “relative charge” field θ ρ − in L C2S2 (see Appendix B).Nonetheless, the transition is still Kosterlitz-Thouless-like [50] [(1 + 1)D XY] and represents a direct, nontrivialtwo-leg analog of the (2 + 1)D scenario recently proposedby Senthil [16, 17].We now present our numerical results, giving strong ev-idence that the above scenario is actually realized. To nu-merically characterize the system, we focus on four mainquantities: the density structure factor (cid:104) δn q δn − q (cid:105) , thespin structure factor (cid:104) S q · S − q (cid:105) , the dimer structure factor (cid:104)B q B − q (cid:105) , and the electron momentum distribution func-tion (cid:104) c † qα c qα (cid:105) , where δn q , S q , B q , and c qα are the Fouriertransforms of the local operators δn ( x ) ≡ n ( x ) − (cid:104) n ( x ) (cid:105) , S ( x ) ≡ (cid:80) α,β c † α ( x ) σ αβ c β ( x ), B ( x ) ≡ S ( x ) · S ( x + 1),and c α ( x ), respectively. In the data presented here, weconsider systems up to L = 96 sites with periodic bound-ary conditions. (See Appendix A for all details, includingdiscussion of the chosen boundary conditions.)We focus first on the density (charge) structure fac-tor (cid:104) δn q δn − q (cid:105) . A crucial aspect of (cid:104) δn q δn − q (cid:105) lies in itsability to distinguish metallic from insulating behaviorat small wavevectors q . For a metallic state, we expect (cid:104) δn q δn − q (cid:105) ∼ | q | for q ∼
0. Specifically, for the two-bandC2S2 metal, the slope of (cid:104) δn q δn − q (cid:105) at q = 0 is relatedto the “Luttinger parameter” g ρ + for the overall chargemode θ ρ + : (cid:104) δn q δn − q (cid:105) = 2 g ρ + | q | /π as q → . (4)Importantly, the quantity g ρ + as determined from Eq. (4)gives a direct measure of the scaling dimension of H :∆[ H ] = 4 g ρ + (see Appendix B 3). Once ∆[ H ] < g ρ + < / g ρ + → (cid:104) δn q δn − q (cid:105) becomes quadratic at small q : (cid:104) δn q δn − q (cid:105) ∼ q in the Mott insulator.In Fig. 3, we show a series of density structure fac-tor measurements ranging from the noninteracting limitat U/t = 0 to deep in the Mott insulating phase at
U/t = 7 .
0. In the inset, we show estimates of g ρ + byplotting (cid:104) δn q δn − q (cid:105) / (2 | q | /π ) [see Eq. (4)]. Based on theabove arguments, we see that the Mott transition occursnear a critical value of U/t = 1 . g ρ + drops below1/2. Note, however, that for these system sizes (cid:104) δn q δn − q (cid:105) still appears linear in q until much larger overall repul-sion, i.e., U/t (cid:39) .
0. Still, we argue that the systembecomes insulating at
U/t = 1 .
6, as this is where H is q/π ! δ n q δ n − q " q/π ! δ n q δ n − q " / ( q / π ) U/t k F k F k F k F − k F π FIG. 3.
Density structure factor: Locating theMott transition and power-law Friedel oscillations ina Mott insulator.
Measurements of the density struc-ture factor, (cid:104) δn q δn − q (cid:105) , allow us to locate the Mott tran-sition near U/t = 1 . ∗ symbols). Theonset of the Mott transition occurs when the overall chargeLuttinger parameter g ρ + drops below 1/2. We measure g ρ + via the slope of (cid:104) δn q δn − q (cid:105) at q = 0, as shown inthe inset [see Eq. (4)]. For U/t > .
6, the system is in-sulating, yet displays power-law singularities in (cid:104) δn q δn − q (cid:105) at finite wavevectors [51] (see black (cid:63) and hexagram sym-bols). Data correspond to a system of length L = 96 with U/t = 0 , . , . , . , . , . , . , . , . , . , . determined to be relevant based on the measurement of g ρ + . That is, we, rather remarkably, have an insulatingstate with a charge correlation length comparable to oursystem size ( L = 96) for 1 . (cid:46) U/t (cid:46) .
0. Indeed, suchlarge correlation lengths are expected in the weak Mottinsulating spinon metal, which we now argue is preciselythe phase realized immediately on the insulating side ofour model. (For more discussion on the finite-size behav-ior of g ρ + , we refer the reader to Appendix B 3.)To this end, we now turn to the spin structure fac-tor (cid:104) S q · S − q (cid:105) in Fig. 4. In the noninteracting limit U/t = 0, we have familiar singularities at wavevec-tors q = 2 k F , k F , π/ , k F − k F originating fromvarious “2 k F ” processes in our two-band system (seeFig. 2). These singularities are simple slope discontinu-ities, i.e., the scaling dimension for the spin operator ateach wavevector is unity as guaranteed by Wick’s theo-rem. As we enter the putative interacting C2S2 metal byturning on finite U/t , the scaling dimensions at wavevec-tors 2 k F , k F , π/ , k F − k F are renormalized slightlybut remain near unity.Near the Mott transition value U/t = 1 . (cid:104) δn q δn − q (cid:105) above, we observe the remark-able result that the singular features in (cid:104) S q · S − q (cid:105) allsurvive, and those at q = 2 k F , k F , π/ enhanced upon entering the insulating phase. In-deed, these are characteristic signatures of the spinonmetal. (In Figs. 3-5, we display characteristic C1S2 q/π ! S q · S − q " q/π ! S q · S − q " / ( q / π ) k F − k F k F U/t k F π FIG. 4.
Spin structure factor: Watching electronsevolve into spinons.
Measurements of the spin structurefactor, (cid:104) S q · S − q (cid:105) , strongly point toward the presence of gap-less spin excitations in both the metal and putative spinonmetal immediately after the Mott transition at U/t = 1 . ∗ symbols). Gapless spin excitations arecharacterized by (cid:104) S q · S − q (cid:105) ∼ | q | as q →
0, and, as shownin the top inset, the opening of a spin gap occurs only for
U/t (cid:38) .
0, at which point the system dimerizes. The “2 k F ”features of the two electron bands in the metallic phase areinherited by the two spinon bands in the spinon metal, and,as highlighted in the bottom inset for q = 2 k F , they are ac-tually enhanced. Data correspond to the same U/t values andcolor scheme as in Fig. 3. spinon metal data at
U/t = 4 . (cid:104) S q · S − q (cid:105) still correspond to the same “2 k F ” processesas in the metallic phase, but with the charge gappedthey now correspond to spinon transfers across the Fermisea. Second, in the spinon metal, we indeed expect thescaling dimensions of the spin operator at wavevectors2 k F , k F , π/ Q = 2 k F , k F , π/ θ ρ + , i.e., S Q ∼ e ± iθ ρ + ( · · · )—see Appendix B 3 andRef. [13]. Thus, pinning of θ ρ + at the Mott transition re-duces the fluctuating content of the spin operator at thesewavevectors, which in turn reduces the scaling dimen-sions and, ultimately, enhances the structure factor sin-gularities. This enhancement is actually a (1+1)D real-ization of “Amperean” attraction between a spinon “par-ticle” and “hole” moving in opposite directions [10, 13].In the density structure factor measurements of Fig. 3,we also have singular features at the “2 k F ” wavevectors q = 2 k F , k F , π/ , k F − k F within the metallic phase,and in fact in the noninteracting U/t = 0 limit, the den-sity and spin structure factors as defined are identical: (cid:104) δn q δn − q (cid:105) = (cid:104) S q · S − q (cid:105) . In the interacting C2S2 metal,the features at q = 2 k F , π/ , k F − k F are still clearly (cid:239) (cid:239) (cid:239) q/π ! B q B − q " q/π ! B q B − q " k F π k F − k F k F k F U/t
FIG. 5.
Dimer structure factor: Period-2 valencebond solid order in the strong Mott insulator.
Mea-surements of the dimer structure factor, (cid:104)B q B − q (cid:105) , show theemergence of a C0S0 period-2 valence bond solid for U/t (cid:38) .
0. Its long-range order is very clearly demonstrated by theprominent Bragg peaks at q = π , as shown in the inset. Datacorrespond to the same U/t values and color scheme as inFigs. 3 and 4. In the main panel (inset), we show data onlyfor the metal and spinon metal (valence bond solid) corre-sponding to values
U/t < . U/t ≥ . visible. In fact, some of these features survive even uponentering the putative insulating spinon metal and remainuntil U/t (cid:39) . (cid:63) symbols in Fig. 3). Thatis, we have power-law density correlations at finite 2 k F wavevectors—a manifestation of which are the famousFriedel oscillations common in metals—even in a Mottinsulator!Indeed, this remarkable result is expected in the two-band spinon metal theory, where, as with the spin op-erator, the slowly varying part of the density opera-tor at wavevectors Q = 2 k F , k F , π/ θ ρ + (but not the wildly fluctuating conjugate field ϕ ρ + ),i.e., δn Q ∼ e ± iθ ρ + ( · · · ). Thus, we should even expectthe scaling dimension of the density operator at thesewavevectors to be reduced due to the same Amperean at-traction mechanism responsible for enhancement of spincorrelations in Fig. 4. However, there are overridingnonuniversal amplitudes that are expected to be smallin a Mott insulator thus preventing observation of thisenhancement—this is likely the case in our data. Fur-thermore, we see development of a feature, though appar-ently weak or with very small amplitude, as anticipated,at a wavevector q = 4 k F = − k F (see black hexagramsymbols in Fig. 3). This feature is again expected fromtheory and is actually a four-fermion contribution to thedensity operator [13] (and thus is extremely weak at weakcoupling). Interestingly, all these power-law density cor-relations in our electronic two-band spinon metal are adirect two-leg analog [52] of the charge Friedel oscillationsexpected on the insulating side of the continuous Motttransition in higher dimensions, as recently stressed byMross and Senthil [51]. q/π U / t ⟨ c † qα c qα ⟩ FIG. 6.
Electronic momentum distribution function:Disappearance of the Fermi surface.
A dense scan of theelectron momentum distribution function, (cid:104) c † qα c qα (cid:105) , over U/t shows the gradual disappearance of the Fermi surface withincreasing interactions, as we move from a two-band C2S2metal (
U/t < .
6) across the insulating C1S2 spinon metal(SM) (1 . < U/t (cid:46) .
0) to the C0S0 valence bond solid in-sulator (
U/t (cid:38) . L = 96 sitesystem as shown in Figs. 3-5. Returning to the spin sector, we can use the small q behavior of (cid:104) S q · S − q (cid:105) to assess whether or not the spinsector is gapless in the realized phases. In analogy withEq. (4), for a spin gapless state we have (cid:104) S q · S − q (cid:105) = 3 g σ + | q | / π as q → , (5)where g σ + is the “Luttinger parameter” associated withthe overall spin mode θ σ + , which for a gapless SU(2) in-variant fixed point is necessarily unity: g σ + = 1 (see Ap-pendix B 3 and also, e.g., Refs. [29, 53]). In the top insetof Fig. 4, we show (cid:104) S q · S − q (cid:105) / (3 | q | / π ), where we see thatfor free electrons g σ + = 1, while increasing U/t pushesthe L = 96 estimate of g σ + above unity—this increas-ing trend continues until U/t (cid:39) .
0, i.e., well beyond theMott critical value of
U/t = 1 .
6. This robust increasingmeasurement of g σ + > g σ + → L → ∞ )well into the insulator is a strong indicator that the spinis gapless on both the metallic and insulating sides ofthe Mott transition, lending strong credence that we areindeed observing the sought-after C2S2 → C1S2 scenariodescribed above. In Appendix B, we discuss these resultsin more depth and make comparisons to how g σ + behavesin the on-site t - t (cid:48) - U Hubbard model at κ = 0.Eventually, above U/t (cid:39) . g σ + dropsbelow unity and (cid:104) S q · S − q (cid:105) ∼ q for small q , indicat-ing the opening of a spin gap. We identify this strongMott insulating phase as a fully gapped (C0S0) period-2 valence bond solid, which is continuously connectedto the dimerized phase realized by the J - J Heisenbergmodel [54] (and also the on-site t - t (cid:48) - U Hubbard modelat large
U/t [55]). To this end, we turn to the dimerstructure factor in Fig. 5. In the inset, we indeed seeclear Bragg peaks developing in (cid:104)B q B − q (cid:105) at q = π for U/t (cid:38) .
0, hence strongly indicative of period-2 valencebond solid order. Furthermore, the operator content ofthe density, δn ( x ), and bond energy, B ( x ), are identical atall wavevectors except π (see Ref. [13] and Appendix B 3).Thus, in the gapless phases (C2S2 and C1S2) we expectsingularities in (cid:104)B q B − q (cid:105) at the same “2 k F ” wavevectorsfor which we find singularities in (cid:104) δn q δn − q (cid:105) (see Fig. 3).Indeed, in the main plot of Fig. 5 we clearly see featuresin (cid:104)B q B − q (cid:105) at q = 2 k F , k F , k F − k F , and 4 k F . Oncein the putative C1S2 insulator, these features are moreapparent in (cid:104)B q B − q (cid:105) than in (cid:104) δn q δn − q (cid:105) since the latterare expected to have small amplitudes in a Mott insula-tor. In our data, this is especially true at wavevectors2 k F and 4 k F , the latter of which is the very nontrivialfour-fermion contribution discussed above.Finally, we discuss the behavior of the electron mo-mentum distribution function (cid:104) c † qα c qα (cid:105) as shown for adense scan of U/t values in Fig. 6. Beyond the Motttransition, when the field θ ρ + gets pinned, we expectthe electron Green’s function to decay exponentially sothat the power-law singularities in (cid:104) c † qα c qα (cid:105) at the fourFermi points q = ± k F , ± k F become gapped. Whileit is not exceedingly apparent that finite correlationlengths emerge at the Fermi points when we cross theMott transition at U/t = 1 . g ρ + measurements—see Fig. 3), we believe this is anothermanifestation of the large charge correlation lengthspresent in the exotic C1S2 insulator. Deep into the pu-tative C1S2 phase though, e.g., for U/t (cid:39) .
0, finite cor-relation lengths are more apparent.
IV. DISCUSSION AND OUTLOOK
In this paper, we have explored the Mott tran-sition between a metal and a quantum spin liquid,presenting strong evidence through large-scale DMRGsimulations in quasi-1D that such a continuous tran-sition can be realized in reasonable electronic mod-els. Our study is strongly motivated by recent exper-iments on the quasi-two-dimensional organic materials κ -(BEDT-TTF) Cu (CN) and EtMe Sb[Pd(dmit) ] ,each of which is a quantum spin liquid that can be driventhrough a Mott transition to a Fermi liquid under pres-sure. We believe our simulations of an extended Hubbardmodel—a model well-motivated by recent ab initio calcu-lations [35, 36] on κ -(BEDT-TTF) Cu (CN) —representan important first step toward numerically characterizingthis transition. While our study is restricted to the two-leg triangular strip, it does show the universal physicsof a clear and direct quasi-1D analog of the continu-ous Mott metal-to-spin liquid transition in two dimen-sions [16]. It is important to point out that the physicsrealized above is markedly distinct from the well-knownstrictly one-dimensional case where a single nested pair ofFermi points gaps out at infinitesimal U/t = 0 + . In our case, we have two unnested pairs of Fermi points whichgap out simultaneously at some finite and intermediatevalue of U/t . Thus, qualitatively speaking, our resultsare remarkably reminiscent of what would happen in fulltwo dimensions where the entire Fermi surface gaps outat the transition [16].Just as importantly, our calculations also elucidate theremarkable properties of the spin-liquid state stabilizedon the insulating side. In many ways, this electronic“spinon metal” weak Mott insulator, as realized in ourmodel, behaves very much like a metal on length scalesshorter than the charge correlation length, and indeedexhibits long-distance density and spin correlations rem-iniscent of the nearby metallic phase (see Figs. 3 and 4).It is precisely this striking similarity between the metallicand insulating states—in basically all properties exceptthe finite charge correlation length in the latter—whichmakes a continuous Mott metal-insulator transition plau-sible, perhaps even likely.Going forward, it would clearly be desirable to movetowards two dimensions and explore the Mott transitionin models such as Eq. (1) on wider ladders and eventu-ally in full 2D, with the goal to make real connectionswith the actual experiments [20, 26, 27, 34]. In the end,the transition may turn out to not be continuous butinstead be weakly first order, as is perhaps realized in κ -(BEDT-TTF) Cu (CN) . Still, our numerical calcula-tions presented here, as well as the recent field theoreticwork of Senthil et al. , suggest that a continuous Motttransition in the ( d + 1)D XY universality class betweena metal and quantum spin liquid is a very real, excitingpossibility. ACKNOWLEDGMENTS
We would like to thank Adrian del Maestro, ArunParamekanti, Federico Becca, Hsin-Hua Lai, DavidMross, William Witczak-Krempa, T. Senthil, and K.Kanoda for useful discussions. R.V.M. is especially grate-ful to Max Metlitski for help on the RG calculationdiscussed in Appendix B 2. This work was supportedby the NSF under grants DMR-1101912 (R.V.M. andM.P.A.F.), PHY11-25915 (R.G.M.), and DMR-1206096(O.I.M.); MICINN through grant FIS2009–13520 (I.G.).R.G.M acknowledges support from NSERC, the CanadaResearch Chair program, the John Templeton Founda-tion, and the Perimeter Institute (PI) for TheoreticalPhysics. Research at PI is supported by the Govern-ment of Canada through Industry Canada and by theProvince of Ontario through the Ministry of EconomicDevelopment & Innovation. We also acknowledge sup-port by the Caltech Institute of Quantum Informationand Matter, an NSF Physics Frontiers Center with thesupport of the Gordon and Betty Moore Foundation(O.I.M. and M.P.A.F.). This work was made possible bythe computing facilities of the Center for Scientific Com-puting from the CNSI, MRL: an NSF MRSEC award(DMR-1121053), and an NSF grant (CNS-0960316); andCESGA. I.G. acknowledges hospitality from the Univer-sity of California, Santa Barbara, during a research staywhen much of this work was done.
Appendix A: Details of DMRG calculations andobservables
We use large-scale DMRG calculations to calculateground state properties of our model Hamiltonian,Eqs. (1)-(2), on finite-size chains of length L sites. Whilewe have performed simulations with both open and peri-odic boundary conditions, we find the latter to be prefer-able for our model in spite of the well-known more chal-lenging convergence properties with periodic boundariesin DMRG calculations. The long-ranged nature of our in-teraction potential [Eq. (2)], however, makes open bound-aries problematic. The issue is that all interactions upto fourth neighbor are chosen to scale with the overallHubbard strength U , so that, at least for the parameterschosen in our study, it is energetically favorable for theend sites of an open chain to become doubly occupied atlarge U/t . That is, even though the system then has topay very large on-site U on the end sites, it gains signif-icant energy by not having to pay as substantial V to V . Therefore, for the calculations on the extended Hub-bard model presented in the main text, we have employedperiodic boundary conditions.To numerically characterize the ground state proper-ties of the system with the DMRG, we calculate the den-sity structure factor (cid:104) δn q δn − q (cid:105) , the spin structure fac-tor (cid:104) S q · S − q (cid:105) , the dimer structure factor (cid:104)B q B − q (cid:105) , andthe electron momentum distribution function (cid:104) c † qα c qα (cid:105) (where α = ↑ , ↓ with no implied summation). In eachcase, the structure factor is defined as the Fourier trans-form of the associated two-point function. Specifically,we have (cid:104) δn q δn − q (cid:105) = 1 L (cid:88) x,x (cid:48) e − iq ( x − x (cid:48) ) (cid:104) δn ( x ) δn ( x (cid:48) ) (cid:105) , (A1) (cid:104) S q · S − q (cid:105) = 1 L (cid:88) x,x (cid:48) e − iq ( x − x (cid:48) ) (cid:104) S ( x ) · S ( x (cid:48) ) (cid:105) , (A2) (cid:104)B q B − q (cid:105) = 1 L (cid:88) x,x (cid:48) e − iq ( x − x (cid:48) ) (cid:104)B ( x ) B ( x (cid:48) ) (cid:105) , (A3) (cid:104) c † qα c qα (cid:105) = 1 L (cid:88) x,x (cid:48) e − iq ( x − x (cid:48) ) (cid:104) c † α ( x ) c α ( x (cid:48) ) (cid:105) , (A4)where n ( x ) ≡ (cid:80) α = ↑ , ↓ c † α ( x ) c α ( x ) is the numberoperator [with δn ( x ) ≡ n ( x ) − (cid:104) n ( x ) (cid:105) ], S ( x ) ≡ (cid:80) α,β c † α ( x ) σ αβ c β ( x ) is the spin operator, and B ( x ) ≡ S ( x ) · S ( x + 1) is the bond energy operator. For sim-plicity, we set (cid:104)B ( x ) B ( x (cid:48) ) (cid:105) = 0 if B ( x ) and B ( x (cid:48) ) sharecommon sites [13]. When presenting all structure factormeasurements, we only show data for q ≥ q = 0.For the dimer structure factor in Eq. (A3), we do notsubtract a product of local averages from the (cid:104)B ( x ) B ( x (cid:48) ) (cid:105) correlations as we do, e.g., for the density structure factorin Eq. (A1). The main reason for this choice is that atlarge U/t (cid:38) . (cid:104)B ( x ) (cid:105) = (cid:104) S ( x ) · S ( x + 1) (cid:105) . This islikely due to the somewhat awkward way in which peri-odic boundaries are implemented in a traditional DMRGsetup which treats the end sites on a different footing.Fourier transforming (cid:104)B ( x ) B ( x (cid:48) ) (cid:105) − (cid:104)B ( x ) (cid:105)(cid:104)B ( x (cid:48) ) (cid:105) thenwashes out the Bragg peaks preasent at q = π . Hence,we just use (cid:104)B ( x ) B ( x (cid:48) ) (cid:105) as the real-space two-point func-tion and exclude plotting (cid:104)B q B − q (cid:105) at q = 0. This bothwell captures the obvious Bragg peaks at q = π in theC0S0 and also gives very clear power-law singularities atthe various “2 k F ” wavevectors as expected in the C1S2insulator (see Fig. 5, Appendix B 3, and Ref. [13]).More generally, we find that the averaging done in ourFourier transforms when summing over both x and x (cid:48) in Eqs. (A1)-(A3) does an effective job of representingthe structure factors in cases where, due to slight lack ofconvergence in the DMRG ground state, the two-pointfunctions depend on both the separation distance x − x (cid:48) and the “origin” x (cid:48) . (Of course, for a perfectly transla-tionally invariant state the two-point functions dependonly on x − x (cid:48) .)In our DMRG calculations, we keep up to m = 6000states and perform at least 6 finite-size sweeps which re-sults in a density matrix truncation error of on the orderof 10 − or smaller. All measurements are well-convergedto the extent necessary to establish the statements madein the main text. To get a feel for the difficulty encoun-tered in obtaining highly accurate data on the stiffnessparameters g ρ + and g σ + (see the main text and Ap-pendix B below), one can observe the data in the insetsof Figs. 3 and 4 at the free electron point U/t = 0—basically the most challenging point for the DMRG. Forfree electrons, we should have g ρ + = g σ + = 1. We seethat there is a rather severe error at the first allowed mo-mentum q = 2 π/L , yet the error for momenta q > π/L is very acceptable, on the order of 1% or less. Appendix B: Luttinger liquid description andsolution by bosonization
In this section, we spell out the effective low-energydescription of the C2S2 metal and C1S2 spinon metaland intervening Kosterlitz-Thouless (KT)-like Mott tran-sition, focusing on those aspects of the theory mostrelevant to the DMRG results presented in the maintext. Some aspects of our presentation follow that ofRefs. [13, 28].
1. Long-wavelength description of C2S2 metal andC1S2 spinon metal
Consider noninteracting electrons at half-filling on thetwo-leg triangular strip (see Fig. 1). When viewed as a 1Dchain with first-neighbor and second-neighbor hopping, t and t (cid:48) , the electron dispersion is given by (see also Fig. 2) (cid:15) ( q ) = − t cos( q ) − t (cid:48) cos(2 q ) − µ. (B1)For t (cid:48) /t > .
5, which is the case of interest here, theground state consists of two disconnected Fermi seas(bands) which we label by a = 1 ,
2. We take the conven-tion that the Fermi velocities v F a are positive (negative)for electrons moving near k F a ( − k F a ), corresponding toright and left movers, respectively. Furthermore, tak-ing the system to be at half-filling gives the sum rule k F + k F = − π/ π .As usual [29], we take the low-energy continuum limitand expand the electron operator in terms of slowly vary-ing continuum fields at the four Fermi points: c α ( x ) = (cid:88) a,P e iP k Fa x c P aα , (B2)where α = ↑ , ↓ denotes the electron spin, and the sumruns over a = 1 , P = R/L =+ / − for the right and left moving electrons at the Fermipoints of each Fermi sea. Although not written explic-itly, the continuum fields of course depend on position x : c P aα = c P aα ( x ).Next, we bosonize [29] the continuum fields accordingto c P aα = η aα e i ( ϕ aα + P θ aα ) , (B3)where ϕ aα and θ aα are the canonically conjugate bosonicphase and phonon fields, respectively. Specifically, wehave [ ϕ aα ( x ) , ϕ bβ ( x (cid:48) )] = [ θ aα ( x ) , θ bβ ( x (cid:48) )] = 0 , (B4)[ ϕ aα ( x ) , θ bβ ( x (cid:48) )] = iπδ ab δ αβ Θ( x − x (cid:48) ) . (B5)The fields η aα are the Klein factors, i.e., Majoranafermions { η aα , η bβ } = 2 δ ab δ αβ , which are necessary to en-sure the correct anticommutation relations among differ-ent fermionic species aα . Finally, the slowly varying com-ponent of the electron density is given by the derivative ofthe θ aα fields: ρ aα = (cid:80) P = ± c † P aα c P aα = ∂ x θ aα /π , where c † P aα c P aα = ∂ x ( θ aα + P ϕ aα ) / (2 π ). Hence, Eq. (B5) isessentially a statement of the density-phase uncertaintyrelation: [ ρ ( x ) , ϕ ( x (cid:48) )] = iδ ( x − x (cid:48) ).Next, we linearize about the Fermi points and expressthe problem in terms of the bosonized fields introduced above. Working in the Euclidean path integral formal-ism, the low-energy continuum Lagrangian density forthe two-band noninteracting electron gas then reads: L free = H free + (cid:88) a,α iπ ( ∂ x θ aα )( ∂ τ ϕ aα ) , (B6)where H free = (cid:88) a,α v F a π (cid:2) ( ∂ x θ aα ) + ( ∂ x ϕ aα ) (cid:3) . (B7)We now introduce the “charge” and “spin” modes foreach band: θ aρ/σ ≡ √ θ a ↑ ± θ a ↓ ) , (B8)and the “overall” and “relative” combinations with re-spect to the two bands: θ µ ± ≡ √ θ µ ± θ µ ) , (B9)where µ = ρ, σ . Fields analogous to Eqs. (B8) and (B9)are also defined for the ϕ ’s. These newly defined fieldssatisfy the same canonical commutation relations as theoriginal fields [Eqs. (B4)-(B5)]. The free-electron La-grangian L free then as usual decouples into charge andspin sectors: L free = L ρ free + L σ free , (B10)where L µ free = H µ free + (cid:88) a iπ ( ∂ x θ aµ )( ∂ τ ϕ aµ ) , (B11) H µ free = (cid:88) a v F a π (cid:2) ( ∂ x θ aµ ) + ( ∂ x ϕ aµ ) (cid:3) . (B12)We are finally in position to discuss interactions. In theinteracting C2S2 Luttinger liquid, the fixed-point theoryis similar to Eq. (B10), i.e., L C2S2 = L ρ C2S2 + L σ C2S2 , (B13)except we have general mode velocities and, in the chargesector, nontrivial Luttinger parameters. For conveniencein the discussion that follows, in the charge sector wework in the ρ ± basis of Eq. (B9) and write the mostgeneral charge sector Lagrangian as L ρ C2S2 = H ρ C2S2 + iπ ∂ x Θ T · ∂ τ Φ , (B14) H ρ C2S2 = 12 π (cid:2) ∂ x Θ T · A · ∂ x Θ + ∂ x Φ T · B · ∂ x Φ (cid:3) , (B15)where Θ T ≡ ( θ ρ + , θ ρ − ) and Φ T ≡ ( ϕ ρ + , ϕ ρ − ); A and B are symmetric, positive definite 2x2 matrices whichencode interactions. Note that even for free electrons, if v F (cid:54) = v F , the charge sector is not diagonal in the ρ ± basis, i.e., A = A (cid:54) = 0, B = B (cid:54) = 0, and in generalthe interacting C2S2 metal will have coupled ρ + and ρ − modes [28].For the spin sector, we stay in the band basis a = 1 , L σ C2S2 = H σ C2S2 + (cid:88) a iπ ( ∂ x θ aσ )( ∂ τ ϕ aσ ) , (B16) H σ C2S2 = (cid:88) a v aσ π (cid:20) g aσ ( ∂ x θ aσ ) + g aσ ( ∂ x ϕ aσ ) (cid:21) . (B17)SU(2) invariance dictates only trivial Luttinger param-eters in the spin sector, i.e., g σ = g σ = 1 (see Ap-pendix B 3), but we keep them general in Eq. (B17)for further analysis below. Our representation of thespin sector here is somewhat schematic in that allowedstrictly marginal chiral interactions will couple the barespin modes [Eq. (B8)] in the quadratic part of the C2S2action. However, the resulting H σ C2S2 is symmetric underinterchanging θ aσ ↔ ϕ aσ and so can easily be broughtback to diagonal form via a simple orthogonal transfor-mation which acts identically on the θ aσ and ϕ aσ fields,hence keeping the Luttinger parameters at their trivialvalues. Thus, for the quadratic part of the C2S2 fixed-point theory, Eq. (B17) is completely general for ourpurposes. Interestingly, the full C2S2 fixed-point theoryalso contains a strictly marginal chiral interband scatter-ing term of the form ( H σ chiral ) ⊥ ∼ cos(2 ϕ σ − ) cos(2 θ σ − ),which is nonharmonic [53]. However, we expect thatthe presence of this, presumably exactly marginal, non-harmonic chiral interaction will not quantitatively alterthe spin sector at the C2S2 (and C1S2; see below) fixedpoint—at least with respect to the Luttinger parametersand contributions to the scaling dimensions of variousoperators (see Appendix B 3). In fact, assuming that( H σ chiral ) ⊥ is exactly marginal already implies trivial spinsector Luttinger parameters, g σ = g σ = 1, which is en-couraging.In addition to such strictly marginal interactions, thereare many nonchiral interactions allowed by symmetrywhich may be added to Eq. (B13) and potentially desta-bilize the C2S2 theory described above. To connect toa given microscopic Hamiltonian, a common approach isto employ a weak-coupling renormalization group (RG)scheme. That is, one can project the microscopic in-teractions onto all continuum symmetry-allowed inter-actions and read off initial conditions for all such cou-plings; these initial conditions can then be subsequentlyused in a controlled RG analysis valid for weak micro-scopic coupling U/t (cid:28)
1. Then, bosonizing the four-fermion interactions—particularly those that may flowto strong coupling, hence destabilizing the “mother”C2S2—emits a direct physical interpretation of the re-sulting phase. This is the approach pioneered many yearsago in Ref. [43], where it was shown (see also Ref. [44])that for the on-site t - t (cid:48) - U Hubbard model, the C2S2metal is generally unstable at weak repulsive interactions to the opening of a spin gap. The basic idea is that theRG flow equations—which are indeed rather complicatedfor the two-band system and in general require a detailednumerical analysis—have a tendency to eventually drive attractive divergent couplings in the spin sector (e.g., theterms denoted g aσ in Ref. [43] or, equivalently, λ σaa inRef. [28]). These divergent couplings conspire to gap outall modes except the overall conducting charge mode θ ρ + ,leaving a one-mode C1S0 conducting Luttinger liquid, es-sentially the quasi-1D analog of a superconductor.However, this spin-gap tendency is not unavoidable.For example, one can fight such pairing tendencies byadding longer-ranged repulsion to the model Hamilto-nian. This approach was recently explored systematicallyin Ref. [28], where it was shown that the C2S2 metal oc-cupies a substantial portion of the weak-coupling phasediagram for the model considered in our work: Eqs. (1)-(2). Stability of the C2S2 metal at weak coupling in-deed seems to be a necessary component for realizingthe C2S2 → C1S2 Mott transition presented numericallyin the main text, and we buttress off the weak-couplingphase diagram presented in Ref. [28] when selecting thespecific parameters of our model Hamiltonian.Finally, as stressed in the main text, our Mott transi-tion is driven at strong interactions by an eight-fermion umklapp term wherein both spin-up and spin-down elec-trons are scattered across each Fermi sea (see Fig. 2): H = u ( c † R ↑ c † R ↓ c † R ↑ c † R ↓ c L ↑ c L ↓ c L ↑ c L ↓ + H . c . ) , (B18)which when written in terms of the bosonized fields sim-ply becomes a cosine of the overall charge field θ ρ + : H = 2 u cos(4 θ ρ + ) . (B19)The C1S2 spinon metal spin liquid corresponds to rele-vance of H so that u flows to strong coupling. That is,the field content of the C1S2 fixed-point theory looksidentical to that of C2S2 but with a massive overallcharge mode θ ρ + . Specifically, we have L C1S2 = L ρ C1S2 + L σ C1S2 , (B20)where the “charge sector” now only contains the ρ − mode: L ρ C1S2 = H ρ C1S2 + iπ ∂ x θ ρ − ∂ τ ϕ ρ − , (B21) H ρ C1S2 = v ρ − π (cid:20) g ρ − ( ∂ x θ ρ − ) + g ρ − ( ∂ τ ϕ ρ − ) (cid:21) , (B22)which physically represents gapless local current fluctu-ations, and the spin sector formally reads the same asbefore: L σ C1S2 = L σ C2S2 , (B23)still with trivial Luttinger parameters, g σ = g σ = 1.For an extensive discussion of the C1S2 phase with re-spect to its features and stability, we refer the reader toRef. [13].0
2. Renormalization group analysis of theC2S2 → C1S2 Mott transition
We now present the details of the critical theory de-scribing our Mott transition. The theory is KT-like witha complication arising because the field θ ρ + , which is be-ing gapped out, is coupled to the field θ ρ − in the Gaussianfixed-point action for the C2S2 [see Eq. (B15)], and θ ρ − is massless on both sides of the transition.From the above considerations, the charge sector La-grangian describing the transition between the C2S2metal and C1S2 spinon metal reads L = L + L cos , (B24)where L = 12 π (cid:2) ∂ x Θ T · C · ∂ x Θ + ∂ τ Θ T · D · ∂ τ Θ (cid:3) (B25)is just L ρ C2S2 from Eq. (B14) with the ϕ ’s integrated out, Θ T ≡ ( θ ρ + , θ ρ − ), and L cos = 2 u cos( nθ ρ + ) (B26)with n = 4 is our eight-fermion umklapp term. It isconvenient to diagonalize the quadratic part of the theory L in a fashion similar to that described in Ref. [28], thusobtaining for the full theory L = 12 π (cid:88) i =1 , (cid:20) v i ( ∂ τ θ i ) + v i ( ∂ x θ i ) (cid:21) , (B27) L cos = 2 u cos( n θ + n θ ) , (B28)where we have absorbed the nontrivial Luttinger param-eters of the two normal modes, θ and θ , into the realcoefficients n and n via a rescaling of the fields. While θ and θ are specific linear combinations of θ ρ + and θ ρ − ,e.g., nθ ρ + = n ( c θ + c θ ) = n θ + n θ , we do notspell out the details here, but instead refer the readerto the Appendix of Ref. [28] for a similar calculation.Ultimately, this linear combination, as well as the veloci-ties and Luttinger parameters of the normal modes in thediagonalized system, are rather complicated, but still an-alytic, functions of the original parameters C and D ofthe coupled system.We have performed a renormalization group (RG)analysis of the above two-mode system, obtaining the fol-lowing leading-order KT-like (see below) flow equationsfor all couplings: dC d(cid:96) = πn Λ v I (cid:18) v v , n (cid:19) u , (B29) dD d(cid:96) = πn Λ v (cid:18) v v (cid:19) − n / I (cid:18) v v , n (cid:19) u , (B30) dud(cid:96) = (cid:20) − (cid:18) n n (cid:19)(cid:21) u, (B31) where I ( α, β ) ≡ (cid:90) π dθ cos θ (cos θ + α sin θ ) β ≥ . (B32)As with ordinary KT, the coupling u renormalizes ac-cording to the scaling dimension of the cosine with re-spect to the quadratic action,∆[cos( nθ ρ + )] = ∆[cos( n θ + n θ )] = n n , (B33)and obtaining its beta function, Eq. (B31), can proceedin a textbook Wilsonian fashion [29]. However, renormal-izing the parameters in L is significantly more involvedand depends on the specific regularization scheme em-ployed. First, note that since L cos contains only the field θ ρ + , it cannot possibly renormalize any terms containing θ ρ − to any order in perturbation theory; hence, the onlynonzero beta functions are those for the couplings C and D . The respective beta functions, Eqs. (B29) and(B30), were obtained using a field-theoretic approach [49]in which we consider insertions into correlation func-tions of the form (cid:104) ∂ x θ i ( x ) ∂ x θ j ( y ) (cid:105) , where x and y arepoints in our (1+1)D space-time. At O ( u ), one hasto integrate over two 2D points from two u insertions,say z and z (cid:48) . Indeed, as z − z (cid:48) becomes small, the in-tegral diverges logarithmically, and we cut it off at ashort-distance scale Λ − . We then compute correctionsto (cid:104) ∂ x θ i ( x ) ∂ x θ j ( y ) (cid:105) from posited “counterterms” in L which are chosen to exactly cancel the aforementionedlogarithmic divergence. This allows us, after an alto-gether somewhat lengthy calculation, to arrive at theabove RG flow equations for C and D .The case of vanishing θ ρ + - θ ρ − coupling in Eq. (B25)corresponds to the limit n →
0, so that θ ∝ θ ρ + and C and D renormalize at the same rate (up to an over-all scale of v ). This of course corresponds to ordinaryKosterlitz-Thouless RG wherein only one parameter in L renormalizes: d ( g − ) d(cid:96) ∼ u , with g the single-modeLuttinger parameter [29].In the general case, the beta functions for C and D involve highly nonuniversal content, and thus we havenot attempted a detailed study of the flows. Still, thetransition is KT-like in nature except that two parame-ters (as opposed to one) in L are renormalized by thesingle cosine, and the transition occurs when the scalingdimension of the cosine equals the space-time dimension:∆[cos( nθ ρ + )] = n + n = 2, where n and n are func-tions of the parameters C and D .We can formally argue for the KT-like nature as fol-lows. From the start, we focus only on the flowing param-eters C , D , and u . Let us denote the (non-negative)factors multiplying u in the beta functions for C and D as A ( C , D ) and B ( C , D ), respectively, andalso denote the coefficient of u in the beta function for u as Γ( C , D ). We emphasize that A , B , and Γ arefunctions of C and D , which, while perhaps compli-cated functions, are analytical and not special. As we1vary in the ( C , D ) plane, we generically expect tofind a line where Γ = 0 separating regions where a small u perturbation is relevant or irrelevant. Let us considerone point on this line, ( C (0)11 , D (0)11 ), and study small de-viations ( δC , δD ) from this point. The RG equationsare, to leading order, d δC d(cid:96) = A (0) u , (B34) d δD d(cid:96) = B (0) u , (B35) dud(cid:96) = (cid:16) α (0) δC + β (0) δD (cid:17) u, (B36)where A (0) and B (0) are the A and B functions evalu-ated at ( C (0)11 , D (0)11 ), while α (0) and β (0) are derivatives ∂ Γ /∂C and ∂ Γ /∂D evaluated at the same point. De-viations satisfying α (0) δC + β (0) δD = 0 correspond tomoving along the Γ = 0 line, while generic deviations willcut across this line. Formally, we can change variables to r = α (0) δC + β (0) δD , s = − β (0) δC + α (0) δD ,which flow as drd(cid:96) = (cid:16) α (0) A (0) + β (0) B (0) (cid:17) u , (B37) dsd(cid:96) = (cid:16) − β (0) A (0) + α (0) B (0) (cid:17) u , (B38) dud(cid:96) = ru. (B39)Thus, the flow equations for the r and u variables havefamiliar KT-like form and subsequent standard analysiscan kick in. On the other hand, the flow of the s variableis simply slaved to u and does not affect the KT analysis.In principle, one should be able to confirm the KTuniversality class from the numerical DMRG data, forexample, by performing Weber-Minnhagen [50] style fitsto finite-size estimates of the scaling dimension of the co-sine in the metallic phase (essentially the stiffness in theXY model context; see also Appendix B 3 below). How-ever, this requires highly accurate data on large systemsizes in the scaling regime, which is currently prohibitivefor our multimode electronic system (see Appendix A).Also, it is not unreasonable to expect that the presenceof two renormalizing parameters in L , instead of one,might make the finite-size effects more severe. In the endthough, this is a rather nonuniversal matter which we donot pursue further analytically.
3. Observables and stiffness parameters
To characterize the system, we have focused on thedensity structure factor, the spin structure factor, thedimer structure factor, and the electron momentum dis-tribution function as presented in the main text and asdefined in Appendix A. In this section, we lay out thedetails of the bosonization treatment which allows us to use these measurements, both at finite and zero wavevec-tors, to probe the nature of the Luttinger liquid phasesrealized by our model Hamiltonian. a. Establishing the result ∆[ H ] = 4 g ρ + As stressed in the main text, we can directly measurethe scaling dimension of the eight-fermion umklapp term[see Eqs. (3) and (B19)] responsible for driving our Motttransition by measuring the slope of the density structurefactor at q = 0 momentum [see Eq. (4)]. We now spell outhow these two quantities, ∆[ H ] and g ρ + , are formallyrelated.The former is defined through the corresponding two-point function: (cid:68) e i θ ρ + ( x ) e − i θ ρ + (0) (cid:69) ∼ | x | H ] , (B40)where, for simplicity, we work at equal (imaginary) timesuch that x is a spatial coordinate only. Assuming thatthe system is in the C2S2 phase so that the charge sec-tor is described by the quadratic Lagrangian L ρ C2S2 ofEq. (B14), we can use a standard identity [29] and write (cid:68) e i θ ρ + ( x ) e − i θ ρ + (0) (cid:69) = e − (cid:104) [ θ ρ + ( x ) − θ ρ + (0)] (cid:105) . (B41)Now, the slowly varying component of the total elec-tron density (measured relative to the average density) isgiven by δn ( x ) = 2 ∂ x θ ρ + /π , so that the long-wavelengthcontribution to the density-density correlation functionin real space is given by (cid:104) δn ( x ) δn (0) (cid:105) = 4 π ∂ x ∂ x (cid:48) (cid:104) θ ρ + ( x ) θ ρ + ( x (cid:48) ) (cid:105)| x (cid:48) =0 + · · · . (B42)The right-hand side can be obtained from Eq. (B41)via straightforward manipulations, which after invokingEq. (B40) gives (cid:104) δn ( x ) δn (0) (cid:105) = − ∆[ H ]2 π x + · · · . (B43)On the other hand, we define the slope of themomentum-space density structure factor as q → (cid:104) δn q δn − q (cid:105) = 2 g ρ + π | q | , (B44)such that g ρ + = 1 corresponds to a two-band nonin-teracting electron gas. After Fourier transformation,Eqs. (B43) and (B44) imply that∆[ H ] = 4 g ρ + , (B45)which is the desired result. Note that g ρ + is not gen-erally a genuine Luttinger parameter due to the cou-pling between the ρ + and ρ − sectors in the C2S2 phase,2 U/t g ρ + ( L , n = ) L = 36 L = 48 L = 60 L = 72 L = 96 metal insulator FIG. 7. Finite-size estimates of g ρ + [see Eq. (B46)] versus U/t . Our bosonized theory predicts that measured values of g ρ + < / g ρ + → L → ∞ ). The somewhat irregularfinite-size behavior in the gapless regions ( U/t (cid:46) .
0) is likelydue to “shell filling” effects, i.e., the thermodynamic phase ismore readily accommodated by some sizes and less by others. but should instead be viewed as a direct measurement of∆[ H ] through the density structure factor.In the main text, we relied heavily upon Eq. (B45)to distinguish between metallic and insulating behavior,where measured g ρ + > / g ρ + < /
2) implies that H is irrelevant (relevant) so that the system is metal-lic (insulating). Of course, if ∆[ H ] <
2, the system isnecessarily insulating and Eq. (B44) no longer applies;instead we have (cid:104) δn q δn − q (cid:105) ∼ q as q →
0. That is, mea-sured g ρ + < / g ρ + →
0. InFig. 3, even well into the insulating phase of our modelas determined by the above arguments, we see on our L = 96 site system that apparently (cid:104) δn q δn − q (cid:105) ∼ | q | ;however, with H relevant, this must be a finite-size ef-fect due to the large charge correlation length present inour weak Mott insulating C1S2.In Fig. 7, we show finite-size estimates of the quan-tity g ρ + obtained with DMRG for the same parametersof the extended Hubbard model used in the main text.Specifically, we define g ρ + ( L, n ) ≡ L n (cid:104) δn q δn − q (cid:105) (cid:12)(cid:12) q = n πL , (B46)and monitor g ρ + ( L, n = 2) while varying
U/t . This datalooks rather far removed from ordinary KT behavior po-tentially indicating strong finite-size effects (see also dis-cussion at the end of the previous section). Still, basedon the above analysis, we must have a Mott transitionnear
U/t = 1 .
6. Eventual g ρ + → U/t (cid:38) .
6, although that is not apparent on these sizesuntil deep in the insulating phase, say
U/t (cid:38) .
0. Wenote that the fully gapped C0S0 state (
U/t (cid:38) .
0) doesshow clear (cid:104) δn q δn − q (cid:105) ∼ q behavior, which is not surpris-ing given the short charge correlation length expected in that state. b. Bosonized representation of operators at finitewavevectors We now give the bosonized expressions for the spin anddensity operators at finite “2 k F ” wavevectors and math-ematically establish the Amperean enhancement mech-anism summarized in the main text. Expanding thespin operator as S ( x ) = (cid:80) Q S Q e iQx , we can easily writethe slowly varying part of the spin operator at variouswavevectors, i.e., S Q = S Q ( x ), in terms of the right andleft moving electron operators defined in Appendix B 1: S k Fa = 12 c † Laα σ αβ c Raβ , (B47) S π/ = 12 c † R α σ αβ c L β + 12 c † R α σ αβ c L β , (B48) S k F − k F = 12 c † R α σ αβ c R β + 12 c † L α σ αβ c L β . (B49)Similarly, for the density operator, we have δn k Fa = c † Laα c Raα , (B50) δn π/ = c † R α c L α + c † R α c L α , (B51) δn k F − k F = c † R α c R α + c † L α c L α . (B52)In each case, summations over spin indices are implied,and S − Q = S † Q and δn − Q = δn † Q . Throughout, our useof denoting wavevectors with either Q or q is an attemptto distinguish the long-wavelength component of an op-erator, O Q , from the actual exact operator used in theDMRG, O q .Bosonizing the above electron bilinears using Eq. (B3)results in the following expressions for the spin: S x k Fa = − iη a ↑ η a ↓ e iθ ρ + e ± iθ ρ − sin( √ ϕ aσ ) , (B53) S y k Fa = − iη a ↑ η a ↓ e iθ ρ + e ± iθ ρ − cos( √ ϕ aσ ) , (B54) S z k Fa = − e iθ ρ + e ± iθ ρ − sin( √ θ aσ ) , (B55) S xπ/ = e − iθ ρ + (cid:104) − iη ↑ η ↓ e − iθ σ − sin( ϕ ρ − + ϕ σ + ) − iη ↓ η ↑ e iθ σ − sin( ϕ ρ − − ϕ σ + ) (cid:105) , (B56) S yπ/ = e − iθ ρ + (cid:104) − iη ↑ η ↓ e − iθ σ − cos( ϕ ρ − + ϕ σ + )+ iη ↓ η ↑ e iθ σ − cos( ϕ ρ − − ϕ σ + ) (cid:105) , (B57) S zπ/ = e − iθ ρ + (cid:104) − iη ↑ η ↑ e − iθ σ + sin( ϕ ρ − + ϕ σ − )+ iη ↓ η ↓ e iθ σ + sin( ϕ ρ − − ϕ σ − ) (cid:105) , (B58)3 S xk F − k F = e − iθ ρ − (cid:104) − iη ↑ η ↓ e − iθ σ + sin( ϕ ρ − + ϕ σ + ) − iη ↓ η ↑ e iθ σ + sin( ϕ ρ − − ϕ σ + ) (cid:105) , (B59) S yk F − k F = e − iθ ρ − (cid:104) − iη ↑ η ↓ e − iθ σ + cos( ϕ ρ − + ϕ σ + )+ iη ↓ η ↑ e iθ σ + cos( ϕ ρ − − ϕ σ + ) (cid:105) , (B60) S zk F − k F = e − iθ ρ − (cid:104) − iη ↑ η ↑ e − iθ σ − sin( ϕ ρ − + ϕ σ − )+ iη ↓ η ↓ e iθ σ − sin( ϕ ρ − − ϕ σ − ) (cid:105) , (B61)and for the density: δn k Fa = 2 ie iθ ρ + e ± iθ ρ − cos( √ θ aσ ) , (B62) δn π/ = 2 e − iθ ρ + (cid:104) − iη ↑ η ↑ e − iθ σ + sin( ϕ ρ − + ϕ σ − ) − iη ↓ η ↓ e iθ σ + sin( ϕ ρ − − ϕ σ − ) (cid:105) , (B63) δn k F − k F = 2 e − iθ ρ − (cid:104) − iη ↑ η ↑ e − iθ σ − sin( ϕ ρ − + ϕ σ − ) − iη ↓ η ↓ e iθ σ − sin( ϕ ρ − − ϕ σ − ) (cid:105) , (B64)where for expressions with ± in the exponent, + refersto band a = 1, while − refers to band a = 2.Perhaps the most important point to take away is thatall operators at Q = 2 k F a , π/ e ± iθ ρ + .Therefore, the fluctuating field content of these operatorsis reduced upon gapping out (pinning of) θ ρ + when cross-ing the Mott transition from the C2S2 metal to C1S2 in-sulator. This leads to lowering of the associated scalingdimensions and subsequent enhancement of the structurefactor singularities. To illustrate this concretely, assumefor the moment that the ρ + and ρ − sectors are decou-pled in the charge sector Lagrangian for the C2S2, i.e., A = A = B = B = 0 in Eq. (B15), with cor-responding Luttinger parameters g ρ + and g ρ − . We thenhave the following for the scaling dimensions of the aboveoperators:∆[ S k Fa ] = ∆[ δn k Fa ] = 12 + g ρ − g ρ + , (B65)∆[ S π/ ] = ∆[ δn π/ ] = 12 + 14 g ρ − + g ρ + , (B66)∆[ S k F − k F ] = ∆[ δn k F − k F ] = 12 + 14 g ρ − + g ρ − , (B67)where we have assumed SU(2) invariance, g σ = g σ = 1(see the next section). Right at the Mott transition g ρ + = 1 /
2, while immediately on the insulating side g ρ + →
0. Therefore, the dimensions in Eqs. (B65)and (B66) corresponding to operators at Q = 2 k F a , π/ decrease at the transition (by an amount of1/8 in the decoupled approximation). Such an enhance-ment of the associated spin structure factor singularitieson the insulating side of the Mott transition is in factdramatically seen in the DMRG data of Fig. 4.Furthermore, stability of the C1S2 insulator requires g ρ − < S π/ ] > ∆[ S k Fa ](and similarly for δn Q ). Thus, for the structure factorsin the C1S2 phase, the features at q = 2 k F a should bemore pronounced than those at q = π/
2. Indeed, this isobserved in the spin structure factor data of Fig. 4 onthe insulating side of the Mott transition in our model.More generally, the presence of clear power-law singular-ities in (cid:104) S q · S − q (cid:105) at finite wavevectors in both the metaland weak Mott insulator points strongly towards to pres-ence of gapless spin excitations in both phases (see alsoAppendix B 3 c).Note that the density operator at Q = 2 k F a , π/ , k F − k F still remains power law when θ ρ + gets pinned, i.e., δn Q does not contain the wildly fluctuating field ϕ ρ + .In fact, for Q = 2 k F a , π/ θ ρ + [see Eqs. (B62), (B63)] and has the same scal-ing dimension as the spin operator: ∆[ δn Q ] = ∆[ S Q ]!Therefore, such Friedel oscillations should actually be en-hanced in the Mott insulator [51]. This enhancement isdifficult to see in the density structure factor DMRG dataof Fig. 3, but that is likely due to the small amplitudesof the features. The power-law nature, however, is stillapparent, at least around q = 2 k F , k F − k F .The bilinears that get enhanced, i.e., those at Q =2 k F a , π/
2, can be predicted by simple “Amperean rules”.Specifically, in the (1+1)D U(1) gauge theory formulationof the C1S2 spinon metal phase [13], θ ρ + corresponds tothe mode that is pinned upon inclusion of gauge fluc-tuations which implements at long wavelengths the con-straint of one spinon per site (in this language, the upand down spinons carry the same gauge charge). We thenexpect that the bilinears that get enhanced upon intro-ducing the gauge fluctuations are those composed fromoperators that produce parallel gauge currents, so-calledAmperean attraction [10, 13]. This is indeed the case forthe spin and density operators at Q = 2 k F a , π/ Q = k F − k F involve op-erators with antiparallel gauge currents and are there-fore not enhanced; indeed these operators do not contain θ ρ + at all. We remark that in our electronic model, theabove “gauge constraint” is implemented dynamically byelectron repulsion upon pinning of the overall conductingcharge mode θ ρ + .In the main text, we have also used the dimer corre-lations, as defined and detailed in Appendix A, to char-acterize the ground state. Following Ref. [13], we canapproximate the bond energy as the electron hoppingenergy, i.e., B ( x ) ∼ − t (cid:80) α (cid:2) c † α ( x ) c α ( x + 1) + H . c . (cid:3) . Infact, in our DMRG measurements it would have been4reasonable to use this as the definition of B ( x ), but weinstead implemented the full B ( x ) = S ( x ) · S ( x + 1),which makes the two-point function (cid:104)B ( x ) B ( x (cid:48) ) (cid:105) a four-spin (eight-electron) measurement. In any case, expan-sion in continuum fields reveals B Q ∼ e iQ/ δn Q , (B68)which holds for all Q (cid:54) = π . Hence, we expect features atthe same wavevectors in measurements of both (cid:104) δn q δn − q (cid:105) and (cid:104)B q B − q (cid:105) . This is indeed observed in Figs. 3 and 5,where in the putative C1S2 insulator the power-law na-ture of the features is, as expected, much more apparentin the dimer correlations than in the density correlations.We further note that (cid:104)B q B − q (cid:105) very clearly picks up afeature at q = 4 k F = − k F , while this feature is muchweaker, though still present, in (cid:104) δn q δn − q (cid:105) . As mentionedin the main text, the wavevector 4 k F = − k F is a four-fermion contribution to the density/bond energy. Specif-ically, δn k F : c † L ↑ c † L ↓ c R ↑ c R ↓ ∼ e i θ ρ + e i θ ρ − , (B69) δn − k F : c † R ↑ c † R ↓ c L ↑ c L ↓ ∼ e − i θ ρ + e i θ ρ − , (B70)both contribute with independent numerical prefactors,and have scaling dimensions in the decoupled ρ ± approx-imation of∆[ δn k F ] = ∆[ B k F ] = g ρ + + g ρ − . (B71)In the C1S2, g ρ + → B k F ] = g ρ − .Gaplessness of the spin sector requires g ρ − < q = 4 k F in (cid:104)B q B − q (cid:105) should be stronger than a slope discontinuity(unit scaling dimension of the associated operator)—thisindeed appears to be the case in our dimer structure fac-tor data of Fig. 5.There is yet another important four-fermion contribu-tion to the spin and density/bond energy at wavevec-tor Q = π . We here focus on the latter, wherefor the bond energy we get contributions such as [13] B π : iδn k F δn k F + H . c . , which when bosonized gives B π ∼ [cos(2 θ σ + ) + cos(2 θ σ − )] sin(2 θ ρ + ) + · · · . (B72)This operator has unit scaling dimension at the C1S2fixed point (∆[ B π ] = 1) and should thus correspond toa slope discontinuity in (cid:104)B q B − q (cid:105) at q = π . Remark-ably, this appears to be consistent with e.g. our charac-teristic C1S2 data point at U/t = 4 . θ ρ + due to relevance of H = 2 u cos(4 θ ρ + ) is such thatsin(2 θ ρ + ) (cid:54) = 0. This is precisely what we would expectif the pinned value of θ ρ + occurs at 4 θ ρ + = π mod 2 π ,which corresponds to the minima of cos(4 θ ρ + ). We thusconclude that u > u < θ ρ + such that 4 θ ρ + = 0 mod 2 π ,i.e., sin(2 θ ρ + ) = 0, thus killing the feature in (cid:104)B q B − q (cid:105) at q = π .At wavevector Q = π , the bond-centered density B π is odd under mirror symmetry ( x → − x ), while the site-centered density δn π is even. Contributions to the latterinclude δn π : δn k F δn k F + H . c . , which in terms of thebosonized fields reads δn π ∼ [cos(2 θ σ + ) + cos(2 θ σ − )] cos(2 θ ρ + ) + · · · . (B73)Hence, the pinning condition 4 θ ρ + = π mod 2 π inferredabove implies cos(2 θ ρ + ) = 0. Indeed, the DMRG datashows no feature in (cid:104) δn q δn − q (cid:105) at q = π within the pu-tative C1S2 phase (see Fig. 3). Again, we conclude thatfor our system with repulsively interacting electrons, wemust have u > H .Finally, presence of a feature at q = π in (cid:104) δn q δn − q (cid:105) in the C1S2 weak Mott insulator would lead to long-range period-2 (site-centered) charge density wave orderin the C0S0 strong Mott insulator at very large U/t .This is indeed very unnatural in our model where the on-site U term is the largest interaction energy scale in theHamiltonian. Instead, the strong Mott insulator realizedin our model develops period-2 long-range order in the bond-centered density, as evidenced by the Bragg peak in (cid:104)B q B − q (cid:105) at q = π . The power-law feature at the samewavevector in the weak Mott insulator [see Eq. (B72)] isthe precursor of this eventual long-range VBS order atlarge U/t .We finally discuss the electron operator itself[Eq. (B3)], which is of course the most primitive oper-ator of all. When written in terms of “ ρ ± ” and “ aσ ”modes, we have c P aα = η aα exp (cid:26) i √ (cid:20) √ ϕ ρ + ± ϕ ρ − ) ± ϕ aσ (cid:21) + iP √ (cid:20) √ θ ρ + ± θ ρ − ) ± θ aσ (cid:21)(cid:27) , (B74)where the first ± on each line refers to a = 1 ,
2, while thesecond refers to α = ↑ , ↓ . Of course, once the θ ρ + field ispinned, the electron Green’s function (cid:104) c † α ( x ) c α (0) (cid:105) is ex-pected to decay exponentially at all wavevectors. Mathe-matically, this is due to its conjugate field ϕ ρ + also beingpresent in the bosonized representation of the electronoperator: By the uncertainty principle, pinning of θ ρ + will cause ϕ ρ + to fluctuate wildly leading to exponentialdecay of the Green’s function. While it is somewhat diffi-cult to ascertain this exponential decay within the puta-tive C1S2 phase for the electron momentum distributionfunction DMRG data of Fig. 6, we again believe this isdue to the excessively large charge correlation lengthspresent in our electronic spinon metal.From Eq. (B74), we also see that gapping of a spinmode will cause the associated electron Fermi point togap out, and thus the electron Green’s function canin principle detect spin-gap behavior. However, this is5rather difficult in practice [45], and in the following sec-tion we discuss a better approach as employed in themain text. c. Assessing gaplessness of the spin sector through g σ + Inspection of the bosonized expressions for the dif-ferent components of the spin operator at wavevectors Q = 2 k F a in Eqs. (B53)-(B55), reveals that in the fixed-point theory for either the C2S2 metal or C1S2 insulatorwe must have only trivial Luttinger parameters in thespin sector: g σ = g σ = 1. Specifically, for arbitrary g aσ as in Eq. (B17) and decoupled ρ + and ρ − modes as inthe illustrative discussion in Appendix B 3 b above, wehave∆[ S x k Fa ] = ∆[ S y k Fa ] = g ρ + g ρ − g aσ , (B75)∆[ S z k Fa ] = g ρ + g ρ − g aσ , (B76)where in the C1S2 insulator we have g ρ + →
0. There-fore, SU(2) spin invariance manifest through isotropicspin-spin correlations functions at wavevectors 2 k F a , i.e.,∆[ S x k Fa ] = ∆[ S y k Fa ] = ∆[ S z k Fa ], indeed dictates that g σ = g σ = 1 , (B77)which constitutes a simple generalization of the well-known one-mode case [29] (see also Ref. [53]).We now show how measurement of the spin structurefactor at zero momentum can assess the condition inEq. (B77). The slowly varying part of the spin densityis S z ( x ) = ∂ x θ σ + /π , hence the long-wavelength part ofthe real-space spin-spin correlation function evaluated inthe fixed-point theory for either the C2S2 or C1S2 [seeEq. (B17)] reads (cid:104) S z ( x ) S z (0) (cid:105) = − g σ + π x + · · · , (B78)where we have defined g σ + ≡ g σ + g σ . (B79)Equation (B78) gives for the spin structure factor as q → (cid:104) S zq S z − q (cid:105) = g σ + π | q | , (B80)which we use in the main text to estimate the parameter g σ + [see Eq. (5) and the inset of Fig. 4]. Clearly thenwithin the fixed-point theory we should have g σ + = 1,while in the presence of a spin gap (cid:104) S zq S z − q (cid:105) ∼ q , so that g σ + →
0. Note that, as with g ρ + above, g σ + is not agenuine Luttinger parameter as even free electrons arenot generally diagonal in the σ ± basis. The above considerations are valid for the fixed pointin the thermodynamic limit. However, there are severalmarginal interactions that need to be irrelevant for thespin sector to remain gapless and the C2S2 and C1S2 tobe stable phases. Thus, the presence of such marginallyirrelevant interactions will affect measurement of g σ + onfinite-size systems. In the case of our C2S2 and C1S2,the residual interactions in the spin sector that mix rightand left movers read H σRL = − (cid:88) a,b ( w σab J Rab · J Lab + λ σab J Raa · J Lbb ) , (B81)where J P ab ≡ c † P aα σ αβ c P bβ . In the C2S2 and C1S2,the w σab terms are strictly irrelevant, while the λ σab termsare only marginally irrelevant [13, 28]. Bosonizing thelatter interactions gives˜ H σRL = V z + V ⊥ , (B82) V z = (cid:88) a λ σaa π (cid:2) ( ∂ x ϕ aσ ) − ( ∂ x θ aσ ) (cid:3) (B83)+ λ σ π [( ∂ x ϕ σ )( ∂ x ϕ σ ) − ( ∂ x θ σ )( ∂ x θ σ )] , (B84) V ⊥ = (cid:88) a λ σaa cos(2 √ θ aσ ) (B85)+ 2 λ σ ˆΓ cos(2 θ σ + ) cos(2 ϕ σ − ) , (B86)where ˆΓ ≡ η ↑ η ↓ η ↑ η ↓ .A necessary condition for the spin to be gapless is thatthe couplings λ σab be initially positive, corresponding tothe system being overall repulsive in the spin sector. Ul-timate stability of the C2S2 and C1S2 corresponds to λ σab renormalizing to zero via slow marginal flows. It shouldin principle be possible to calculate precise flows (andfinite-size scaling behavior) of our effective g σ + param-eter by analyzing the behavior of the zero-momentumpiece of the spin structure factor perturbatively in the λ σab . We do not pursue this here, but instead to geta rough, initial feel for the trends within our Abelianbosonization, imagine for the moment naively ignoringthe V ⊥ cosines and λ σ cross terms. Then, the quadratic V z terms effectively feed into renormalizing the g aσ Lut-tinger parameters above (below) unity for λ σaa positive(negative), hence effectively corresponding to g σ + > g σ + <
1) on a finite-size system. This is indeed theexpected trend for overall repulsion in the spin sector.On the other hand, the flows for the C1S0 supercon-ductor (the main instability of the C2S2) correspond to λ σaa eventually becoming negative (attraction in the spinsector) and then diverging to −∞ . All modes then even-tually get gapped out except the overall conducting ρ +mode [43, 44], so that for the spin structure factor wehave (cid:104) S zq S z − q (cid:105) ∼ q as q →
0, i.e., g σ + →
0. On a finite-size system, we thus expect the spin gap to be manifestas a measured g σ + <
1. Note, though, that due to initial repulsion in the spin sector [ λ σab ( (cid:96) = 0) > U/t g σ + ( L , n = ) L = 36 L = 48 L = 72 L = 96 C1S2 insulator,spin liquid C0S0insulator,VBSC2S2 metal
FIG. 8. Finite-size estimates of g σ + [see Eq. (B87)] versus U/t for the same parameters in the extended Hubbard modelat the focus of the main text. The putative realized phases(see text) are labeled with separating vertical dashed-dottedlines. At
U/t = 0, our DMRG calculations give g σ + ( L, n =2) = 1 to within 1% for all sizes; this serves as a very usefulcheck on our convergence since free electrons are, ironically,very challenging to converge in the DMRG. on relatively short length scales, i.e., measured g σ + > t - t (cid:48) - U Hub-bard model [45]. We stress, however, that in our modelwith longer-ranged repulsion—a model which is known tobe spin gapless at weak coupling (
U/t (cid:28)
1) for our chosenparameters [28]—measurements of g σ + still strongly in-dicate spin gaplessness all the way up to U/t (cid:39) .
0, wellpast the Mott critical value of
U/t = 1 .
6. In the nextsection, we contrast this with the behavior of the on-site t - t (cid:48) - U Hubbard model at κ = 0 in which the metal andinsulator are presumably both spin gapped.Finally, we again mention that the observed power-law singularities in the spin structure factor at the var-ious “2 k F ” wavevectors (see the main text and Ap-pendix B 3 b) provide complementary evidence that thespin sector is gapless in both the metal (C2S2) and weakMott insulator (C1S2) of our model.
4. Further analysis of g σ + DMRG data
Here we present more data of our DMRG measure-ments of the parameter g σ + discussed in the previous sec-tion. Specifically, we define a finite-size estimate of g σ + via Eq. (5) by evaluating the slope of the spin structurefactor at a momentum q = n πL with n a small integer: g σ + ( L, n ) ≡ L n (cid:104) S q · S − q (cid:105) (cid:12)(cid:12) q = n πL , (B87)where in what follows we choose n = 2.In Figs. 8 and 9, we show g σ + ( L, n = 2) versus
U/t onseveral system sizes L for the extended Hubbard model aspresented in the main text [Eqs. (1)-(2) with t (cid:48) /t = 0 . U/t g σ + ( L , n = ) L = 36 L = 48 L = 72 L = 96 L = 128 C1S0 metal C0S0 insulator,VBS
FIG. 9. Finite-size estimates of g σ + [see Eq. (B87)] ver-sus U/t for the on-site t - t (cid:48) - U Hubbard model at t (cid:48) /t = 0 . U/t = 3 . g ρ + measurements (notshown; see Ref. [56]). This value is in good agreement withearlier studies of the half-filled t - t (cid:48) - U Hubbard model [45, 55].Here, we use open boundary conditions which gives very goodconvergence, though at the expense of some small systematicerror in determining g σ + from the momentum-space struc-ture factor; e.g., g σ + ( L, n = 2) is slightly less than one at
U/t = 0 which is due entirely to the usage of open boundaryconditions. κ = 0 . γ = 0 .
2] and the on-site t - t (cid:48) - U Hubbard model[Eqs. (1)-(2) with t (cid:48) /t = 0 . κ = 0], respectively. In theformer case, we use periodic boundary conditions due tothe reasons discussed in Appendix A, while in the lattercase we use standard open boundary conditions. Notethat the L = 96 data in Fig. 8 corresponds to the second( q = 2 π ) data points in the inset of Fig. 4.We first focus on the extended Hubbard model dataas shown in Fig. 8. Here, g σ + ( L ) increases above unityas we turn on U/t and continues to do so well past theputative Mott transition from the C2S2 metal to C1S2insulator at
U/t = 1 .
6. Rather remarkably, the data doesnot start renormalizing visibly downwards until
U/t (cid:38) .
0. Around
U/t (cid:39) .
0, the system starts showing signsof spin-gap behavior (e.g., a Bragg peak in the dimerstructure factor; see Fig. 5) near which g σ + ( L ) finallystarts bending downward. While the data points on thelarge sizes are still not fully converged due to the periodicboundary conditions and inherent difficulty involved inconverging such a quantity at small momenta, we believethat as L → ∞ we would find g σ + = 1 for U/t (cid:46) . g σ + = 0 for U/t (cid:38) . S tot in the ground state (we work only in the S z tot = 0sector in the DMRG) by evaluating the computed spinstructure factor at q = 0: (cid:104) S q · S − q (cid:105) (cid:12)(cid:12) q =0 = 1 L (cid:104) S (cid:105) = 1 L S tot ( S tot + 1) . (B88)7In simulations of Eqs. (1)-(2) with periodic boundaryconditions, we often find for S tot some small nonintegervalue on the order of unity. For example, on L = 96 siteswith m = 6000 states, at the free electron point U/t = 0,we find S tot = 0 .
60, and at the characteristic C1S2 spinonmetal point
U/t = 4 .
0, we find S tot = 0 .
46. However, webelieve this is just a benign effect of our inability to fullyconverge the DMRG and the eventual ground state at m → ∞ will be a spin-singlet with S tot = 0. We knowthis to be true at U/t = 0, while all indications pointtoward a spin-singlet C1S2 for 1 . < U/t (cid:46) .
0, e.g.,the features at 2 k F and 2 k F are symmetrically locatedabout q = π/ (cid:104) S q · S − q (cid:105) (see Fig. 4).In fact, this convergence difficulty is to be expected inour parameter regime of t (cid:48) /t = 0 .
8, as realization of thetwo-band spinon metal in a pure spin model with ringexchanges (Ref. [13]) found similar DMRG convergenceproblems in the corresponding parameter regime of thatmodel.Also, these difficulties are likely responsible for thesmall “jumps” in the data in Fig. 8, since measured fi-nite total spin will have a small, somewhat unpredictable,quantitative effect on our g σ + ( L, n ) values. For instance,we are able to converge to a singlet for all
U/t on the L = 36 site system, and hence its curve is smooth. Onthe other hand, on the L = 48 site system, the mea-sured total spin starts abruptly dropping toward zeronear U/t = 4 .
4, and we believe this behavior is respon-sible for the corresponding feature in the L = 48 curveof Fig. 8. Ultimately, however, these convergence prob-lems will almost certainly have no qualitative effect onour conclusions being drawn from the g σ + data.In Fig. 9, we show analogous g σ + ( L, n = 2) measure-ments for the ordinary on-site t - t (cid:48) - U Hubbard model at t (cid:48) /t = 0 .
8. This model has a spin gap at weak coupling
U/t (cid:28) U , thesystem is initially repulsive (stable) in the spin sector,while the eventual gapping out of both the spin modesand the “ ρ − ” mode happens due to a delicate interplayof all channels (see Fig. 3 of Ref. [28]). We believe thisinitial repulsion in the spin sector is responsible for mea-sured g σ + > g σ + < J S i · S j that favors a spin-gapped(Luther-Emery) liquid (see Ref. [56]).The Mott transition in the t - t (cid:48) - U Hubbard model willalso be driven by the same eight-fermion umklapp termdiscussed above. By measuring its scaling dimension inthe same fashion as we have done for the extended model(see Fig. 3 and Appendix B 3), we have determined thatfor the U -only Hubbard model at t (cid:48) /t = 0 . U/t = 3 .
5, after which period-2VBS order sets in immediately (see Ref. [56] for moredetails). We see, however, that g σ + ( L ) already startsbending downward well before then. We stress that thisis in sharp contrast to the data of Fig. 8 in which ourmodel with longer-ranged repulsion shows no signs of aspin gap until well past the Mott transition. In that case,the intervening phase is the spin gapless C1S2 spin liquidinsulator. [1] N. F. Mott, Metal-Insulator Transitions (Taylor & Fran-cis Inc., USA, 1990).[2] M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys. , 1039 (1998).[3] A. Zylbersztejn and N. F. Mott, Phys. Rev. B , 4383(1975).[4] P. Limelette, A. Georges, D. Jrome, P. Wzietek, P. Met-calf, and J. M. Honig, Science , 89 (2003).[5] P. Limelette et al. , Phys. Rev. Lett. , 016401 (2003).[6] F. Kagawa, T. Itou, K. Miyagawa, and K. Kanoda, Phys.Rev. B , 064511 (2004).[7] F. Kagawa, K. Miyagawa, and K. Kanoda, Nature ,534 (2005).[8] P. Anderson, Materials Research Bulletin , 153 (1973).[9] P. W. Anderson, Science , 1196 (1987).[10] P. A. Lee, N. Nagaosa, and X.-G. Wen, Rev. Mod. Phys. , 17 (2006).[11] L. Balents, Nature , 199 (2010).[12] S.-S. Lee and P. Lee, Phys. Rev. Lett. , 036403 (2005).[13] D. N. Sheng, O. I. Motrunich, and M. P. A. Fisher, Phys.Rev. B , 205112 (2009).[14] O. I. Motrunich, Phys. Rev. B , 045105 (2005). [15] S. Florens and A. Georges, Phys. Rev. B , 035114(2004).[16] T. Senthil, Phys. Rev. B , 045109 (2008).[17] D. Mross and T. Senthil, Phys. Rev. B , 165126 (2011).[18] M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S.Fisher, Phys. Rev. B , 546 (1989).[19] Y. Shimizu, K. Miyagawa, K. Kanoda, M. Maesato, andG. Saito, Phys. Rev. Lett. , 107001 (2003).[20] Y. Kurosaki, Y. Shimizu, K. Miyagawa, K. Kanoda, andG. Saito, Phys. Rev. Lett. , 177001 (2005).[21] T. Itou, A. Oyamada, S. Maegawa, M. Tamura, and R.Kato, Phys. Rev. B , 104413 (2008).[22] S. Yamashita et al. , Nature Phys. , 459 (2008).[23] M. Yamashita et al. , Science , 1246 (2010).[24] M. Yamashita et al. , Nature Phys. , 44 (2009).[25] W. Witczak-Krempa, P. Ghaemi, T. Senthil, and Y. B.Kim, Phys. Rev. B , 245102 (2012).[26] K. Kanoda, talk at “KITP Conference on Exotic Phasesof Frustrated Magnets”, October 2012, http://online.kitp.ucsb.edu/online/fragnets_c12/kanoda/ .[27] T. Furukawa, K. Miyagawa, H. Taniguchi, R. Kato, andK. Kanoda, Nature Phys. , 221 (2015). [28] H.-H. Lai and O. I. Motrunich, Phys. Rev. B , 045105(2010).[29] T. Giamarchi, Quantum Physics in One Dimension , In-ternational Series of Monographs on Physics (OxfordUniversity Press, New York, 2003).[30] T. Yoshioka, A. Koga, and N. Kawakami, Phys. Rev.Lett. , 036401 (2009).[31] H.-Y. Yang, A. M. L¨auchli, F. Mila, and K. P. Schmidt,Phys. Rev. Lett. , 267204 (2010).[32] J. Kokalj and R. H. McKenzie, Phys. Rev. Lett. ,206402 (2013).[33] Z. Y. Meng, T. C. Lang, S. Wessel, F. F. Assaad, and A.Muramatsu, Nature , 847 (2010).[34] R. Kato, A. Tajima, A. Nakao, N. Tajima, and M.Tamura, in
Multifunctional Conducting Molecular Mate-rials , Vol. 306 of
RSC Special Publication Series , editedby G. Saito, F. Wudl, R. C. Haddon, and K. Tanigaki(The Royal Society of Chemistry, Cambridge, 2006), pp.31–36.[35] K. Nakamura, Y. Yoshimoto, T. Kosugi, R. Arita, andM. Imada, Journal of the Physical Society of Japan ,083710 (2009).[36] T. Koretsune and C. Hotta, Phys. Rev. B , 045102(2014).[37] S. R. White, Phys. Rev. Lett. , 2863 (1992).[38] S. R. White, Phys. Rev. B , 10345 (1993).[39] U. Schollw¨ock, Rev. Mod. Phys. , 259 (2005).[40] M. S. Block, D. N. Sheng, O. I. Motrunich, and M. P. A.Fisher, Phys. Rev. Lett. , 157202 (2011).[41] R. V. Mishmash, M. S. Block, R. K. Kaul, D. N. Sheng,O. I. Motrunich, and M. P. A. Fisher, Phys. Rev. B , 245127 (2011).[42] H.-C. Jiang, M. S. Block, R. V. Mishmash, J. R. Garri-son, D. N. Sheng, O. I. Motrunich, and M. P. A. Fisher,Nature , 39 (2013).[43] L. Balents and M. P. A. Fisher, Phys. Rev. B , 12133(1996).[44] K. Louis, J. V. Alvarez, and C. Gros, Phys. Rev. B ,113106 (2001).[45] G. I. Japaridze, R. M. Noack, D. Baeriswyl, and L. Tin-cani, Phys. Rev. B , 115118 (2007).[46] F. D. M. Haldane, Phys. Rev. Lett. , 1840 (1981).[47] H.-H. Lin, L. Balents, and M. P. A. Fisher, Phys. Rev. B , 1794 (1998).[48] M. P. A. Fisher, in Topological aspects of low dimensionalsystems , Vol. 69 of
Les Houches Lecture Series , editedby A. Comtet, T. Jolicoeur, S. Ouvry, and F. David(Springer, Berlin, 1999), pp. 575–641.[49] D. J. Amit, Y. Y. Goldschmidt, and S. Grinstein, Journalof Physics A: Mathematical and General , 585 (1980).[50] H. Weber and P. Minnhagen, Phys. Rev. B , 5986(1988).[51] D. F. Mross and T. Senthil, Phys. Rev. B , 041102(2011).[52] H.-H. Lai and O. I. Motrunich, Phys. Rev. B , 235120(2009).[53] N. Sedlmayr, P. Korell, and J. Sirker, Phys. Rev. B ,195113 (2013).[54] S. R. White and I. Affleck, Phys. Rev. B , 9862 (1996).[55] L. F. Tocchio, F. Becca, and C. Gros, Phys. Rev. B81