# Magnetic phases of the anisotropic triangular lattice Hubbard model

Yang Yu  Department of Physics, University of Michigan, Ann Arbor, MI 48109,
USA  Shaozhi Li  Materials Science and Technology Division, Oak Ridge National
Laboratory, Oak Ridge, TN 37831, USA  Sergei Iskakov  Department of Physics,
University of Michigan, Ann Arbor, MI 48109, USA  Emanuel Gull  Department of
Physics, University of Michigan, Ann Arbor, MI 48109, USA

###### Abstract

The Hubbard model on an anisotropic triangular lattice in two dimensions, a
fundamental model for frustrated electron physics, displays a wide variety of
phases and phase transitions. This work investigates the model using the
ladder dual fermion approximation which captures local correlations non-
perturbatively but approximates non-local correlations. We find metallic, one-
dimensional antiferromagnetic, non-collinear antiferromagnetic, square-lattice
antiferromagnetic, and spiral phases but no evidence of collinear
antiferromagnetic order in different parts of the phase diagram. Analyzing the
spin susceptibility in detail, we see both regions of agreement and of
discrepancy with previous work. The case of Cs2CuCl4 is discussed in detail.

##  I Introduction

Correlation physics on triangular lattices is a subject of intense interest.
Many triangular lattice magnets including Cs2CuCl4 [1], $\kappa$-(BEDT-
TTF)2Cu2(CN)3 [2], EtMe3Sb[Pd(dmit)2]2 [3], NiGa2S4 [4], Ba3CuSb2O9 [5],
Ba3CoSb2O9 [6], YbMgGaO4 [7], and NaYbO2 [8] have been experimentally
investigated. Transition metal dichalcogenide moiré superlattices [9, 10, 11]
and triangular optical lattices [12, 13] have also recently been proposed as
platforms to study correlation physics on the triangular lattice.
Theoretically, both triangular lattice spin models and the triangular lattice
Hubbard model have been investigated to describe frustrated quantum magnets.
Nevertheless, the phase diagrams of triangular lattice systems are still under
debate both experimentally [14, 15, 16] and theoretically [17, 18]. Quantum
spin liquid (QSL) states [19, 20, 21, 22, 23] and Mott transitions [24, 25,
26, 27] proposed in triangular lattice systems therefore remain elusive.

The case of spatially anisotropic hopping/exchange interaction is particularly
interesting. Previous studies show inconsistencies for both the anisotropic
triangular lattice Heisenberg and Hubbard model. Among those is the existence
of collinear antiferromagnetic order in the limit of weakly coupled one-
dimensional chains [28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41,
42, 43, 44, 45, 46, 47, 48].

Here we present a dual fermion (DF) study of the anisotropic triangular
lattice Hubbard model. The DF method, as one of the diagrammatic extensions of
the dynamical mean field theory (DMFT), is unbiased (in the sense that it does
not favor certain types of magnetic orders) and provides high-resolution spin
susceptibilities, which renders it a valuable tool in discussing the above-
mentioned problem. We survey the phase diagram, connect it to ground state
simulations, and discuss in detail the parameter believed to be relevant for
Cs2CuCl4.

The remainder of this paper proceeds as follows. In Sec. II, we give a brief
introduction to the methodology by describing the model, the dual fermion
method, the observables of interest, and the possible magnetic orders. In Sec.
III, we present our phase diagram, discuss the phases and phase transitions,
and apply the method to the concrete example of Cs2CuCl4. Conclusions are
drawn in Sec. IV.

##  II Method

###  II.1 Hubbard model

We study the Hubbard model [49, 50]

$\hat{H}=-\sum_{\langle
ij\rangle\sigma}t_{ij}\hat{c}_{i\sigma}^{\dagger}\hat{c}_{j\sigma}+\text{H.c.}+U\sum_{i}\hat{n}_{i\uparrow}\hat{n}_{i\downarrow}-\mu\sum_{i\sigma}\hat{n}_{i\sigma}$
(1)

on an anisotropic triangular lattice in two dimensions at half filling. Here
$\hat{c}_{i\sigma}$ ($\hat{c}_{i\sigma}^{\dagger}$) are fermion annihilation
(creation) operators on site $i$ with spin $\sigma=\uparrow,\downarrow$.
$\hat{n}_{i\sigma}=\hat{c}_{i\sigma}^{\dagger}\hat{c}_{i\sigma}$ is the
particle number operator. $\langle ij\rangle$ denotes nearest-neighbor pairs
on the triangular lattice. $U$ is the on-site interaction. The anisotropic
hopping is set to be $t_{ij}=t$ for two of the three directions and
$t_{ij}=t^{\prime}$ for the third direction as illustrated in Fig. 1(c). When
$t^{\prime}\gg t$, the system approaches the limit of decoupled one-
dimensional chains, and when $t^{\prime}\ll t$, the system becomes a square
lattice. The chemical potential $\mu$ is chosen such that the system is close
to half filling (we adjust $\mu$ such that
$|\sum_{\sigma}\langle\hat{n}_{i\sigma}\rangle-1|<0.01$). All the data in this
paper are calculated at fixed temperature $T=0.1t$ ($\beta=1/T=10/t$) unless
stated otherwise. In the remainder of the paper, the hopping $t$ is set to $1$
and used as an overall energy scale.

###  II.2 Dual fermion expansion and ladder dual fermion approximation

The DF expansion [51] is a diagrammatic technique developed to describe non-
local correlation effects in models with local interactions. It can be
understood as an expansion around a local limit in terms of ‘dual’ diagrams.
At the lowest order, only local correlations are captured, but the method
converges to the exact limit when all ‘dual’ diagrams to all orders are
summed. In contrast to other diagrammatic techniques, it provides both a non-
perturbative strongly correlated local self-energy and non-local contributions
on a continuous momentum grid, as well as direct access to one- and two-
particle correlation functions. A detailed discussion of the technique is
presented in Refs. [51, 52, 53, 54] and in a review [55].

In this paper, we employ the ladder dual fermion approximation (LDFA) with a
dynamical mean field theory [56] starting point (see Refs. [57, 58, 59, 60,
61, 62] for previous DF work on the triangular lattice Hubbard model). This
approximation only contains dual particle-hole ladder diagrams, which are
summed to infinite order. All non-ladder diagrams and all higher-order
vertices are neglected. The quality of the approximation is assessed in Refs.
[63, 64, 65, 66, 67] (see also [68]). We use the implementation provided by
the open-source OpenDF code [69] based on ALPSCore [70, 71].

As a starting point, the self-consistent vertices and self-energies of a DMFT
impurity problem are used as the input of the LDFA [72]. These quantities are
obtained using numerically exact continuous-time auxiliary field quantum Monte
Carlo [73, 74, 75, 76] (for analytic expressions of the non-interacting
triangular lattice local density of states see Appendix A). The LDFA method
proceeds by solving the Bethe-Salpeter equation (BSE) to obtain dual vertex
functions. Once the dual vertex functions are obtained, the Schwinger Dyson
equation (SDE) and the Dyson equation are used to obtain dual self-energies,
which are used to update the dual Green’s functions. The BSE, SDE, and the
Dyson equation are iterated until convergence [69]. This iterative procedure
method may encounter convergence issues where the DMFT starting point is far
from the exact solution, in particular in the vicinity of phase transitions
and at low temperatures.

LDFA results, such as one-body Green’s function $G$ and spin susceptibility
$\chi_{\text{s}}$, are obtained as a function of Matsubara frequencies.
Dynamical quantities, such as spectral functions and dynamical spin
susceptibilities, require analytic continuation from the Matsubara to the real
axis. This procedure is ill-conditioned and uncontrolled [77]. In this paper,
where possible, we therefore present quantities that can be directly extracted
from the Matsubara results. Where we show analytically continued quantities,
we use the maximum entropy method implemented in the open-source MaxEnt code
[78] based on ALPSCore [70, 71].

Despite being approximate, the LDFA has advantages that make it a unique tool
for studying correlated phases. First, since it contains strong local
correlations non-perturbatively, it can be used to study strongly correlated
systems outside the perturbative regime. Second, since it gives access to
generalized susceptibilities with almost continuous momentum resolution (we
use systems of size $24\times 24$ and $48\times 48$ in this work), it can be
used to detect the emergence of phases with long-range and incommensurate
orders. Finally, since the formalism is based on detecting large fluctuations,
rather than entering an ordered state, all possible orders are treated on
equal footing without an a priori assumption of an expected phase.

###  II.3 Spin structure factor

The main quantities of interest in this paper are the dynamical spin
susceptibilities and the dynamical spin structure factor. The magnetic
scattering cross section is directly related to the dynamical spin structure
factor,

$\mathcal{S}^{\alpha\alpha}(\bm{q},\omega)=\frac{1}{N}\int_{-\infty}^{\infty}\frac{\mathrm{d}t}{2\pi}\mathrm{e}^{\mathrm{i}\omega
t}\left\langle\hat{S}_{\bm{q}}^{\alpha}(t)\hat{S}_{-\bm{q}}^{\alpha}(0)\right\rangle,$
(2)

where
$\hat{S}^{\alpha}_{\bm{q}}=\sum_{i}\mathrm{e}^{-\mathrm{i}\bm{q}\cdot\bm{r}_{i}}\hat{S}^{\alpha}_{i}$
denotes the Fourier transform of the spin operator with $\alpha=x,y,z$, and
$N$ is the number of lattice sites. The total scattering cross section that
integrates over all frequencies is related to the equal-time spin structure
factor,

$\mathcal{S}^{\alpha\alpha}_{0}(\bm{q})=\frac{1}{N}\left\langle\hat{S}_{\bm{q}}^{\alpha}(0)\hat{S}_{-\bm{q}}^{\alpha}(0)\right\rangle.$
(3)

Our calculation is performed in the normal state, where SU(2) symmetry is
preserved, which allow us to focus on the $z$ component. The frequency- and
momentum-dependent dynamical spin susceptibility,

$\chi_{\text{s}}(\bm{q},\omega)=-\frac{\mathrm{i}}{N}\int_{0}^{\infty}\mathrm{d}t\mathrm{e}^{\mathrm{i}(\omega+\mathrm{i}0^{+})t}\left\langle\left[\hat{S}^{z}_{\bm{q}}(t),\hat{S}^{z}_{-\bm{q}}(0)\right]\right\rangle,$
(4)

which describes the response to a weak externally applied magnetic
perturbation, can be extracted from the spin susceptibility on the Matsubara
axis,

$\chi_{\text{s}}(\bm{q},\omega_{n})=-\frac{1}{N}\int_{0}^{\beta}\mathrm{d}\tau\mathrm{e}^{\mathrm{i}\omega_{n}\tau}\langle
T_{\tau}\hat{S}^{z}_{\bm{q}}(\tau)\hat{S}^{z}_{-\bm{q}}(0)\rangle,$ (5)

via analytic continuation. (Note a minus sign difference with respective to
the physical spin susceptibilities.) The dynamical spin structure factor can
then be written as

$\mathcal{S}^{zz}(\bm{q},\omega)=-\frac{1}{\pi}\left[1+n_{B}(\omega)\right]\operatorname{Im}\chi_{\text{s}}(\bm{q},\omega)$
(6)

via the fluctuation-dissipation relation [79, 80], where
$n_{B}(\omega)=1/\left(\mathrm{e}^{\beta\omega}-1\right)$ denotes the Bose-
Einstein distribution. The equal-time spin structure factor can be calculated
from a sum over the Matsubara susceptibility as

$\mathcal{S}^{zz}_{0}(\bm{q})=-\frac{1}{\beta}\sum_{\omega_{n}}\chi_{\text{s}}(\bm{q},\omega_{n}).$
(7)

###  II.4 Magnetic orders

Figure 1: Ordering wave vectors for magnetic orders on the triangular lattice
are represented by the different symbols or lines in (a) where the $\Gamma$
point is set to be the origin and the hexagon indicates the first Brillouin
zone with side length equaling $4\pi/3$. The corresponding spin configurations
in the real space are illustrated in: (c) square-lattice antiferromagnetic
(SAFM), (d) and (e) collinear antiferromagnetic (CAFM1, CAFM2), (f) spiral,
(g) non-collinear antiferromagnetic (NCAFM), and (h) one-dimensional
antiferromagnetic (1DAFM) orders. The ground state of the classical Heisenberg
model is illustrated in (b) where the relative angle between the spin on each
site and the spin on the left-bottom site is indicated.

Upon approaching a magnetically ordered state, a divergence in the equal-time
spin structure factor with the ordering wave vectors indicates the onset of
the corresponding magnetic order. Since the LDFA calculation performed here
always remains in the normal state, large well-defined peaks in the equal-time
spin structure factor indicate strong magnetic correlations. These are viewed
as a precursor of magnetic order. To compare our results with previous studies
of related models at zero temperature, we identify strong magnetic
correlations at the lowest temperature studied ($T=0.1$) with the
corresponding magnetic orders at zero temperature in this work.

For the symmetry-broken magnetically ordered states, we can directly evaluate
the ordering wave vectors from the definition in Eq. (3). Figs. 1(c)- 1(h)
show the real-space spin configurations for the magnetic orders discussed in
this paper. The corresponding ordering wave vectors and equivalent
points/lines are indicated by the different symbols in Fig. 1(a). Square-
lattice antiferromagnetic (SAFM) order (Fig. 1(c)) with wave vector at
M′$(0,2\pi/\sqrt{3})$ and two degenerate collinear antiferromagnetic (CAFM1,
CAFM2) orders (Figs. 1(d) and 1(e)) with wave vectors at M$(\pi,\pi/\sqrt{3})$
and M′′$(\pi,-\pi/\sqrt{3})$, have no difference in spin configurations except
for the direction of the ferromagnetically ordered chains. The spiral order
with wave vector at K$(4\pi/3,0)$ and the non-collinear antiferromagnetic
(NCAFM) order with wave vector at X$(\pi,0)$ are shown in Figs. 1(f) and 1(g)
respectively. Additionally, the one-dimensional antiferromagnetic (1DAFM)
order is shown in Fig. 1(h). In the 1DAFM order, no correlation exists between
the horizontal antiferromagnetic chains, resulting in ordering wave vectors
uniformly spread over the line connecting M and M′. The values of equal-time
spin structure factor on this line for the 1DAFM order are expected to be much
smaller than the peak values for the 2D magnetically ordered states mentioned
above [32].

##  III Results

###  III.1 Phase diagram

Figure 2: Phase diagram of the anisotropic triangular lattice Hubbard model.
Green, cyan, red, orange, and yellow regions correspond to metal, one-
dimensional antiferromagnetic (1DAFM), non-collinear antiferromagnetic
(NCAFM), spiral, and square-lattice antiferromagnetic (SAFM) orders
respectively. Plus markers: places with converged results. Square markers:
LDFA fails to converge. Star markers: data points shown in Sec III.3. All the
points are assigned to the closest commensurate orders or the metallic phase.

The phase diagram of the anisotropic triangular lattice Hubbard model, as
identified by LDFA, is shown in Fig. 2. At $4\leq U\leq 6$, the metallic phase
dominates around the isotropic lattice limit and transforms into the SAFM
order when $t^{\prime}$ is decreased and into the 1DAFM order when
$t^{\prime}$ is increased. At $6\leq U\leq 10$, the metal-insulator
transitions from the metallic phase to SAFM, spiral, and NCAFM orders emerge
at different $t^{\prime}$ due to the increase of $U$. A transition from 1DAFM
to CAFM also begins in this region with larger $t^{\prime}$. At $10\leq U\leq
12$, different magnetic phases appear in the order of SAFM, spiral, NCAFM, and
1DAFM orders upon the increase of $t^{\prime}$ going from the square lattice
limit to the isotropic lattice limit to the limit of decoupled chains.

A $24\times 24$ lattice is used to scan this phase diagram. The physical
quantities obtained from the $24\times 24$ lattice are found to have no
qualitative difference from the ones obtained from a $48\times 48$ lattice
which have higher resolution and are therefore shown in the remainder of the
paper. To identify the different phases in Fig. 2, we use two quantities that
reflect the properties of the system in the charge sector and in the spin
sector respectively: the spectral weight at the Fermi energy,
$A(\nu=0)=-\operatorname{Im}G(\mathrm{i}0^{+})/\pi$, and the static spin
susceptibility, $\chi_{\text{s}}(\bm{q},\omega_{n}=0)$.

We first identify the 2D magnetically ordered states according to the absolute
values of the static spin susceptibility at the possible magnetic wave vectors
mentioned in Sec II.4, i.e., the K, M, M′, M′′, X points shown in Fig. 1(a).
If the largest absolute value at those points is larger than $1$, then we
identify the system as belonging to the corresponding magnetically ordered
state. As we will show in Sec. III.3, incommensurate magnetic orders where the
ordering wave vectors do not belong to those points, and a coexisting phase
where the features for the 1DAFM and NCAFM orders coexist, are also observed.
For simplicity, we attribute them as the closest commensurate magnetic orders
in Fig. 2.

To further classify the remainder of the phase diagram, we define that a point
where the spectral weight at the Fermi energy is larger than 0.1 lies in a
metallic phase, and otherwise lies in a phase for the 1DAFM order. This way of
identification of the metallic phase is consistent with the above method used
for the 2D magnetically ordered states. Note that using other threshold values
for the static spin susceptibility and the spectral weight instead of $1$ and
$0.1$ will generate slightly different phase diagrams but the overall shape of
the phase diagram remains similar for those values.

Although our phase diagram does not give precise phase boundaries, our results
can reflect the rich phases of the anisotropic triangular lattice systems. A
precise determination of the location of the boundaries of phases in the
ground state would need additional careful studies of the temperature
dependence and of corrections to LDFA which is beyond our current focus. Note
that we use the terminology “1DAFM order/state” in this work to refer to the
state that has a straight line feature connecting M and M′′ in the equal-time
spin structure factor, such as the case of the uncorrelated
antiferromagnetically ordered chains illustrated in Sec. II D. In fact, the
quantum state with a such 1D feature found by the DF method is the previously
identified 1D QSL which is not a magnetically ordered state [28, 29, 30, 31,
32, 33, 34, 35, 36, 37, 44, 45, 46].

###  III.2 Discussion of previous studies

Before discussing the details of our results, we review previous numerical
studies on related systems. Both the Hubbard model and the Heisenberg model
are known to have a rich ground state phase diagram on the anisotropic
triangular lattice. In the discussion of the Heisenberg model, we will focus
on $J-J^{\prime}$ model where $J$ is used to denote the nearest-neighbor
exchange interaction along two directions and $J^{\prime}$ the one along the
remaining direction. We will not discuss other possible interactions in the
spin model such as in-plane isotropy, next-nearest-neighbor exchange, and ring
exchange that may or may not be related to the Hubbard model at intermediate
$U$ [81].

In the classical ground state solution of the anisotropic triangular lattice
Heisenberg model, the next nearest neighbor chains in the vertical direction
are ferromagnetically correlated. The ordering wave vector is given by
$(2\theta,0)$ where $\theta=\cos^{-1}[-J/(2J^{\prime})]$ for $0\leq
J/J^{\prime}\leq 2$ and $\theta=\pi$ for $J/J^{\prime}>2$ [82, 83]. The real
space spin configuration is illustrated in Fig. 1(b), where $\theta$ is the
angle between two nearest neighbor spins in the two diagonal directions while
$2\theta$ is the angle between two nearest neighbor spins in the horizontal
direction. When $J/J^{\prime}$ increases from $0$ to $1$ to infinity, the
system undergoes the 1D to 2D transition from NCAFM to spiral to SAFM orders
with the ordering wave vector changing smoothly from X to K to M′. (Note that
one of the equivalent points of M′ is located at $(2\pi,0)$ so there is no
jump in the ordering wave vector during the transition.) For the NCAFM order
in this classical solution, the spins on the nearest neighbor chains in the
vertical direction are orthogonal to each other since $\theta=\pi/2$. Note
that a non-orthogonal configuration gives the same energy in the classical
solution. Further, the peak at X in the equal-time spin structure factor is
only related to the fact that the next nearest neighbor antiferromagnetic
chains are ferromagnetically correlated in the vertical direction and are not
related to the relative configuration between the nearest neighbor chains.
Therefore, when the real-space illustration for the NCAFM order is presented
in Fig. 1(g), we intentionally show a non-orthogonal configuration to leave
some ambiguity.

Unlike the classical case, the quantum solution of the Heisenberg model is
complicated and is still under debate. In the square lattice limit
($J^{\prime}/J=0$) and the isotropic lattice ($J^{\prime}/J=1$), the SAFM
order [84, 85, 86] and the spiral order [87, 88, 89, 90, 91] are still the
ground states, as in the classical problem. The incommensurate orders between
the spiral order and SAFM order are also observed [82, 32, 92, 35, 37].
However, the transition point between the SAFM order and the
spiral/incommensurate orders are found to be shifted above the classical value
$J^{\prime}/J=0.5$ [82, 93, 28, 94, 32, 30, 31, 35, 36, 37, 95]. Furthermore,
a dimerized phase is proposed by series expansion (SE) [82] and a putative QSL
phase is proposed by modified spin wave theory (MSWT) [30, 31], variational
Monte Carlo (VMC) [35], exact diagonalization (ED) [36], and high-temperature
SE [95] in the region between the SAFM order and the spiral order. Going from
the isotropic lattice ($J^{\prime}/J=1$) to the limit of decoupled chains
($J^{\prime}/J=\infty$), the phase diagram becomes more involved, especially
between the limit of decoupled chains and the region where the incommensurate
orders are found [82, 96, 40, 34, 42, 32, 43, 35, 37]. While the Tomonaga-
Luttinger physics in the 1D antiferromagnetic Heisenberg chain has been clear
for some time [97, 98, 99, 100], there is some controversy in the 1D-2D
crossover region. On one hand, a gapless QSL, which corresponds to the 1DAFM
order in our notation, is expected to occur as a remnant of the Tomonaga-
Luttinger liquid in the crossover region and is expected to occupy a large
portion of the phase diagram due to dimensional reduction by frustration and
quantum fluctuation [101]. This gapless QSL, along with a strong 1D feature,
is confirmed by density matrix renormalization group (DMRG) [28], resonating
valence bond mean field theory [29], MSWT [30, 31], functional renormalization
group (FRG) [32], VMC [33, 34, 35], ED [36], and Schwinger boson theory [37].
On the other hand, an unexpected CAFM order, which does not have a classical
counterpart, is purposed to be stabilized in the crossover region by
renormalization group (RG) studies [38, 39] and is later supported by SE [40],
coupled cluster method [41], RG [42], and ED [43]. Besides, DMRG [94] and ED
[43] studies point out that incommensurate correlations may exist in the
entire region going from the square lattice limit to the limit of decoupled
chains. In addition, other QSL states besides the 1D QSL are also proposed in
the crossover region [102, 33, 34, 30, 36].

For the Hubbard model at half-filling, a metallic phase at small $U$ and
intermediate anisotropy is expected to be connected with the isotropic non-
interacting triangular lattice. The possibility of a superconducting phase in
the nearby region is also proposed [103, 104, 105] but remains controversial
[106, 107, 108, 61]. Recent DMRG studies find that a Luther-Emery liquid (LEL)
[109] could be the true ground state in this region [110, 111, 46]. Another
recent DMRG study also points out that the superconducting phase only exists
with large doping at intermediate $U$ [112]. Approaching the square lattice
limit at small $U$, a metal-insulator transition (MIT) is expected since the
SAFM order is known to be the ground state of the square lattice Hubbard model
for all $U$ [50]. Approaching the limit of decoupled chains at small $U$, a
MIT is reported from variational cluster approximation (VCA) [44] and VMC [45]
studies while recent DMRG work also shows the possibility of transiting from
metal or LEL to 1D metal [46]. Besides the two MITs just mentioned, another
interaction-driven MIT around the isotropic lattice limit upon the increasing
of $U$ is heavily investigated [113, 103, 114, 115, 106, 116, 117, 57, 107,
118, 119, 120, 121, 108, 59, 44, 47, 45, 122, 123, 61, 124, 48, 125, 126, 111,
62, 127, 46, 128]. Many studies further propose a QSL or a non-magnetic
insulating (NMI) phase existing between the metallic region at small $U$ and
the ordered phase at large $U$, which makes a Mott transition between the
metallic phase and the QSL or NMI phase possible [113, 103, 115, 106, 116,
119, 58, 59, 61, 123, 125, 126, 111, 46, 127]. Recent studies show evidence of
chiral order in the intermediate $U$ region, leading to a conclusion of chiral
QSL [111, 129, 130, 46, 127, 112, 131, 132]. At large $U$, similar to the
Heisenberg model, the phases at $t^{\prime}/t<1$ are more or less clear with
an observation that the transition line between the SAFM order and the
spiral/incommensurate orders tend to move above the classical value
$t^{\prime}/t=\sqrt{2}/2$ [113, 104, 107, 108, 47, 61, 124, 46]. The proposal
of a QSL or a NMI phase existing between the SAFM and spiral orders is also
aroused by VMC [118, 108], and VCA [47]. For $t^{\prime}/t>1$, the phase
diagram is as unclear as the corresponding part of the Heisenberg model. VCA
[47, 44], VMC [45], cluster dynamical mean-field theory [48], and DMRG [46]
studies find a CAFM order above the spiral/incommensurate orders. The 1DAFM
order, i.e., the 1D QSL phase, is also found in a similar region by VCA [44],
VMC [45], and DMRG [46]. Note that the CAFM order and 1DAFM order are found to
belong to different regions of the phase diagram in Ref. [44, 45] but they
were found in different DMRG simulation configurations in Ref. [46]. The
additional QSL phases in the 1D-2D crossover region besides the 1D QSL are
also found in some configurations in Ref. [46].

###  III.3 Magnetic phases

Our results for the equal-time spin structure factor as a function of
$t^{\prime}$ at small $U$ and large $U$ are given in Fig. 3.

When $t^{\prime}<1$, the SAFM order is observed in a large portion of the
parameter space admitting the competition with the spiral order when
approaching the isotropic lattice limit. As shown in Figs. 3(l) and 3(n), at
relatively large $U$, a peak at M′ is established in equal-time spin structure
factor at $t^{\prime}=0.50$ and is then split into two peaks when an
incommensurate order is formed at $t^{\prime}=0.75$. A similar transition with
less obvious peaks is also observed in the metallic region as shown in Figs.
3(k) and 3(m).

Figure 3: Equal-time spin structure factor $\mathcal{S}^{zz}_{0}(\bm{q})$ at
small $U$ (left column) and large $U$ (right column) with different
$t^{\prime}$ ranging from $0.5$ (bottom row) to $5$ (top row). $U$ and
$t^{\prime}$ are as indicated. The Brillouin zone is shown by the white dashed
lines.

For the isotropic lattice, the peaks at the six corners of the first Brillouin
zone are formed when $U$ is gradually increased, as shown in Figs. 3(i) and
3(j). No obvious Mott insulating phase is observed at the temperature studied,
which agrees with the previous DF results [59] and a recent multimethod study
[127]. For the potential spin liquid in the middle of the MIT, the multimethod
study shows that related roton-like excitations are developed when the
temperature decreases from $T=0.1$ [127], which can be inferred from the
equal-time spin structure factor at the M point. In a similar region, we find
convergence difficulties at half-filling below $T=0.1$. Therefore, the non-
magnetic insulator phase is absent on our phase diagram. A different starting
point for the DF method may alleviate the convergence problem at low
temperature [58].

Approaching the limit of decoupled chains with large $U$, incommensurate order
is observed between $t^{\prime}=1$ and $t^{\prime}\approx 2$ with the peak
value in the equal-time spin structure factor evolving from the K point to the
X point as shown in Figs. 3(j), 3(h), and 3(f). At even larger $t^{\prime}$,
as shown in Figs. 3(b) and 3(d), the peak at X decreases while a background
broadening ranging from M to M′ emerges and finally becomes comparable to the
peak at X. A similar trend is also observed at small $U$ as shown in Figs.
3(i), 3 (g), 3(e), 3(c), and 3(a), where the broadening 1D feature
corresponding to the 1DAFM order appears much early with the increase of
$t^{\prime}$. Despite the occurrence of the 1DAFM order, we note that the peak
at X seems to be stable and may extend to large $t^{\prime}$ for both small
and large $U$.

Our equal-time spin structure factor results at large $U$ clearly show a
smooth transition from the NCAFM order to the 1DAFM order, which is remarkably
similar to the susceptibility or structure factor shown in the FRG study of
the Heisenberg model [32] and the ED study of the XX model [133]. The
unexpected CAFM order proposed in the Heisenberg model [38, 39, 40, 41, 42,
43] and the Hubbard model [47, 44, 45, 48, 46] is not found in our
calculation. Furthermore, the coexisting state, where both the peak at X and
the line going from M to M′ are visible in the equal-time spin structure
factor such as the one shown in Fig. 3(d), is found to connect the NCAFM order
and the 1DAFM order despite that we classify them to the nearest ordered phase
in Fig. 2. To our best knowledge, this coexisting state has not yet been
discussed and could be directly related to the controversy in the large
$t^{\prime}/t$ or large $J^{\prime}/J$ region. We find some literature indeed
providing indirect evidence for this coexisting state. For example, both DMRG
[94] and ED [43] studies for the Heisenberg model support that the
incommensurate correlation exists in the entire region going from the
isotropic lattice to the limit of decoupled chains which agrees with the
observation here that the peak at X tend to exist at very large $t^{\prime}$.
Although the ED study [43] observes a phase transition from the incommensurate
order to the CAFM order instead of the NCAFM order as $J^{\prime}$ increases,
the energy of the NCAFM order is found to be larger than the CAFM order only
by a margin of $\sim 10^{-7}$ at large $J^{\prime}$, and this small energy
difference even decreases with system size. Besides, a recent DMRG study on
the Hubbard model also finds a similar phase like the coexisting phase found
here in the YC6 cylinder [46]. The supplemental materials of this work also
point out that the 1D QSL phase found by the work may be non-trivial because
of the isolated gapless points (probably at $\Gamma$ and X). We note that this
is consistent with the picture of the NCAFM order coexisting with the 1DAFM
order presented here. Furthermore, a recent VMC work on the Hubbard model find
a similar coexisting state of the CAFM order and the 1DAFM order where
additional anisotropy is introduced for the hopping terms in the model such
that CAFM can be stabilized in the nearby regions [134]. In the VMC work, this
coexisting phase is further identified as a new type of QSL. Lastly, a smooth
transition from the 1DAFM order to the SAFM order is also indicated in the
auxiliary field quantum Monte Carlo (AFQMC) study of the anisotropic square
lattice Hubbard model with additional frustration from the next-nearest
neighbor hopping [135].

Besides the 1DAFM order and the coexisting state discussed above, the phases
slightly above and below the isotropic lattice limit at large $U$ may also not
have counterparts in the classical solution of the Heisenberg model. As shown
in Figs. 3(h) and 3(l), the amplitude of the equal-time spin structure factor
decreases in these two regions and the peaks get blurred. We note that it is
possible to relate the characteristics of the equal-time spin structure factor
for the coexisting state and the phases shown in Figs. 3(h) and 3(l) to QSL
phases [34, 133, 134]. However, we do not have other reliable methods to
classify QSL within the framework of LDFA and we cannot obtain convergence for
even larger $U$ to discuss the stability of these phases at the infinite $U$
limit. Therefore, we conservatively classify them as the nearest ordered
phases in the phase diagram, Fig. 2.

###  III.4 Lifshitz transitions

Figure 4: Electron occupation $n(\bm{k})$ at $U=4$ (left column) and $U=6$
(right column) with different $t^{\prime}$ ranging from $0.5$ (bottom row) to
$1.5$ (top row). $U$ and $t^{\prime}$ are as indicated. Black lines show the
non-interacting Fermi surfaces. The Brillouin zone is shown by the white
dashed lines.

For the small $U$ region of the phase diagram, we identify the 1DAFM order
from the metallic phase according to the spectral weight at the Fermi energy
in Sec. III.1. In fact, the same boundary between these two phases is also
related to the Lifshitz transitions which in our case is the Fermi surface
opening at the boundary of the Brillouin zone. Fig. 4 shows the electron
occupation
$n(\bm{k})=\sum_{\sigma,\omega_{n}}\mathrm{e}^{-\mathrm{i}\omega_{n}0^{-}}G_{\sigma}(\mathrm{i}\omega_{n},\bm{k})/\beta$
at $U=4$ and $U=6$ with different $t^{\prime}$ ranging from $0.5$ to $1.5$.
The black line indicates the non-interacting Fermi surface
$\varepsilon(\bm{k})-\mu=0$. We see that at $U=4$ the Fermi surface (white
regions) roughly matches the non-interacting Fermi surface at all $t^{\prime}$
and two Lifshitz transitions inherited from the non-interacting limit occur at
$t^{\prime}$ slightly above $1.5$ and below $0.5$. At $U=6$, with the changing
of $t^{\prime}$ the Fermi surface may have already opened at $t^{\prime}=0.5$
and $t^{\prime}=1.5$, as indicated by the white regions in Figs. 4(b) and
4(f). In the recent DMRG study of the same system, the Fermi surface opening
is also observed at very close regions and the trend that the Fermi surface at
larger $U$ opens at $t^{\prime}$ closer to $1$ also agrees well [46]. Besides,
as shown in Fig. 4, in the presence of interactions the saddle points of
Lifshitz transitions for the non-interacting case become broadened white
regions and at larger $U$ the Fermi surface encloses effectively a larger
region. This phenomenon is well understood in a Landau-type theory [136] as
flat bands near the conventional Lifshitz transition. Similar phenomena are
also found in the isotropic triangular lattice Hubbard model at van Hove
filling [60] and the doped Kondo lattice model [137].

We note that in Ref. [46], although a similar Lifshitz transition is found,
the small $U$ and large $t^{\prime}$ region (left top corner of Fig. 2) is
recognized as 1D metal instead of 1DAFM. In our calculation, lower temperature
data are needed to distinguish the insulating phase and a correlated metal.
This is beyond the scope of the method. We note that the insulating phase is
found in VCA [44] and VMC [45] in similar regions, which agrees with the 1DAFM
order at large $t^{\prime}$ we found. Besides, a similar transition from a 1D
insulating phase to a higher-dimensional metallic phase with antiferromagnetic
spin correlations is reported from the AFQMC study of the anisotropic square
lattice Hubbard model [138]. Therefore, the metallic or insulating nature at
the small $U$ and large $t^{\prime}$ region may need further study.

Figure 5: Equal-time spin structure factor $\mathcal{S}^{zz}_{0}(\bm{q})$ at
$U=11$ and $t^{\prime}=1.715$ with different temperatures $T$ ranging from
$0.1$ (left) to $10$ (right). $T$ is as indicated. The Brillouin zone is shown
by the white dashed lines. Figure 6: Dynamical spin structure factor
$\mathcal{S}^{zz}(\bm{q},\omega)$ at $U=11$ and $t^{\prime}=1.715$ with
different temperatures $T$ ranging from $0.1$ (top row) to $10$ (bottom row).
$T$ is as indicated. Left column, middle column, and right column show
$\mathcal{S}^{zz}(\bm{q},\omega)$ at $q_{y}=0$, $q_{y}=\sqrt{3}\pi$, and
$q_{x}=-\pi/2$, respectively. White lines shows the locations of peaks for
each fixed $\bm{q}$. White open squares are peak values at $T<0.1$K reported
from the neutron scattering experiment with a rescaling by a factor of two
[1]. Values of peaks at $q_{y}=2\pi/\sqrt{3}$ and $q_{y}=4\pi/\sqrt{3}$ are
shown in the right column with white solid dots and numbers.

###  III.5 Spin excitation of $\text{Cs}_{2}\text{CuCl}_{4}$

The anisotrpic triangular lattice Hubbard model at $t^{\prime}=1.715$
(estimated from $J/J^{\prime}=0.34$) is believed to describe the physics of
Cs2CuCl4. We present the spin structure factors and compare them with the
inelastic neutron scattering experiment of Cs2CuCl4 [1]. Our results overall
show good agreement with Ref. [1]. Besides, we find a region where a distorted
1D feature and the feature of an incommensurate order coexist in the equal-
time spin structure factor. This is similar to the coexisting states of the
NCAFM and 1DAFM orders discussed in Sec. III.3.

Fig. 5 shows the equal-time spin structure factor at $U=11$ and
$t^{\prime}=1.715$ with different temperatures. Upon increasing the
temperature from $T=0.1$ to $T=1$, the peaks at incommensurate ordering wave
vector shown in Fig. 5(a) gradually evolve into broadening 1D curves as shown
in Fig. 5(c). At intermediate temperature $T=0.5$ as shown in Fig. 5(b), both
the blurred peaks at the incommensurate ordering wave vector and a distorted
1D feature are observed. Such coexistence actually also appears at $T=0.1$ and
$T=1$ with either the 1D feature or the peak at incommensurate ordering wave
vector much less obvious. At high temperature $T=10$, as shown in Fig. 5(d),
the structure factor becomes almost featureless considering the minor
difference between the maximum and the minimum. The finite temperature
transition presented above is consistent with the transition from an
incommensurate magnetic order to a spinon excitation continuum to a
paramagnetic state observed in the inelastic neutron scattering experiment
[1]. The coexisting states also appear to be consistent with the dispersion
relation given in the experiment.

To further illustrate the spin excitations, we also present the dynamical spin
structure factor in Fig. 6 with momentum cuts at $q_{y}=0$,
$q_{y}=\sqrt{3}\pi$, and $q_{x}=-\pi/2$, following the ones used in Fig. 3 of
Ref. [1]. The comparison of the spin excitations from the experiment and
calculations based on the Heisenberg model are well known. The peak dispersion
of neutron scattering is well described by series expansion and $1/S$
expansion [139, 140, 141, 142]. A perturbative calculation based on 1D spin
chains further reveals the excitation continuum and the peak dispersion
originate from lower dimension spinons [143].

Due to the ill-conditioned nature of the analytic continuation problem, the
details of the dynamical spin structure factor depend on the default model
used in maximum entropy calculation and are also sensitive to the error
estimation of the calculation which is difficult to obtain in the LDFA
calculation. This limits the reliability of the dynamical quantities and thus
makes the direct comparison with experimental data difficult. Since fine-
tuning the parameters in the analytic continuation process may introduce
misleading agreement, we conservatively choose to use the same model and same
error estimation for all momentum points $\bm{q}$ and all temperatures $T$
studied. This results in anomalies that could be removed by tuning the
parameters such as high energy tails around $q_{x}=\pi,3\pi$ in Fig. 6(a) and
around $q_{x}=\pi$ in Fig. 6(b) where the Matsubara data are significantly
larger than at other $\bm{q}$ points.

As shown in the top row of Fig. 6, the peak dispersion (white line) at low
temperature $T=0.1$ qualitatively matches the peak dispersion (white open
square) below $0.1$K given in the neutron scattering experiment [1]. More
precisely, the minimum close to $q_{x}=\pi$ and $q_{x}=3\pi$ shown in Fig.
6(a), the asymmetric dispersion with respect to $q_{x}=\pi$ and $q_{x}=3\pi$
shown in Fig. 6(a), the symmetric dispersion with respect to $q_{x}=\pi$ shown
in Fig. 6(b), and the asymmetric dispersion with respect to
$q_{y}=\sqrt{3}\pi$ shown in Fig. 6(c) all agree with the dispersion relation
obtained in the experiment. Note that to highlight this agreement, the
experimental results have been scaled by a factor of two. Using the estimate
$J=4t^{2}/U$ and $J=0.128(5)$meV [1], one can show $T=0.1$ used in the
calculation corresponds to a physical temperature above $0.1$K but below
$T_{N}=0.62$K, the transition temperature from an incommensurate order to
excitation continuum reported in Ref. [1]. Adjusting temperature $T$ (and/or
$U,t^{\prime}$) therefore may yield a more accurate absolute agreement here.
Besides the agreement for dispersion, the peak value at
$\bm{q}=(-\pi/2,2\pi/\sqrt{3})$ is found to be larger than the peak value at
$\bm{q}=(-\pi/2,4\pi/\sqrt{3})$ as indicated by the white solid dots and
numbers in Fig. 6(c). The ratio of these two peak values is found to be about
$4$ in the experiment (see scan E, F in Fig. 5 of Ref. [1]) and is well
explained by the perturbative study [143]. Although the ratio found in our
calculation is much smaller than the experiment value, we find the feature of
a larger peak value occurring at $\bm{q}=(-\pi/2,2\pi/\sqrt{3})$ is robust
upon different default models and error estimation.

When temperature $T$ is further increased, the dispersion tends to become flat
as shown in the second to the fourth row of Fig. 6 which signatures the
features of the incommensurate order at lower temperature gradually disappear.
The ratio of peak values at $\bm{q}=(-\pi/2,2\pi/\sqrt{3})$ and
$\bm{q}=(-\pi/2,4\pi/\sqrt{3})$ becomes larger at $T=0.5$ and then decreases,
which agrees with the picture that spinon excitations are enhanced at
intermediate temperature and the paramagnetic state appears at higher
temperature [1].

##  IV Conclusion

We have presented a comprehensive study of the phase diagram of the
anisotropic triangular lattice Hubbard model as a function of $t^{\prime}$ and
$U$ at half-filling.

The ladder dual fermion approximation used here captures all local
correlations while treat non-local correlations perturbatively. It provides a
fine k-space resolution of the susceptibility and is able to resolve spin
fluctuations without an a priori assumption of possible orders.

We found a rich phase diagram with metallic, SAFM, spiral, NCAFM, and 1DAFM
phases. While many aspects of the phase diagram have been found before, we
resolve the discrepancies between previously published results. More
precisely, the disputed CAFM order (see Sec. III.2) near the limit of
decoupled chains is not found by our method. Instead, the NCAFM order and the
coexisting phase of the 1DAFM and NCAFM orders are found to be dominant in the
corresponding region. Some indirect support in previous studies related to the
newly found coexisting phase is discussed in Sec. III.3. We also investigated
the Lifshitz transition on the anisotropic triangular lattice in the presence
of interaction.

Our study investigated the physics of the model for parameters relevant to
Cs2CuCl4 in more detail. We found a transition from an incommensurate magnetic
order to a spinon excitation continuum to a paramagnetic state upon the
increase of temperature, which is in agreement with the inelastic neutron
scattering experiment [1]. Despite the limitation of analytic continuation,
the dynamical spin structure factor at $U=11$ and $t^{\prime}=1.715$ shows a
dispersion relation consistent with results obtained in Ref. [1].

## Acknowledgement

We thank Wei-Ting Lin and André Erpenbeck for helpful discussions. This work
was supported by NSF DMR 2001465. This work used Expanse at SDSC through
allocation DMR130036 from the Extreme Science and Engineering Discovery
Environment (XSEDE), which was supported by National Science Foundation grant
number #1548562.

##  Appendix A Local density of states for the non-interacting anisotropic
triangular lattice

The local density of states,
$\rho(\nu)=-\operatorname{Im}G^{0}(\nu+\mathrm{i}0^{+})/\pi$, of the non-
interacting anisotropic triangular lattice is used as the input of DMFT (see
Ref [144] for the expression of the non-interacting Green’s function $G^{0}$).
It can be expressed analytically in terms of the complete elliptic integral of
the first kind $K(m)=\int_{0}^{\pi/2}d\phi\left[1-m\sin^{2}\phi\right]^{-1/2}$
as:

$\rho(\nu;t,t^{\prime})=\frac{1}{\pi^{2}t^{\prime}\sqrt{z_{0}}}K\left(\frac{z_{1}}{z_{0}}\right),$
(8) $z_{0}=\left\\{\begin{aligned} &q\qquad&0<p\leq q\\\ &p\qquad&0<q<p\\\
&p-q\qquad&q<0,\end{aligned}\right.$ (9) $z_{1}=\left\\{\begin{aligned}
&q-p\qquad&0<p\leq q\\\ &p-q\qquad&0<q<p\\\ &p\qquad&q<0,\end{aligned}\right.$
(10)

where we define two dimensionless parameters $u=t/t^{\prime}\geq 0$ and
$E=\nu/t$ along with the following notation:

$\displaystyle r$ $\displaystyle\equiv u\sqrt{u^{2}-Eu+2},$ (11)
$\displaystyle p$ $\displaystyle\equiv 4r,$ $\displaystyle q$
$\displaystyle\equiv\frac{(r-u^{2})^{2}(r^{2}-4u^{2}+2ru^{2}+u^{4})}{4u^{4}}.$

The range of $\nu$ is indicated by

$\displaystyle-4-\frac{2}{u}\leq E\leq u+\frac{2}{u}\qquad$ $\displaystyle
0<u\leq 2,$ (12) $\displaystyle-4-\frac{2}{u}\leq E\leq 4-\frac{2}{u}\qquad$
$\displaystyle u>2.$

## References

  * Coldea  _et al._ [2003] R. Coldea, D. A. Tennant, and Z. Tylczynski, Physical Review B 68, 134424 (2003).
  * Shimizu  _et al._ [2003] Y. Shimizu, K. Miyagawa, K. Kanoda, M. Maesato, and G. Saito, Physical Review Letters 91, 107001 (2003).
  * Itou  _et al._ [2008] T. Itou, A. Oyamada, S. Maegawa, M. Tamura, and R. Kato, Physical Review B 77, 104413 (2008).
  * Nakatsuji  _et al._ [2005] S. Nakatsuji, Y. Nambu, H. Tonomura, O. Sakai, S. Jonas, C. Broholm, H. Tsunetsugu, Y. Qiu, and Y. Maeno, Science 309, 1697 (2005).
  * Zhou  _et al._ [2011] H. D. Zhou, E. S. Choi, G. Li, L. Balicas, C. R. Wiebe, Y. Qiu, J. R. D. Copley, and J. S. Gardner, Physical Review Letters 106, 147204 (2011).
  * Shirata  _et al._ [2012] Y. Shirata, H. Tanaka, A. Matsuo, and K. Kindo, Physical Review Letters 108, 057205 (2012).
  * Li  _et al._ [2015] Y. Li, G. Chen, W. Tong, L. Pi, J. Liu, Z. Yang, X. Wang, and Q. Zhang, Physical Review Letters 115, 167203 (2015).
  * Bordelon  _et al._ [2019] M. M. Bordelon, E. Kenney, C. Liu, T. Hogan, L. Posthuma, M. Kavand, Y. Lyu, M. Sherwin, N. P. Butch, C. Brown, M. J. Graf, L. Balents, and S. D. Wilson, Nature Physics 15, 1058 (2019).
  * Wu  _et al._ [2018] F. Wu, T. Lovorn, E. Tutuc, and A. H. MacDonald, Physical Review Letters 121, 026402 (2018).
  * Zhang  _et al._ [2021] Y.-H. Zhang, D. N. Sheng, and A. Vishwanath, Physical Review Letters 127, 247701 (2021).
  * Regan  _et al._ [2022] E. C. Regan, D. Wang, E. Y. Paik, Y. Zeng, L. Zhang, J. Zhu, A. H. MacDonald, H. Deng, and F. Wang, Nature Reviews Materials 7, 778 (2022).
  * Struck  _et al._ [2011] J. Struck, C. Ölschläger, R. Le Targat, P. Soltan-Panahi, A. Eckardt, M. Lewenstein, P. Windpassinger, and K. Sengstock, Science 333, 996 (2011).
  * Yang  _et al._ [2021] J. Yang, L. Liu, J. Mongkolkiattichai, and P. Schauss, PRX Quantum 2, 020344 (2021).
  * Xu  _et al._ [2016] Y. Xu, J. Zhang, Y. S. Li, Y. J. Yu, X. C. Hong, Q. M. Zhang, and S. Y. Li, Physical Review Letters 117, 267202 (2016).
  * Ni  _et al._ [2019] J. M. Ni, B. L. Pan, B. Q. Song, Y. Y. Huang, J. Y. Zeng, Y. J. Yu, E. J. Cheng, L. S. Wang, D. Z. Dai, R. Kato, and S. Y. Li, Physical Review Letters 123, 247204 (2019).
  * Bourgeois-Hope  _et al._ [2019] P. Bourgeois-Hope, F. Laliberté, E. Lefrançois, G. Grissonnanche, S. R. de Cotret, R. Gordon, S. Kitou, H. Sawa, H. Cui, R. Kato, L. Taillefer, and N. Doiron-Leyraud, Physical Review X 9, 041051 (2019).
  * Savary and Balents [2016] L. Savary and L. Balents, Reports on Progress in Physics 80, 016502 (2016).
  * Zhou  _et al._ [2017] Y. Zhou, K. Kanoda, and T.-K. Ng, Reviews of Modern Physics 89, 025003 (2017).
  * Anderson [1973] P. W. Anderson, Materials Research Bulletin 8, 153 (1973).
  * Anderson [1987] P. W. Anderson, Science 235, 1196 (1987).
  * Balents [2010] L. Balents, Nature 464, 199 (2010).
  * Knolle and Moessner [2019] J. Knolle and R. Moessner, Annual Review of Condensed Matter Physics 10, 451 (2019).
  * Broholm  _et al._ [2020] C. Broholm, R. J. Cava, S. A. Kivelson, D. G. Nocera, M. R. Norman, and T. Senthil, Science 367, eaay0668 (2020).
  * Lee  _et al._ [2006] P. A. Lee, N. Nagaosa, and X.-G. Wen, Reviews of Modern Physics 78, 17 (2006).
  * Kanoda and Kato [2011] K. Kanoda and R. Kato, Annual Review of Condensed Matter Physics 2, 167 (2011).
  * Powell and McKenzie [2011] B. J. Powell and R. H. McKenzie, Reports on Progress in Physics 74, 056501 (2011).
  * Furukawa  _et al._ [2015] T. Furukawa, K. Miyagawa, H. Taniguchi, R. Kato, and K. Kanoda, Nature Physics 11, 221 (2015).
  * Weng  _et al._ [2006] M. Q. Weng, D. N. Sheng, Z. Y. Weng, and R. J. Bursill, Physical Review B 74, 012407 (2006).
  * Hayashi and Ogata [2007] Y. Hayashi and M. Ogata, Journal of the Physical Society of Japan 76, 053705 (2007).
  * Hauke  _et al._ [2011] P. Hauke, T. Roscilde, V. Murg, J. I. Cirac, and R. Schmied, New Journal of Physics 13, 075017 (2011).
  * Hauke [2013] P. Hauke, Physical Review B 87, 014415 (2013).
  * Reuther and Thomale [2011] J. Reuther and R. Thomale, Physical Review B 83, 024402 (2011).
  * Yunoki and Sorella [2006] S. Yunoki and S. Sorella, Physical Review B 74, 014408 (2006).
  * Heidarian  _et al._ [2009] D. Heidarian, S. Sorella, and F. Becca, Physical Review B 80, 012404 (2009).
  * Ghorbani  _et al._ [2016] E. Ghorbani, L. F. Tocchio, and F. Becca, Physical Review B 93, 085111 (2016).
  * Castells-Graells  _et al._ [2019] D. Castells-Graells, A. Yuste, and A. Sanpera, Physical Review B 100, 155119 (2019).
  * Gonzalez  _et al._ [2020] M. G. Gonzalez, E. A. Ghioldi, C. J. Gazza, L. O. Manuel, and A. E. Trumper, Physical Review B 102, 224410 (2020).
  * Starykh and Balents [2007] O. A. Starykh and L. Balents, Physical Review Letters 98, 077205 (2007).
  * Starykh  _et al._ [2010] O. A. Starykh, H. Katsura, and L. Balents, Physical Review B 82, 014421 (2010).
  * Pardini and Singh [2008] T. Pardini and R. R. P. Singh, Physical Review B 77, 214433 (2008).
  * Bishop  _et al._ [2009] R. F. Bishop, P. H. Y. Li, D. J. J. Farnell, and C. E. Campbell, Physical Review B 79, 174405 (2009).
  * Ghamari  _et al._ [2011] S. Ghamari, C. Kallin, S.-S. Lee, and E. S. Sørensen, Physical Review B 84, 174415 (2011).
  * Thesberg and Sørensen [2014] M. Thesberg and E. S. Sørensen, Physical Review B 90, 115117 (2014).
  * Yamada [2014a] A. Yamada, Physical Review B 90, 235138 (2014a).
  * Tocchio  _et al._ [2014] L. F. Tocchio, C. Gros, R. Valentí, and F. Becca, Physical Review B 89, 235107 (2014).
  * Szasz and Motruk [2021] A. Szasz and J. Motruk, Physical Review B 103, 235132 (2021).
  * Yamada [2014b] A. Yamada, Physical Review B 89, 195108 (2014b).
  * Acheche  _et al._ [2016] S. Acheche, A. Reymbaut, M. Charlebois, D. Sénéchal, and A. M. S. Tremblay, Physical Review B 94, 245133 (2016).
  * Arovas  _et al._ [2022] D. P. Arovas, E. Berg, S. A. Kivelson, and S. Raghu, Annual Review of Condensed Matter Physics 13, 239 (2022).
  * Qin  _et al._ [2022] M. Qin, T. Schäfer, S. Andergassen, P. Corboz, and E. Gull, Annual Review of Condensed Matter Physics 13, 275 (2022).
  * Rubtsov  _et al._ [2008] A. N. Rubtsov, M. I. Katsnelson, and A. I. Lichtenstein, Physical Review B 77, 033101 (2008).
  * Brener  _et al._ [2008] S. Brener, H. Hafermann, A. N. Rubtsov, M. I. Katsnelson, and A. I. Lichtenstein, Physical Review B 77, 195105 (2008).
  * Rubtsov  _et al._ [2009] A. N. Rubtsov, M. I. Katsnelson, A. I. Lichtenstein, and A. Georges, Physical Review B 79, 045133 (2009).
  * Hafermann  _et al._ [2009] H. Hafermann, G. Li, A. N. Rubtsov, M. I. Katsnelson, A. I. Lichtenstein, and H. Monien, Physical Review Letters 102, 206401 (2009).
  * Rohringer  _et al._ [2018] G. Rohringer, H. Hafermann, A. Toschi, A. A. Katanin, A. E. Antipov, M. I. Katsnelson, A. I. Lichtenstein, A. N. Rubtsov, and K. Held, Reviews of Modern Physics 90, 025003 (2018).
  * Georges  _et al._ [1996] A. Georges, G. Kotliar, W. Krauth, and M. J. Rozenberg, Reviews of Modern Physics 68, 13 (1996).
  * Lee  _et al._ [2008] H. Lee, G. Li, and H. Monien, Physical Review B 78, 205117 (2008).
  * Antipov  _et al._ [2011] A. E. Antipov, A. N. Rubtsov, M. I. Katsnelson, and A. I. Lichtenstein, Physical Review B 83, 115126 (2011).
  * Li  _et al._ [2014] G. Li, A. E. Antipov, A. N. Rubtsov, S. Kirchner, and W. Hanke, Physical Review B 89, 161118 (2014).
  * Yudin  _et al._ [2014] D. Yudin, D. Hirschmeier, H. Hafermann, O. Eriksson, A. I. Lichtenstein, and M. I. Katsnelson, Physical Review Letters 112, 070403 (2014).
  * Laubach  _et al._ [2015] M. Laubach, R. Thomale, C. Platt, W. Hanke, and G. Li, Physical Review B 91, 245125 (2015).
  * Li and Gull [2020] S. Li and E. Gull, Physical Review Research 2, 013295 (2020).
  * Katanin [2013] A. A. Katanin, Journal of Physics A: Mathematical and Theoretical 46, 045002 (2013).
  * Simons Collaboration on the Many-Electron Problem  _et al._ [2015] Simons Collaboration on the Many-Electron Problem, J. P. F. LeBlanc, A. E. Antipov, F. Becca, I. W. Bulik, G. K.-L. Chan, C.-M. Chung, Y. Deng, M. Ferrero, T. M. Henderson, C. A. Jiménez-Hoyos, E. Kozik, X.-W. Liu, A. J. Millis, N. V. Prokof’ev, M. Qin, G. E. Scuseria, H. Shi, B. V. Svistunov, L. F. Tocchio, I. S. Tupitsyn, S. R. White, S. Zhang, B.-X. Zheng, Z. Zhu, and E. Gull, Physical Review X 5, 041041 (2015).
  * Iskakov  _et al._ [2016] S. Iskakov, A. E. Antipov, and E. Gull, Physical Review B 94, 035102 (2016).
  * Ribic  _et al._ [2017] T. Ribic, P. Gunacker, S. Iskakov, M. Wallerberger, G. Rohringer, A. N. Rubtsov, E. Gull, and K. Held, Physical Review B 96, 235127 (2017).
  * LeBlanc  _et al._ [2019] J. P. F. LeBlanc, S. Li, X. Chen, R. Levy, A. E. Antipov, A. J. Millis, and E. Gull, Physical Review B 100, 075123 (2019).
  * Gukelberger  _et al._ [2017] J. Gukelberger, E. Kozik, and H. Hafermann, Physical Review B 96, 035152 (2017).
  * Antipov  _et al._ [2015] A. E. Antipov, J. P. F. LeBlanc, and E. Gull, Physics Procedia Proceedings of the 28th Workshop on Computer Simulation Studies in Condensed Matter Physics (CSP2015), 68, 43 (2015).
  * [70] M. Wallerberger, S. Iskakov, A. Gaenko, J. Kleinhenz, I. Krivenko, R. Levy, J. Li, H. Shinaoka, S. Todo, T. Chen, X. Chen, J. P. F. LeBlanc, J. E. Paki, H. Terletska, M. Troyer, and E. Gull, arXiv:1811.08331 [cond-mat, physics:physics] .
  * Gaenko  _et al._ [2017] A. Gaenko, A. E. Antipov, G. Carcassi, T. Chen, X. Chen, Q. Dong, L. Gamper, J. Gukelberger, R. Igarashi, S. Iskakov, M. Könz, J. P. F. LeBlanc, R. Levy, P. N. Ma, J. E. Paki, H. Shinaoka, S. Todo, M. Troyer, and E. Gull, Computer Physics Communications 213, 235 (2017).
  * Rohringer  _et al._ [2012] G. Rohringer, A. Valli, and A. Toschi, Physical Review B 86, 125114 (2012).
  * Gull  _et al._ [2008] E. Gull, P. Werner, O. Parcollet, and M. Troyer, EPL (Europhysics Letters) 82, 57003 (2008).
  * Gull  _et al._ [2011a] E. Gull, A. J. Millis, A. I. Lichtenstein, A. N. Rubtsov, M. Troyer, and P. Werner, Reviews of Modern Physics 83, 349 (2011a).
  * Gull  _et al._ [2011b] E. Gull, P. Staar, S. Fuchs, P. Nukala, M. S. Summers, T. Pruschke, T. C. Schulthess, and T. Maier, Physical Review B 83, 075122 (2011b).
  * Chen  _et al._ [2015] X. Chen, J. P. F. LeBlanc, and E. Gull, Physical Review Letters 115, 116402 (2015).
  * Jarrell and Gubernatis [1996] M. Jarrell and J. E. Gubernatis, Physics Reports 269, 133 (1996).
  * Levy  _et al._ [2017] R. Levy, J. P. F. LeBlanc, and E. Gull, Computer Physics Communications 215, 149 (2017).
  * Mahan [2000] G. D. Mahan,  _Many-Particle Physics_ (Springer US, Boston, MA, 2000).
  * Sachdev [2011] S. Sachdev,  _Quantum Phase Transitions_, 2nd ed. (Cambridge University Press, 2011).
  * Yang  _et al._ [2010] H.-Y. Yang, A. M. Läuchli, F. Mila, and K. P. Schmidt, Physical Review Letters 105, 267204 (2010).
  * Zheng  _et al._ [1999] W. Zheng, R. H. McKenzie, and R. R. P. Singh, Physical Review B 59, 14367 (1999).
  * Merino  _et al._ [1999] J. Merino, R. H. McKenzie, J. B. Marston, and C. H. Chung, Journal of Physics: Condensed Matter 11, 2965 (1999).
  * Anderson [1952] P. W. Anderson, Physical Review 86, 694 (1952).
  * Chakravarty  _et al._ [1989] S. Chakravarty, B. I. Halperin, and D. R. Nelson, Physical Review B 39, 2344 (1989).
  * Manousakis [1991] E. Manousakis, Reviews of Modern Physics 63, 1 (1991).
  * Huse and Elser [1988] D. A. Huse and V. Elser, Physical Review Letters 60, 2531 (1988).
  * Jolicoeur and Le Guillou [1989] T. Jolicoeur and J. C. Le Guillou, Physical Review B 40, 2727 (1989).
  * Bernu  _et al._ [1992] B. Bernu, C. Lhuillier, and L. Pierre, Physical Review Letters 69, 2590 (1992).
  * Bernu  _et al._ [1994] B. Bernu, P. Lecheminant, C. Lhuillier, and L. Pierre, Physical Review B 50, 10048 (1994).
  * Capriotti  _et al._ [1999] L. Capriotti, A. E. Trumper, and S. Sorella, Physical Review Letters 82, 3899 (1999).
  * Harada [2012] K. Harada, Physical Review B 86, 184421 (2012).
  * Chung  _et al._ [2001] C. H. Chung, J. B. Marston, and R. H. McKenzie, Journal of Physics: Condensed Matter 13, 5159 (2001).
  * Weichselbaum and White [2011] A. Weichselbaum and S. R. White, Physical Review B 84, 245130 (2011).
  * Gonzalez  _et al._ [2022] M. Gonzalez, B. Bernu, L. Pierre, and L. Messio, SciPost Physics 12, 112 (2022).
  * Bocquet  _et al._ [2001] M. Bocquet, F. H. L. Essler, A. M. Tsvelik, and A. O. Gogolin, Physical Review B 64, 094425 (2001).
  * Bethe [1931] H. Bethe, Zeitschrift für Physik 71, 205 (1931).
  * des Cloizeaux and Pearson [1962] J. des Cloizeaux and J. J. Pearson, Physical Review 128, 2131 (1962).
  * Endoh  _et al._ [1974] Y. Endoh, G. Shirane, R. J. Birgeneau, P. M. Richards, and S. L. Holt, Physical Review Letters 32, 170 (1974).
  * Tennant  _et al._ [1993] D. A. Tennant, T. G. Perring, R. A. Cowley, and S. E. Nagler, Physical Review Letters 70, 4003 (1993).
  * Powell  _et al._ [2017] B. J. Powell, E. P. Kenny, and J. Merino, Physical Review Letters 119, 087204 (2017).
  * Alicea  _et al._ [2005] J. Alicea, O. I. Motrunich, and M. P. A. Fisher, Physical Review Letters 95, 247203 (2005).
  * Kyung and Tremblay [2006] B. Kyung and A. M. S. Tremblay, Physical Review Letters 97, 046402 (2006).
  * Sahebsara and Sénéchal [2006] P. Sahebsara and D. Sénéchal, Physical Review Letters 97, 257004 (2006).
  * Watanabe  _et al._ [2006] T. Watanabe, H. Yokoyama, Y. Tanaka, and J.-i. Inoue, Journal of the Physical Society of Japan 75, 074707 (2006).
  * Clay  _et al._ [2008] R. T. Clay, H. Li, and S. Mazumdar, Physical Review Letters 101, 166403 (2008).
  * Watanabe  _et al._ [2008] T. Watanabe, H. Yokoyama, Y. Tanaka, and J. Inoue, Physical Review B 77, 214505 (2008).
  * Tocchio  _et al._ [2013] L. F. Tocchio, H. Feldner, F. Becca, R. Valentí, and C. Gros, Physical Review B 87, 035143 (2013).
  * Luther and Emery [1974] A. Luther and V. J. Emery, Physical Review Letters 33, 589 (1974).
  * Gannot  _et al._ [2020] Y. Gannot, Y.-F. Jiang, and S. A. Kivelson, Physical Review B 102, 115136 (2020).
  * Szasz  _et al._ [2020] A. Szasz, J. Motruk, M. P. Zaletel, and J. E. Moore, Physical Review X 10, 021042 (2020).
  * Zhu  _et al._ [2022] Z. Zhu, D. N. Sheng, and A. Vishwanath, Physical Review B 105, 205110 (2022).
  * Morita  _et al._ [2002] H. Morita, S. Watanabe, and M. Imada, Journal of the Physical Society of Japan 71, 2109 (2002).
  * Aryanpour  _et al._ [2006] K. Aryanpour, W. E. Pickett, and R. T. Scalettar, Physical Review B 74, 085117 (2006).
  * Koretsune  _et al._ [2007] T. Koretsune, Y. Motome, and A. Furusaki, Journal of the Physical Society of Japan 76, 074719 (2007).
  * Sahebsara and Sénéchal [2008] P. Sahebsara and D. Sénéchal, Physical Review Letters 100, 136402 (2008).
  * Ohashi  _et al._ [2008] T. Ohashi, T. Momoi, H. Tsunetsugu, and N. Kawakami, Physical Review Letters 100, 076402 (2008).
  * Tocchio  _et al._ [2009] L. F. Tocchio, A. Parola, C. Gros, and F. Becca, Physical Review B 80, 064419 (2009).
  * Yoshioka  _et al._ [2009] T. Yoshioka, A. Koga, and N. Kawakami, Physical Review Letters 103, 036401 (2009).
  * Galanakis  _et al._ [2009] D. Galanakis, T. D. Stanescu, and P. Phillips, Physical Review B 79, 115116 (2009).
  * Kokalj and McKenzie [2013] J. Kokalj and R. H. McKenzie, Physical Review Letters 110, 206402 (2013).
  * Dang  _et al._ [2015] H. T. Dang, X. Y. Xu, K.-S. Chen, Z. Y. Meng, and S. Wessel, Physical Review B 91, 155101 (2015).
  * Mishmash  _et al._ [2015] R. V. Mishmash, I. González, R. G. Melko, O. I. Motrunich, and M. P. A. Fisher, Physical Review B 91, 235140 (2015).
  * Goto  _et al._ [2016] S. Goto, S. Kurihara, and D. Yamamoto, Physical Review B 94, 245145 (2016).
  * Shirakawa  _et al._ [2017] T. Shirakawa, T. Tohyama, J. Kokalj, S. Sota, and S. Yunoki, Physical Review B 96, 205130 (2017).
  * Misumi  _et al._ [2017] K. Misumi, T. Kaneko, and Y. Ohta, Physical Review B 95, 075124 (2017).
  * Wietek  _et al._ [2021] A. Wietek, R. Rossi, F. Šimkovic, M. Klett, P. Hansmann, M. Ferrero, E. M. Stoudenmire, T. Schäfer, and A. Georges, Physical Review X 11, 041013 (2021).
  * [128] P.-O. Downey, O. Gingras, M. Charlebois, C.-D. Hébert, and A. M. S. Tremblay, arXiv:2207.08008 [cond-mat, physics:physics] .
  * Chen  _et al._ [2022] B.-B. Chen, Z. Chen, S.-S. Gong, D. N. Sheng, W. Li, and A. Weichselbaum, Physical Review B 106, 094420 (2022).
  * Cookmeyer  _et al._ [2021] T. Cookmeyer, J. Motruk, and J. E. Moore, Physical Review Letters 127, 087201 (2021).
  * Kadow  _et al._ [2022] W. Kadow, L. Vanderstraeten, and M. Knap, Physical Review B 106, 094417 (2022).
  * [132] V. Zampronio and T. Macrì, arXiv:2210.13551 [cond-mat, physics:quant-ph] .
  * Yuste  _et al._ [2017] A. Yuste, M. Moreno-Cardoner, and A. Sanpera, Physical Review B 95, 195167 (2017).
  * Ido  _et al._ [2022] K. Ido, K. Yoshimi, T. Misawa, and M. Imada, npj Quantum Materials 7, 1 (2022).
  * Raczkowski and Assaad [2013] M. Raczkowski and F. F. Assaad, Physical Review B 88, 085120 (2013).
  * Volovik [2017] G. E. Volovik, Low Temperature Physics 43, 47 (2017).
  * Peters and Kawakami [2015] R. Peters and N. Kawakami, Physical Review B 92, 075103 (2015).
  * Raczkowski  _et al._ [2015] M. Raczkowski, F. F. Assaad, and L. Pollet, Physical Review B 91, 045137 (2015).
  * Veillette  _et al._ [2005] M. Y. Veillette, A. J. A. James, and F. H. L. Essler, Physical Review B 72, 134429 (2005).
  * Dalidovich  _et al._ [2006] D. Dalidovich, R. Sknepnek, A. J. Berlinsky, J. Zhang, and C. Kallin, Physical Review B 73, 184403 (2006).
  * Zheng  _et al._ [2006] W. Zheng, J. O. Fjærestad, R. R. P. Singh, R. H. McKenzie, and R. Coldea, Physical Review Letters 96, 057201 (2006).
  * Fjærestad  _et al._ [2007] J. O. Fjærestad, W. Zheng, R. R. P. Singh, R. H. McKenzie, and R. Coldea, Physical Review B 75, 174447 (2007).
  * Kohno  _et al._ [2007] M. Kohno, O. A. Starykh, and L. Balents, Nature Physics 3, 790 (2007).
  * Horiguchi [1991] T. Horiguchi, Physica A: Statistical Mechanics and its Applications 178, 351 (1991).