Chiral spin liquid with spinon Fermi surfaces in spin-$1/2$ triangular Heisenberg model

TL;DR

This paper uncovers a gapless chiral spin liquid phase with spinon Fermi surfaces in a triangular Heisenberg model, revealing a novel quantum state driven by competing interactions and characterized by emergent fractionalized excitations.

Contribution

It introduces the first example of a gapless CSL with coexisting spinon Fermi surfaces and chiral edge states, supported by large-scale DMRG and a parton mean-field theory.

Findings

Discovery of a gapless chiral spin liquid phase in the model.

Identification of spinon Fermi surfaces and chiral edge modes.

Numerical evidence of emergent fractionalized excitations.

Abstract

We study the interplay of competing interactions in spin-$1/2$ triangular Heisenberg model through tuning the first- ($J_1$), second- ($J_2$), and third-neighbor ($J_3$) couplings. Based on large-scale density matrix renormalization group calculation, we identify a quantum phase diagram of the system and discover a new {\it gapless} chiral spin liquid (CSL) phase in the intermediate $J_2$ and $J_3$ regime. This CSL state spontaneously breaks time-reversal symmetry with finite scalar chiral order, and it has gapless excitations implied by a vanishing spin triplet gap and a finite central charge on the cylinder. Moreover, the central charge grows rapidly with the cylinder circumference, indicating emergent spinon Fermi surfaces. To understand the numerical results we propose a parton mean-field spin liquid state, the $U(1)$ staggered flux state, which breaks time-reversal symmetry with chiral edge modes by adding a Chern insulator mass to Dirac spinons in the $U(1)$ Dirac spin liquid. This state also breaks lattice rotational symmetries and possesses two spinon Fermi surfaces driven by nonzero $J_2$ and $J_3$, which naturally explains the numerical results. To our knowledge, this is the first example of a gapless CSL state with coexisting spinon Fermi surfaces and chiral edge states, demonstrating the rich family of novel phases emergent from competing interactions in triangular-lattice magnets.

Equations

$$H=J_1\sum _{⟨i,j⟩}{\mathbf{S}}_i\cdot {\mathbf{S}}_j+J_2\sum _{⟨⟨i,j⟩⟩}{\mathbf{S}}_i\cdot {\mathbf{S}}_j+J_3\sum _{⟨⟨⟨i,j⟩⟩⟩}{\mathbf{S}}_i\cdot {\mathbf{S}}_j,$$

$${\mathbf{S}}_i=\frac{1}{4}\text{Tr}\left({\psi }_i^{†}{\psi }_i{\mathbf{\sigma }}^T\right),{\psi }_i=\left(\begin{array}{cc}{f}_{i,\uparrow }& f_{i,\downarrow }\\ f_{i,\downarrow }^{†}& -f_{i,\uparrow }^{†}\end{array}\right)$$

$$H_{\mathbf{MF}}=\frac{1}{8}\sum _{ij}{J}_{ij}\text{Tr}\left({\psi }_i^{†}{u}_{ij}{\psi }_j+h.c.\right)+\frac{1}{8}\sum _{ij}{J}_{ij}\text{Tr}\left(u_{ij}^{†}{u}_{ij}\right)$$

$$u_{ij}=\left(\begin{array}{cc}-{\bar{\chi }}_{ij}& 0\\ 0& {\chi }_{ij}\end{array}\right),$$

$$H_{\mathbf{MF}}=\frac{J}{4}\sum _{⟨ij⟩}\sum _{\alpha }\left(-{\bar{\chi }}_{ij}{f}_{i,\alpha }^{†}{f}_{j,\alpha }+h.c.\right)+\frac{J}{4}\sum _{⟨ij⟩}\left(|{\chi }_{ij}{|}^2\right)$$

$$\sum _{\alpha =\uparrow ,\downarrow }{f}_{i\alpha }^{†}{f}_{i\alpha }=1,\forall i.$$

$$f_{\mathbf{r},\alpha }\stackrel{T_2}{\to }(-{)}^{r_1}{f}_{\mathbf{r}+{\mathbf{a}}_2,\alpha }^{†},f_{\mathbf{r},\alpha }\stackrel{T_1}{\to }{f}_{\mathbf{r}+{\mathbf{a}}_1,\alpha }.$$

$$S(\mathbf{k})=\frac{1}{N}\sum _{i,j}⟨{\mathbf{S}}_i\cdot {\mathbf{S}}_j⟩e^{i\mathbf{k}\cdot ({\mathbf{r}}_i-{\mathbf{r}}_j)}.$$

$$D_l(\mathbf{k})=\frac{1}{N}\sum _{i,j}(⟨{\mathcal{B}}_{i,i+{\delta }_l}{\mathcal{B}}_{j,j+{\delta }_l}⟩-⟨{\mathcal{B}}_{i,i+{\delta }_l}⟩⟨{\mathcal{B}}_{j,j+{\delta }_l}⟩)e^{i\mathbf{k}\cdot ({\mathbf{r}}_i-{\mathbf{r}}_j)},$$

$${\delta }_1$$

$${\delta }_2$$

$${\delta }_3$$

$$H_1$$

$$+e^{\text{i}{\varphi }_2}{f}_{\mathbf{r}}^{†}{f}_{\mathbf{r}+{\delta }_2}-e^{-\text{i}{\varphi }_2}{f}_{\mathbf{r}+{\delta }_2}^{†}{f}_{\mathbf{r}+2{\delta }_2}$$

$$+f_{\mathbf{r}}^{†}{f}_{\mathbf{r}+{\delta }_3}+f_{\mathbf{r}+{\delta }_2}^{†}{f}_{\mathbf{r}+{\delta }_3+{\delta }_2}+h.c.),$$

$$f_{\mathbf{r}p}=\frac{1}{\sqrt{N}}\sum _{\mathbf{k}}{e}^{\text{i}\mathbf{k}\cdot \mathbf{r}}{f}_{\mathbf{k}p},$$

$$H_1$$

$$h^{(1)}(\mathbf{k})$$

$$h_{00}^{(1)}(\mathbf{k})$$

$$h_{11}^{(1)}(\mathbf{k})$$

$$h_{01}^{(1)}(\mathbf{k})$$

$$h_{10}^{(1)}(\mathbf{k})$$

$$H_2$$

$$+e^{\text{i}{\phi }_2}{f}_{\mathbf{r}}^{†}{f}_{\mathbf{r}+{\delta }_2+{\delta }_3}+e^{-\text{i}{\phi }_2}{f}_{\mathbf{r}+{\delta }_2}^{†}{f}_{\mathbf{r}+2{\delta }_2+{\delta }_3}$$

$$+e^{\text{i}{\phi }_3}{f}_{\mathbf{r}}^{†}{f}_{\mathbf{r}+{\delta }_3-{\delta }_1}-e^{-\text{i}{\phi }_3}{f}_{\mathbf{r}+{\delta }_2}^{†}{f}_{\mathbf{r}+{\delta }_2+{\delta }_3-{\delta }_1}+h.c.).$$

$$H_2$$

$$h^{(2)}(\mathbf{k})$$

$$h_{00}^{(2)}(\mathbf{k})$$

$$h_{11}^{(2)}(\mathbf{k})$$

$$h_{01}^{(2)}(\mathbf{k})$$

$$h_{10}^{(2)}(\mathbf{k})$$

$$H_3$$

$$+e^{\text{i}{\gamma }_2}{f}_{\mathbf{r}}^{†}{f}_{\mathbf{r}+2{\delta }_2}-e^{-\text{i}{\gamma }_2}{f}_{\mathbf{r}+{\delta }_2}^{†}{f}_{\mathbf{r}+3{\delta }_2}$$

$$+e^{\text{i}{\gamma }_3}{f}_{\mathbf{r}}^{†}{f}_{\mathbf{r}+2{\delta }_3}-e^{-\text{i}{\gamma }_3}{f}_{\mathbf{r}+{\delta }_2}^{†}{f}_{\mathbf{r}+{\delta }_2+2{\delta }_3}+h.c.).$$

$$H_3$$

$$h^{(3)}(\mathbf{k})$$

$$h_{00}^{(3)}(\mathbf{k})$$

$$+\cos (2\mathbf{k}\cdot {\delta }_3+{\gamma }_3),$$

$$h_{11}^{(3)}(\mathbf{k})$$

$$-\cos (2\mathbf{k}\cdot {\delta }_3-{\gamma }_3).$$

Full text

Introduction. Quantum spin liquids (QSLs) are novel quantum phases of matter, which do not exhibit any symmetry-breaking orders even at zero temperature Balents (2010); Savary and Balents (2016); Zhou et al. (2017) but feature long-range entanglement and fractionalized excitations Wen (1991); Senthil and Fisher (2000, 2001); Wen (2013). QSLs have been studied extensively in the past few decades, due to their important role in understanding strongly correlated materials and potential application in topological quantum computation Anderson (1973); Read and Sachdev (1991); Sachdev (1992); Lee et al. (2006); Kitaev (2006). While gapped QSLs have been classified and characterized systematically, there is much less understanding on gapless QSLs and how they could be realized in materials. Although a gapless QSL with Dirac cones of spinons has been shown to exist in the exactly soluble Kitaev model Kitaev (2006), so far there is no definitive evidence that a Dirac spin liquid has been realized in any magnetic materials Savary and Balents (2016); Hermanns et al. (2018). A more exotic state is the gapless QSL with spinon Fermi surfaces (SFSs) Ioffe and Larkin (1989); Nagaosa and Lee (1990); Lee and Nagaosa (1992). Such a QSL has an extensive number of low-energy excitations, and was shown to be stabilized by four-spin ring-exchange couplings that can arise from strong charge fluctuations in weak Mott insulators Misguich et al. (1999); Motrunich (2005); Sheng et al. (2009); Block et al. (2011); He et al. (2018).

Experimentally, many QSL candidate materials fall into the family of layered spin-$1/2$ magnets on the triangular lattice, such as the organic salts Shimizu et al. (2003); Kurosaki et al. (2005); Yamashita et al. (2008, 2009, 2010) and the transition metal dichalcogenides Law and Lee (2017); Yu et al. (2017); Ribak et al. (2017); Klanjšek et al. (2017). Specific heat and thermal transport measurements point towards the presence of extensive mobile gapless spin excitations, which appear to be consistent with a gapless QSL with SFSs Yamashita et al. (2008, 2010); Yu et al. (2017). These materials are considered to be weak Mott insulators with strong charge fluctuations, which may induce such gapless QSL behaviors Misguich et al. (1999); Motrunich (2005); Sheng et al. (2009); Block et al. (2011); He et al. (2018). However, a direct study on the triangular Hubbard model suggests a possible gapped chiral spin liquid (CSL) phase in the intermediate $U$ region Szasz et al. (2018). Therefore, a clear theoretical understanding on the mechanism to realize gapless QSLs in these layered quasi-two-dimensional magnets is still lacking.

Another route to QSL is through competing interactions between different neighboring sites, such as the kagome compound kapellasite Fåk et al. (2012) and $J_{1}$-$J_{2}$-$J_{3}$ kagome model Gong et al. (2014, 2015). Recently, competing interactions have also been found essential to understand possible QSLs in the triangular-lattice rare-earth compounds Paddison et al. (2016); Zhu et al. (2018); Zhang et al. (2018); Li et al. (2015a); Shen et al. (2018) and delafossite oxides Kimura et al. (2006); Ye et al. (2007); Kadowaki et al. (1990); Seki et al. (2008). Indeed, a QSL phase has been found in the spin-$1/2$ $J_{1}$-$J_{2}$ triangular Heisenberg model (THM) although its nature has not been established Kaneko et al. (2014); Iqbal et al. (2016); Li et al. (2015b); Zhu and White (2015); Hu et al. (2015); Zheng et al. (2015); Lu (2016); Saadatmand and McCulloch (2016). Therefore, understanding how QSL phases emerge from competing interactions is an important issue in order to discover new QSL materials Ivanov and Richter (1995); Yao et al. (2018); Iaconis et al. (2018).

In this Rapid communication, we systematically study the spin-$1/2$ $J_{1}$-$J_{2}$-$J_{3}$ THM using density-matrix renormalization group (DMRG) method and the parton construction. The model Hamiltonian is given as

[TABLE]

where $J_{1},J_{2},J_{3}$ are the first-, second-, and third NN interactions as shown in the inset of Fig. 1(a). We choose $J_{1}=1.0$ as the energy scale. In the coupling range $0\leq J_{2}/J_{1}\leq 0.7,0\leq J_{3}/J_{1}\leq 0.4$, besides the previously found $J_{1}$-$J_{2}$ spin liquid and different magnetic orders, we identify a new gapless CSL phase as shown in Fig. 1(a). This CSL state spontaneously breaks time-reversal symmetry (TRS) with a finite scalar chiral order. We also observe spin pumping upon flux insertion, similar to the charge pumping in Laughlin-type fractional quantum Hall states, indicative of a chiral edge mode, which is further confirmed by the entanglement spectrum. Finite-size scaling of spin triplet gap on the square-like clusters shows a vanished spin gap. The gapless nature is further supported by the bipartite entanglement entropy, which exhibits a logarithmic correction of the area law versus subsystem length, leading to a finite central charge. The central charge grows with the cylinder circumference consistent with a QSL with emergent SFSs. We propose a staggered flux state in the Abrikosov-fermion representation of spin-$1/2$ operators, which explains the coexistence of chiral edge mode and SFSs observed in this gapless CSL.

We study the system by using DMRG with $SU(2)$ symmetry White (1992); McCulloch and Gulácsi (2002). We use cylinder geometry with periodic boundary conditions along circumference direction and open boundary conditions along extended direction. The lattice vectors are defined as ${\bf a}_{1}=(1,0)$ and ${\bf a}_{2}=(\frac{1}{2},\frac{\sqrt{3}}{2})$. Two geometries named YC and XC cylinders are studied, both having extended direction along ${\bf a}_{1}$. For the YC and XC cylinders, circumference direction is along ${\bf a}_{2}$ and perpendicular to ${\bf a}_{1}$, respectively. The cylinders are denoted as YC$L_{y}-L_{x}$ and XC$L_{y}-L_{x}$ with $L_{y}$ and $L_{x}$ being the numbers of sites along circumference and extended directions. We study the systems with $L_{y}=5-12$ by keeping up to $8000$ $SU(2)$ states (equivalent to about 24000 $U(1)$ states) to obtain accurate results with truncation error less than $10^{-5}$ in most calculations.

Quantum phase diagram. We demonstrate the quantum phase diagram in Fig. 1(a). With growing $J_{2}$ and $J_{3}$, we find different magnetically ordered phases and QSL phases. In Fig. 1(a), the black dotted lines denote the classical phase boundaries of the 1200, stripe, and zigzag orders. We also find an incommensurate (IC) magnetic order in the neighbor of the zigzag order, consistent with previous spin-wave calculations Ivanov and Richter (1995). The incommensurate order might be considered as the zigzag order with an incommensurate modulation (see Supplemental Material sup ). In the presence of quantum fluctuations, we find a new gapless CSL phase near the triple point of the classical orders, which sits at the neighbor of the $J_{1}$-$J_{2}$ SL. By computing spin and dimer correlation functions, we find featureless spin and dimer structure factors that indicate the absence of spin and dimer orders in the CSL state sup . Next, we further characterize the nature of this new CSL state.

Spontaneous time-reversal symmetry breaking. To detect spontaneous TRS breaking, we use complex-valued wavefunction, which has been applied in DMRG to find chiral ground states in different systems Gong et al. (2014); Zhu et al. (2016). If TRS is spontaneously broken, the system is featured by finite scalar chiral order $\langle\chi\rangle=\langle{\bf S}_{1}\cdot({\bf S}_{2}\times{\bf S}_{3})\rangle$, where ${\bf S}_{i}$ $(i=1,2,3)$ label the three spins on each triangle. On the YC cylinder with both even and odd $L_{y}$, we find a nonzero chiral order in the bulk of cylinder with a large circumference, as shown in Fig. 2(a) for $J_{2}=0.3,J_{3}=0.15$. In these states, the chiral orders of all the up- and down-triangles have the same sign, and the chiral order grows more robust as the circumference increases. On the XC cylinder shown in Fig. 2(b), the chiral order vanishes in the bulk for small circumference but becomes stable on the wide XC12 cylinder. Combining these results we conclude a CSL state with spontaneous TRS breaking.

Spin triplet gap and entanglement characterization. We calculate the spin triplet gap by obtaining the ground state (in the $S=0$ sector) on long cylinder and then sweeping the $S=1$ sector for the middle $N_{x}$ columns Yan et al. (2011), which gives the gap of the middle $N_{x}\times L_{y}$ system. We find that the gap versus $1/N_{x}$ shows length dependence sup . To estimate the gap in the 2D limit and avoid 1D physics, we extrapolate the gap data of the square-like clusters as shown in Fig. 3(a). The gap drops fast as a function of $1/L_{y}$ and smoothly scales to zero, suggesting gapless spin-triplet excitations.

Furthermore, we study entanglement entropy versus subsystem length $l_{x}$ by cutting the cylinder into two parts. Since the real-valued wavefunction is a superposition of the two chiral states with opposite chiralities, it has a higher entanglement entropy and is harder to converge to; thus we also use complex-valued wavefunction to compute entanglement entropy. As shown in Fig. 3(b) for $J_{2}=0.3,J_{3}=0.15$ on the YC6 cylinder, the entropy shows a logarithmic correction of the area law and follows the behavior $S(l_{x})=(c/6)\ln[(L_{x}/\pi)\sin(l_{x}\pi/L_{x})]+g$ Calabrese and Cardy (2004), where $S(l_{x})$ is the bipartite entanglement entropy, $c$ is the central charge, and $g$ is a non-universal constant. The YC6 cylinders with different $L_{x}$ give a consistent central charge $c\simeq 1$. For the YC8 cylinder, we choose $L_{x}=16$ (the entropy for larger $L_{x}$ is much harder to converge and we show the results for $L_{x}=24$ in Supplemental Material sup , which are consistent with the fitted central charge $c=5$). As shown in Fig. 3(c), the entropy continues to grow with kept state number and converges very well by keeping $6000$ $SU(2)$ states, giving a large central charge of $c\simeq 5$. The finite central charge supports the gapless nature of the CSL state. Once a 2D quantum state is confined to a 1D cylinder, the finite circumference quantizes the momentum around the cylinder. The central charge of the 1D cylinder needs to sum over contributions from all quantized momenta. Take the $U(1)$ Dirac spin liquid for example, the cylinder central charge $c=2-1=1$ if the quantized momenta for each spin species only cross one Dirac cone, where the extra $-1$ accounts for the $U(1)$ gauge field fluctuations which gaps out the total spinon density fluctuation Geraedts et al. (2016); sup . Similarly $c\leq 3$ if the quantized momenta cross two Dirac cones (for each spin species), which is an upper bound for the central charge on a cylinder of any $L_{y}$. The large central charge we found from DMRG is therefore inconsistent with the $U(1)$ Dirac spin liquid on triangular lattice Kaneko et al. (2014); Iqbal et al. (2016); Lu (2016), but provides a strong evidence supporting emergent SFSs Ioffe and Larkin (1989); Nagaosa and Lee (1990); Motrunich (2005); Sheng et al. (2009). Now that each pair of crossings (one right mover and one left mover) between the quantized momenta and the SFSs contributes a unit of central charge, the total central charge of SFSs generally grows with $L_{y}$, with an upper bound of $c\leq N_{w}-1$ where $2N_{w}$ is the total number of crossings sup .

For gapped CSL states, entanglement spectrum has a one-to-one correspondence with physical edge spectrum Li and Haldane (2008). Interestingly, for this gapless CSL state entanglement spectrum also shows a quasi-degenerate group of levels with the counting $\{1,1,2,3,5,\cdots\}$ agreeing with chiral $SU(2)_{1}$ conformal field theory Francesco et al. (2012), as shown in Fig. 3(d). This may be the first example of such novel states for interacting system, which demonstrates similar edge physics as the non-interacting $p+ip$ chiral superconductor with a gapless bulk spectrum Dubail and Read (2015); Poilblanc (2017).

The staggered flux state. To understand the DMRG results, we propose a staggered flux state, whose mean-field ansatz is constructed in the Abrikosov-fermion representation of spin-$1/2$ operators Abrikosov (1965)

[TABLE]

The Heisenberg Hamiltonian $H=\sum_{\langle{ij}\rangle}J_{ij}\mathbf{S}_{i}\cdot\mathbf{S}_{j}$ is decoupled into the mean-field form as

[TABLE]

with the mean-field amplitude $u_{ij}=\langle\psi_{i}\psi_{j}^{\dagger}\rangle=u_{ij}^{\dagger}$. In the $U(1)$ QSL states all spinon pairing terms will vanish and thus

[TABLE]

where $\chi_{ij}=\sum_{\alpha}\langle f_{i,\alpha}^{\dagger}f_{j,\alpha}\rangle=\bar{\chi}_{ji}$. Then the mean-field ansatz can be simplified as

[TABLE]

where the mean-field ground state is at half-filling due to the single-occupancy constraint on the parton Hilbert space

[TABLE]

We consider a $U(1)$ spin liquid known as the staggered flux state Wen (2002); Li et al. (2017); Bieri et al. (2016), where fermionic spinons transform under translations as follows

[TABLE]

Although the mean-field ansatz doubles the unit cell (along ${\bf a}_{2}$ direction), the projected wavefunction preserves the lattice translation symmetries by ${\bf a}_{1,2}$.

Since the spin model has couplings up to the 3rd NN sites, we consider the symmetry-allowed mean-field ansatz with hopping terms up to the 3rd NN, which are shown in Supplemental Material sup . The NN hopping ansatz reduces to the $\pi$-flux $U(1)$ QSL state Lu (2016) in the case of $\phi_{1}=\phi_{2}=\pi/2$ ($\phi_{1},\phi_{2}$ are the phases of the NN hoppings), with a pair of Dirac spinons at half filling for each spin species. The 2nd and 3rd NN hoppings can open up a direct gap at each Dirac cone, leading to a Chern number $C=\pm 1$ of the lower spinon band and the chiral edge states shown in Fig. 4(b). Meanwhile the 3rd NN hoppings can break the degeneracy of two Dirac cones, giving rise to one particle-like SFS around one Dirac point (blue in Figs. 4(c-d)) and a hole-like SFS around the other Dirac point (pink in Figs. 4(c-d)). Due to single-occupancy constraint Eq. (4), the particle-like SFS and hole-like SFS are perfectly compensated at half filling. Choosing mean-field parameters as $\chi=1.0,\phi_{1}=\phi_{2}=\pi/2,\lambda=1.0,\varphi_{1}=\varphi_{2}=\varphi_{3}=0,\rho=3.0,\gamma_{1}=\gamma_{2}=\gamma_{3}=\pi/2$ sup , the mean-field dispersion and edge spectrum of fermionic spinons are shown in Figs. 4(a-b).

For further comparison with DMRG, we follow the YC cylinder geometry with quantized momentum $k_{2}=2\pi l_{2}/L_{y}$ along $\mathbf{b}_{2}$ direction. In Fig. 4(c) and (d) we depict how the 1D channels with quantized momenta $k_{2}=2\pi l_{2}/L_{y}$ intersect with the two SFSs in the reduced Brillouin zone of the staggered flux state. On the YC6 cylinder, as shown in Fig. 4(c), there are $N_{w}=2\times 2=4$ pairs of gapless 1D modes crossing the SFSs (counting both spin species), constraining the central charge to be $c\leq N_{w}-1=3$. On the YC8 cylinder, as shown in Fig. 4(d), there are $N_{w}=2\times 4=8$ pairs of gapless 1D modes crossing the SFSs, restricting the central charge as $c\leq N_{w}-1=7$. This is consistent with the observed $c\approx 1$ on YC6-24 cylinder (Fig. 3(b)) and $c\approx 5$ on YC8-16 cylinder (Fig. 3(c)). Note that the number $N_{w}-1$ only bound the actual central charge from above, since symmetric backscatterings between these gapless 1D channels can further reduce the total central charge from $N_{w}-1$ sup .

Discussion. The spin structure factor of the gapless CSL phase in Fig. 1(b) resembles that of the $U(1)$ Dirac spin liquid Iqbal et al. (2016). Specifically, it exhibits high intensities on the edge and at the corner of the hexagonal Brillouin zone, which are associated with fermion bilinears and monopoles respectively in the $U(1)$ Dirac spin liquid Song et al. (2018a, b). This suggests the proximity of the gapless CSL to the $U(1)$ Dirac state, which is indeed the case for the proposed staggered flux state. We have also studied the phase transition from the $J_{1}$-$J_{2}$ SL phase to the gapless CSL phase. The ground-state energy versus couplings is very smooth, suggesting a possible continuous phase transition sup . Interestingly, in the $J_{1}$-$J_{2}$ SL entanglement entropy also shows a logarithmic correction of the area law, which leads to a finite central charge sup . A new insight for its ground state could be a gapless spin liquid with SFSs but preserving TRS, which we leave for future work.

Summary. We have studied the spin-$1/2$ $J_{1}$-$J_{2}$-$J_{3}$ THM by extensive DMRG calculations. We identify a CSL state spontaneously breaking TRS, featuring a chiral edge mode and spin pumping upon flux insertion. The vanishing spin triplet gap and finite central charge reveal the gapless nature of this state. The central charge which grows with system circumference further indicates emergent SFSs. While the competing $J_{2},J_{3}$ couplings lead to a gapped CSL on kagome lattice He et al. (2014); Gong et al. (2015), they induce a gapless CSL on triangular lattice. On the mean-field level we propose a staggered flux state driven by the $J_{2},J_{3}$ couplings, which breaks TRS and forms SFS, providing a theoretical understanding for such a novel gapless phase. The discovery of this gapless CSL reveals the novel possibility for the coexistence of chiral edge modes and SFSs in a gapless QSL, emergent from competing interactions in a frustrated two-dimensional magnet.

Note added. After completion of this work, we became aware of a work by Shijie Hu et al. Hu et al. (2019), who studied the $J_{1}$-$J_{2}$ spin liquid. Compared to their work, our work focuses on the $J_{1}$-$J_{2}$-$J_{3}$ model and found a staggered flux state driven by further-neighbor interactions.

**Supplemental Material **

I Ground-state energy

We first demonstrate the bulk energy on the YC cylinder in the gapless chiral spin liquid phase by choosing the parameter point at $J_{2}=0.3,J_{3}=0.15$ as a representative. We show the energy versus cylinder circumference for $L_{y}=5-11$ in Fig. 1(a). On the YC5 cylinder, the ground state is a dimer state that breaks lattice translational symmetry. On the other YC cylinders, the bond energy is translational invariant in the bulk of cylinder. One may find that the energy shows an even-odd oscillation versus $L_{y}$ but changes small, approaching to the same value with increasing circumference. By a rough estimation, we propose the ground-state energy in the thermodynamic limit near $-0.4685$ for this parameter point.

II Zigzag order and incommensurate order

We have shown the quantum phase diagram of the model in the main text, which has different magnetic order phases. Besides the $120^{0}$ order and the stripe order which have been found in the $J_{1}-J_{2}$ triangular Heisenberg model, we find another two new orders. One order is called as zigzag order in our paper, as shown in Fig. 2(a). We have confirmed this order on both YC and XC cylinders. Here we just show the results on the XC cylinder. The spins along the zigzag path in each column are parallel to each other.

Another order seems to be an incommensurate magnetic order. As shown in Fig. 2(b), the magnitudes of the magnetic moments have oscillation between large and small values. The overall pattern of spin correlations seems compatible with the zigzag order but with oscillation, which suggests this order as the zigzag order with incommensurate modulation. On the small XC cylinders such as the XC6 and XC8 cylinders, the spin correlations decay fast, and the order becomes stable on the larger XC10 cylinder. On the YC cylinder, this incommensurate order is also found.

III Absent spin and dimer order in the gapless chiral spin liquid state

In Fig. 3, we show spin correlation function and structure factor of the gapless chiral spin liquid state. We show the results for $J_{2}=0.3,J_{3}=0.15$. As shown in Fig. 3(a), one can find that spin correlation functions decay with distance like an exponential way, showing the absent magnetic order. We define spin structure factor as

[TABLE]

In Fig. 3(b), we show spin structure factor $S({\bf k})$ on the YC9-18 cylinder, which shows broad peaks near the corners of the Brillouin zone. To explicitly show spin structure factor, we plot $S({\bf k})$ along the momentum lines in the first Brillouin zone in Fig. 3(d). The momentum lines are denoted by the dashed lines in Fig. 3(c). The $x$ label of Fig. 3(d) $l_{2}L_{x}+l_{1}$ denotes that the momentum points follow the sequence from left to right and from bottom to top in the first Brillouin zone in Fig. 3(c). The momentum positions of the small peaks are marked as blue dots in Fig. 3(c). One can find that there is no sharp peak of $S({\bf k})$, consistent with absent magnetic order. While the dominant peaks appear at the corners of the Brillouin zone, some smaller peaks appear along one edge of the Brillouin zone. On the other YC cylinders such as YC8 and YC10, spin structure factors are similar to that of the YC9 cylinder. The dominant peaks also appear at the corners of the Brillouin zone, but the smaller peaks may not exactly locate along the edge of Brillouin zone but have small deviation duo to the geometric size.

To study dimer order in the spin liquid phase, we calculate dimer correlation functions and their structure factors. Since the nearest-neighbor bonds on the triangular lattice are along three different directions as shown in the inset of Fig. 4, we calculate their dimer structure factors separately. The dimer structure factors $D_{l}({\bf k})$ $(l=1,2,3)$ can be defined as

[TABLE]

where $\mathcal{B}_{i,i+\delta_{l}}={\bf S}_{i}\cdot{\bf S}_{i+\delta_{l}}$ with $\delta_{l}$ denoting the three bond directions $\delta_{1}=(1,0)$, $\delta_{2}=(\frac{1}{2},\frac{\sqrt{3}}{2})$, $\delta_{3}=(-\frac{1}{2},\frac{\sqrt{3}}{2})$. The obtained dimer structure factors on different YC cylinders such as YC8, YC9, and YC10 are similar. Here we show the results on the YC8 cylinder in Fig. 4 as an example. The dimer structure factors along the three directions are all featureless, showing the vanished dimer order in the spin liquid phase.

III.1 Size scaling of spin triplet gap

In the main text, we have shown the $1/L_{y}$ size scaling of the spin triplet gap obtained on the square-like clusters for $J_{2}=0.3,J_{3}=0.15$. Here we demonstrate the $1/N_{x}$ size scaling of the gap on the YC6 and YC8 cylinders. We calculate the ground state on a long cylinder and then sweep the middle $N_{x}\times L_{y}$ sites in the bulk of cylinder by targeting the total spin $S=1$ sector. The obtained spin gap versus $1/N_{x}$ is shown in Fig. 5. We find that the gap depends on the system length and seems to decrease with $N_{x}$, which is consistent with the 1D gapless feature. To avoid the 1D physics, we use the square-like clusters for $1/L_{y}$ size scaling to estimate the spin triplet gap in the thermodynamic limit.

III.2 Entanglement spectrum and chiral edge mode

Next, we study the entanglement spectrum in the gapless chiral spin liquid state. In chiral topological states, it was conjectured that there is a one-to-one correspondence between entanglement spectrum and physical edge spectrum Li and Haldane (2008). In Fig. 6, we show the entanglement spectrum on the YC8 cylinder, which is labeled by the momentum quantum number along the ${\bf a}_{2}$ direction $\Delta k_{y}$ in different total $S^{z}$ sectors. The spectrum shows quasi-degenerate groups of levels with counting $\{1,1,2,3,5,\cdots\}$, which agrees with the chiral SU(2)1 Wess-Zumino-Witten conformal field theory Francesco et al. (2012).

The spectrum result seems to suggest a gapless chiral spin liquid state with chiral edge mode. Such chiral states with chiral edge mode and gapless bulk properties have been studied recently for non-interacting systems Dubail and Read (2015), such as the non-interacting $p+ip$ topological chiral state. However, for interacting systems, such example is rare. This gapless chiral spin liquid state may be the first example to realize such a novel quantum state for interacting systems, which should be further studied in future work.

III.3 Quantum phase transition from the $J_{1}-J_{2}$ spin liquid phase to the gapless chiral spin liquid phase

Since this gapless chiral spin liquid phase emerges at the neighbor of the $J_{1}-J_{2}$ spin liquid phase, it would be very interesting to study how the system evolves from one state to the other one. To study this phase transition, we choose a parameter line that crosses the typical parameter points in the two phases that are $J_{2}=0.125,J_{3}=0$ and $J_{2}=0.3,J_{3}=0.15$. Thus the parameter line is defined as $J_{3}/J_{1}=(6/7)(J_{2}/J_{1})-3/28$, which is shown in Fig. 7(a) as the blue solid line. First of all, we calculate the scalar chiral order in the bulk of the YC8-36 and the YC10-36 cylinders along this parameter line as shown in Fig. 7(b). While the chiral order decreases with circumference for $J_{2}<0.3$, it becomes stable with growing circumference at $J_{2}=0.3,J_{3}=0.15$, showing the absent and emergent scalar chiral order in either of the two phases. Furthermore, we calculate the ground-state energies in the bulk of cylinder along the parameter line as shown in Fig. 7(c), which change very small from the YC8 to the YC10 cylinder. With growing couplings, the energy seems very smooth, suggesting a possible continuous phase transition between the two phases although a weak first-order transition cannot be excluded.

Next we reexamine the entanglement entropy in the $J_{1}-J_{2}$ spin liquid state. A typical point $J_{2}=0.125,J_{3}=0$ is shown in Fig. 7(d) on the YC8-16 cylinder. With growing number of kept states, the entropy keeps to increase. We can see that the entropy also has a logarithmic correction of the area law. The data by keeping $6000$ $SU(2)$ states could be fitted quite well by the formula $S(l_{x})=(c/6)\ln[(L_{x}/\pi)\sin(l_{x}\pi/L_{x})]+g$ Calabrese and Cardy (2004), where $S(l_{x})$ is the bipartite entanglement entropy, $c$ is the central charge, and $g$ is a nonuniversal constant. The fitting gives a finite central charge $c\sim 5$, which is consistent with the central charge of the point $J_{2}=0.3,J_{3}=0.15$ (in the gapless chiral spin liquid phase) on the same system size. The similar behavior of entropy suggests possible spinon Fermi surface in the $J_{1}-J_{2}$ spin liquid state. A new insight for the ground-state of the $J_{1}-J_{2}$ spin liquid could be a gapless spin liquid with small spinon Fermi surface, which however preserves time-reversal symmetry. In the presence of the third-neighbor $J_{3}$ coupling, time-reversal symmetry is broken and the system evolves to the gapless chiral spin liquid phase. This new insight needs more studies in future work.

IV Details in the mean-field analysis for the staggered flux state

The ground state is at the half-filling of fermionic spinons. We know that the U(1) Dirac state does not break translational symmetry after projection. However, on mean-field level, the unit cell has to be doubled due to magnetic translation PhysRevLett.98.117205. Now therefore, we choose to label the lattice as ${\bf a}_{1}^{\prime}=(1,0),{\bf a}_{2}^{\prime}=(1,\sqrt{3}),{\bf r}=r_{1}{\bf a}_{1}^{\prime}+r_{2}{\bf a}_{2}^{\prime},r_{1},r_{2}\in{\bf Z}$. The dual reciprocal space are labeled as $\mathbf{b}_{1}^{\prime}=(1,-\sqrt{3}/3)$, $\mathbf{b}_{2}^{\prime}=(0,\sqrt{3}/3),\mathbf{k}=k_{1}\mathbf{b}_{1}^{\prime}+k_{2}\mathbf{b}_{2}^{\prime}$, $k_{1}=2\pi l_{1}/N_{1}$, $k_{2}=2\pi l_{2}/N_{2}$, $l_{1},l_{2}\in{\bf Z}$, where $N_{1},N_{2}$ are the numbers of the unit cells along ${\bf a}_{1}^{\prime}$- and ${\bf a}_{2}^{\prime}$-direction, respectively. The first Brillouin zone (FBZ) can be identified as $k_{1},k_{2}\in[-\pi,\pi)$. We also define the three NN bridge vectors in Descartes coordinate system:

[TABLE]

In the next we only consider one of the two branches of U(1) spin liquid fermionic spinons, which are spinless. First of all, we consider the mean-field Hamiltonian with the NN hopping terms according to Eq. (5) of the main text, which is given as

[TABLE]

We set the NN hopping amplitude $\chi=1.0$ as unit. In reciprocal lattice with the Fourier transformation on the torus

[TABLE]

we have

[TABLE]

where $\eta_{\mathbf{k}}=(f_{\mathbf{k}+},f_{\mathbf{k}-})^{T}$ and

[TABLE]

We may assume real and equal amplitude $\lambda$ in the 2nd NN hopping, thus

[TABLE]

In reciprocal space similarly,

[TABLE]

with

[TABLE]

We may also assume real and equal amplitude $\rho$ for the 3rd NN couplings, thus

[TABLE]

In the reciprocal space we have

[TABLE]

with

[TABLE]

Finally, $H=J_{1}H_{1}+J_{2}H_{2}+J_{3}H_{3}$, where $J_{1,2,3}$ follows as same as the parameters in DMRG simulation. The 2nd and 3rd NN hoppings can open up a direct gap at each Dirac cone, leading to a Chern number $C=\pm 1$ of the lower spinon band and the chiral edge states. Meanwhile the 3rd NN hoppings can break the degeneracy of two Dirac cones, giving rise to one particle-like spinon Fermi surface (SFS) around one Dirac point and a hole-like SFS around the other Dirac point. Due to single-occupancy constraint, the particle-like SFS and hole-like SFS are perfectly compensated at half filling.

V Central charge on the YC cylinders

Below we discuss the central charge of the staggered flux state on YC cylinders, as illustrated in Fig. 4 (c)-(d) of the main text. Consider the 1D channels of quantized momenta (dashed orange lines in Fig. 4 of the main text) cross each fermi pockets $2n$ times for each spin species, so that in total there are $N_{w}=4n$ pairs of right and left movers among the gapless 1d modes. They can be labeled by their spin index $\sigma=\uparrow,\downarrow$, particle/hole pocket index $\mu=\pm 1$, chirality index $\chi=\pm 1$ for right and left movers, and flavor index $1\leq\rho\leq n$. Each gapless 1d mode of fermionic spinons on the YC cylinder can be bosonized as

[TABLE]

where $\kappa$ represents the Klein factors. Due to “staggered” magnetic translation symmetry (Eq. (5) of the main text), the momenta of the 1D modes on particle and hole pockets are related by

[TABLE]

where ${\bf b}_{1,2}$ are reciprocal lattice vectors dual to Bravais lattice vectors ${\bf a}_{1,2}$. And the boson fields transform under translations as

[TABLE]

Due to the single-occupancy constraint of the spinons on each site, the total density fluctuation of the 1D gapless modes is forbidden:

[TABLE]

This is a consequence of the $U(1)$ gauge field coupled to the fermionic spinons. Therefore the low effective theory of the mean-field spinon Fermi surfaces coupled to a $U(1)$ gauge field, when placed on a finite YC cylinder, can be written as

[TABLE]

where $v_{\chi,\rho}>0$ is the fermi velocity of each 1D mode $\phi_{\chi,\mu,\rho,\sigma}$. The cosine term $\propto g_{0}$ plays the role of $U(1)$ gauge field, removing the total density fluctuation Eq. (22) from the low-energy physics. Clearly the effective action preserves both translations $T_{1,2}$ in Eq. (21) and $SU(2)$ spin rotational symmetry.

Notice that the $U(1)$ gauge fluctuations (i.e. the $g_{0}$ term in effective action) reduces the total central charge $c$ of the YC cylinder from $N_{w}=4n$ to $N_{w}-1$, gapping out the total density mode. Also, we have not included any interaction between fermionic spinons in the above effective 1d action, which can further reduce the central charge to an integer smaller than $N_{w}-1$. Therefore we have shown that the central charge for the mean-field spinon Fermi surfaces in Fig. 4 (c)-(d) of the main text, is upper bounded by $c\leq N_{w}-1=4n-1$. This means $c\leq 3$ on YC6 cylinder ($n=1,N_{w}=4$) in Fig. 4(c), and $c\leq 7$ on YC8 cylinder ($n=2,N_{w}=8$) in Fig. 4(d) of the main text.