# Floquet Fractional Chern Insulators and Competing Phases in Twisted Bilayer
Graphene

Peng-Sheng Hu  Zhao Liu  Zhejiang Institute of Modern Physics, Zhejiang
University, Hangzhou 310027, China

###### Abstract

We study the many-body physics in twisted bilayer graphene coupled to periodic
driving of a circularly polarized light. In the limit of high driving
frequency $\Omega$, we use Floquet theory to formulate the system by an
effective static Hamiltonian truncated to the order of $\Omega^{-2}$, which
consists of a single-electron part and the screened Coulomb interaction. We
numerically simulate this effective Hamiltonian by extensive exact
diagonalization. Remarkably, under appropriate twist angle and driving
strength, we identify Floquet fractional Chern insulator states in the
partially filled Floquet valence band. We characterize these topologically
ordered states by ground-state degeneracy, spectral flow, and entanglement
spectrum. We also find topologically trivial charge density waves and
interaction-induced Fermi liquid which strongly compete with fractional Chern
insulator states.

##  I Introduction

Van-der-Waals heterostructures with moiré patterns [1, 2] have attracted
tremendous attention over the last few years. When two atomic layers are
stacked with each other, the mismatch between two crystals due to different
lattice constants and/or a twist angle generates a large-scale superlattice
and affects the interlayer coupling. Such moiré systems, whose band structures
are highly controllable, are promising to host numerous exotic phenomena. A
representative example is the twisted bilayer graphene (TBG), consisting of
two sheets of monolayer graphene with a twist angle in between [3]. One of the
most striking features of TBG is that, at certain special twist angles (called
magic angles), the low-energy bands near charge neutrally point (CNP) can be
tuned to be flat, thus providing an ideal platform to investigate correlated
physics. Indeed, unconventional superconductivity and ferromagnetism have been
discovered in TBG [4, 5, 6, 7, 8, 9]. Another salient progress in this
direction is the observation of the quantum anomalous Hall effect at zero
magnetic field (also called Chern insulators) in TBG aligned with a hexagonal
boron nitride (hBN) [10]. In this case, hBN gaps out the protected Dirac
points of TBG, such that the flat bands around CNP is isolated and acquire
nonzero Chern numbers [11, 12, 13, 6, 14, 10]. Motivated by this topological
band structure, a lot theoretical efforts have been made to demonstrate the
possibility of realizing the fractional Chern insulators (FCIs) [15, 16, 17]
in TBG-hBN when the flat bands near CNP are partially filled by interacting
electrons [18, 19, 20, 21]. Excitingly, evidence of FCIs in TBG-hBN has been
reported in a recent experiment at weak magnetic fields [22].

External driving fields of light provide an alternative avenue to obtain
topological band structures even if the original bands in the absence of
driving is topologically trivial [23, 24, 25, 26, 27, 28, 29, 30]. Unlike
nondriven systems, systems periodically driven by light do not have well-
defined static band structures. However, the effect of light enters an
effective static Hamiltonian that captures the dynamics of the system on time
scales much longer than the driving period [31, 24, 32, 33]. This effective
static Hamiltonian, which describes photo-dressed band structures, has been
extensively used to predict the out-of-equilibrium topological properties of
various periodically driven systems [23, 24, 34, 27, 28]. In particular,
accompanying the rapid development of moiré materials, exploring the effects
of light on moiré band structures has become an intriguing direction recently
[35, 36, 37, 38, 39, 40, 41, 42]. Just like in monolayer graphene [23, 24,
43], a circularly polarized light can open a band gap at the Dirac points of
TBG and give rise to Floquet topological flat bands [35, 36, 37]. In this
sense, for TBG the light driving plays the role of hBN in the static case.

While exciting progress has been made on Floquet moiré materials, the
interactions between electrons are not yet taken into account in previous
works [35, 36, 37, 38, 39, 40, 41, 42]. Hence whether correlated topological
phases can be induced by light driving in these systems remains a crucial open
question. In the context of photon-dressed topological band structure, a
natural candidate of correlated topological states is the Floquet FCI which
was first proposed in monolayer graphene [44, 45]. The Floquet FCI is formed
when interacting particles partially fill in a topological flat band of the
effective static Hamiltonian of the driven system.

In this work, to investigate the possibility of stabilizing Floquet FCIs in
moiré materials, we consider TBG driven by a circularly polarized laser light.
We focus on the limit of high driving frequency $\Omega$, in which case we
truncate the effective static Hamiltonian to the order of $\Omega^{-2}$. Then
we use exact diagonalization to search for the evidence of Floquet FCIs in the
parameter space spanned by the twist angle and the driving strength.
Strikingly, for electrons interacting via the screened Coulomb potential, we
find a wide region in which FCIs are promising to exist when the Floquet
valence band is occupied at filling $\nu=1/3$. We also observe charge density
waves (CDWs) and Fermi liquid (FL)-like states in the neighboring regions of
FCIs. The overall phase diagram is similar for different interlayer coupling
strength of TBG. Because the effective Hamiltonian of Floquet TBG in the high-
frequency limit is analogous to that of static TBG-hBN, our results could also
apply to the many-body physics in the later system.

##  II Model

We consider the low-energy dynamics of TBG at small twist angles ($\theta\sim
1^{\circ}$) by following Bistritzer and MacDonald’s continuum model [46].
After twisting, the original first Brillouin zone of single-layer graphene is
folded into moiré Brillouin zones (MBZ). We focus on the MBZ near the valley
$\mathbf{K}_{+}$ of the single-layer graphene (Fig. 1) and assume spin
polarization. In this MBZ, the $\mathbf{K}_{+}$ points of the top and bottom
graphene layers are located at ${\bf K}_{+}^{t}=R_{\theta/2}{\bf K}_{+}$ and
${\bf K}_{+}^{b}=R_{-\theta/2}{\bf K}_{+}$, respectively, with $R_{\theta}$ a
counter-clockwise rotation around the $z$-axis in the momentum space. ${\bf
G}_{1}$ and ${\bf G}_{2}$ are the reciprocal lattice vectors of the MBZ.

Figure 1:  Moiré Brillouin zone of TBG. Two large hexagons are original first
Brillouin zones of top and bottom graphene layers. The small hexagon
represents the MBZ resulting from twist. The $\mathbf{K}_{+}^{t,b}$ points,
the vectors ${\bf q}_{0,1,2}$ in Eq. (1) and the MBZ reciprocal lattice
vectors ${\bf G}_{1}$ and ${\bf G}_{2}$ are given.

###  II.1 Static system

In the absence of driving, the single-electron Hamiltonian of TBG in the
momentum space [47, 18] takes the form of

$\displaystyle H_{\text{kin}}=$
$\displaystyle\sum_{\mathbf{k}}\left(\psi_{t}^{\dagger}(\mathbf{k})h_{-\theta/2}\left(\mathbf{k}-\mathbf{K}_{+}^{t}\right)\psi_{t}(\mathbf{k})\right.$
(1)
$\displaystyle\left.+\psi_{b}^{\dagger}(\mathbf{k})h_{\theta/2}\left(\mathbf{k}-\mathbf{K}_{+}^{b}\right)\psi_{b}(\mathbf{k})\right)$
$\displaystyle+\sum_{\mathbf{k}}\sum_{j=0}^{2}\left(\psi_{t}^{\dagger}\left(\mathbf{k}-\mathbf{q}_{0}+\mathbf{q}_{j}\right)T_{j}\psi_{b}(\mathbf{k})+\text{H.c.}\right),$

where $\psi_{t}({\bf k})=\begin{pmatrix}\psi_{tA}(\mathbf{k})\\\
\psi_{tB}(\mathbf{k})\end{pmatrix}$ and $\psi_{b}({\bf
k})=\begin{pmatrix}\psi_{bA}(\mathbf{k})\\\
\psi_{bB}(\mathbf{k})\end{pmatrix}$ are spinors of annihilation operators for
electrons in top ($t$) and bottom ($b$) graphene layers, respectively, and $A$
and $B$ correspond to the two sublattices in single-layer graphene. The first
two terms in $H_{\text{kin}}$ are the Dirac Hamiltonians of top and bottom
graphene layers, for which $h_{\theta}(\mathbf{k})=h(R_{\theta}\mathbf{k})$
with $h(\mathbf{k})=-\hbar v_{F}(k_{x}\sigma_{x}+k_{y}\sigma_{y})$. Here
$\hbar v_{F}=\sqrt{3}at_{0}/2$, with $t_{0}=2.62\ {\rm eV}$ and $a=0.246\
\text{nm}$ the nearest-neighbor hopping amplitude and lattice constant of
single-layer graphene, respectively. The third term describes the moiré
tunneling between the two graphene layers. Such tunneling is encoded in the
matrix

$T_{j}=w_{0}-w_{1}e^{i(2\pi/3)j\sigma_{z}}\sigma_{x}e^{-i(2\pi/3)j\sigma_{z}}$
(2)

with the momenta
$\mathbf{q}_{0}=R_{-\theta/2}\mathbf{K}_{+}-R_{\theta/2}\mathbf{K}_{+}$,
$\mathbf{q}_{1}=R_{2\pi/3}\mathbf{q}_{0}$ and
$\mathbf{q}_{2}=R_{-2\pi/3}\mathbf{q}_{0}$ (Fig. 1). $w_{0}$ and $w_{1}$ in
$T_{j}$ are the tunneling strengths between $AA$ and $AB$ sites, respectively.
Ab initio numerics gives $w_{1}\approx 110\ {\rm meV}$ [48], however, some
works suggested smaller values [49, 50, 18, 21]. In this paper, we consider
two situations with $w_{1}=90\ \text{meV}$ and $w_{1}=110\ \text{meV}$ to
account for variations of different theoretical models and realistic samples.
Furthermore, we fix $w_{0}=0.7w_{1}$ to include the effects of lattice
relaxation [51, 52] and corrugation [53, 54, 50]. To calculate the band
structure at momentum ${\bf k}_{0}$ in the MBZ, we write ${\bf k}$ in
$H_{\text{kin}}$ as ${\bf k}_{0}+m{\bf G}_{1}+n{\bf G}_{2}$ with integers
$m,n=-d,...,d$ then diagonalize $H_{\text{kin}}$, where $d$ is a suitably
chosen cutoff. As there is no alignment with the hexagonal boron nitride (hBN)
substrate, the Dirac band touching exists at the corners of the MBZ.

We simulate the interaction between electrons via the screened Coulomb
potential

$H_{\text{int}}=\frac{1}{2}\sum_{\mathbf{q}}V(\mathbf{q}):\rho(\mathbf{q})\rho(-\mathbf{q}):,$
(3)

where $\rho(\mathbf{q})$ is the valley and spin projected density operator and
$:\ :$ means the normal order. We choose the Yukawa potential $V({\bf
q})=\frac{e^{2}}{4\pi\epsilon_{r}\epsilon_{0}S}\frac{2\pi}{\sqrt{|{\bf
q}|^{2}+\kappa^{2}}}$ to describe the screening, where $e$ is the electron
charge, $\epsilon_{r}$ is the relative dielectric constant of the material,
$\epsilon_{0}$ is the dielectric constant of vacuum, $S$ is the area of the
moiré superlattice, and $\kappa$ measures the screening strength. Throughout
this work, we fix $\epsilon_{r}=4$ [55, 56] and $\kappa=1/a_{M}$, with
$a_{M}=a/(2\sin(\theta/2))$ the lattice constant of TBG.

###  II.2 Floquet system

Now we consider the coupling of TBG with light in the scenario of periodic
driving. We suppose the system is driven by a circularly polarized light which
shines vertically and uniformly across the surface of TBG. The light field is
represented by an electric field rotating in-plane as
$\bm{\mathcal{E}}=\mathcal{E}_{0}\left(\sin\Omega t,\cos\Omega t\right)$,
where $\Omega$ is the driving frequency. The corresponding vector potential
${\bf A}=A_{0}\left(\cos\Omega t,-\sin\Omega t\right)$, satisfying
$\bm{\mathcal{E}}=-\frac{\partial{\bf A}}{\partial t}$. For the single-
electron Hamiltonian, the light field only affects the intralayer hopping.
This is because the interlayer tunneling is dominated by hopping between atoms
that are exactly on top of each other, thus mostly contributed by
$z$-component of the vector potential which is absent in our setup [36, 57,
35]. We include the effect of light using a Peierls substitution
$\mathbf{k}\rightarrow\mathbf{k}+e\mathbf{A}(t)/\hbar$ in the intralayer terms
of Eq. (1), resulting in a time-dependent single-particle Hamiltonian
$H_{\text{kin}}(t)$. The interaction Hamiltonian Eq. (3) remains as in the
static case since it is expressed in terms of the density operator [44, 58].
Combining both terms, we get a new time-periodic Hamiltonian

$H(t)=H_{\text{kin}}(t)+H_{\text{int}}.$ (4)

to describe our system irradiated by the circularly polarized light, where
$H(t)=H(t+2\pi/\Omega)$.

According to Floquet theory, the stroboscopic evolution of the system, upon a
unitary transformation, can be captured by an effective static Hamiltonian
$H_{\rm eff}$ that does not depend on initial conditions [31, 32, 33]. While
in general it is complicated to evaluate $H_{\rm eff}$, we consider the limit
where the driving frequency $\Omega$ is large compared to other characteristic
energy scales in the system. In this case, $H_{\rm eff}$ can be represented by
a series expansion of $1/\Omega$ [32, 45, 31, 58, 33]:

$H_{\rm eff}=H_{\rm eff}^{(0)}+H_{\rm eff}^{(1)}+H_{\rm eff}^{(2)}+\dots,$ (5)

with

$\displaystyle H_{\rm eff}^{(0)}=$ $\displaystyle H_{0},$ (6a) $\displaystyle
H_{\rm eff}^{(1)}=$
$\displaystyle\frac{1}{\hbar\Omega}\sum_{m=1}^{\infty}\frac{\left[H_{m},H_{-m}\right]}{m},$
(6c) $\displaystyle H_{\rm eff}^{(2)}=$
$\displaystyle\frac{1}{(\hbar\Omega)^{2}}\Bigg{\\{}\sum_{m=1}^{\infty}\frac{\left[H_{m},\left[H_{0},H_{-m}\right]\right]}{2m^{2}}+\sum^{\infty}_{m,m^{\prime}=1\atop
m\neq m^{\prime}}$
$\displaystyle\frac{\left[H_{-m^{\prime}},\left[H_{m^{\prime}-m},H_{m}\right]\right]-\left[H_{m^{\prime}},\left[H_{-m^{\prime}-m},H_{m}\right]\right]}{3mm^{\prime}}$
$\displaystyle+\text{H.c.}\Bigg{\\}}.$

Here $H_{m}$ is the Fourier transform of $H(t)$, i.e.,
$H(t)=\sum_{m}H_{m}e^{im\Omega t}$. For our model, $H_{m}$ is nonzero only
when $m=0,\pm 1$.

The zeroth order term $H_{\rm eff}^{(0)}=H_{0}$ is just the static Hamiltonian
$H_{\rm kin}+H_{\rm int}$. The first order term is

$\displaystyle H_{\rm
eff}^{(1)}=\frac{(eA_{0}v_{F})^{2}}{\hbar\Omega}\sum_{\mathbf{k}}\left(\psi_{t}^{\dagger}(\mathbf{k})\sigma_{z}\psi_{t}(\mathbf{k})+\psi_{b}^{\dagger}(\mathbf{k})\sigma_{z}\psi_{b}(\mathbf{k})\right),$

which is the same as that derived for non-interacting TBG [37, 35, 36, 59]
because the interaction is time-independent in our model. $H_{\rm eff}^{(1)}$
is a single-particle term which introduces a staggered potential of strength
$P=(eA_{0}v_{F})^{2}/(\hbar\Omega)=(3a^{2}t_{0}^{2}e^{2}\mathcal{E}_{0}^{2})/(4\hbar^{3}\Omega^{3})$
in both graphene layers. As this stagger potential breaks the $C_{2}$
sublattice symmetry, it gaps the Dirac cones of the static TBG. Then the two
flat bands near the CNP may be isolated and carry non-zero Chern number in
specific range of $P$ and twist angle $\theta$ [37, 35]. Note that $H_{\rm
eff}^{(1)}$ plays a similar role to the alignment with an hBN substrate,
which, to the lowest order, also introduces a staggered potential of about
$15\ {\rm meV}$ [14]. The advantage of light driving compared to hBN is that
we can easily control the strength of this staggered potential by tuning the
electric field $\mathcal{E}_{0}$ or frequency $\Omega$.

Previous works studying high-frequency driving in non-interacting lattice
models often neglect $H_{\rm eff}^{(2)}$ and other higher-order terms in
$H_{\rm eff}$, because their corrections to the single-particle Hamiltonian is
quite small. However, once interactions are considered, one should be very
careful when dealing with these high-order terms, because they can include
effective interactions even though the original interaction $H_{\rm int}$ is
time-independent. To the leading order, these effective interactions are
present in $H_{\rm eff}^{(2)}$ if $\left[H_{m},\left[H_{\rm
int},H_{-m}\right]\right]\neq 0$. While being much weaker than $H_{\rm
eff}^{(0)}$ and $H_{\rm eff}^{(1)}$, these effective interactions may still be
comparable to the many-body gap protecting the ground state of $H_{\rm eff}$,
thus having essential influences on the low-energy stroboscopic physics.
Indeed, it was found in some Floquet topological lattice models that the
effective interactions in $H_{\rm eff}^{(2)}$ led by original density-density
repulsions destabilize topologically ordered FCIs [58]. In our model, we
carefully evaluate $\left[H_{1},\left[H_{\rm int},H_{-1}\right]\right]$.
Remarkably, we find it is zero due to the special forms of our $H_{\pm 1}$ and
$H_{\rm int}$ (see Appendix). Therefore, by contrast to Ref. [58], $H_{\rm
eff}^{(2)}$ in our model is still a single-particle correction without
effective interactions:

$\displaystyle H_{\rm eff}^{(2)}$ $\displaystyle=$
$\displaystyle\frac{(eA_{0}v_{F})^{2}}{(\hbar\Omega)^{2}}\Bigg{[}-H_{\rm
kin}+$
$\displaystyle\sum_{\mathbf{k}}\sum_{j}\left(\psi_{t}^{\dagger}\left(\mathbf{k}-\mathbf{q}_{0}+\mathbf{q}_{j}\right)\mathbb{W}_{\theta}\psi_{b}(\mathbf{k})+\text{H.c.}\right)\Bigg{]},$

where

$\mathbb{W}_{\theta}=w_{0}\begin{pmatrix}e^{i\theta}&0\\\
0&e^{-i\theta}\end{pmatrix}.$ (9)

In the following, we choose a high frequency $\hbar\Omega=1.5\ {\rm eV}$,
which significantly exceeds the energy scale of low-energy bands in static
TBG. When $\Omega$ is fixed, $P$ quantifies the driving strength. We consider
$P$ up to $60\ {\rm meV}$, corresponding to a strong electric field
$\mathcal{E}_{0}\approx 8\ \rm{MV/cm}$. Then the prefactor
$\frac{(eA_{0}v_{F})^{2}}{(\hbar\Omega)^{2}}$ of $H_{\rm eff}^{(2)}$ is only
$\sim 4\%$ at most. In this sense, our Floquet model could also simulate the
static TBG aligned with hBN both on the top and at the bottom.

##  III Floquet band structure

Now we have obtained the effective static Hamiltonian

$\displaystyle H_{\rm eff}=H_{\rm kin}+H_{\rm int}+H_{\rm eff}^{(1)}+H_{\rm
eff}^{(2)}$ (10)

describing the stroboscopic nature of our Floquet system, which includes the
original interaction and a single-particle part corrected by driving. Before
we dive into the interaction induced many-body physics, let us first analyze
the properties of the Floquet bands. In the following, we focus on the valence
band below the CNP.

In Fig. 2, we present the indirect band gap $\Delta_{\rm SP}$ and bandwidth
$W_{\rm SP}$ of the valence band as functions of $\theta$ and $P$. For most
parameters that we consider, we find positive $\Delta_{\rm SP}$, meaning that
the Floquet valence band is isolated from other bands at these parameters.
However, there are some lines in the $(\theta,P)$ space along which the band
gap vanishes. After calculating the Chern number $C$ of the Floquet valance
band, we find these gap-vanishing lines correspond to the transition between
$C=-1$ and $C=0$ (Fig. 3). In the parameter range that we consider, the
largest band gaps appear in the topological $C=-1$ case.

Figure 2:  The indirect band gap $\Delta_{\text{SP}}$ and the bandwidth
$W_{\rm SP}$ of the Floquet valence band for $w_{1}=90\ {\rm meV}$ [(a),(c)]
and $w_{1}=110\ {\rm meV}$ [(b),(d)].  Figure 3:  Chern number of the Floquet
valence band for (a) $w_{1}=90\ {\rm meV}$ and (b) $w_{1}=110\ {\rm meV}$.
Solid dots represent $C=-1$, and circles represent $C=0$.

##  IV Many-body physics

Because the Floquet valence band is well isolated and carries Chern number
$C=-1$ in a wide range of parameters $\theta$ and $P$ (Figs. 2 and 3), it is
promising to host Floquet FCIs. Now we consider the situation in which this
band is partially filled by $N$ interacting electrons to confirm this
possibility. We choose a finite periodic system with $N_{1}$ and $N_{2}$ moiré
unit cells in the directions of the two primitive moiré lattice vectors, such
that the band filling factor $\nu$ is defined as $N/(N_{1}N_{2})$. Due to the
periodic boundary condition, each energy level of this finite system can be
labeled by the total two-dimensional (2D) momentum $(K_{1},K_{2})$, with
integers $K_{1}=0,\cdots,N_{1}-1$ and $K_{2}=0,\cdots,N_{2}-1$. Motivated by
the observations of robust FCIs at $\nu=1/3$ (lattice analogs of the
celebrated Laughlin state [60]) in various static $|C|=1$ topological flat
bands [61, 62], we also fix $\nu=1/3$ in our Floquet system.

We then focus on the $C=-1$ region in the $(\theta,P)$ space to numerically
study the effective static Hamiltonian Eq. (10) to search for the evidence of
FCIs. Since the Floquet valence band is well isolated in the $C=-1$ region, it
is fair to project Eq. (10) to this active band, leading to

$\displaystyle H_{\rm eff}^{\text{proj}}=\sum_{{\bf k}}E({\bf k})c_{\bf
k}^{\dagger}c_{\bf
k}+\sum_{\\{\mathbf{k}_{i}\\}}V_{\mathbf{k}_{1}\mathbf{k}_{2}\mathbf{k}_{3}\mathbf{k}_{4}}c^{\dagger}_{\mathbf{k}_{1}}c^{\dagger}_{\mathbf{k}_{2}}c_{\mathbf{k}_{3}}c_{\mathbf{k}_{4}},$

where $c_{\mathbf{k}}^{\dagger}$ ($c_{\mathbf{k}}$) is the operator creating
(annihilating) an electron with momentum $\mathbf{k}$ in the Floquet valence
band, $E(\mathbf{k})$ is the corresponding band dispersion, and the matrix
element $V_{\mathbf{k}_{1}\mathbf{k}_{2}\mathbf{k}_{3}\mathbf{k}_{4}}$ is
given by [18]

$\displaystyle
V_{\mathbf{k}_{1}\mathbf{k}_{2}\mathbf{k}_{3}\mathbf{k}_{4}}=\frac{1}{2}\delta^{\prime}_{\mathbf{k}_{1}+\mathbf{k}_{2},\mathbf{k}_{3}+\mathbf{k}_{4}}\sum_{\bf
G}V(\mathbf{k}_{1}-\mathbf{k}_{4}+\mathbf{G})$ $\displaystyle\times\langle
u({\bf k}_{1})|u({\bf k}_{4}-\mathbf{G})\rangle\langle u({\bf k}_{2})|u({\bf
k}_{3}+\mathbf{G}+\delta\mathbf{G})\rangle.$ (12)

Here $\delta^{\prime}_{{\bf k},{\bf k}^{\prime}}$ is the 2D periodic Kronecker
delta function with period of primitive reciprocal lattice vectors, $|u({\bf
k})\rangle$ is the Floquet valence band eigenvector, $\delta\mathbf{G}={\bf
k}_{1}+{\bf k}_{2}-{\bf k}_{3}-{\bf k}_{4}$, and the sum of ${\bf G}$ is over
the entire reciprocal space rather than only in the MBZ. $E(\mathbf{k})$ and
$|u({\bf k})\rangle$ can be obtained by diagonalizing the single-particle part
of $H_{\rm eff}$. We then use exact diagonalization to extract the low-energy
physics of this projected effective Hamiltonian Eq. (IV) for $w_{1}=90\ {\rm
meV}$ and $w_{1}=110\ {\rm meV}$, respectively. We have examined that the
ground-state energy gap obtained from diagonalizing Eq. (IV) is always smaller
than the single-electron band gap of the Floquet valence band, thus justifying
the validity of band projection.

Figure 4:  The FCI gap $\Delta_{\rm FCI}$, FCI splitting $W_{\rm FCI}$, and
their ratio $\Delta_{\rm FCI}/W_{\rm FCI}$ for $w_{1}=90\ \text{meV}$
[(a)-(c)] and $w_{1}=110\ \text{meV}$ [(d)-(f)]. We compute the data for
$N=10$ and $N_{1}\times N_{2}=5\times 6$. In (a) and (d), we also give the
tentative phase diagram. As shown later, the CDW phase is identified by the
structure factor, and the Fermi liquid (FL)-like phase is characterized by the
step structure in the $n({\bf k})-E_{h}({\bf k})$ curve. The markers indicate
the representative parameter points that we choose in Figs. 5, 6, 7, 8 and 9.
Figure 5: Evidence for the $\nu=1/3$ Laughlin FCIs in the Floquet valence band
at $\theta=1.05^{\circ}$, $P=10\ \text{meV}$ when $w_{1}=90\text{meV}$. (a)
The low-lying energy spectra for $N=8,10,12$ electrons on the $N/2\times 6$
lattice. (b) The finite-size scaling of the energy gap $E_{4}-E_{1}$ and the
ground-state splitting $E_{3}-E_{1}$ for $N=4,5,6,8,10,12$ electrons. (c) The
spectral flow for $N=10$, $N_{1}\times N_{2}=5\times 6$, where $\Phi_{y}$ is
the magnetic flux insertion in the $\mathbf{a}_{2}$ direction. (d) The
particle entanglement spectrum for $N=10$, $N_{1}\times N_{2}=5\times 6$ and
$N_{A}=5$, with $23256$ levels below the entanglement gap (the dashed line).

###  IV.1 $w_{1}=90\ \text{meV}$

Let us first choose $w_{1}=90\ \text{meV}$. On the torus geometry, an
essential feature of the $\nu=1/3$ FCIs is the robust three-fold ground-state
degeneracy in the momentum sectors determined by the Haldane statistics of the
$\nu=1/3$ Laughlin state [61, 63, 64]. Therefore, we compute the FCI gap
$\Delta_{\rm FCI}$ in the $C=-1$ region as the energy difference between the
fourth and the first eigenvalues of the projected effective Hamiltonian Eq.
(IV), where all eigenvalues are sorted in ascending order. If the lowest three
eigenstates are not in the momentum sectors predicted by the Haldane
statistics, we simply set the FCI gap to be zero. Meanwhile, we also measure
the FCI splitting $W_{\rm FCI}$, quantified by the energy difference between
the third and the first eigenvalues when they are in the FCI momentum sectors.
The result for $N=10$ electrons is demonstrated in Figs. 4(a) and 4(b).
Strikingly, we can identify a wide range of parameters for which the lowest
three eigenstates are located in FCI momentum sectors, protected by a
significant gap, and approximately degenerate [i.e., very large $\Delta_{\rm
FCI}/W_{\rm FCI}$, as shown in Fig. 4(c)].

To examine the robustness of such ground-state degeneracies, we examine their
dependence on the system size and the boundary condition. In Fig. 5(a), we
demonstrate the low-energy spectra of the projected effective Hamiltonian for
$N=8$, $10$, and $12$ electrons at a representative parameter point
$(\theta,P)=(1.05^{\circ},10\ \text{meV})$ [labelled by the star in Figs.
4(a)-4(c)]. For each system size, we observe excellent three-fold ground-state
degeneracy, i.e., the splitting of the three ground states are much smaller
than their separation from higher-energy levels. A finite-size scaling of the
ground-state splitting and the energy gap suggests that both the three-fold
degeneracy and the ground-state gap are very likely to survive in the
thermodynamic limit [Fig. 5(b)]. By inserting magnetic flux through the
handles of the toroidal system, we find the three-fold ground-state degeneracy
persists [Fig. 5(c)]. All these data further confirm the robustness of the
three-fold topological degeneracy. Remarkably, the observed gap corresponds to
a temperature of about $20$ Kelvin, which is an order of magnitude higher than
required by the conventional fractional quantum Hall states in two-dimensional
electron gases.

To further corroborate the nontrivial topological properties of the ground
state, we also investigate the particle entanglement spectrum (PES), which
encodes the information of quasihole excitations of the system and
distinguishes FCIs from other competing phases [61, 65, 66]. In Fig. 5(d), we
divide the whole system into $N_{A}$ and $N-N_{A}$ electrons and label each
PES level by the total momentum $(K_{1}^{A},K_{2}^{A})$ of those $N_{A}$
electrons. A clear entanglement gap appears separating the low-lying PES
levels from higher ones, and the number of levels below the gap exactly
matches the pertinent counting of quasihole excitations in the $\nu=1/3$
Laughlin state [61, 63, 64]. This entanglement spectroscopy, together with the
low-energy spectrum, strongly suggests that the most robust $\nu=1/3$ Floquet
FCIs exist in the region with $\theta\approx 1.0^{\circ}-1.1^{\circ}$ and
$P\approx 5\ {\rm meV}-30\ \text{meV}$ when $w_{1}=90\ {\rm meV}$.

Figure 6:  (a) The structure factor at $\theta=0.95^{\circ}$ as a function of
$P$ for $N=6$, $N_{1}\times N_{2}=3\times 6$. The values of $S({\bf q})$ at
${\bf q}={\bf K}_{+}^{t,b}$ are significantly larger than those at ${\bf
q}\neq{\bf K}_{+}^{t,b}$. (b) Distribution of $S(\mathbf{q})$ in the MBZ at
$(\theta,P)=(0.95^{\circ},20\ \textrm{meV})$ for $N=12$, $N_{1}\times
N_{2}=6\times 6$. (c) The many-body energy spectrum at the same parameters
with (b). The lowest three states are lying at
$(K_{1},K_{2})=(0,0),(2,2),(4,4)$, separated by $\Delta{\bf K}=(2,2)$ and
$\Delta{\bf K}=(4,4)$.

There are also regions in Fig. 4(a) in which the three-fold topological
degeneracy of the ground states becomes poor and eventually collapses. On the
left side of the FCI phase (with smaller $\theta$), we find a pronounced
sensitivity of the energy spectrum to the lattice size. For $N=6,N_{1}\times
N_{2}=3\times 6$ and $N=12,N_{1}\times N_{2}=6\times 6$, we observe a new kind
of three-fold ground-state degeneracies in different momentum sectors from the
$\nu=1/3$ Laughlin FCIs. These momentum sectors are separated exactly by the
moiré Dirac point momenta $\mathbf{K}_{+}^{b}$ and $\mathbf{K}_{+}^{t}$. For
example, while the three Laughlin FCIs of $N=12,N_{1}\times N_{2}=6\times 6$
all carry $(K_{1},K_{2})=(0,0)$ [Fig. 5(a)], the three ground states at the
parameter point $(\theta,P)=(0.95^{\circ},20\ \textrm{meV})$ [labelled by the
dot in Figs. 4(a)-4(c)] of the same system size are located in the
$(K_{1},K_{2})=(0,0),(2,2)$ and $(4,4)$ sectors [Fig. 6(c)], which are
separated by momentum $\Delta{\bf
K}=(2,2)=(\mathbf{G}_{1}+\mathbf{G}_{2})/3\sim\mathbf{K}_{+}^{b}$ and
$\Delta{\bf K}=(4,4)=2(\mathbf{G}_{1}+\mathbf{G}_{2})/3\sim\mathbf{K}_{+}^{t}$
(Fig. 1). Here $\sim$ means “equal to” up to a MBZ reciprocal lattice vector.
The distribution of degenerate ground states over equally spaced momenta is a
signal of charge density waves. To further confirm this, we compute the
structure factor $S(\mathbf{q})$ which can reveal the CDW order. We define
$S(\mathbf{q})$ as

$\displaystyle S(\mathbf{q})=\frac{1}{N_{1}N_{2}}\left(\langle\bar{\rho}({\bf
q})\bar{\rho}(-{\bf q})\rangle-N^{2}\delta_{{\bf q},{\bf 0}}\right),$ (13)

where

$\displaystyle\bar{\rho}(\mathbf{q})=\sum_{{\bf k}\in{\rm MBZ}}\langle
u(\mathbf{k})|u(\mathbf{k}-{\bf q})\rangle c^{\dagger}_{\mathbf{k}}c_{{\bf
k}-{\bf q}}.$ (14)

Remarkably, we find pronounced peaks at the corners of the MBZ [Figs. 6(a) and
6(b)], revealing the underlying phase is the CDW with the order momentum
$\mathbf{K}_{+}^{t,b}$. Such a $\mathbf{K}$-CDW phase was also identified in
Refs. [21, 67] for static TBG-hBN as a competing phase against FCIs. These
$\mathbf{K}$-CDW states corresponds to a Wigner crystal, whose unit cell is
tripled compared to the original moiré lattice [21]. On the other hand, we do
not see the CDW degeneracy among the lowest states for $N=8,N_{1}\times
N_{2}=4\times 6$ and $N=10,N_{1}\times N_{2}=5\times 6$. This is because the
$\mathbf{K}_{+}^{t,b}$ points are absent in the MBZ of these finite lattices.
Remember that the single-electron momentum ${\bf k}$ can only take ${\bf
k}=\frac{m_{1}}{N_{1}}{\bf G}_{1}+\frac{m_{2}}{N_{2}}{\bf G}_{2}$ for a finite
periodic system. Therefore, both $N_{1}$ and $N_{2}$ must be divisible by
three if $\mathbf{K}_{+}^{t,b}$ belong to this set of allowed ${\bf k}$.
Otherwise, the $\mathbf{K}$-CDW cannot develop in the finite system.

Figure 7:  The ground-state occupation $\braket{n(\mathbf{k})}$ at momentum
${\bf k}$ as a function of $-E_{h}(\mathbf{k})$ for three representative
parameter points. (a) and (b) correspond to the already identified CDW and FCI
phase, respectively.

Finally, we go to the region on the right side of the FCI phase (larger
$\theta$), where we see neither the FCI topological degeneracy nor the CDW
degeneracy in the low-energy spectra. To explore the nature of the ground
states in this region, we switch from the picture of electrons to that of
holes, i.e., changing $c(\mathbf{k})$ to $d^{\dagger}(\mathbf{k})$ in the
Hamiltonian Eq. (IV). Here $d^{\dagger}(\mathbf{k})$ creates a hole with
momentum ${\bf k}$ in the Floquet valence band. As shown in Refs. [18, 68],
this particle-hole transformation induces an effective hole dispersion

$\displaystyle E_{h}(\mathbf{k})$ $\displaystyle=$
$\displaystyle\sum_{\mathbf{k}^{\prime}}\left(V_{\mathbf{k}^{\prime}\mathbf{k}\mathbf{k}^{\prime}\mathbf{k}}+V_{\mathbf{k}\mathbf{k}^{\prime}\mathbf{k}\mathbf{k}^{\prime}}-V_{\mathbf{k}\mathbf{k}^{\prime}\mathbf{k}^{\prime}\mathbf{k}}-V_{\mathbf{k}^{\prime}\mathbf{k}\mathbf{k}\mathbf{k}^{\prime}}\right)$
(15) $\displaystyle-E(\mathbf{k})$

in the dual Hamiltonian of holes. As $E_{h}(\mathbf{k})$ in general is non-
constant, it breaks the particle-hole symmetry of Eq. (IV). If such a
dispersion dominates over the hole-hole interaction, the system can form a
Fermi liquid-like state, that is, electrons (holes) tend to occupy ${\bf k}$
points with large (small) $E_{h}({\bf k})$. This can be manifest by a sharp
step in the electronic occupation number $\braket{n(\mathbf{k})}$ at an
electron’s “Fermi energy” $-E_{h}(\mathbf{k})$. In Fig. 7, we display the
ground-state occupation $\braket{n(\mathbf{k})}$ at momentum ${\bf k}$ as a
function of $-E_{h}(\mathbf{k})$ for the three representative parameter points
[indicated by the dot, star, and triangle in Figs. 4(a)-4(c), respectively].
The $\braket{n(\mathbf{k})}$ curves show obviously distinct features in these
three cases. In the FCI phase, the correlation between
$\braket{n(\mathbf{k})}$ and $-E_{h}(\mathbf{k})$ is quite weak, being
consistent with the expectation of $\braket{n(\mathbf{k})}\approx\nu$ for the
FCIs [Fig. 7(b)]. The dependence of $\braket{n(\mathbf{k})}$ on
$-E_{h}(\mathbf{k})$ becomes stronger in the CDW phase [Fig. 7(a)], as
electrons tend to occupy ${\bf K}_{+}^{t,b}$ points. Eventually, for ground
states at parameter points on the right side of the FCI phase, we find
striking Fermi surface-like structures, that is,
$\braket{n(\mathbf{k})}\approx 1$ for small $-E_{h}(\mathbf{k})$ and drops to
very small values at a critical $-E_{h}(\mathbf{k})$ [Fig. 7(c)]. This feature
strongly suggests that the region on the right side of the FCI phase is a
Fermi-liquid like phase.

Based on numerical results above, we present a tentative phase diagram for
$w_{1}=90\ {\rm meV}$ in Fig. 4(a). Because of the limited system sizes which
we can reach by exact diagonalization, the phase boundaries are only roughly
determined by the topological degeneracy of FCIs for $10$ electrons, so they
should not be thought as being precise.

Figure 8: Evidence for the $\nu=1/3$ Laughlin FCIs in the Floquet valence band
at $\theta=1.27^{\circ}$, $P=10\ \text{meV}$ when $w_{1}=110\ \text{meV}$. (a)
The low-lying energy spectra for $N=8,10,12$ electrons on the $N/2\times 6$
lattice. (b) The finite-size scaling of the energy gap $E_{4}-E_{1}$ and the
ground-state splitting $E_{3}-E_{1}$ for $N=4,5,6,8,10,12$ electrons. (c) The
spectral flow for $N=10$, $N_{1}\times N_{2}=5\times 6$, where $\Phi_{y}$ is
the magnetic flux insertion in the $\mathbf{a}_{2}$ direction. (d) The
particle entanglement spectrum for $N=10$, $N_{1}\times N_{2}=5\times 6$ and
$N_{A}=5$, with $23256$ levels below the entanglement gap (the dashed line).

###  IV.2 $w_{1}=110\ \text{meV}$

We have repeated the investigations in Sec. IV.1 for stronger interlayer
tunneling $w_{1}=110\ \text{meV}$. We find the overall phase diagram [Figs.
4(d)-4(f)] is very similar to that for $w_{1}=90\ \text{meV}$. However, the
FCI phase exists at larger twist angles, which are still small but obviously
deviate from the magic angle. This can be explained by the distribution of
Berry curvature in the MBZ – we find that the twist angle has to be increased
to achieve flatter Berry curvature when $w_{1}$ is strengthened from $90\
\text{meV}$ to $110\ \text{meV}$. The best three-fold degeneracies of ten
electrons appear at $\theta\approx 1.20^{\circ}-1.35^{\circ}$ and $P\approx 5\
{\rm meV}-35\ {\rm meV}$ [Fig. 4(f)]. We demonstrate the energy spectrum and
the PES at a representative parameter point [labelled by the star in Figs.
4(d)-4(f)] in this region (Fig. 8), where the energy gap protecting FCI ground
states is even larger than that at $w_{1}=90\ \text{meV}$. Meanwhile, the
competition between FCI and CDW remains for $w_{1}=110\ \text{meV}$. The
signals of the CDW phase and the Fermi liquid-like phase, namely, the
structure factor and the $\braket{n(\mathbf{k})}-E_{h}(\mathbf{k})$ curve, are
shown in Fig. 9 for two representative parameter points [labelled by the dot
and triangle in Figs. 4(d)-4(f), respectively].

Figure 9:  (a) Distribution of $S(\mathbf{q})$ in the MBZ for $N=12$,
$N_{1}\times N_{2}=6\times 6$ at $(\theta,P)=(1.20^{\circ},55\ \textrm{meV})$.
(b) The ground-state occupation $\braket{n(\mathbf{k})}$ at momentum ${\bf k}$
as a function of $-E_{h}(\mathbf{k})$ for $N=12$, $N_{1}\times N_{2}=6\times
6$ at $(\theta,P)=(1.50^{\circ},10\ \textrm{meV})$.

##  V Discussion

In this work, we investigate the interaction effect in twisted bilayer
graphene irradiated with monochromatic circularly polarized light. We work in
the regime of high driving frequency. When the Floquet valence band obtained
from the effective static Hamiltonian is partially occupied by electrons at
$\nu=1/3$ filling, we find compelling numerical evidence that Floquet
fractional Chern insulators exist for a wide range of twist angle and driving
strength, as characterized by the robust ground-state topological degeneracy
and the counting of quasihole excitations. By calculating the structure factor
and the electron occupation numbers, the intriguing interplay of Floquet FCIs
with charge density waves and Fermi liquid-like states are also identified in
a single phase diagram. Our results demonstrate the rich many-body physics in
Floquet TBG. As the effective Hamiltonian in our setup is analogous to the
static TBG-hBN system with hBN alignment on both top and bottom sides of TBG,
the phase diagram obtained by us may also be helpful to the experimental
realization of zero-field FCI in the static TBG-hBN system. Our results imply
that the zero-field FCI could be stabilized in static TBG-hBN at twist angles
larger than the magic value.

There are several interesting theoretical future directions following our
present work. First, it would be interesting to take into account more model
parameters, such as the ratio $w_{0}/w_{1}$ and the dielectric constant, to
explore the many-body phase diagram in a larger parameter space. Second, while
we only consider $\nu=1/3$ filling in this work, it remains unclear whether
FCIs at other fillings, especially non-Abelian ones, can be stabilized in this
Floquet system. Finally, inspired by the recent development in Floquet band
structures of other moiré materials beyond TBG [38, 39, 40, 41], it is natural
to study the Floquet many-body physics thereof.

Considering that the Floquet Chern insulator has been realized in experiments
of monolayer graphene [43], it is possible to observe it also in TBG driven by
the light field, which is the first step towards realizing the Floquet FCIs
predicted in our work. In fact, the parameters we choose here should be within
experimentally realizable parameters. For example, both the driving frequency
$\hbar\Omega=1.5\ {\rm eV}$ and the driving strength $\mathcal{E}_{0}=8\
\rm{MV/cm}$ (corresponding to $P\approx 60\ {\rm meV}$) are accessible by
current laser technology. However, obvious challenges exist. For instance, the
Floquet heating out of the long-lived prethermal regime results a featureless
infinite-temperature system. Moreover, one should keep the driving off-
resonant, otherwise the direct absorption of photons by electrons makes the
effective static Hamiltonian insufficient to capture the out-of-equilibrium
properties (such as the transport) of the system [24].

###### Acknowledgements.

This project is supported by the National Key Research and Development Program
of China through Grant No. 2020YFA0309200 and the National Natural Science
Foundation of China through Grant No. 11974014.

##  Appendix A Effective Floquet Hamiltonian

In this Appendix, we give a detailed derivation of the static effective
Hamiltonian $H_{\text{eff}}$. Note that we do this derivation before the band
projection. Within each layer of TBG, the Dirac Hamiltonian of monolayer
graphene under driving becomes

$\displaystyle\begin{aligned}
h(\mathbf{k})\rightarrow&h\left(\mathbf{k}+\frac{e}{\hbar}\mathbf{A}(t)\right)&\\\
=&-\hbar v_{F}\begin{pmatrix}0&k_{-}+\frac{eA_{0}}{\hbar}e^{i\Omega t}\\\
k_{+}+\frac{eA_{0}}{\hbar}e^{-i\Omega t}&0\end{pmatrix},&\end{aligned}$ (16)

where $k_{\pm}=k_{x}\pm ik_{y}$. The interlayer tunneling remains unchanged
since we have neglected its transverse components. Then the only non-zero
Fourier components $H_{m}=\frac{1}{T}\int_{0}^{T}H(t)e^{-im\Omega t}dt$ are

$\displaystyle H_{0}$ $\displaystyle=$ $\displaystyle
H_{\text{kin}}+H_{\text{int}},$ $\displaystyle H_{1}$ $\displaystyle=$
$\displaystyle\left(-eA_{0}v_{F}\right)$
$\displaystyle\times\sum_{\mathbf{k}}\left(e^{i\theta/2}\psi_{tA}^{\dagger}(\mathbf{k})\psi_{tB}(\mathbf{k})+e^{-i\theta/2}\psi_{bA}^{\dagger}(\mathbf{k})\psi_{bB}(\mathbf{k})\right),$
$\displaystyle H_{-1}$ $\displaystyle=$ $\displaystyle H_{1}^{\dagger},$ (17)

where $H_{\text{kin}}$ is static TBG single-particle Hamiltonian defined in
main text. Using the relation
$\left[c^{\dagger}_{m}c_{n},c^{\dagger}_{k}c_{l}\right]=\delta_{n,k}c^{\dagger}_{m}c_{l}-\delta_{m,l}c^{\dagger}_{k}c_{n}$,
we have the first-order term of $H_{\text{eff}}$ as

$\displaystyle H_{\text{eff}}^{(1)}$ $\displaystyle=$
$\displaystyle\frac{\left[H_{1},H_{-1}\right]}{\hbar\Omega}$ (18)
$\displaystyle=$ $\displaystyle
P\sum_{\mathbf{k}}\left(\psi^{\dagger}_{t}(\mathbf{k})\sigma_{z}\psi_{t}(\mathbf{k})+\psi^{\dagger}_{b}(\mathbf{k})\sigma_{z}\psi_{b}(\mathbf{k})\right),$

with $P=\frac{(eA_{0}v_{F})^{2}}{\hbar\Omega}$. The second-order term of
$H_{\text{eff}}$ is

$\displaystyle H^{(2)}_{\text{eff}}$ $\displaystyle=$
$\displaystyle\frac{1}{2(\hbar\Omega)^{2}}\left[H_{1},\left[H_{0},H_{-1}\right]\right]+\textrm{H.c.}$
(19) $\displaystyle=$
$\displaystyle\frac{1}{2(\hbar\Omega)^{2}}\left(\left[H_{1},\left[H_{\text{kin}},H_{-1}\right]\right]+\left[H_{1},\left[H_{\text{int}},H_{-1}\right]\right]\right).$
$\displaystyle+\textrm{H.c.}$

A straightforward calculation gives

$\displaystyle\left[H_{1},\left[H_{\rm
kin},H_{-1}\right]\right]+\textrm{H.c.}=2(eA_{0}v_{F})^{2}\Bigg{[}-H_{\rm
kin}$
$\displaystyle+\sum_{\mathbf{k}}\sum_{j=1}^{2}\left(\psi_{t}^{\dagger}(\mathbf{k}-\mathbf{q}_{0}+\mathbf{q}_{j})\mathbb{W}_{\theta}\psi_{b}(\mathbf{k})+\textrm{H.c.}\right)\Bigg{]}$
(20)

with

$\displaystyle\mathbb{W}_{\theta}=w_{0}\begin{pmatrix}e^{i\theta}&0\\\
0&e^{-i\theta}\end{pmatrix}.$ (21)

The calculation of $\left[H_{1},\left[H_{\text{int}},H_{-1}\right]\right]$ is
more tedious. We write the interaction Eq. (3) in the second-quantized form

$\displaystyle\begin{aligned}
H_{\text{int}}=&\sum_{\begin{subarray}{c}\\{\mathbf{k}_{i}\\},\mathbf{q}\end{subarray}}\sum_{\alpha,\beta}\delta_{\mathbf{k}_{1}-\mathbf{k}_{2},\mathbf{q}}\delta_{\mathbf{k}_{1}+\mathbf{k}_{2},\mathbf{k}_{3}+\mathbf{k}_{4}}&\\\
&\times
V(\mathbf{q})\psi^{\dagger}_{\alpha}(\mathbf{k}_{1})\psi^{\dagger}_{\beta}(\mathbf{k}_{2})\psi_{\beta}(\mathbf{k}_{3})\psi_{\alpha}(\mathbf{k}_{4}),&\end{aligned}$
(22)

where $\alpha,\beta\in(tA,tB,bA,bB)$ are layer and sublattice indices. For
fermionic creation and annihilation operators, we have

$\displaystyle\left[c^{\dagger}_{\alpha}(\mathbf{k}_{1})c^{\dagger}_{\beta}(\mathbf{k}_{2})c_{\gamma}(\mathbf{k}_{3})c_{\delta}(\mathbf{k}_{4}),\sum_{\mathbf{k}}h(\mathbf{k})c^{\dagger}_{\mu}(\mathbf{k})c_{\nu}(\mathbf{k})\right]$
(23) $\displaystyle=$
$\displaystyle\delta_{\mu,\delta}h(\mathbf{k}_{3})c^{\dagger}_{\alpha}(\mathbf{k}_{1})c^{\dagger}_{\beta}(\mathbf{k}_{2})c_{\gamma}(\mathbf{k}_{3})c_{\nu}(\mathbf{k}_{4})$
$\displaystyle+$
$\displaystyle\delta_{\mu,\gamma}h(\mathbf{k}_{4})c^{\dagger}_{\alpha}(\mathbf{k}_{1})c^{\dagger}_{\beta}(\mathbf{k}_{2})c_{\nu}(\mathbf{k}_{3})c_{\delta}(\mathbf{k}_{4})$
$\displaystyle-$
$\displaystyle\delta_{\nu,\alpha}h(\mathbf{k}_{1})c^{\dagger}_{\mu}(\mathbf{k}_{1})c^{\dagger}_{\beta}(\mathbf{k}_{2})c_{\gamma}(\mathbf{k}_{3})c_{\delta}(\mathbf{k}_{4})$
$\displaystyle-$
$\displaystyle\delta_{\nu,\beta}h(\mathbf{k}_{2})c^{\dagger}_{\alpha}(\mathbf{k}_{1})c^{\dagger}_{\mu}(\mathbf{k}_{2})c_{\gamma}(\mathbf{k}_{3})c_{\delta}(\mathbf{k}_{4}).$

We then apply Eq. (23) to $\left[H_{\text{int}},H_{-1}\right]$, and get

$\displaystyle\left[H_{\text{int}},H_{-1}\right]$ $\displaystyle=$
$\displaystyle\left(-eA_{0}v_{F}\right)\sum_{\begin{subarray}{c}\\{\mathbf{k}_{i}\\},\mathbf{q}\end{subarray}}\sum_{\mathbf{k}}\sum_{\alpha,\beta}\delta_{\mathbf{k}_{1}-\mathbf{k}_{2},\mathbf{q}}\delta_{\mathbf{k}_{1}+\mathbf{k}_{2},\mathbf{k}_{3}+\mathbf{k}_{4}}$
(24) $\displaystyle\times$ $\displaystyle
V(\mathbf{q})\big{\\{}e^{-i\theta/2}\psi^{\dagger}_{tB}(\mathbf{k}_{1})\psi^{\dagger}_{\alpha}(\mathbf{k}_{2})\psi_{\alpha}(\mathbf{k}_{3})\psi_{tA}(\mathbf{k}_{4})$
$\displaystyle+e^{-i\theta/2}\psi^{\dagger}_{\alpha}(\mathbf{k}_{1})\psi^{\dagger}_{tB}(\mathbf{k}_{2})\psi_{tA}(\mathbf{k}_{3})\psi_{\alpha}(\mathbf{k}_{4})$
$\displaystyle-e^{-i\theta/2}\psi^{\dagger}_{tB}(\mathbf{k}_{1})\psi^{\dagger}_{\alpha}(\mathbf{k}_{2})\psi_{\alpha}(\mathbf{k}_{3})\psi_{tA}(\mathbf{k}_{4})$
$\displaystyle-e^{-i\theta/2}\psi^{\dagger}_{\alpha}(\mathbf{k}_{1})\psi^{\dagger}_{tB}(\mathbf{k}_{2})\psi_{tA}(\mathbf{k}_{3})\psi_{\alpha}(\mathbf{k}_{4})$
$\displaystyle+e^{i\theta/2}\psi^{\dagger}_{bB}(\mathbf{k}_{1})\psi^{\dagger}_{\alpha}(\mathbf{k}_{2})\psi_{\alpha}(\mathbf{k}_{3})\psi_{bA}(\mathbf{k}_{4})$
$\displaystyle+e^{i\theta/2}\psi^{\dagger}_{\alpha}(\mathbf{k}_{1})\psi^{\dagger}_{bB}(\mathbf{k}_{2})\psi_{bA}(\mathbf{k}_{3})\psi_{\alpha}(\mathbf{k}_{4})$
$\displaystyle-e^{i\theta/2}\psi^{\dagger}_{bB}(\mathbf{k}_{1})\psi^{\dagger}_{\alpha}(\mathbf{k}_{2})\psi_{\alpha}(\mathbf{k}_{3})\psi_{bA}(\mathbf{k}_{4})$
$\displaystyle-e^{i\theta/2}\psi^{\dagger}_{\alpha}(\mathbf{k}_{1})\psi^{\dagger}_{bB}(\mathbf{k}_{2})\psi_{bA}(\mathbf{k}_{3})\psi_{\alpha}(\mathbf{k}_{4})\big{\\}}$
$\displaystyle=$ $\displaystyle 0.$

Therefore, only $H_{\rm kin}$ contributes to $H_{\text{eff}}^{(2)}$ in our
model, making $H_{\text{eff}}^{(2)}$ a non-interacting part in
$H_{\text{eff}}$. This is different from the result in Ref. [45], whose
authors obtained effective interactions in $H_{\text{eff}}^{(2)}$ for another
Floquet lattice model. The reason for this difference is that, in our model
the factor in front of $\psi^{\dagger}(\mathbf{k})\psi(\mathbf{k})$ in
$H_{\text{eff}}^{(1)}$ does not depend on $\mathbf{k}$, and $\alpha$ and
$\beta$ in the interaction can take the same value in Eq. (22). On the
contrary, in the model of Ref. [45] the corresponding prefactor in
$H_{\text{eff}}^{(1)}$ does depend on ${\bf k}$ and $\alpha$ and $\beta$ in
the interaction cannot take the same value. Hence,
$\left[H_{\text{int}},H_{-1}\right]$ is different between our model and the
model in Ref. [45].

## References

  * Novoselov  _et al._ [2016] K. S. Novoselov, A. Mishchenko, A. Carvalho, and A. H. Castro Neto, 2D materials and van der Waals heterostructures, Science 353, 10.1126/science.aac9439 (2016).
  * Geim and Grigorieva [2013] A. K. Geim and I. V. Grigorieva, Van der Waals heterostructures, Nature 499, 419 (2013), 1307.6718 .
  * Andrei and MacDonald [2020] E. Y. Andrei and A. H. MacDonald, Graphene bilayers with a twist, Nature Materials 19, 1265 (2020).
  * Cao  _et al._ [2018a] Y. Cao, V. Fatemi, A. Demir, S. Fang, S. L. Tomarken, J. Y. Luo, J. D. Sanchez-yamagishi, K. Watanabe, T. Taniguchi, E. Kaxiras, R. C. Ashoori, and P. Jarillo-herrero, Correlated insulator behaviour at half-filling in magic-angle graphene superlattices, Nature Publishing Group 10.1038/nature26154 (2018a).
  * Cao  _et al._ [2018b] Y. Cao, V. Fatemi, S. Fang, K. Watanabe, T. Taniguchi, E. Kaxiras, and P. Jarillo-herrero, magic-angle graphene superlattices 10.1038/nature26160 (2018b).
  * Sharpe  _et al._ [2019] A. L. Sharpe, E. J. Fox, A. W. Barnard, J. Finney, K. Watanabe, T. Taniguchi, M. A. Kastner, and D. Goldhaber-Gordon, Emergent ferromagnetism near three-quarters filling in twisted bilayer graphene, Science 365, 605 (2019), 1901.03520 .
  * Lu  _et al._ [2019] X. Lu, P. Stepanov, W. Yang, M. Xie, M. A. Aamir, I. Das, C. Urgell, K. Watanabe, T. Taniguchi, G. Zhang, A. Bachtold, A. H. MacDonald, and D. K. Efetov, Superconductors, orbital magnets and correlated states in magic-angle bilayer graphene, Nature 574, 653 (2019), 1903.06513 .
  * Xu  _et al._ [2020] Y. Xu, S. Liu, D. A. Rhodes, K. Watanabe, T. Taniguchi, J. Hone, V. Elser, K. F. Mak, and J. Shan, Correlated insulating states at fractional fillings of moiré superlattices, Nature 587, 214 (2020).
  * Yankowitz  _et al._ [2019] M. Yankowitz, S. Chen, H. Polshyn, Y. Zhang, K. Watanabe, T. Taniguchi, D. Graf, A. F. Young, and C. R. Dean, Tuning superconductivity in twisted bilayer graphene, Science 363, 1059 (2019), 1808.07865 .
  * Serlin  _et al._ [2020] M. Serlin, C. L. Tschirhart, H. Polshyn, Y. Zhang, J. Zhu, K. Watanabe, T. Taniguchi, L. Balents, and A. F. Young, Intrinsic quantized anomalous Hall effect in a moiré heterostructure, Science 367, 900 (2020), 1907.00261 .
  * Hunt  _et al._ [2013] B. Hunt, J. D. Sanchez-Yamagishi, A. F. Young, M. Yankowitz, B. J. LeRoy, K. Watanabe, T. Taniguchi, P. Moon, M. Koshino, P. Jarillo-Herrero, and R. C. Ashoori, Massive dirac fermions and hofstadter butterfly in a van der waals heterostructure, Science 340, 1427 (2013), https://www.science.org/doi/pdf/10.1126/science.1237240 .
  * San-Jose  _et al._ [2014] P. San-Jose, A. Gutiérrez-Rubio, M. Sturla, and F. Guinea, Spontaneous strains and gap in graphene on boron nitride, Phys. Rev. B 90, 075428 (2014).
  * Jung  _et al._ [2015] J. Jung, A. M. DaSilva, A. H. MacDonald, and S. Adam, Origin of band gaps in graphene on hexagonal boron nitride, Nature Communications 6, 6308 (2015).
  * Zhang  _et al._ [2019a] Y.-H. Zhang, D. Mao, and T. Senthil, Twisted bilayer graphene aligned with hexagonal boron nitride: Anomalous hall effect and a lattice model, Phys. Rev. Research 1, 033126 (2019a).
  * Regnault and Bernevig [2011] N. Regnault and B. A. Bernevig, Fractional chern insulator, Phys. Rev. X 1, 021014 (2011).
  * Parameswaran  _et al._ [2013] S. A. Parameswaran, R. Roy, and S. L. Sondhi, Fractional quantum hall physics in topological flat bands, Comptes Rendus Physique 14, 816 (2013), topological insulators / Isolants topologiques.
  * BERGHOLTZ and LIU [2013] E. J. BERGHOLTZ and Z. LIU, Topological flat band models and fractional chern insulators, International Journal of Modern Physics B 27, 1330017 (2013), https://doi.org/10.1142/S021797921330017X .
  * Abouelkomsan  _et al._ [2020] A. Abouelkomsan, Z. Liu, and E. J. Bergholtz, Phys. Rev. Lett. 124, 106803 (2020).
  * Ledwith  _et al._ [2020] P. J. Ledwith, G. Tarnopolsky, E. Khalaf, and A. Vishwanath, Fractional Chern insulator states in twisted bilayer graphene: An analytical approach, Physical Review Research 2, 10.1103/PhysRevResearch.2.023237 (2020), 1912.09634 .
  * Repellin and Senthil [2020] C. Repellin and T. Senthil, Chern bands of twisted bilayer graphene: Fractional Chern insulators and spin phase transition, Physical Review Research 2, 1 (2020), 1912.11469 .
  * Wilhelm  _et al._ [2021] P. Wilhelm, T. C. Lang, and A. M. Läuchli, Phys. Rev. B 103, 125406 (2021).
  * Xie  _et al._ [2021] Y. Xie, A. T. Pierce, J. M. Park, D. E. Parker, E. Khalaf, P. Ledwith, Y. Cao, S. H. Lee, S. Chen, P. R. Forrester, K. Watanabe, T. Taniguchi, A. Vishwanath, P. Jarillo-Herrero, and A. Yacoby, Fractional chern insulators in magic-angle twisted bilayer graphene, Nature 600, 439 (2021).
  * Oka and Aoki [2009] T. Oka and H. Aoki, Photovoltaic hall effect in graphene, Phys. Rev. B 79, 081406 (2009).
  * Kitagawa  _et al._ [2011] T. Kitagawa, T. Oka, A. Brataas, L. Fu, and E. Demler, Transport properties of nonequilibrium systems under the application of light: Photoinduced quantum hall insulators without landau levels, Phys. Rev. B 84, 235108 (2011).
  * Lindner  _et al._ [2011] N. H. Lindner, G. Refael, and V. Galitski, Floquet topological insulator in semiconductor quantum wells, Nature Physics 7, 490 (2011).
  * Wang  _et al._ [2013] Y. H. Wang, H. Steinberg, P. Jarillo-Herrero, and N. Gedik, Observation of floquet-bloch states on the surface of a topological insulator, Science 342, 453 (2013), https://www.science.org/doi/pdf/10.1126/science.1239834 .
  * Cayssol  _et al._ [2013] J. Cayssol, B. Dóra, F. Simon, and R. Moessner, Floquet topological insulators, physica status solidi (RRL) - Rapid Research Letters 7, 101 (2013).
  * Kitagawa  _et al._ [2010] T. Kitagawa, E. Berg, M. Rudner, and E. Demler, Topological characterization of periodically driven quantum systems, Phys. Rev. B 82, 235114 (2010).
  * Fläschner  _et al._ [2016] N. Fläschner, B. S. Rem, M. Tarnowski, D. Vogel, D.-S. Lühmann, K. Sengstock, and C. Weitenberg, Experimental reconstruction of the berry curvature in a floquet bloch band, Science 352, 1091 (2016), https://www.science.org/doi/pdf/10.1126/science.aad4568 .
  * Bao  _et al._ [2022] C. Bao, P. Tang, D. Sun, and S. Zhou, Light-induced emergent phenomena in 2d materials and topological materials, Nature Reviews Physics 4, 33 (2022).
  * Rahav  _et al._ [2003] S. Rahav, I. Gilary, and S. Fishman, Effective Hamiltonians for periodically driven systems, Physical Review A 68, 18 (2003), 0301033 [nlin] .
  * Goldman and Dalibard [2014] N. Goldman and J. Dalibard, Periodically driven quantum systems: Effective Hamiltonians and engineered gauge fields, Physical Review X 4, 1 (2014), 1404.4373 .
  * Bukov  _et al._ [2015] M. Bukov, L. D’Alessio, and A. Polkovnikov, Universal high-frequency behavior of periodically driven systems: from dynamical stabilization to floquet engineering, Advances in Physics 64, 139 (2015), https://doi.org/10.1080/00018732.2015.1055918 .
  * Jiang  _et al._ [2011] L. Jiang, T. Kitagawa, J. Alicea, A. R. Akhmerov, D. Pekker, G. Refael, J. I. Cirac, E. Demler, M. D. Lukin, and P. Zoller, Majorana fermions in equilibrium and in driven cold-atom quantum wires, Phys. Rev. Lett. 106, 220402 (2011).
  * Katz  _et al._ [2020] O. Katz, G. Refael, and N. H. Lindner, Optically induced flat bands in twisted bilayer graphene, Physical Review B 102, 155123 (2020), 1910.13510 .
  * Vogl  _et al._ [2020a] M. Vogl, M. Rodriguez-Vega, and G. A. Fiete, Effective floquet hamiltonians for periodically driven twisted bilayer graphene, Phys. Rev. B 101, 235411 (2020a).
  * Li  _et al._ [2020] Y. Li, H. A. Fertig, and B. Seradjeh, Floquet-engineered topological flat bands in irradiated twisted bilayer graphene, Phys. Rev. Research 2, 043275 (2020).
  * Rodriguez-Vega  _et al._ [2020] M. Rodriguez-Vega, M. Vogl, and G. A. Fiete, Floquet engineering of twisted double bilayer graphene, Phys. Rev. Research 2, 033494 (2020).
  * Assi  _et al._ [2021] I. A. Assi, J. P. F. LeBlanc, M. Rodriguez-Vega, H. Bahlouli, and M. Vogl, Floquet engineering and nonequilibrium topological maps in twisted trilayer graphene, Phys. Rev. B 104, 195429 (2021).
  * Topp  _et al._ [2021] G. E. Topp, C. J. Eckhardt, D. M. Kennes, M. A. Sentef, and P. Törmä, Light-matter coupling and quantum geometry in moiré materials, Phys. Rev. B 104, 064306 (2021).
  * Vogl  _et al._ [2021] M. Vogl, M. Rodriguez-Vega, B. Flebus, A. H. MacDonald, and G. A. Fiete, Floquet engineering of topological transitions in a twisted transition metal dichalcogenide homobilayer, Phys. Rev. B 103, 014310 (2021).
  * Benlakhouy  _et al._ [2022] N. Benlakhouy, A. Jellal, H. Bahlouli, and M. Vogl, Chiral limits and effect of light on the hofstadter butterfly in twisted bilayer graphene, Phys. Rev. B 105, 125423 (2022).
  * McIver  _et al._ [2020] J. W. McIver, B. Schulte, F. U. Stein, T. Matsuyama, G. Jotzu, G. Meier, and A. Cavalleri, Light-induced anomalous Hall effect in graphene, Nature Physics 16, 38 (2020), arXiv:1811.03522 .
  * Grushin  _et al._ [2014] A. G. Grushin, Á. Gómez-León, and T. Neupert, Floquet fractional chern insulators, Physical Review Letters 112, 1 (2014).
  * Anisimovas  _et al._ [2015] E. Anisimovas, G. Žlabys, B. M. Anderson, G. Juzeliūnas, and A. Eckardt, Phys. Rev. B 91, 245135 (2015).
  * Bistritzer and MacDonald [2011] R. Bistritzer and A. H. MacDonald, Moiré bands in twisted double-layer graphene, Proc. Natl. Acad. Sci. USA 108, 12233 (2011), 1009.4203 .
  * Zhang  _et al._ [2019b] Y.-H. Zhang, D. Mao, Y. Cao, P. Jarillo-Herrero, and T. Senthil, Nearly flat chern bands in moiré superlattices, Phys. Rev. B 99, 075127 (2019b).
  * Jung  _et al._ [2014] J. Jung, A. Raoux, Z. Qiao, and A. H. MacDonald, Ab initio theory of moiré superlattice bands in layered two-dimensional materials, Phys. Rev. B 89, 205414 (2014).
  * Moon and Koshino [2013] P. Moon and M. Koshino, Optical absorption in twisted bilayer graphene, Physical Review B 87, 1 (2013), 1302.5218 .
  * Koshino  _et al._ [2018] M. Koshino, N. F. Yuan, T. Koretsune, M. Ochi, K. Kuroki, and L. Fu, Maximally Localized Wannier Orbitals and the Extended Hubbard Model for Twisted Bilayer Graphene, Physical Review X 8, 31087 (2018), 1805.06819 .
  * Popov  _et al._ [2011] A. M. Popov, I. V. Lebedeva, A. A. Knizhnik, Y. E. Lozovik, and B. V. Potapkin, Commensurate-incommensurate phase transition in bilayer graphene, Physical Review B 84, 1 (2011).
  * Nam and Koshino [2017] N. N. Nam and M. Koshino, Lattice relaxation and energy band modulation in twisted bilayer graphene, Physical Review B 96, 1 (2017), 1706.03908 .
  * Uchida  _et al._ [2014] K. Uchida, S. Furuya, J.-I. Iwata, and A. Oshiyama, Atomic corrugation and electron localization due to moiré patterns in twisted bilayer graphenes, Phys. Rev. B 90, 155451 (2014).
  * van Wijk  _et al._ [2015] M. M. van Wijk, A. Schuring, M. I. Katsnelson, and A. Fasolino, Relaxation of moiré patterns for slightly misaligned identical lattices: graphene on graphite, 2D Materials 2, 034010 (2015).
  * Bessler  _et al._ [2019] R. Bessler, U. Duerig, and E. Koren, The dielectric constant of a bilayer graphene interface, Nanoscale Adv. 1, 1702 (2019).
  * Guinea and Walet [2018] F. Guinea and N. R. Walet, Electrostatic effects, band distortions, and superconductivity in twisted graphene bilayers, Proc. Natl. Acad. Sci 115, 13174 (2018).
  * Vogl  _et al._ [2020b] M. Vogl, M. Rodriguez-Vega, and G. A. Fiete, Floquet engineering of interlayer couplings: Tuning the magic angle of twisted bilayer graphene at the exit of a waveguide, Physical Review B 101, 241408 (2020b).
  * Eckardt and Anisimovas [2015] A. Eckardt and E. Anisimovas, High-frequency approximation for periodically driven quantum systems from a Floquet-space perspective, New Journal of Physics 17, 93039 (2015), 1502.06477 .
  * Lu  _et al._ [2021] M. Lu, J. Zeng, H. Liu, J.-H. Gao, and X. C. Xie, Valley-selective floquet chern flat bands in twisted multilayer graphene, Phys. Rev. B 103, 195146 (2021).
  * Laughlin [1983] R. B. Laughlin, Anomalous quantum hall effect: An incompressible quantum fluid with fractionally charged excitations, Phys. Rev. Lett. 50, 1395 (1983).
  * Regnault and Andrei Bernevig [2011] N. Regnault and B. Andrei Bernevig, Fractional Chern Insulator, Physical Review X 1, 1 (2011), 1105.4867 .
  * Wu  _et al._ [2012] Y.-L. Wu, B. A. Bernevig, and N. Regnault, Zoology of fractional chern insulators, Phys. Rev. B 85, 075116 (2012).
  * Bernevig and Regnault [2012] B. A. Bernevig and N. Regnault, Emergent many-body translational symmetries of Abelian and non-Abelian fractionally filled topological insulators, Physical Review B 85, 1 (2012), 1110.4488 .
  * Wu  _et al._ [2014] Y. L. Wu, N. Regnault, and B. A. Bernevig, Haldane statistics for fractional Chern insulators with an arbitrary Chern number, Physical Review B 89, 1 (2014).
  * Sterdyniak  _et al._ [2011] A. Sterdyniak, N. Regnault, and B. A. Bernevig, Extracting excitations from model state entanglement, Physical Review Letters 106, 1 (2011).
  * Li and Haldane [2008] H. Li and F. D. Haldane, Entanglement spectrum as a generalization of entanglement entropy: Identification of topological order in non-Abelian fractional quantum hall effect states, Physical Review Letters 101, 1 (2008), arXiv:0805.0332 .
  * Abouelkomsan  _et al._ [2022] A. Abouelkomsan, K. Yang, and E. J. Bergholtz, Quantum metric induced phases in moiré materials 10.48550/ARXIV.2202.10467 (2022).
  * Liu  _et al._ [2013] Z. Liu, E. J. Bergholtz, R. Moessner, and A. M. La, Hierarchy of Fractional Chern Insulators and Competing Compressible States,  126802, 1 (2013).