Ferromagnetism and Spin-Valley liquid states in Moir\'{e} Correlated Insulators

TL;DR

This paper investigates ferromagnetism and exotic quantum phases in Moiré superlattice insulators, revealing a ferromagnetic region, potential metal-insulator transitions, and topological superconductivity upon doping, using analytical, numerical, and theoretical methods.

Contribution

It introduces a two-orbital model with SU(2) symmetries to explain ferromagnetism and spin-liquid phases in Moiré systems, and predicts topological superconductivity after doping.

Findings

Existence of ferromagnetic order despite antiferromagnetic couplings

Zeeman field induces a metal-insulator transition as observed experimentally

Doped system likely becomes a spin-triplet/valley-singlet topological superconductor

Abstract

Motivated by the recent observation of evidences of ferromagnetism in correlated insulating states in systems with Moir\'{e} superlattices, we study a two-orbital quantum antiferromagnetic model on the triangular lattice, where the two orbitals physically correspond to the two valleys of the original graphene sheet. For simplicity this model has a SU(2)$^s$$\otimes$SU(2)$^v$ symmetry, where the two SU(2) symmetries correspond to the rotation within the spin and valley space respectively. Through analytical argument, Schwinger boson analysis and also DMRG simulation, we find that even though all the couplings in the Hamiltonian are antiferromagnetic, there is still a region in the phase diagram with fully polarized ferromagnetic order. We argue that a Zeeman field can drive a metal-insulator transition in our picture, as was observed experimentally. We also construct spin liquids and topological ordered phases at various limits of this model. Then after doping this model with extra charge carriers, the system most likely becomes spin-triplet/valley-singlet $d+id$ topological superconductor as was predicted previously.

Figures

Equations

$$H=\sum _{<i,j>}\sum _{a,b=1}^{3}J{T}_i^{ab}{T}_j^{ab}+J^s{\sigma }_i^a{\sigma }_j^a+J^v{\tau }_i^b{\tau }_j^b,$$

$$H$$

$$H_t$$

$$H_U$$

$$H_V$$

$$\mathcal{P}=\prod _j(-1)\frac{1}{6}{n}_j\left(n_j-2\right)\left(n_j-3\right)\left(n_j-4\right).$$

$$\mathcal{P}={\mathcal{P}}_{ss}+{\mathcal{P}}_{st}+{\mathcal{P}}_{ts}+{\mathcal{P}}_{tt},$$

$$H_{\text{eff}}=\mathcal{P}{H}_t\frac{1}{E_0-H_U-H_V}{H}_t\mathcal{P},$$

$$H_{\text{eff}}=-\frac{2t^2}{U+4V}{\mathcal{P}}_{st}-\frac{2t^2}{U-4V}{\mathcal{P}}_{ts}.$$

$$J^s$$

$$J^v$$

$$J$$

$$|{\Psi }_{\mathrm{FM}}⟩=\left(\prod _i|{\sigma }_i^z=+1⟩\right)\otimes |\mathrm{AF}\mathrm{of}{\vec{\tau }}_{\{i\}}⟩.$$

$$J_{\mathrm{eff}}=J^s+J⟨{\vec{\tau }}_i\cdot {\vec{\tau }}_j⟩.$$

$$\sum _{\alpha =1}^4{b}_{j,\alpha }^{†}{b}_{j,\alpha }=\kappa .$$

$$H$$

$$+$$

$$\left({\vec{\Delta }}_{ij}^{ts},{\vec{\Delta }}_{ij}^{st}\right)=b_i^t\left(\mathrm{i}{\sigma }^{32},{\sigma }^{02},\mathrm{i}{\sigma }^{12},{\sigma }^{23},\mathrm{i}{\sigma }^{20},{\sigma }^{21}\right)b_j,$$

$$J_{ts}=-\frac{1}{4}(3J_v+3J-J_s),J_{st}=-\frac{1}{4}(3J+3J_s-J_v),$$

$$J_{ss}=\frac{3}{4}(3J-J_s-J_v),J_{tt}=\frac{1}{4}(J+J_s+J_v).$$

$$H_{\mathrm{MF}}$$

$$+$$

$$\vec{M}\sim \mathrm{i}\vec{\varphi }\times {\vec{\varphi }}^{\ast }\sim {\vec{\varphi }}_1\times {\vec{\varphi }}_2.$$

$$H_{\mathrm{MF}}$$

$$+$$

$$H_{\mathrm{MF}}$$

$$+$$

$$\mathcal{T}:b_j\to \mathrm{i}{\sigma }^{21}{b}_j,{\vec{\Delta }}_{ij}^{ts}\to {\vec{\Delta }}_{ij}^{ts},$$

$$\sum _{\alpha =1,2}{b}_{j,\alpha }^{s,†}{b}_{j,\alpha }^s=b_{j,\alpha }^{v,†}{b}_{j,\alpha }^v=1.$$

$$\sum _{\alpha =1,2}{f}_{j,\alpha }^{s,†}{f}_{j,\alpha }^s=f_{j,\alpha }^{v,†}{f}_{j,\alpha }^v=1.$$

$$H=\sum _{ij}J\left(-\frac{5}{4}({\vec{\Delta }}_{ij}^{†})\cdot ({\vec{\Delta }}_{ij})+\cdots \right),$$

$$H_{\mathrm{MF}}=\sum _{ij}J\left(-\frac{5}{4}\vec{\varphi }\cdot {\vec{\Delta }}_{ij}+H.c.\right)+\frac{5}{4}J|\vec{\varphi }{|}^2.$$

Full text

I Introduction

Recently a series of surprising correlated physics such as superconductivity and insulators at commensurate fractional charge fillings have been discovered in multiple systems with Moiré superlattices Chen et al. (2018); Cao et al. (2018a, b); Yankowitz et al. (2019); Chen et al. (2019); Sharpe et al. (2019); Kim (2019); Shen et al. (2019); Liu et al. (2019). These discoveries have motivated very active theoretical studies Xu and Balents (2018); Po et al. (2018a); Dodaro et al. (2018); Yuan and Fu (2018); Kang and Vafek (2018a); Padhi et al. (2018); Padhi and Phillips (2018); Baskaran (2018); Liu et al. (2018a); Rademaker and Mellado (2018); Isobe et al. (2018); Koshino et al. (2018); You and Vishwanath (2018); Wu et al. (2018); Zhu et al. (2018); Lian et al. (2018); Guo et al. (2018); Zhang and Senthil (2018a); Po et al. (2018b); Kang and Vafek (2018b); Tang et al. (2018); Thomson et al. (2018); Hejazi et al. (2019); Liu et al. (2018b); Wu et al. (2019); Seo et al. (2018); Bultinck et al. (2019); Lee et al. (2019); Liu and Dai (2019); Wu and Sarma (2019). For different reasons, these systems can have narrow electron bandwidth near charge neutrality Bistritzer and MacDonald (2011); Suárez Morell et al. (2010); Fang and Kaxiras (2016); Trambly de Laissardière et al. (2012); Shen et al. (2019); Liu et al. (2019), hence interaction effects are significantly enhanced. A consensus of the mechanism for the observed insulator and superconductor has not yet been reached. A minimal two-orbital extended Hubbard model on the triangular lattice was proposed in Ref. Xu and Balents, 2018 (physically the orbital space corresponds to the space of two Dirac valleys in the original Brillouin zone of graphene), which at least describes the trilayer graphene and hexagonal BN heterostructure (TLG/h-BN) Chen et al. (2018, 2019), as well as the twisted double bilayer graphene (TDBG) Kim (2019); Shen et al. (2019); Liu et al. (2019) with certain twisted angle and out-of-plane electric field (displacement field), since in these cases there is no symmetry protected band touching below the fermi energy, and the isolated narrow band has trivial quantum valley topological number Po et al. (2018a); Zhang et al. (2019); Zhang and Senthil (2018b); Chittari et al. (2019); Lee et al. (2019); Liu and Dai (2019). This minimal model would then naturally predict either a spin-triplet Xu and Balents (2018) or spin-singlet Dodaro et al. (2018) $d+id$ topological superconductor, depending on the sign of the on-site Hund’s coupling.

Evidence of spin-triplet pairing predicted previously Xu and Balents (2018) was recently found in TDBG Kim (2019); Liu et al. (2019). Evidence of ferromagnetic correlated insulator at half filling away from charge neutrality was discovered in the same system Kim (2019); Shen et al. (2019); Liu et al. (2019). In TDBG, besides clear ferromagnetism signature observed at the $1/2$-filling insulator Kim (2019); Shen et al. (2019); Liu et al. (2019), it was also observed that correlated insulators at $1/4$ and $3/4$ fillings emerge under inplane magnetic field Liu et al. (2019), whose main effect is likely just a spin-polarizing Zeeman effect. This observation implies that the TDBG at $1/4$ and $3/4$ filling is rather close to a ferromangetic correlated insulator, and a Zeeman field would drive a metal-insulator transition. Evidence of ferromagnetism at the insulating phase with 3/4 filling away from the charge neutrality was also found in another Moiré system Sharpe et al. (2019).

Motivated by these experiments, in this work we investigate a quantum spin-valley model on the triangular lattice with one fermion per site, which corresponds to either $1/4$ filling or $3/4$ filling on the Moiré superlattice. The Hamiltonian of this model reads

[TABLE]

where $\sigma^{a}$ and $\tau^{b}$ are Pauli operators in the spin and valley spaces, and $T^{ab}=\sigma^{a}\otimes\tau^{b}$. When $J^{s}=J^{v}=J$, this model becomes the ${\rm SU}(4)$ quantum antiferromagnetic model with fundamental representation on each site. The ${\rm SU}(4)$ symmetry is broken by the Hund’s coupling Xu and Balents (2018), which in general makes $J^{v}>J>J^{s}$, if we choose the standard sign of the Hund’s coupling which favors large spin on each site. But we assume that the ${\rm SU}(4)$ breaking effect is not strong enough to change the sign of $J^{s}$, $J^{v}$ and $J$, namely we keep all three coupling constants positive, $i.e.$ antiferromagnetic. Indeed, since the Hund’s coupling originates from the exchange coupling which involves overlap between wave functions at the two valleys, the Hund’s coupling should be a relatively weak effect since the inter-valley wave function overlap is expected to be small because large momentum transfer between the two valleys is suppressed by the long wavelength modulation of the Moiré superlattice. For simplicity we ignore other mechanisms that break the ${\rm SU}(4)$ symmetry, such as valley-dependent hopping Po et al. (2018a), hence in the spin-valley model Eq. 1 the valley space has its own ${\rm SU}(2)^{v}$ symmetry.

II Derivation of the Spin-Valley model

The model Eq. 1 can be derived in the standard perturbation theory starting with a Hubbard model plus an on-site Hund’s coupling on the triangular lattice. As we mentioned in the introduction, this model at least applies to TLG/h-BN Chen et al. (2018, 2019) and TDBG Kim (2019); Shen et al. (2019); Liu et al. (2019) with certain twisted angle and displacement field, since in these cases there is no symmetry protected band touching below the fermi energy, and the narrow band with correlated physics has trivial quantum valley topological number (although the exact valley topological number computed depends on the models used in the literature) Po et al. (2018a); Zhang et al. (2019); Zhang and Senthil (2018b); Chittari et al. (2019); Lee et al. (2019); Liu and Dai (2019):

[TABLE]

where $c_{j,\alpha}^{\dagger},c_{j,\alpha}$ are electron creation annihilation operators which have four flavors including both the spin and valley indices, $\hat{n}_{j}=\sum_{\alpha}c_{j,\alpha}^{\dagger}c_{j,\alpha}$ is the total particle number per site, and $\frac{1}{2}\hat{\sigma}_{j}^{a}=\frac{1}{2}c_{j}^{\dagger}\sigma^{a}c_{j}$ is the total on-site spin operator. The $H_{U}$ term is the on-site repulsive interaction, and we assume the most natural sign of the Hund’s coupling with $V>0$, as a result of exchange interaction. We treat the kinetic term $H_{t}$ perturbatively. This amounts to integrating out the charge degree of freedom to obtain an effective spin-valley model in the correlated insulator phase.

We follow the standard approach of degenerate perturbation theory. At quarter-filling, the ground state of $H_{U}+H_{V}$ has precisely one electron per site, and the projection operator to the ground state manifold reads

[TABLE]

Considering any pair of nearest neighbor sites on the Moiré superlattice, the ground state manifold can be further divided into four sectors which correspond to spin-singlet/triplet and valley-singlet/triplet states. We can write

[TABLE]

where (for example) $\mathcal{P}_{st}$ means the projection to spin-single/valley-triplet states. The effective Hamiltonian can be calculated as

[TABLE]

where $E_{0}$ is the ground state energy for the two-site problem. A detailed analysis of the intermediate states considering the virtual hopping process can be found in Ref. Xu and Balents, 2018. We find that only $\mathcal{P}_{st}$ and $\mathcal{P}_{ts}$ contribute to the effective Hamiltonian which takes a diagonal form in this basis

[TABLE]

Rewritten in terms of the ${\rm SU}(4)$ generators on the nearest neighbor sites, the effective Hamiltonian is equivalent to the spin-valley model Eq. 1 with the coupling constants given by

[TABLE]

There is a $Z_{2}$ symmetry regarding the sign of the Hund’s coupling. The coupling constants transform as $J^{s}/J\leftrightarrow J^{v}/J$ when we change $V\leftrightarrow-V$, as is naturally expected from the form of the Hund’s coupling Eq. 5.

III The $\mathrm{FM}\otimes 120^{\circ}$ state

At least in certain limit, $i.e.$ $J^{v}\gg J\gg J^{s}>0$, it is fairly easy to see why ferromagnetism would emerge in model Eq. 1 with all antiferromagnetic coupling constants. First of all, the following state will always be an eigenstate of the Hamlitonian:

[TABLE]

This state is a direct product of two parts: the first part is a fully-polarized ferromagnetic state of the spin $\vec{\sigma}_{i}$ space; the second part is the ground state of the nearest-neighbor antiferromagnetic quantum Heisenberg model on the triangular lattice in the $\vec{\tau}_{i}$ space. Although we cannot write down the explicit form of the exact microscopic wave-function $|\mathrm{AF}\ \mathrm{of}\ \vec{\tau}_{\{i\}}\ \rangle$, we do know that this state has a $120^{\circ}$ antiferromagnetic order with reduced moment due to quantum fluctuation and geometric frustration. This state Eq. 15 is always the eigenstate of Eq. 1 because a fully polarized ferromagnetic spin state is the eigenstate of operator $\vec{\sigma}_{i}\cdot\vec{\sigma}_{j}$ on every link $<i,j>$. Then in the limit of $J^{v}\gg J\gg J^{s}>0$, this eigenstate $|\Psi_{\mathrm{FM}}\rangle$ is also the ground state, because intuitively on every link the spin $\vec{\sigma}_{i}$ will see a background “effective” ferromagnetic coupling

[TABLE]

Because $\langle\vec{\tau}_{i}\cdot\vec{\tau}_{j}\rangle<0$ for the $120^{\circ}$ state of $\vec{\tau}_{i}$, for large enough $J$ the spins will see an effective ferromagnetic coupling, even though in the original model Eq. 1 all the couplings are antiferromagnetic.

With a fixed large $J^{v}$, while increasing $J^{s}$, eventually $J\langle\vec{\tau}_{i}\cdot\vec{\tau}_{j}\rangle$ will not be strong enough to overcome the antiferromagnetic coupling $J_{s}$, hence we expect to see a transition from the “$\mathrm{FM}\otimes 120^{\circ}$” state to another state without ferromagnetic order. Numerically Bernu et al. (1993), $\langle\vec{\tau}_{i}\cdot\vec{\tau}_{j}\rangle$ is found to be $\sim-0.73$ for the triangular lattice quantum antiferromagnet. If we evaluate the energy of $|\Psi_{\mathrm{FM}}\rangle$ in Eq. 15, while increasing $J_{s}/J$, this state is no longer the ground state when $J_{s}/J>0.73$. Hence the intuitive argument gives an upper bound for the transition point of $J^{s}/J$.

—DMRG simulation of the spin-valley model

We now provide numerical evidence for the fact that the $\mathrm{FM}\otimes 120^{\circ}$ state is the ground state in the $J^{v}\gg J\gg J^{s}$ limit and obtain the phase boundary of this state. To this end we use the density matrix renormalization group (DMRG) method White (1992); Schollwöck (2005). We note that in finite systems, boundaries can introduce strong oscillations in the expectation value of $\langle\vec{\tau}_{i}\cdot\vec{\tau}_{j}\rangle$ on nearest-neighbor bonds (which are expected to be uniform when the valley degree of freedom is in the $120^{\circ}$ state), and thus affect the effective coupling seen by the spins. To avoid such boundary effects we use infinite DMRG McCulloch (2008). To observe uniform bond expectation values when the valley degree of freedom is in the $120^{\circ}$ state wide enough cylinders have to be considered. We perform our analysis on cylinders of circumference $L_{y}=6$, for which we obtain mean bond expectation value $\langle\vec{\tau}_{i}\cdot\vec{\tau}_{j}\rangle\approx-0.74$, consistent with Ref. Bernu et al., 1993, with spatial variations below one percent.

For our numerical simulations we use the ITensor library ITe . We assume a 3-site unit cell along the cylinder to allow for the formation of a $120^{\circ}$ state in the valley and/or spin degrees of freedom. The valley degree of freedom is initialized in the total $\tau^{z}=0$ sector and $\tau^{z}$ quantum number conservation is used. The spin degree of freedom is initialized either in the fully polarized state, or a classical anti-ferromagnetic state with total $\sigma^{z}=0$. The maximal bond dimension in our simulations is $M=1000$.

We find that at large $J^{v}$ and small $J^{s}$, the system indeed converges to a fully polarized spin-FM and a $120^{\circ}$-valley ordered state, independent of the initial conditions. At larger $J^{s}$ we observe a state without any net magnetization. At this stage we cannot conclude as to the nature of the entire region with zero magnetization (with larger $J^{s}$ and $J^{s}<J$, for the total $\sigma^{z}=0$ sector and going to bond dimensions of up to 4000 we have not identified a clear order for the spin degree of freedom), but in the next section we will propose some possible interesting liquid states and topological orders for this region of phase diagram. The ground state energy obtained using iDMRG for a fixed $J^{v}/J=7$ as function of $J^{s}$, for the two initial states, is shown in the inset of Fig. 1. The solid blue line on the same plot indicates the energy expected for the $\rm{FM}\otimes 120^{\circ}$ state, calculated using a uniform bond expectation value $\langle\vec{\tau}_{i}\cdot\vec{\tau}_{j}\rangle\approx-0.74$ that we obtain for the $120^{\circ}$ state on infinite cylinders of circumference $L_{y}=6$ as mentioned above. We estimate the position of the phase boundary for each $J^{v}/J$ to be at the $J^{s}/J$ for which the lowest energy obtained using iDMRG drops below the one expected for the $\rm{FM}\otimes 120^{\circ}$ state. The resulting phase boundary as function of $J^{v}/J$ and $J^{s}/J$ is shown in the main Fig. 1. Our results for the ground state energy and the magnetization across the phase boundary both suggest that the transition between the ferromagnetic order and the paramagnet is a first order level-crossing.

—Schwinger boson analysis

We can also construct the $\mathrm{FM}\otimes 120^{\circ}$ state using the Schwinger boson formalism. We first define a four component Schwinger boson $b_{j,\alpha}$ on every site which forms fundamental representation under both the spin and valley ${\rm SU}(2)$ symmetry, and also a fundamental representation of the enlarged ${\rm SU}(4)$ symmetry. The Schwinger boson Hilbert space is subject to a local constraint

[TABLE]

Physically $\kappa=1$, but in the Schwinger boson mean field calculation $\kappa$ is often treated as a tuning parameter. The Hamiltonian Eq. 1 can be reorganized into the following form:

[TABLE]

The operator $\vec{\Delta}^{ts}_{ij}$ and $\vec{\Delta}^{st}_{ij}$ are the spin-triplet/valley-singlet, and spin-singlet/valley-triplet pairing operator between Schwinger boson $b_{\alpha}$ on site $i,j$:

[TABLE]

where $\sigma^{ab}=\sigma^{a}\otimes\tau^{b}$, $\sigma^{0}=\tau^{0}=\mathbf{1}_{2\times 2}$. Then $\Delta^{ss}$ and $\Delta^{tt}$ are the singlet/singlet, and triplet/triplet pairing respectively, for example $\Delta^{ss}_{ij}=b^{t}_{i}\sigma^{22}b_{j}$.

In Eq. 20,

[TABLE]

With the most natural parameter region $J_{v}>J>J_{s}>0$, $J_{ts}$ is always negative, and it corresponds to the strongest mean field channel, while none of the other parameters are guaranteed to be negative (for example $J_{tt}$ is always positive). Thus for the purpose of mean field analysis, we will just keep the first term of Eq. 20, and ignore all the rest three terms of Eq. 20. The mean field Hamiltonian reads

[TABLE]

where $\vec{\phi}$ is a complex vector under ${\rm SU}(2)^{s}$. Here we choose a uniform ansatz of $\vec{\phi}$ on the entire lattice, and all the links $i,j$ are included in the sum with the convention $j=i+\hat{e}$ with $\hat{e}=(1,0)$, $(-1/2,\pm\sqrt{3}/2)$, $i.e.$ the mean field ansatz explicitly preserves the translation and rotation by $2\pi/3$ symmetry of the triangular lattice, while all the crystal symmetries are preserved as projected symmetry group (PSG). $\mu$ is another variational parameter of the mean field calculation which guarantees that the filling of Schwinger boson is fixed at $\kappa$ on every site.

The $\mathrm{FM}\otimes 120^{\circ}$ ordered state corresponds to the mean field ansatz with $\vec{\phi}=\vec{\phi}_{1}+\mathrm{i}\vec{\phi}_{2}$, and the real vectors $\vec{\phi}_{1}$ and $\vec{\phi}_{2}$ orthogonal with each other. For example, when $\vec{\phi}=\phi(1,\mathrm{i},0)$, only the spin-up ($\sigma^{3}=+1$) Schwinger bosons participate in this mean field analysis. The ferromagnetic order parameter corresponds to the following gauge invariant quantity:

[TABLE]

With only spin-up Schwinger bosons, the mean field calculation reduces precisely to the SU(2) spin-1/2 Heisenberg model on the triangular lattice Sachdev (1992) with Heisenberg coupling $J_{ij}=-4J_{ts}$ (The $\mathrm{Sp}(N)$ Heisenberg model defined in Ref. Sachdev, 1992 has Hamiltonian $H=\sum_{i,j}-\frac{1}{2N}J_{ij}\Delta^{\dagger}_{ij}\Delta_{ij}$ where $\Delta_{ij}$ is the $\mathrm{Sp}(N)$ singlet pairing between Schwinger bosons on sites $i,j$):

[TABLE]

Because the spin-down Schwinger bosons do not contribute to the mean field decomposition when $\vec{\phi}\sim(1,\mathrm{i},0)$, we replace the constraint in Eq. 17 by $\sum_{\alpha=1}^{2}b^{\dagger}_{j,\uparrow,\alpha}b_{j,\uparrow,\alpha}=\kappa$ in Eq. 31. Now technically the mean field theory Eq. 31 corresponds to the “zero-flux state” in Ref. Wang and Vishwanath, 2006, which has lower mean field energy than other mean field ansatz Wang and Vishwanath (2006) for this nearest neighbor model, and it makes the minima of the Schwinger boson band structure locate at the corner of the Brillouin zone $\vec{Q}=(\pm 4\pi/3,0)$. The mean field solution gives $\mu>0$, which is consistent with the fact that we set $\sum_{\alpha=1}^{2}b^{\dagger}_{j,\downarrow,\alpha}b_{j,\downarrow,\alpha}=0$. And at the mean field level, when the filling of the Schwinger boson $\kappa$ is greater than 0.34 Sachdev (1992), $b_{\alpha}$ condenses, which leads to a fully polarized FM in the spin space, and also $120^{\circ}$ state in the valley space.

If the mean field value of $\vec{\phi}$ is real (or equivalently if $\vec{\phi}_{1}$ is parallel to $\vec{\phi}_{2}$, for example, $\vec{\phi}\sim\phi(0,0,1)$, both spin-up and spin-down Schwinger bosons participate in the mean field analysis, and the mean field analysis is technically equivalent to the calculations in Ref. Sachdev, 1992 for the $\mathrm{Sp}(2)\sim\mathrm{SO}(5)$ antiferromagnet on the triangular lattice also with Heisenberg coupling $J_{ij}=-4J_{ts}$, because $\left(\vec{\Delta}^{ts}_{ij},\vec{\Delta}^{st}_{ij}\right)$ together form a $\mathrm{SO}(6)$ vector, and condensing each component of the vector breaks the $\mathrm{SO}(6)$ down to $\mathrm{SO}(5)\sim\mathrm{Sp}(2)$. Each component of $\vec{\Delta}_{ts}$ can be viewed as the $\mathrm{Sp}(2)$ singlet introduced in Ref. Sachdev, 1992:

[TABLE]

Quoting the results in Ref. Sachdev, 1992, the $\mathrm{FM}\otimes 120^{\circ}$ state with the previous mean field ansatz with $\phi_{1}\perp\phi_{2}$ has a lower mean field ground state energy density, which is consistent with our analytical observation and also numerical simulation.

—Zeeman field driven Metal-Insulator transition

Since the insulator has a fully polarized ferromagnetic order, its energy can be tuned by an external Zeeman field. An inplane magnetic field, whose main effect is the Zeeman coupling can drive a first order metal-insulator transition (a level-crossing) between the unpolarized metal and the fully polarized ferromagnetic insulator, as was observed experimentally at the 1/4 and 3/4 filling of TDBG Kim (2019); Liu et al. (2019).

There is another possible mechanism of metal-insulator transition driven by a Zeeman field. At the metalic side at the transition, the system is likely described by a $t-J$ model with a similar $J,J^{s},J^{v}$ terms as Eq. 1. The Zeeman field tends to polarize the spin, which effectively increases the antiferromagnetic coupling in the valley space $J^{v}_{eff}=J^{v}+J\langle\vec{\sigma}_{i}\cdot\vec{\sigma}_{j}\rangle$. Thus at certain temperature, the magnitude of the $120^{\circ}$ order in the valley space is tunable and enhanced by an external Zeeman field. If the insulating behavior of the system is a consequence of the finite momentum valley order which folds the Brillouin zone and partially gaps out the Fermi surface, an increasing magnitude of the $120^{\circ}$ order in the valley space can gap out larger portion of the Fermi surface, decrease the charge carrier density, and hence eventually drive a continuous metal insulator transition.

IV Liquids and topological phases

When $J^{v}\sim J^{s}\sim J$, it would be rather difficult for the system to form any semiclassical order due to “double frustration”: the $J^{s}$ and $J^{v}$ term of Eq. 1 are both already frustrated due to the geometry of the triangular lattice, while the $J$ term further frustrates/disfavors the simultaneous $120^{\circ}$ semiclassical order of $\vec{\sigma}_{i}$ and $\vec{\tau}_{i}$. Since there is an obvious Lieb-Shultz-Matthis theorem which forbids a completely trivial disordered phase, we expect this “double frustration” effect to lead to either a completely disordered spin-valley liquid state, or a partially ordered state with certain topological order. In this section we explore several possible spin-valley liquids or topological orders in the region $J^{v}\sim J^{s}\sim J$.

—Spin nematic $Z_{2}$ topological phase

Let us get back to the mean field Hamiltonian Eq. 27. As we discussed before, if the mean field value of $\vec{\phi}$ is real (or equivalently if $\vec{\phi}_{1}$ is parallel to $\vec{\phi}_{2}$, for example, $\vec{\phi}\sim(0,0,1)$, both spin-up and spin-down Schwinger bosons participate in the mean field analysis, and the mean field analysis is technically equivalent to the calculations in Ref. Sachdev, 1992 for the $\mathrm{Sp}(4)$ antiferromagnet on the triangular lattice. And with large spin symmetry, the quantum fluctuation makes it more difficult for $b_{\alpha}$ to condense. If $b_{\alpha}$ is not condensed, the mean field order parameter $\vec{\phi}$ already breaks the ${\rm SU}(2)^{s}$, and also break the ${\rm U}(1)$ gauge symmetry down to $Z_{2}$ gauge degree of freedom.

The nature of the state with condensed $\vec{\phi}$ but uncondensed $b_{\alpha}$ depends on the nature of $\vec{\phi}$ under time-reversal. The transformation of $b_{\alpha}$ under time-reversal can be inferred by the fact that $\vec{\sigma}\rightarrow-\vec{\sigma}$, $(\tau^{1},\tau^{2},\tau^{3})\rightarrow(\tau^{1},\tau^{2},-\tau^{3})$:

[TABLE]

as long as $\vec{\phi}$ is a real vector (or $\vec{\phi}_{1}$ parallel with $\vec{\phi}_{2}$), time-reversal is preserved, and this state is a spin nematic $Z_{2}$ topological order. By contrast, if $\vec{\phi}_{1}\perp\vec{\phi}_{2}$ the time-reversal is broken.

—$Z_{2}\times Z_{2}$ spin-valley liquid

More states can be constructed by introducing two flavors of Schwinger bosons $b^{s}_{j,\alpha}$ and $b^{v}_{j,\alpha}$ for the spin and valley space on each site respectively, which are subject to the constraint

[TABLE]

These two constraints introduces two ${\rm U}(1)$ gauge symmetries. It is fairly straightforward to construct the $\mathrm{FM}\otimes 120^{\circ}$ state using this type of Schwinger bosons: $b^{s}_{\alpha}$ condenses at zero momentum, while simultaneously $b^{v}_{\alpha}$ condenses at the corner of the Brillouin zone.

In fact, due to the “double frustration” effect, both the spin and valley space can form a $Z_{2}$ topological order (overall speaking the system is in a $Z_{2}\times Z_{2}$ spin-valley liquid state), whose $e$ particles carry the fundamental representation of ${\rm SU}(2)^{s}$ and ${\rm SU}(2)^{v}$ respectively, as long as neither $b^{s}_{j,\alpha}$ nor $b^{v}_{j,\alpha}$ introduced in Eq. 36 condenses when the mean field parameters break both ${\rm U}(1)$ gauge symmetries down to $Z_{2}$.

Starting from the $Z_{2}\times Z_{2}$ spin-valley liquid state, one can also construct a spin-valley liquid with only one $Z_{2}$ topological order. This can be formally obtained by forming bound state of the “visons” (the $m$ excitations) of both $Z_{2}$ topological orders, and condense the bound state. This condensate will confine $b^{s}_{\alpha}$ and $b^{v}_{\alpha}$ separately, but their bound state is still deconfined, and becomes the $e$ particle of the new $Z_{2}$ topological order. This final $Z_{2}$ topological order preserves all the symmetries of the system, and it can also be constructed using the same mean field formalism as Eq. 20, as long as one condenses the spin-singlet/valley-singlet pairing operator $\Delta^{ss}_{ij}$ in Eq. 20.

—${\rm U}(1)\times{\rm U}(1)$ Dirac spin-valley liquid

More exotic spin-valley liquid states can be constructed by introducing fermionic slave particles $f^{s}_{j,\alpha}$ and $f^{v}_{j,\alpha}$ which are subject to the constraints

[TABLE]

In Ref. Lu, 2016, a Dirac spin liquid with ${\rm U}(1)$ gauge field and $N_{f}=4$ flavors of Dirac fermions was constructed for spin-1/2 systems on the triangular lattice. And this Dirac spin liquid is the parent state of both the $120^{\circ}$ ordered state and the valence bond solid state Lu (2016); Song et al. (2018a, b), and it could be a deconfined quantum critical point between these two different ordered states Jian et al. (2018).

In our case, both spin and valley space can form the Dirac liquid phase mentioned above, due to the double frustration effect. Thus in total there are eight flavors of Dirac fermions and two ${\rm U}(1)$ gauge fields.

— The SU(4) point

At the point $J^{v}=J^{s}=J$, this model has a ${\rm SU}(4)\sim\mathrm{SO}(6)$ symmetry. Although semiclasical approach such as nonlinear sigma model were studied before for ${\rm SU}(N)$ antiferromagnet with other representations Read and Sachdev (1989), with a fundamental representation on every site, this model has no obvious semiclassical limit to start with, and it is expected to be a nontrivial spin liquid or topological order. At this point, it is most convenient to define a four component Schwinger boson $b_{j,\alpha}$ on every site which forms fundamental representation under both the spin and valley ${\rm SU}(2)$ symmetry, and there is a constraint $\sum_{\alpha=1}^{4}b^{\dagger}_{j,\alpha}b_{j,\alpha}=1$.

Unlike a ${\rm SU}(2)$ spin system, one can prove that at the ${\rm SU}(4)$ point there cannot be a fully symmetric $Z_{2}$ spin liquid whose $e$ particle is the $b_{\alpha}$ slave particle. The reason is that all the local spin excitations can be written as $b^{\dagger}_{j,\alpha}b_{j,\beta}$ with different $\alpha,\beta=1\cdots 4$, hence all the local spin excitations are invariant under the $Z_{4}$ center of the ${\rm SU}(4)$ group. In a $Z_{2}$ topological order, two of the $e$ particles should merge into a local excitations, while two $b_{\alpha}$ slave particle cannot fuse into a representation that is invariant under the $Z_{4}$ center. This argument also shows that a $Z_{2}$ topological order whose $e$ particle is a $\mathrm{SO}(6)$ vector is allowed.

On the other hand, using the slave particle $b_{\alpha}$ one can construct a $Z_{2}$ topological order with certain spontaneous ${\rm SU}(4)$ symmetry breaking. At the ${\rm SU}(4)$ point, the model Eq. 1 can be written as

[TABLE]

where $\vec{\Delta}_{ij}$ is a six component vector pairing between $b_{\alpha}$. One can introduce a six component complex $\mathrm{SO}(6)$ vector mean field parameter $\vec{\phi}$:

[TABLE]

The complex vector $\vec{\phi}=\vec{\phi}_{1}+\mathrm{i}\vec{\phi}_{2}$, where its real and imaginary parts $\vec{\phi}_{1}$ and $\vec{\phi}_{2}$ can be either parallel or orthogonal to each other. If the Schwinger boson does not condense, both mean field theories would lead to a $Z_{2}$ topological order on top of the spontaneous ${\rm SU}(4)$ symmetry breaking.

Our DMRG simulation actually suggest that the ${\rm SU}(4)$ point of the spin-valley model is a spin-valley liquid state with a Fermi surface of fermionic slave particles, which will be presented in detail in another work.

V Conclusion

In this work we demonstrated both analytically and numerically that a quantum spin-valley model with all antiferromagnetic interaction can have a fully polarized ferromagnetic order in its phase diagram. We propose possible mechanism for an inplane Zeeman field to drive a metal-insulator transition, as was observed experimentally at the 1/4 and 3/4 filling of TDBG. We also discussed various possible nontrivial spin-valley liquid state and topological order of this model.

We would like to acknowledge several previous theoretical works that studied the ferromagnetism in Moiré systems using different approaches and different models Kang and Vafek (2018b); Seo et al. (2018); Bultinck et al. (2019); Wu and Sarma (2019). For example, in Ref. Seo et al., 2018, a spin-valley model with ferromagnetic couplings on an effective honeycomb Moiré lattice was derived for the twisted bilayer graphene system. While our work (which aims to understand a different Moiré system, $i.e.$ the twisted double bilayer graphene) demonstrated that ferromagnetism can emerge from the spin-valley model on a triangular lattice with fully antiferromagnetic interaction.

Within our framework, under doping, again the system is likely described by the $t-J$ model with the similar $J,J^{s},J^{v}$ terms as Eq. 1. Then the analysis in Ref. Xu and Balents, 2018 still applies: the spin-triplet/valley-singlet pairing channel between electrons would become the strongest pairing channel. Due to a strong on-site Hubbard interaction, the system would still prefer to become a $d+id$ topological superconductor with spin triplet pairing.

Anna Keselman and Chao-Ming Jian are supported by the Gordon and Betty Moore Foundations EPiQS Initiative through Grant GBMF4304. Cenke Xu is supported by the David and Lucile Packard Foundation. Use was made of the computational facilities administered by the Center for Scientific Computing at the CNSI and MRL (an NSF MRSEC; DMR-1720256) and purchased through NSF CNS-1725797.