Quantum state preparation of topological chiral spin liquids
via Floquet engineering

Matthieu Mambrini1, Didier Poilblanc1*

1 Laboratoire de Physique Théorique UMR5152, C.N.R.S. and Université de
Toulouse,

118 rte de Narbonne, 31062 Toulouse, FRANCE

*<EMAIL_ADDRESS>
## Abstract

In condensed matter, Chiral Spin Liquids (CSL) are quantum spin analogs of
electronic Fractional Quantum Hall states (in the continuum) or Fractional
Chern Insulators (on the lattice). As the latter, CSL are remarquable states
of matter, exhibiting topological order and chiral edge modes. Preparing CSL
on quantum simulators like cold atom platforms is still an open challenge.
Here we propose a simple setup on a finite cluster of spin-1/2 located at the
sites of a square lattice. Using a Resonating Valence Bond (RVB) non-chiral
spin liquid as initial state on which fast time-modulations of strong nearest-
neighbor Heisenberg couplings are applied, following different protocols (out-
of-equilibrium quench or semi-adiabatic ramping of the drive), we show the
slow emergence of such a CSL phase. An effective Floquet dynamics, obtained
from a high-frequency Magnus expansion of the drive Hamiltonian, provides a
very accurate and simple framework fully capturing the out-of-equilibrium
dynamics. An analysis of the resulting prepared states in term of Projected
Entangled Pair states gives further insights on the topological nature of the
chiral phase. Finally, we discuss possible applications to quantum computing.

###### Contents

  1. 1 Introduction
  2. 2 Setups and methods
    1. 2.1 General considerations on time evolution
    2. 2.2 Initial and targeted spin liquid states
    3. 2.3 Fast chiral periodic drive
    4. 2.4 Effective Floquet Hamiltonian
    5. 2.5 Quench and adiabatic evolution
  3. 3 Results
    1. 3.1 Validity of the Magnus high-frequency expansion
    2. 3.2 Quench of the fast drive
      1. 3.2.1 Long-time average and diagonal ensembles
      2. 3.2.2 Time evolution at small and intermediate times
    3. 3.3 Adiabatic ramping of the fast drive
    4. 3.4 Fidelities of the final CSL state
      1. 3.4.1 Comparison with the adiabatic limit
      2. 3.4.2 Comparison with chiral PEPS
  4. 4 Discussions and conclusions
    1. 4.1 Realizing the drive Hamiltonian
    2. 4.2 Preparation of the initial state
    3. 4.3 Final remarks
  5. A PEPS representation of spin liquids
  6. B Derivation of the Floquet effective Hamiltonian
    1. B.1 $H_{0}=0$ case
    2. B.2 Finite constant $H_{0}$
  7. C Derivation of the kick operator

##  1 Introduction

The search for spin liquids in condensed matter materials is a very rapidly
developing area of two-dimensional (2D) quantum magnetism [1]. Chiral spin
liquids (CSLs) are the quantum spin analogs of the topological Fractional
Quantum Hall States (FQHS) realized in the strongly interacting two-
dimensional (2D) electron gas [2]. A zoo of 2D topological Abelian and non-
Abelian CSLs, breaking time-reversal symmetry (T), have already been
identified for SU(2) spins on various (frustrated) lattices [3, 4, 5, 6], and
for spins of higher ${\text{SU}}(N)$ internal symmetry, $N>2$ [7, 8]. While
T-breaking can also occur spontaneously [9, 10, 11], most CSLs are stabilized
by an explicit T-breaking chiral perturbation.

In 1982, Richard Feynman first proposed that one quantum system could be used
to simulate the dynamics of other quantum systems of interest [12]. Nowadays,
quantum simulators based on cold atoms on two-dimensional (2D) optical
lattices are being realized experimentally with the realistic perspective of
emulating simple models of condensed matter physics [13, 14, 15, 16, 17] or
studying topological quantum matter [18]. Rydberg atoms platforms are also
used to realize 2D quantum phases of matter [19] like dimer liquids or spin
liquids [20, 21].

Designing atomic platforms of quantum spins with a T-breaking term to emulate
CSLs is still an open challenge. Very recently, using Floquet engineering, the
emulation of the non-Abelian CSL of the Kitaev honeycomb model in the presence
of a chiral term has been proposed theoretically using cold atoms or Rydberg
simulators [22, 23]. Preparing CSL with full spin-rotational symmetry is also
very challenging and of great interest. Proposals based on synthetic gauge
fields in N-colors SU($N$) Hubbard models on the square lattice with 1
fermion/site were shown to be effective for $N\geq 3$ [24]. Here we rather
focus on the preparation of SU(2)-symmetric spin-1/2 CSLs.

Quench dynamics is an active area of study in recent years [25, 26, 27, 28,
29, 30], relevant to a wide range of physical systems, including condensed
matter systems, ultracold atomic gases, and quantum field theories. While
simple quenches provide a way to probe the non-equilibrium dynamics of quantum
systems, they cannot serve to prepare quantum phases of matter under specific
conditions (low temperature, etc…). In contrast, Floquet engineering
consisting of applying a fast and strong periodic drive to the system –
described by Floquet theory [31] – is an appealing route towards phase/state
preparation [32], and even towards non-equilibrium phases without an
equivalent counterpart [33]. For example, such a scheme has been proposed to
trigger non-trivial energy bands via engineered gauge fields [34]. Our aim is
here to design some simple protocols based on Floquet theory to prepare
genuine spin-1/2 CSLs. For that purpose, we have developed new theoretical
tools to compute faithfully the time evolution of our 2D quantum spin systems
in order to address Floquet dynamics in various experimental setups, e.g. a
sudden quench of the drive or a semi-adiabatic ramping of the drive. Starting
from a paradigmatic spin liquid, the Resonating Valence Bond (RVB) state [35],
the finite system is slowly driven out of its equilibrium state, and the time
evolution of the system becomes very complex, showing rapid oscillations
(micromotion). However, using a high-frequency Magnus expansion [36, 37, 38,
39], we analytically construct an effective Floquet Hamiltonian which
describes accurately the (stroboscopic) time evolution and makes it much
easier to compute numerically.

Under the application of the fast drive, preserving the spin-rotational
invariance (singlet character) of the state, we observe a steady increase of
the plaquette chirality (more rapid for a quench than for a semi-adiabatic
ramping) signaling the emergence of a CSL. Using tensor network techniques in
2D [40], more specifically a chiral Projected Entangled Pair State (PEPS)
ansatz [41] (see Appendix A for details), the properties of the final state
are uncovered, showing topological features and chiral edge physics
characteristic of the Abelian SU(2)1 CSL [42].

##  2 Setups and methods

###  2.1 General considerations on time evolution

The considerations below are fairly general and apply to a closed system (i.e.
with no external bath) defined by some initial many-body state
$|\Psi_{0}\big{>}$ and its unitary evolution
$|\Psi(t)\big{>}=U(t,t_{0})|\Psi_{0}\big{>}$. In practice, we assume the
system is finite. Let us first consider a time-independent Hamiltonian $H$ –
typically a Floquet effective Hamiltonian – applied (abruptly) at time
$t_{0}$. The unitary time evolution operator $U(t_{0}+t,t_{0})$ is just given
by $\exp{(-iHt)}$ (setting $\hbar=1$ here) and its action on any state
$|\Phi\big{>}$ can be obtained by constructing the Krylov space
$\\{H^{m}|\Phi\big{>}\\}$ recursively. Then, in principle, one can numerically
obtain the exact evolution of the quantum state on a finite system, within
machine precision (provided that the computer memory is large enough to
accommodate the exponentially large Hilbert space).

In the case of a (more microscopic) time-dependent Hamiltonian $H(t)$, the
unitary time evolution operator involves time-ordering, i.e.
$U(t_{0}+t,t_{0})={\cal
T}_{t^{\prime}}\exp{(-i\int_{t_{0}}^{t_{0}+t}H(t^{\prime})dt^{\prime})}$. In
practice, one needs to discretize the time interval by regular time steps
$\tau=t/M$. Then, a Trotter-Suzuki (TS) decomposition [43, 44] of the unitary
time evolution can be performed in term of elementary gates,

$U(t_{0}+t,t_{0})\simeq\exp{(-iH(t_{M})\tau)}\cdots\exp{(-iH(t_{2})\tau)}\exp{(-iH(t_{1})\tau)},$
(1)

where $t_{n}=t_{0}+(n-\frac{1}{2})\tau$. The time-ordering operator ${\cal
T}_{t^{\prime}}$ can then be seen as the $\tau\rightarrow 0$ limit of Eq. 1.
In practice, $\tau=t/M$ is a small time step, significantly smaller that the
smallest relevant time scale in the problem, typically the drive period,
($\tau\ll T_{\rm drive}$). The TS decomposition involves Trotter errors of
order $\tau^{2}$ which can be controlled by choosing $\tau$ small enough. The
states $|\Psi(\tau)\big{>}$, $|\Psi(2\tau)\big{>}$, … $|\Psi(M\tau)\big{>}$
are obtained recursively using
$|\Psi(n\tau)\big{>}=\exp{(-iH(t_{n})\tau)}|\Psi((n-1)\tau)\big{>}$. After
Taylor expanding $\exp{(-iH(t_{n})\tau)}$ in powers of $H(t_{n})$, the
calculation is done in the Krylov’s basis $\\{|\Phi_{m}^{n}\big{>}\\}$,

$\displaystyle|\Psi(n\tau)\big{>}$ $\displaystyle=$
$\displaystyle\sum_{0}^{m_{\rm
max}}\frac{(-i)^{m}}{m!}|\Phi_{m}^{n}\big{>}\,,$ (2) $\displaystyle{\rm
where}\,\,\,\,\,|\Phi_{m}^{n}\big{>}$ $\displaystyle=$ $\displaystyle
H(t_{n})^{m}|\Psi((n-1)\tau)\big{>},$

using the normalisation condition $|\big{<}\Psi(n\tau)|\Psi(n\tau)\big{>}|=1$
as a test of convergence of the power-expansion up to order $m_{\rm max}$. The
Krylov’s basis is computed recursively using
$|\Phi_{m}^{n}\big{>}=H(t_{n})|\Phi_{(m-1)}^{n}\big{>}$ starting from
$|\Phi_{0}^{n}\big{>}=|\Psi((n-1)\tau)\big{>}$, for each $n$.

###  2.2 Initial and targeted spin liquid states

We now move more specifically to our system constituted by an assembly of
spin-1/2 on a square lattice. The initial prepared state considered in this
work is a simple spin-1/2 nearest-neighbor (NN) Resonating Valence Bond (RVB)
state [35]. It can be view as an equal-weight superposition of NN singlet
(hard-core) covering of the lattice. Conveniently, this spin liquid has a
simple Projected Entangled Pair State (PEPS) representation [45, 46] (see
Appendix A) which enables to construct it easily on finite systems or on the
infinite plane. The NN RVB state has very short-range spin-spin correlations
(with a spin correlation length of the order of the lattice spacing) but
critical dimer-dimer correlations (on any bipartite lattice). However, any
admixture of longer singlets generates a finite dimer correlation length [47,
48]. In that case, the gapped RVB spin liquid acquires ${\mathbb{Z}}_{2}$
topological order, similarly to the Kitaev toric code [49].

Most interestingly, possible (closely related) experimental realizations have
been proposed recently using Rydberg platforms [21]. Hereafter, we shall
consider a $4\times 4$ finite torus which conveniently offers lattice
translation symmetry. However, we believe our setup proposal and our findings
are also relevant to lattice systems with open boundaries, a geometry
encountered in experiments.

After applying a (fast) drive over a sufficiently long time interval we wish
to prepare the system in an (Abelian) Chiral Spin Liquid (CSL) phase, possibly
at low effective temperature (if not in the ground state). It is known that
such a phase is realized in a simple chiral spin-1/2 Heisenberg Hamiltonian on
the square lattice, which we view as a “target” Hamiltonian,

$H_{\rm target}=H_{0}+\lambda_{\rm
chiral}\,i\sum_{\square}(P_{ijkl}-P_{ijkl}^{-1}),$ (3)

where the last term involves cyclic (let’s say anti-clockwise) permutations
$P_{ijkl}$ on all plaquettes $\square$ and breaks time-reversal (T) and parity
(P), and $H_{0}$ is the spin-1/2 quantum antiferromagnet with NN
antiferromagnetic coupling $J_{1}$ and, optionally, next-NN (NNN)
antiferromagnetic coupling $J_{2}<J_{1}$,

$H_{0}=J_{1}\sum_{(ij)}{\bf S}_{i}\cdot{\bf S}_{j}+J_{2}\sum_{((ij))}{\bf
S}_{i}\cdot{\bf S}_{j}\,.$ (4)

In fact, the phase diagram of this model hosts the desired (gapped)
topological CSL phase over a significant region of its parameter space [50,
51, 52, 53], even for a vanishing frustrating coupling $J_{2}=0$, but for
finite $\lambda_{\rm chiral}$. Interestingly, the corresponding CSL phase can
be represented accurately by PEPS [41, 42, 51, 54], which conveniently give
crucial insights on the topological features and on the nature of the chiral
edge modes, via the bulk-edge correspondence PEPS theorem [55].

###  2.3 Fast chiral periodic drive

Figure 1: Left panel – Sketch of the $H_{\rm drive}$ Hamiltonian given by Eq.
5. Dephased terms $H_{a},H_{b},H_{c}$ and $H_{d}$ are represented as blue,
yellow, red and green bonds. Note that the bond modulations around the sites A
(black dots) and sites B (white dots) of the (bipartite) square lattice are
dephased by $\pi$. Right panel – Floquet effective Hamiltonian $H_{F}$ given
by Eq. 7 acting on square plaquettes as circular anticlockwise (resp.
clockwise) permutations with a factor $i$ (resp. $-i$).

The (ideal) monochromatic drive we consider involves dephased sinusoidal
modulations of the spin-spin interaction on the NN $a$, $b$, $c$ and $d$ bonds
around all the sites of one of the two sublattices (let’s say $A$). More
precisely, we assume,

$\displaystyle H_{\rm drive}(t)$ $\displaystyle=$ $\displaystyle
H_{a}(t)+H_{b}(t)+H_{c}(t)+H_{d}(t),$ (5) $\displaystyle H_{a}(t)$
$\displaystyle=$ $\displaystyle J\cos{(\omega t)}\,\sum_{i\in A}{\bf
S}_{i}\cdot{\bf S}_{i+e_{x}},$ $\displaystyle H_{b}(t)$ $\displaystyle=$
$\displaystyle J\cos{(\omega t-\pi/2)}\,\sum_{i\in A}{\bf S}_{i}\cdot{\bf
S}_{i+e_{y}},$ $\displaystyle H_{c}(t)$ $\displaystyle=$ $\displaystyle
J\cos{(\omega t-\pi)}\,\sum_{i\in A}{\bf S}_{i}\cdot{\bf S}_{i-e_{x}},$
$\displaystyle H_{d}(t)$ $\displaystyle=$ $\displaystyle J\cos{(\omega
t-3\pi/2)}\,\sum_{i\in A}{\bf S}_{i}\cdot{\bf S}_{i-e_{y}}.$

where $J$ is the coupling constant and $e_{x}$, $e_{y}$ are the unit vectors
along the crystal axis. As shown in Fig. 1 the bond pattern is staggered
leading to two types of plaquettes $adcb$ and $cbad$. Interestingly, $H_{\rm
drive}(t)$ is chiral, breaking parity $P$ and time-reversal $T$ but not the
product PT. Indeed, $P$ (reflection w.r.t. a crystal axis) changes the
plaquette bond ordering, e.g. $adcb\rightarrow abcd$, i.e. the sign of the
phases in Eq. 5. Time-reversal $T:t\rightarrow-t$ has the same effect. Note
that the sign of the chirality in the two types of plaquettes is the same so
that the overall chirality is uniform. Note also that the initial state is a
spin-singlet and the Hamiltonian in Eq. 5 is invariant under SU($2$) spin
rotations. The time evolution therefore occurs in the spin-singlet manifold.

It can be noticed that the sequence we are considering has similarities to the
driving protocol introduced in Ref. [56] for a tight-binding model of non-
interacting particles on the square lattice, in which hopping amplitudes are
varied in a spatially homogeneous but chiral, time-periodic way resulting in
chiral edge states.

Here the drive is assumed to be very fast compared to the natural time scale
$\hbar/J$, i.e. $\hbar\omega\gg J$. Hereafter we set $\hbar=1$ (unless useful
to keep it) and take $\omega=2\pi f$ with a frequency $f=10$. Hence, the drive
time-periodicity $T_{\rm drive}$ is 0.1 (in units of $1/J$ taken as the unit
of time throughout), which is the smallest time scale in the problem. Hence,
to compute the time evolution using a TS decomposition we have to divide the
drive period into a sufficiently large number $N_{\tau}$ of time steps
multiple of $\tau=T_{\rm drive}/N_{\tau}$. Note that such a monochromatic
drive may not be easy to implement on experimental cold atom platforms but has
been chosen for simplicity. We shall discuss variations more adapted to
experiment later in the conclusion section.

###  2.4 Effective Floquet Hamiltonian

In the case of fast drives, different types of high-frequency expansions can
be performed in powers of $1/\omega$. For practical reasons (which will become
clear soon), we follow here the Magnus expansion for the stroboscopic Floquet
Hamiltonian [39]. In that case, the effective (time-independent) Floquet
Hamiltonian describes the slow motion at stroboscopic times,
$T_{p}=t_{0}+pT_{\rm drive}$, $p\in{\mathbb{N}}$. The fast micromotion between
stroboscopic times is set by the periodic kick operator (not computed) which
vanishes at $t=nT$ (choosing $t_{0}=0$ i.e. $|\Psi(0)\big{>}=|\Psi_{\rm
RVB}\big{>}$). The drive Hamiltonian has no constant part, $(1/T_{\rm
drive})\int_{t_{0}}^{t_{0}+T_{\rm drive}}H_{\rm drive}(t)dt=0$, and just
single-harmonics components,

$\displaystyle H_{\rm drive}(t)$ $\displaystyle=$ $\displaystyle\cos{(\omega
t)}H_{x}+\sin{(\omega t)}H_{y},$ (6) $\displaystyle{\rm where}\hskip
8.5359ptH_{x}$ $\displaystyle=$ $\displaystyle J\sum_{i\in A}({\bf
S}_{i}\cdot{\bf S}_{i+e_{x}}-{\bf S}_{i}\cdot{\bf S}_{i-e_{x}})$
$\displaystyle{\rm and}\hskip 8.5359ptH_{y}$ $\displaystyle=$ $\displaystyle
J\sum_{i\in A}({\bf S}_{i}\cdot{\bf S}_{i+e_{y}}-{\bf S}_{i}\cdot{\bf
S}_{i-e_{y}}).$

The leading-order effective Floquet Hamiltonian is therefore [39],

$\displaystyle H_{F}^{0}$ $\displaystyle=$
$\displaystyle\frac{i}{2\omega}[H_{x},H_{y}],$ (7) $\displaystyle=$
$\displaystyle J_{F}\,i\sum_{\square}(P_{ijkl}-P_{ijkl}^{-1}),$

where $J_{F}=\frac{J^{2}}{4\omega}$ and $P_{ijkl}$ is the 4-site cyclic
permutation which applies to all plaquettes $\square$. The details of the
derivation is left in the Appendix B. It is interesting that the next order of
the expansion is ${\cal O}(1/\omega^{3})$ so that we expect $H_{F}^{0}$ to
capture accurately the stroboscopic motion. $H_{F}^{0}$ is a chiral
Hamiltonian, acquiring its P and T symmetries from the microscopic model
$H_{\rm drive}$: it transforms into its opposite under $P$ or $T$ and is
invariant under $PT$. Interestingly, it is fully symmetric – invariant under
lattice translations and under all discrete symmetries of the $C_{4v}$ point
group – in contrast to $H_{\rm drive}$ which bears two types of plaquettes (of
same chirality) as seen Fig. 1. The Floquet Hamiltonian is introducing the
longest time scale in the problem
$t_{F}=\hbar/J_{F}=(4\hbar\omega/J)\frac{\hbar}{J}$. Hence we are in the
regime where,

$t_{F}\gg\hbar/J\gg T_{\rm drive}\,,$ (8)

with 3 well-separated time (or energy) scales. Next, we shall also consider
the case where a constant term $\kappa H_{0}$ is added to the drive
Hamiltonian leading to a more generic Floquet Hamiltonian $H_{F}$. For
convenience we shall denote $|\Psi_{\rm Floquet}(t)\big{>}$ the state evolved
according to the Floquet Hamiltonian $\exp{(-iH_{F}t)}|\Psi_{\rm RVB}\big{>}$,
to be distinguished from the state $|\Psi_{\rm drive}(t)\big{>}={\cal
T}_{t^{\prime}}\exp{(-i\int_{0}^{t}(\kappa H_{0}+H_{\rm
drive}(t^{\prime}))dt^{\prime})}|\Psi_{\rm RVB}\big{>}$ following the exact
evolution.

Figure 2: Infidelity $1-|\big{<}\Psi_{\rm RVB}|U(t,0)|\Psi_{\rm RVB}\big{>}|$
w.r.t. the initial state, the NN RVB spin liquid, computed on a $4\times 4$
torus in the case of a pure chiral drive $H(t)=H_{\rm drive}(t)$. (a)
Comparison of the unitary evolutions with the drive Hamiltonian and with the
(time-independent) effective Floquet Hamitonian obtained on a $4\times 4$
cluster. Zooms at small times (b) or around $t=30$ (c) are shown. Trotter
errors are negligible at short times (b) but become sizeable for $N_{\tau}=8$
at intermediate times (c).

###  2.5 Quench and adiabatic evolution

The objective of the setup is to emulate the low-energy physics (possibly the
groundstate) of the (static) chiral antiferromagnetic Hamiltonian $H_{\rm
target}$ on the square lattice defined in Eq. 3. It is known that the phase
diagram of this model bears a wide CSL phase, even for a vanishing frustrating
coupling $J_{2}=0$. While emulating the NN coupling $J_{1}$ between the spins
can be done experimentally, generating the necessary chiral term is a
difficult task. For this aim, we add the Hamiltonian $H_{0}$ to the drive
Hamiltonian with some (small) amplitude $\kappa$ to be adjusted,

$H(t)=\kappa H_{0}+\lambda(t)H_{\rm drive}(t)\,.$ (9)

We also consider the possibility of slowly ramping the amplitude $J$ of the
drive by some time-dependent factor $0\leq\lambda(t)\leq 1$. Two scenarios
will be considered, either (i) a quench setup where $\lambda(t)=1$ from the
very beginning of the drive application at $t=0$ or (ii) a quasi-adiabatic
setup where $\lambda$ is increased smoothly from $0$ at $t=0$ to $1$ at
$t=t_{\rm ramp}$.

The most interesting regime to consider is the so-called strong (and fast)
drive regime $\kappa J_{1}\ll J$, i.e. $\kappa\ll 1$ if $J$ and $J_{1}$ are
both taken as energy reference, $J=J_{1}=1$. If $\lambda(t)$ is not too small
(i.e. excluding the short initial branching of the drive), the new high-
frequency effective Floquet Hamiltonian is simply,

$H_{F}(t)\simeq\kappa
H_{0}+\lambda^{2}(t)J_{F}\,i\sum_{\square}(P_{ijkl}-P_{ijkl}^{-1})\,,$ (10)

which is exactly of the form of the target Hamiltonian (3). Note that in the
case of a slow ramp, the Floquet Hamiltonian itself acquires a slow time
dependence. When a quench is realized, $\lambda(t)=1$, one can choose $\kappa$
so that $H_{F}=\kappa H_{\rm target}$. Using (3) and (10), one gets

$\kappa=\frac{J_{F}}{\lambda_{\rm chiral}}=\frac{J^{2}}{4\omega\lambda_{\rm
chiral}}\ll 1\,.$ (11)

Note that, because of the presence of $H_{0}$ in Eq. 9, there are corrections
to Eq. 10 of order ${\cal O}(1/\omega)$ instead of ${\cal O}(1/\omega^{3})$
(see Appendix B). However, these additional terms have amplitudes $\kappa
J_{1}J/\omega$ i.e, from Eq. 11, of the order of $J_{1}J^{2}/\omega^{2}$ and,
hence, can be negleted compared to the terms in Eq. 10.

##  3 Results

###  3.1 Validity of the Magnus high-frequency expansion

In order to check the relevance of the Floquet Hamiltonian to describe the
effective dynamics, we first consider the simplest possible setup, namely the
application of the fast chiral drive alone on the NN RVB state, i.e. we set
$x=0$ and $\lambda(t)=1$ in Eq. 9. Our initial state is invariant under $P$
and $T$ while $H_{F}$ is odd so we expect
$\big{<}H_{F}\big{>}\equiv\big{<}\Psi_{\rm Floquet}(t)|H_{F}|\Psi_{\rm
Floquet}(t)\big{>}=0$, i.e. the stroboscopic motion occurs in the zero-
chirality manifold. Considering a $4\times 4$ finite cluster with PBC (torus),
we have checked numerically that this property is even exactly realized during
the micromotion under $H_{\rm drive}(t)$ at every time $t$. Indeed, if one
considers any 4-site oriented plaquette $(i,j,k,l)$, $\big{<}\Psi(t)|{\bf
S}_{i}\cdot({\bf S}_{j}\times{\bf S}_{k})|\Psi(t)\big{>}=-\big{<}\Psi(t)|{\bf
S}_{k}\cdot({\bf S}_{l}\times{\bf S}_{i})|\Psi(t)\big{>}$ at all times, i.e
facing triangles in all plaquettes have always opposite chiralities (of small
magnitudes). Hence, the overall plaquette chirality
$\big{<}\Psi(t)|i(P_{ijkl}-P_{ijkl}^{-1})|\Psi(t)\big{>}$ is exactly vanishing
at all times.

To have a more precise comparison between the exact and the effective
dynamics, we have computed the infidelity w.r.t. to the initial state, i.e.
${\cal I}(t)=1-|\big{<}\Psi(t)|\Psi_{\rm RVB}\big{>}|$, where
$|\Psi(t)\big{>}$ is computed using either $H_{\rm drive}(t)$
($|\Psi(t)\big{>}=|\Psi_{\rm drive}(t)\big{>}$) or $H_{F}$
($|\Psi(t)\big{>}=|\Psi_{\rm Floquet}(t)\big{>}$), starting from the same
initial state $|\Psi_{\rm RVB}\big{>}$ at $t_{0}=0$. The results in Fig. 2(a)
show an excellent agreement. A zoom at small time in Fig. 2(b) shows a clear
oscillation within every time period (micromotion), while the motion at
(discrete) stroboscopic times $(p/N_{\tau}+n)T_{\rm drive}$, p fixed, falls on
smooth curves when varying $n$. The result at $p=0$, that is at times exact
multiples of $T_{\rm drive}$, compares very well with the calculation obtained
using $H_{F}$, i.e. approximating $|\Psi(t)\big{>}$ by $|\Psi_{\rm
Floquet}(t)\big{>}=\exp{(-iH_{F}t)}|\Psi_{\rm RVB}\big{>}$.

However, a zoom at a longer time $t\sim 30$ in Fig. 2(c) shows small
deviations between the stroboscopic motion computed with the drive Hamiltonian
using two choices of $N_{\tau}$ and the effective Floquet motion (computed
exactly). These deviations may have two origins: first, the Trotter errors in
the time discretization of the drive oscillations and, secondly, the neglected
higher order terms in the Floquet Hamiltonian. The truncation of the Floquet
Hamiltonian leads to relative errors of order $1/\omega^{2}$, i.e. $\sim
2\times 10^{-4}$ which can be neglected. In contrast, the Trotter errors are
still sizeable at $N_{\tau}=8$. Using $N_{\tau}=16$ Trotter steps per drive
period seems to give already a much better accuracy, typically within $2\%$ of
the (exact) effective Floquet motion. Interestingly, we also observe that the
amplitude of the micromotion oscillations becomes quite small compared to the
average behavior, except at very small time.

###  3.2 Quench of the fast drive

In this section we shall assume the drive Hamiltonian (10) is suddenly switch
on at $t=0$ on the RVB initial state, i.e. we assume $\lambda(t)=1$,
independent of $t$. In that case, the effective motion follows the unitary
evolution given by $U(t)=\exp{(-i\kappa H_{\rm target}t)}$, which depends on
the control parameters $J_{1}$, $J_{2}$ and $\lambda_{\rm chiral}$. The case
$J_{1}=J_{2}=0$ corresponds to the previously studied case $H_{0}=0$ i.e. a
purely oscillatory drive. If the constant term $H_{0}$ is present we set
$J_{1}=1$.

Figure 3:  Diagonal ensembles: comparison between distributions of diagonal
density matrices $\rho_{\mu}$ vs $E_{\mu}$, the eigenvalue spectrum of $H_{\rm
target}=\frac{\lambda_{\rm chiral}}{J_{F}}H_{F}$; (a) pure drive $\kappa=0$
i.e. $H_{\rm target}$ given by $J_{1}=J_{2}=0$ and $\lambda_{\rm chiral}=0.5$;
(b) $\kappa\neq 0$ given by Eq. 11, $H_{\rm target}$ given by $J_{1}=1$,
$J_{2}=0$ and $\lambda_{\rm chiral}=0.7$; (c) Same as (b) but for $J_{2}=0.4$
and $\lambda_{\rm chiral}=0.5$. Energy levels are represented by black dotted
lines and their weight by segments whose 1/2-length is equal to $\rho$. Note
that $\sum_{\mu}\rho_{\mu}=1$. The first singlet excited state, forming
together with the GS the pair of relevant CFT states [50], is shown by a red
line in (b) and (c).

####  3.2.1 Long-time average and diagonal ensembles

Since the effective Floquet dynamics is expected to be very accurate, one can
use it to investigate the physics in the $t\rightarrow\infty$ limit, after a
quench of the drive Hamiltonian, i.e. assuming $\lambda=1$ in Eq. 10.
Decomposing the initial (non-chiral) spin liquid state in terms of the Floquet
eigenstates $|e_{\mu}\big{>}$ of $H_{F}$, $|\Psi_{\rm
RVB}\big{>}=\sum_{\mu}w_{\mu}|e_{\mu}\big{>}$, the infinite-time average of
observables $\cal O$ can be written in terms of the diagonal density matrix
$\rho_{\rm diag}=\sum_{\mu}\rho_{\mu}|e_{\mu}\big{>}\big{<}e_{\mu}|$ [39],

$\big{<}{\cal O}\big{>}_{\rm ave}=\sum_{\mu}\overline{{\cal
O}}_{\mu\mu}\rho_{\mu}\,,$ (12)

where $\rho_{\mu}=|w_{\mu}|^{2}$. Within a time period, the wavefunction
$|\Psi(nT+t)\big{>}$ is subject to a unitary rotation $\exp{(-iK(t))}$ where
$K(t)$ is a periodic kick operator (see explicit expression in Appendix C).
Hence, the infinite-time average in Eq. 12 involves the averaged operator,

$\overline{\cal O}=\frac{1}{T}\int_{0}^{T}\exp{(iK(t))}{\cal
O}\exp{(-iK(t))}dt\,.$ (13)

Note that, for the observables used in this work, the effect of the kick
rotations can be neglected so that $\overline{{\cal O}}_{\mu\mu}\simeq{\cal
O}_{\mu\mu}$ in Eq. 12. Note also that, since the RVB state is a singlet state
invariant under all space group operations, its decomposition in terms of
Floquet eigenstates is limited to the most symmetric singlet sector.

When $\kappa=0$ (sinusoidal drive alone) the initial RVB state is a highly
excited state of $H_{F}^{0}$. Indeed, $H_{F}^{0}$ has a symmetric spectrum
with $\pm E_{\mu}$ opposite eigenvalues associated to pairs of complex-
conjugate eigenvectors $|e_{\mu}^{0},-\big{>}=|e_{\mu}^{0},+\big{>}^{*}$. The
RVB state being a real (T-invariant) wavefunction, its decomposition is of the
form

$|\Psi_{\rm
RVB}\big{>}=\sum_{\mu}c_{\mu}^{0}(|e_{\mu}^{0},-\big{>}+|e_{\mu}^{0},+\big{>})\,,$
(14)

where $c_{\mu}^{0}$ are real coefficients. The $E_{\mu}\leftrightarrow-
E_{\mu}$ symmetry of the weight distribution is indeed clear in Fig. 3(a). One
also immediately gets $\big{<}H_{F}^{0}\big{>}_{\rm ave}=0$, in agreement with
the fact that the chirality remains zero under time evolution, as noted
before. Interestingly, we see that the initial RVB state has a very large
overlap with the zero energy states located at the middle of the many-body
spectrum of $H_{F}^{0}$. For larger $L\times L$ systems of increasing
$N=L^{2}$ sites, a broader weight distribution is expected, still centered
around the spectrum center $E=0$, so that the system can be viewed as being at
infinite temperature. A more exotic scenario would be realized if the RVB
state keeps a finite overlap on the $E=0$ state manifold suggesting the
existence of low-entangled ”scar” states [57, 58] (i.e. following the area
law) in the center of the spectrum and, hence, a lack of thermalization.

Weight distributions for $\kappa\neq 0$, i.e. when a constant term is present
in the drive Hamiltonian, are shown in Fig. 3(b,c), in the cases corresponding
to $H_{\rm target}$ given by $J_{1}=1$, $J_{2}=0$ and $\lambda_{\rm
chiral}=0.7$, and by $J_{1}=1$, $J_{2}=0.4$ and $\lambda_{\rm chiral}=0.5$,
respectively. Then, the spectrum of $H_{\rm target}=\frac{1}{\kappa}H_{F}$ is
no longer symmetric and one gets finite values for the time-averaged plaquette
chirality ${\rm Im}\big{<}P_{ijkl}\big{>}_{\rm ave}\sim 0.2496$ and $0.2334$
in the (b) and (c) cases, respectively. The GS and the 1st singlet excited
state (shown as a red line in Figs. 3(b,c)) have been identified from a
Conformal Field Theory (CFT) construction as the two relevant states forming
the CSL doubly-degenerate GS manifold on the infinite-size torus [50]. Note
that the RVB state has no overlap with the second CFT state which possesses
different quantum numbers. However, large weights are found on the GS and on
the next lowest fully symmetric singlet eigenstates of $H_{F}$.

####  3.2.2 Time evolution at small and intermediate times

Figure 4: Plaquette chirality vs time for a quench of the drive computed on a
$4\times 4$ torus. Parameters are $J_{1}=1$, $J_{2}=0$ and $\lambda_{\rm
chiral}=0.7$ or $J_{1}=1$, $J_{2}=0.4$, and $\lambda_{\rm chiral}=0.5$. (a)
Short and intermediate times behavior; all computations involving the fast
drive are done with $N_{\tau}=16$ intervals in each drive period but only data
at times $kT_{\rm drive}$ and $(k+1/2)T_{\rm drive}$, $k\in{\mathbb{N}}$, are
reported since the amplitude of the oscillations within each time period
(micromotion) remains very small. Computations using the effective Floquet
dynamics with $U(t)=\exp{(-i\kappa H_{\rm target}t)}$ shown as dashed-dotted
curves agree very well. Results involving a ramping of the drive are also
shown for comparison. (b,c) Long time behavior obtained using the effective
Floquet dynamics; the chirality (blue curve) is shown to oscillate around the
average value (black dashed line) computed independently using Eq. (12). The
time average $\frac{1}{t}\int_{0}^{t}\\!\text{Im}\langle
P_{\square}\rangle(\tau)\,\mathrm{d}\tau$ is displayed as a red line.

We now turn to the actual computation of the time evolution in the case of a
quench of the drive. Results obtained on the $4\times 4$ torus are shown in
Fig. 4(a). The constant $H_{0}$ of the quenched Hamiltonian (10) involves
either $J_{1}=1$ and $J_{2}=0$ or $J_{1}=1$ and $J_{2}=0.4$. The amplitude
$\kappa$ of $H_{0}$ is set using Eq. 11 i.e.
$\kappa=J^{2}/(4\omega\lambda_{\rm chiral})$, and we choose either
$\lambda_{\rm chiral}=0.7$ when $J_{2}=0$ or $\lambda_{\rm chiral}=0.5$ when
$J_{2}=0.4$, corresponding to the two realistic $H_{\rm target}$ whose spectra
are shown in Fig. 3(b,c). We observe a rapid increase of the plaquette
chirality, initially vanishing in the RVB state. Note that the amplitude of
the micromotion fast oscillations remains quite small, barely visible on the
scale of the plot, so that the full time evolution can be compared directly to
the effective Floquet dynamics, not only at stroboscopic times. We observe
that the effective Floquet approach is able to capture accurately the slow
variations of the chirality expectation value up to long times. Even more, we
believe that the small deviation at intermediate times with the full time-
dependent simulations is primarily due to the Trotter error in the later
rather than to the truncation of the Floquet Hamiltonian at order $1/\omega$.
At times $t>t_{F}$ we observe in Fig. 4(b) some strong oscillations around the
expected long-time average values computed from Eq. 12. This is clearly a
finite size effect and one expects the amplitude of the oscillations to
decrease rapidly with increasing system size. In any case, we believe that
quenching the drive Hamiltonian is not the correct procedure to achieve the
goal of preparing the system in the final CSL state and we now move to a
different setup.

###  3.3 Adiabatic ramping of the fast drive

Figure 5:  Plaquette chirality vs reduced time $t/t_{\rm ramp}$ for adiabatic
ramping of the drive using different values of $t_{\rm ramp}$, computed on a
$4\times 4$ torus. (a) The RVB state is used as initial state, approximating
the GS of the frustrated Heisenberg quantum antiferromagnet $H_{0}$ defined by
$J_{1}=1$ and $J_{2}=0.4$. $\kappa$ is given by Eq. 11 with $\lambda_{\rm
chiral}=0.5$. All time-dependent computations (shown as continuous curves) are
done with $N_{\tau}=16$ intervals in each drive period. (b) The GS of the
Heisenberg quantum antiferromagnet $H_{0}$ ($J_{1}=1$, $J_{2}=0$) is used as
initial state. $\kappa$ is given by Eq. 11 with $\lambda_{\rm chiral}=0.7$.
(a,b) Results using the (weakly time-dependent) effective Floquet Hamiltonian
of Eq. 10 are shown as dashed-dotted curves. The (ideal) adiabatic limits are
shown by green dotted lines.

The new setup we shall follow now consists in branching the fast drive in a
quasi-adiabatic way. For that we have chosen a simple linear ramp,

$\lambda(t)=t/t_{\rm ramp}\,,$ (15)

for $t\in[0,t_{\rm ramp}]$. In a truly infinitely-slow ramping $\lambda(t)$
from 0 to 1 (e.g. taking $t_{\rm ramp}\rightarrow\infty$), the state
$|\Psi(t)\big{>}$ is expected to follow adiabatically the GS of

$H_{F}(t)\simeq\kappa[H_{0}+\lambda^{2}(t)(H_{\rm target}-H_{0})]\,,$ (16)

provided that $|\Psi(0)\big{>}$ is the GS of $H_{0}$. In practice, deviations
from this ideal scenario are present. First, the parent Hamiltonian of the RVB
state in the thermodynamic limit involves longer range interactions [45, 46]
than the NN and NNN interactions considered here. However, we have checked
that, on a finite system, the RVB state is a good approximation of the exact
GS of $H_{0}$ for $J_{2}=0.4J_{1}$, with an overlap of order $0.9960$ on a
$4\times 4$ cluster. Secondly, $t_{\rm ramp}$ is finite leading to a deviation
from the adiabatic behavior. It is then interesting to try to quantify the
effect of a finite ramp time. Note that, although the branching of the drive
is linear in $t$, the effective coupling to $H_{\rm target}$ in (16) is
smoother, increasing only as $t^{2}$ at small time. Note also that the initial
RVB state has a sizeable overlap with the targeted final state, as shown in
Fig. 3(c), which should enable an adiabatic evolution in a shorter time.

We have computed the time evolution (using $N_{\tau}=16$ TS steps per drive
period $T_{\rm drive}=0.1$) and the emergence of the plaquette chirality for
three ramp times of increasing length, $t_{\rm ramp}=80,160$, and $240$ as
plotted in Fig. 4(a), showing a smooth behavior of the observable, in contrast
to the quench setups. As expected, the increase of the chirality is slower and
slower when increasing the ramp time but the final chirality reached at
$t=t_{\rm ramp}$ becomes larger and larger.

Fig. 5(a) shows similar data as a function of the reduced time $t/t_{\rm
ramp}$. The adiabatic groundstate $|\Psi_{\rm adiab}(\lambda)\big{>}$ is
obtained by following the GS of the Floquet Hamiltonian (16), varying
$\lambda$ ”by hand” from $0$ to $1$. It should be reached asymptotically by
taking the limit $t_{\rm ramp}\rightarrow\infty$ keeping the ratio $t/t_{\rm
ramp}=\lambda$ fixed. In fact, the GS chirality of $H_{\rm target}$ is
$\sim\\!0.478$ while the achieved chirality at $t=t_{\rm ramp}=320$ is only
$\sim\\!0.383$. This indicates that longer ramp times are needed to achieve
the goal of approaching more closely the GS of $H_{\rm target}$. However,
running the full drive simulation over longer times is expensive in CPU-time
and, in fact, not necessary. Indeed, we have checked that for $t_{\rm
ramp}=320$, the effective evolution with the slowly varying Floquet
Hamiltonian $H_{F}(t)$ gives very similar results as seen in Fig. 5(a). In
particular, the overlap $|\big{<}\Psi_{\rm Floquet}(t)|\Psi_{\rm
drive}(t)\big{>}|$ between the two final states at $t=t_{\rm ramp}$ is
$0.99908$, very close to 1. The effective Floquet dynamics allows to use much
longer Trotter steps, typically of order 0.5, compared to $T_{\rm
drive}/N_{\tau}=6.25\times 10^{-3}$ for the full dynamics, and is therefore a
much cheaper alternative for longer times. Hereafter, further results will
then be obtained using Floquet dynamics, showing no appreciable difference on
the plots and performed easily up to $t_{\rm ramp}=960$ (i.e. around $10\,000$
drive time periods!). For $t_{\rm ramp}>500$, while approaching the adiabatic
limit, one starts to see oscillations in the time evolution of the chirality
which are attributed to the deviation of the initial RVB state from the true
GS of $H_{0}$.

To confirm the origin of the spurious oscillations found above, we consider
the simple quantum Heisenberg antiferromagnet, setting $J_{1}=1$, $J_{2}=0$ in
$H_{0}$, and use its GS as initial state, instead of the RVB state. Note that,
for this new choice of $H_{0}$, the overlap between the RVB state and the true
GS is only $0.9256$ on the $4\times 4$ torus and, hence, the RVB is no longer
an appropriate choice for the initial state. Results are shown in Fig. 5(b),
showing a smooth evolution towards the adiabatic limit when increasing $t_{\rm
ramp}$, even at large $t_{\rm ramp}$, in contrast to Fig. 5(a).

###  3.4 Fidelities of the final CSL state

The results shown in the previous sections strongly suggest that the ramp
setup can give efficient implementations of the CSL state. Here we try to
analyse more quantitatively the properties of the final state.

Figure 6: Overlap $|\langle\Psi_{t_{\rm ramp}}(t)|\Psi_{\rm
adiab}(\lambda=t/t_{\rm ramp})\rangle|$ between the time-evolved state and the
adiabatic GS computed at $\lambda=t/t_{\rm ramp}$, plotted versus $t/t_{\rm
ramp}$. Here the time-evolved state is obtained using the simplified Floquet
dynamics with $J_{1}=1$, $J_{2}=0$ and $\lambda_{\rm chiral}=0.7$ and various
finite ramp times $t_{\rm ramp}$ up to $t_{\rm ramp}=960$ are considered. The
GS of $H_{0}$ is taken as initial state. The inset is a zoom of the same
quantity for $t/t_{\rm ramp}\leq 0.5$.

####  3.4.1 Comparison with the adiabatic limit

First, we show in Fig. 6, for various values of $t_{\rm ramp}$, the overlap of
the time-evolved state $|\Psi_{t_{\rm ramp}}(t)\big{>}$ with the adiabatic GS
$|\Psi_{\rm adiab}(\lambda)\big{>}$, computed at $\lambda=t/t_{\rm ramp}$, as
a function of the reduced time $t/t_{\rm ramp}$. As expected, we observe that,
for a given ramp time, the fidelity deteriorates as a function of time.
However, for a fixed ratio $\lambda=t/t_{\rm ramp}$ defining a target (chiral)
state $|\Psi_{\rm adiab}(\lambda)\big{>}$ GS of a target Hamiltonian where
$\lambda_{\rm chiral}$ is replaced by $\lambda^{2}\,\lambda_{\rm chiral}$, the
fidelity quickly increases with $t_{\rm ramp}$. For $\lambda=1$, corresponding
to $\lambda_{\rm chiral}=0.7$ in $H_{\rm target}$, a fidelity of
$\sim\\!0.854$ can be obtained using $t_{\rm ramp}\simeq 1000$. For
$\lambda=0.8$, corresponding to $\lambda_{\rm chiral}\simeq 0.45$, a fidelity
of $\sim\\!0.953$ can be reached with the same value of $t_{\rm ramp}$, but
after a shorter time $\lambda\,t_{\rm ramp}\simeq 800$.

Figure 7: Overlap $|\big{<}\Psi_{t_{\rm ramp}}|\Psi_{\rm PEPS}\big{>}|$
between the time-evolved state at $t=t_{\rm ramp}$ and the best optimized PEPS
states, plotted versus $t_{\rm ramp}$. (a) $J_{1}=1$, $J_{2}=0.4$ and
$\lambda_{\rm chiral}=0.5$, and NN RVB state as initial state; (b) $J_{1}=1$,
$J_{2}=0$ and $\lambda_{\rm chiral}=0.7$, and GS of $H_{0}$ as initial state.
Different symbols correspond to different PEPS manifold, as indicated in the
legend. The overlap $|\big{<}\Psi_{t_{\rm ramp}}|\Psi_{\rm init}\big{>}|$ with
the initial state is shown as a reference. The adiabatic limit expected for
$t_{\rm ramp}\rightarrow\infty$ is shown on the right of the picture.

####  3.4.2 Comparison with chiral PEPS

We expect the final state to be in a CSL phase, in the close neighborhood of
the GS of the target Hamiltonian. The later has been shown to be represented
faithfully by chiral PEPS ansatze [42, 54]. It is therefore instructive to
compare our final states prepared at $t=t_{\rm ramp}$ with optimized PEPS. We
shall consider fully symmetric, i.e. SU(2) and translationally invariant, PEPS
with bond virtual spaces $\nu=0\oplus 1/2$ ($D=3$), $\nu=0\oplus 1/2\oplus
1/2$ ($D=5$) and $\nu=0\oplus 1/2\oplus 1$ ($D=6$). Within these symmetric
PEPS manifolds, restricted $A_{1}+iA_{2}$ chiral sub-manifolds (with same bond
dimensions) are designed to represent our CSL (see Appendix A for more
details). While generic manifolds possess $(n_{\rm tensors}-1)$ complex
independent parameters, where $n_{\rm tensors}$ is the total number of $A_{1}$
and $A_{2}$ tensors (see Table 1 for numbers), the chiral PEPS manifolds (with
a fixed relative phase between the $A_{1}$ and $A_{2}$ tensors) have only half
degrees of freedom i.e. $(n_{\rm tensors}-1)$ real parameters to be optimized.

As shown in Figs. 7, as expected the $D=3$ PEPS cannot accommodate very well,
for increasing $t_{\rm ramp}$, the increase of entanglement of $|\Psi_{\rm
t_{\rm ramp}}\big{>}$ signaled by the rapid decrease of the fidelity w.r.t.
the initial state. However, one already sees that the performance of the
chiral PEPS (with only two variational parameters) becomes comparable to the
generic (symmetric) PEPS for $t_{\rm ramp}>250$, suggesting a crossover
towards the CSL phase. The extension of the chiral PEPS to $D=5$ or $D=6$
unfortunately does not improve the overlaps significantly (see Table 1) but a
similar trend is found as shown in Fig. 7(b) for $D=6$.

##  4 Discussions and conclusions

We first discuss possible experimental realisations of our setup, how to
realize the fast drive and how to prepare the initial state.

###  4.1 Realizing the drive Hamiltonian

The Heisenberg interactions in Eq. 5 can be realized with cold atoms loaded on
an optical lattice [59], for which a microscopic Hubbard-like description
applies [60, 61]. The excellent isolation of the optical lattice from its
environment enables to realize ideal unitary dynamics. We assume one atom per
site (on average) and a large on-site repulsion $U$ to be in the Mott
insulating phase [62], and different hopping terms $t_{v}(t)$ on the four NN
bonds $v=a,b,c,d$ of Fig. 1. Our setup is relatively simple as it involves
only the modulation in time and space of the NN (effective) Heisenberg
couplings. A (selective) modulation in time of the Heisenberg couplings
$J_{v}(t)$ of $H_{v}$ can be realized via a (selective) modulation of the
hoppings between sites $t_{v}(t)$, i.e. changing the tunneling barriers of the
optical lattice appropriately. Note that this procedure is conceptually
different from the synthetic gauge field approach [18], in which the chiral
term would appear in the Mott insulating phase only in fourth-order of the
(complex) hopping $t$, i.e. $\lambda_{\rm chiral}\sim t^{4}/U^{3}$. However,
in fermionic Hubbard systems the Heisenberg couplings can only be
antiferromagnetic ($J_{v}=4(t_{v}(t))^{2}/U>0$), not ferromagnetic ($J_{v}<0$)
and, strictly speaking, the single-harmonic modulations of (5) cannot be
realized 111A constant NN antiferromagnetic interaction $J$ cannot be added to
all bonds since it will contribute to the constant $H_{0}$ part of $H_{F}$.. A
more practical alternative would then be to use, on the optical lattice, two-
component bosons [63, 64] which (i) are more easily brought to low
temperatures and (ii) are more versatile than fermions to realize alternating
antiferromagnetic and ferromagnetic couplings, tuning the
$U_{\uparrow\uparrow}=U_{\downarrow\downarrow}$ and $U_{\uparrow\downarrow}$
on-site Hubbard interactions.

It is interesting to mention that a simplified drive sequence consists in
replacing the sinusoidal modulation by a square wave i.e.

$\displaystyle{\tilde{H}}_{\rm drive}(t)$ $\displaystyle=$ $\displaystyle
H_{x}\,\hskip 8.5359pt\text{for}\hskip
8.5359ptt\in[-\frac{T}{4},\frac{T}{4}],$ $\displaystyle{\tilde{H}}_{\rm
drive}(t)$ $\displaystyle=$ $\displaystyle\,H_{y}\hskip
8.5359pt\text{for}\hskip 8.5359ptt\in[0,\frac{T}{2}],$
$\displaystyle{\tilde{H}}_{\rm drive}(t)$ $\displaystyle=$
$\displaystyle-\,H_{x}\hskip 8.5359pt\text{for}\hskip
8.5359ptt\in[\frac{T}{4},\frac{3T}{4}],$ $\displaystyle{\tilde{H}}_{\rm
drive}(t)$ $\displaystyle=$ $\displaystyle-\,H_{y}\hskip
8.5359pt\text{for}\hskip 8.5359ptt\in[\frac{T}{2},T],$ (17)

This new type of drive sequence introduces higher harmonics
$\frac{2}{n\pi}(\exp{(in\omega t)\pm c.c.)}$ in the amplitude modulations of
$H_{x}$ and $iH_{y}$ which modify the prefactor $J_{F}$ in Eq. 7 into

${\tilde{J}}_{F}=\frac{14\zeta(3)}{\pi^{2}}J_{F}\simeq 1.72J_{F},$ (18)

and the chiral term in Eqs. 7 and 10 accordingly. Note also the emergence of
chiral terms in the (targeted) Floquet effective Hamiltonian at order ${\cal
O}(\omega^{-2})$ instead of ${\cal O}(\omega^{-3})$. Nevertheless, we believe
the main ideas/results developed in this work are still valid qualitatively in
the case of (17). Finally, we would like to mention that Rydberg atom arrays
can also emulate quantum spin systems [65], but keeping SU(2) spin rotational
invariance might be challenging. Nevertheless, it would be interesting to
investigate whether our setup can be exported to the case of anisotropic XXZ
dipolar spin interactions.

###  4.2 Preparation of the initial state

Spin liquids of the RVB type can be realized on cold atom platforms i.e. in
Mott insulating ultracold bosons [66]. Using a two-component one-dimensional
Bose–Hubbard system, a highly entangled Heisenberg antiferromagnet with SU(2)
symmetry has been experimentally realized [67]. We believe that this promising
approach could be extended to two dimensions. Note however that, as shown by
our studies, it is important that the initial state in the ramp process is as
close as possible to the GS of the Hamiltonian at $t=0$ i.e. $H(0)=H_{0}$.
Using cold atoms on an optical square lattice, the Mott phase can be realized
on fermionic systems of hundred of sites [68, 69] with magnetic correlation
lengths of several lattice spacings [70]. Similarly, systems of two-component
bosons could also be brought in Mott phases at large Hubbard repulsions
$U_{\sigma\sigma^{\prime}}$ [64]. Owing to rapid experimental progress in
decreasing the temperature, we believe the quantum state of such systems may
soon approach the (finite size) GS of a plain spin-1/2 NN quantum Heisenberg
model (i.e. with $J_{2}=0$) and could serve as an initial state for the second
type of ramping. Tweezer arrays of fermions have also been proposed [71], that
could be used to reach the same goal. We point out that the precise nature of
the initial state (e.g. spin liquid vs Néel state) is not completely settled
on finite size systems which are always quantum-mechanically disordered (i.e.
are in a global singlet state) and possess a small gap in their many-body
spectrum. In other words, the distinction between different phases is not
sharp on finite size systems.

###  4.3 Final remarks

In summary, we have studied different protocols to prepare a (finite size)
system of quantum spin-1/2 on a square lattice into a spin-1/2 chiral spin
liquid phase. We start from a simple low-entangled (SU(2)-symmetric) initial
state like a (non-chiral) spin liquid (e.g. a RVB spin liquid) or the (finite
size) GS of a simple NN Heisenberg quantum antiferromagnet which, we believe,
could (or will soon) be realized experimentally in platforms of cold atoms on
an optical lattice. A simple drive involving dephased periodic modulations of
four types of NN bonds is proposed. We focus on the high-frequency/large
amplitude (strong drive) regime for which an effective Floquet dynamics is
shown to give very accurate results. We argue that our goal cannot be realized
by branching the drive suddenly, but rather by a semi-adiabatic ramping.
Quantitative results are given for a periodic $4\times 4$ system. In that
case, we show that the (finite size) CSL ground state of a relevant target
Hamiltonian can be reached with high fidelity typically within less than
$\\!10$ thousands drive oscillations, still beyond current experimental
possibilities.

Recent experiments are now performed on systems with hundred sites or even
more, and with a box-like confining potential (corresponding theoretically to
open boundary conditions). The use of periodic boundary conditions in our
system was in fact dictated by convenience, translation invariance making
calculations simpler. However, we believe no conceptual difference is to be
expected for open boundaries. The case of a harmonic confining potential is
more tedious and left for future investigations. In any case, the system size
$N$ is an important relevant parameter. Indeed, in the transition region
between the non-chiral initial phase and the CSL, the gap should reach a
minimum value, typically of order $1/N$, or a power of $1/N$. Hence, the
minimum ramp time to reach a good fidelity should qualitatively scales like
$N$, which could be a limitation for $N$ of a few hundred sites. However we
believe that preparation times could be drastically reduced using quantum
control algorithms, an interesting issue left for future studies.

Note that our computation has been performed on the square lattice but similar
setups can also be realized in frustrated 2D lattices like kagome or
triangular lattices which can be obtained experimentally [72]. In that case,
dephased (fast) modulations of 3 or 6 groups of equivalent NN bonds can be
performed leading to chiral permutations on triangles or hexagons. The
resulting topological CSL state and its anyonic excitations may be used for
topological quantum computations following Kitaev’s proposal [49].

Lastly, we note that the implementation of a dimer spin liquid state (of the
toric code universality class) has been proposed on the kagome lattice by
using Rydberg atom arrays [65]. It would then be interesting to investigate
whether a dimer CSL could be realized, starting from such an initial state, by
Floquet engineering involving dimer degrees of freedom (represented by excited
Rydberg atoms in the blockade regime) instead of the original spin-1/2 degrees
of freedom considered in this work.

##  Appendix A PEPS representation of spin liquids

PEPS are defined as a tensor network of local onsite tensors of size $dD^{4}$
where $d=2$ is the local physical Hilbert space dimension (spin 1/2) and $D$
is the bond dimension i.e. the number of virtual states on each tensor leg. By
contracting the tensors w.r.t. the bond indices, one generates an entangled
ansatz (whose entanglement entropy per bond is limited by $\ln{D}$ and, hence,
controlled by $D$). PEPS offer an extremely efficient variational scheme to
address local Hamiltonians [73]. In the infinite-PEPS (iPEPS) algorithm [40]
infinite-size systems can be handled.

Both the initial state and the Floquet Hamiltonian governing the effective
dynamics of our problem are SU(2) symmetric and transform according to the
trivial representation of the square lattice space group222Note that the time-
dependent drive (5) breaks the $C_{4v}$ point group symmetry causing small
(negligible) lattice asymmetry of the time-evolved state at short time
$t<1/J$. . It is therefore relevant to enforce such symmetries at every step
of our scheme. A symmetric spin liquid on the square lattice is simply
obtained by choosing a $C_{4v}$-symmetric site tensor placed uniformly on
every lattice sites. As explained in details in Ref. [74], the continuous
SU(2) of PEPS can be implemented at the level of local tensors by imposing
that the $D$-dimensional subspace $\nu$ of each virtual leg is a direct sum of
SU(2) irreducible representations. Working with SU(2)-symmetric PEPS allows to
greatly reduce the number of parameters describing the PEPS family and is a
major advantage for computations based on optimization.

The NN RVB state [35] which consists of an equal weight superposition of all
NN singlet (valence bond) coverings333On a bipartite lattice, NN singlets are
defined in the bosonic representation oriented from one sublattice to the
other.  has a simple representation using $\nu=0\oplus 1/2$ virtual bonds
$(D=3)$ [45, 46]. In fact, the $D=3$ symmetric PEPS manifold is represented by
two linearly independent fully symmetric ($A_{1}$) real site tensors
(describing the fusion of either one or three virtual spin-1/2 into a physical
spin-1/2) and generically describes gapped ${\mathbb{Z}}_{2}$ topological RVB
spin liquids which include (thanks to “teleportation”) longer range singlet
bonds [47]. The critical NN RVB state is therefore a special point in this
two-parameter manifold [48], with exponentially large correlation lengths in
its neighborhood [47, 48].

The simplest PEPS representation of the (Abelian) spin-1/2 CSL – of the
Kalmeyer-Laughlin type [2] – is in fact possible within the $\nu=0\oplus 1/2$
virtual space ($D=3$) manifold by adding a third site tensor of different
$A_{2}$ orbital symmetry with a pure-imaginary amplitude, resulting into a
$A_{1}+iA_{2}$ site tensor [41, 42]. Such an ansatz breaks time-reversal T
($i\rightarrow-i$) and parity P ($A_{2}\rightarrow-A_{2}$) while being
invariant under their product PT, as expected for a CSL. Also, despite its
simplicity, it gives a faithful description of the chiral edge modes – a key
feature of the CSL – following closely the prediction of a SU($2$)1 Wess-
Zumino-Witten (WZW) conformal field theory (CFT). Moreover, the CSL ground-
state [75] of the target Hamiltonian (3) has been successfully studied using
similar PEPS [51, 54]. Increasing the bond dimension from $D=3$ to $D=5$
($\nu=0\oplus 1/2\oplus 1/2$) or $D=6$ ($\nu=0\oplus 1/2\oplus 1$) provides a
small improvement of the variational energy. We show in Table 1 the optimal
PEPS overlaps on a $4\times 4$ torus with the adiabatic GS considered in this
work.

D | 3 | 5 | 6
---|---|---|---
no of $A_{1}$ tensors | 2 | 10 | 11
no of $A_{2}$ tensors | 1 | 8 | 8
no of variational parameters | 2 | 17 | 18
$|\big{<}\Psi_{\rm adiab}^{(a)}|\Psi_{\rm PEPS}\big{>}|$ | 0.86470 | 0.86648 |
$|\big{<}\Psi_{\rm adiab}^{(b)}|\Psi_{\rm PEPS}\big{>}|$ | 0.90525 | 0.91256 | 0.91319
Table 1: Number of $A_{1}$ and $A_{2}$ tensors and number of variational
parameters (fixing one tensor amplitude) used to construct the chiral PEPS, as
a function of the bond dimension. The last lines show the corresponding
overlaps with the two adiabatic GS considered in this work.

##  Appendix B Derivation of the Floquet effective Hamiltonian

###  B.1 $H_{0}=0$ case

We start from the expression of the Floquet Hamiltonian at order $1/\omega$ –
see Eq. (42) in [39] – decomposing $H_{\rm drive}$ as $H_{1}\exp{(i\omega
t)}+H_{-1}\exp{(-i\omega t)}$,

$\displaystyle H_{F}^{0}$ $\displaystyle=$
$\displaystyle\frac{1}{\omega}[H_{1},H_{-1}]\,,$ (19) $\displaystyle=$
$\displaystyle\frac{i}{2\omega}[H_{x},H_{y}]\,,$ $\displaystyle=$
$\displaystyle 2J_{F}\,i[h_{x},h_{y}]\,,$

where

$\displaystyle h_{x}$ $\displaystyle=$ $\displaystyle H_{x}/J=\sum_{i\in
A}({\bf S}_{i}\cdot{\bf S}_{i+e_{x}}-{\bf S}_{i}\cdot{\bf S}_{i-e_{x}}),$
$\displaystyle h_{y}$ $\displaystyle=$ $\displaystyle H_{y}/J=\sum_{i\in
A}({\bf S}_{i}\cdot{\bf S}_{i+e_{y}}-{\bf S}_{i}\cdot{\bf S}_{i-e_{y}})\,.$
(20)

It is interesting to note that, in this simple case of an harmonic drive
alone, corrections only appear at order $1/\omega^{3}$. One can then expand
the commutator as $[h_{x},h_{y}]=\sum_{i\in A}(C_{i}+C_{i}^{\prime})$ where,

$\displaystyle C_{i}$ $\displaystyle=$
$\displaystyle(S_{i+e_{x}}^{\alpha}-S_{i-e_{x}}^{\alpha})[S_{i}^{\alpha},S_{i}^{\beta}](S_{i+e_{y}}^{\beta}-S_{i-e_{y}}^{\beta})\,,$
$\displaystyle C_{i}^{\prime}$ $\displaystyle=$
$\displaystyle(S_{i}^{\alpha}-S_{i+2e_{x}}^{\alpha})[S_{i+e_{x}}^{\alpha},S_{i+e_{x}}^{\beta}](S_{i+e_{x}-e_{y}}^{\beta}-S_{i+e_{x}+e_{y}}^{\beta})\,.$
(21)

Using
$[S_{i}^{\alpha},S_{i}^{\beta}]=i\epsilon_{\alpha\beta\gamma}S_{i}^{\gamma}$
we get,

$\displaystyle C_{i}$ $\displaystyle=$ $\displaystyle i{\bf S}_{i}\cdot({\bf
S}_{i+e_{x}}\times{\bf S}_{i+e_{y}})+i{\bf S}_{i}\cdot({\bf
S}_{i+e_{y}}\times{\bf S}_{i-e_{x}})$ (22) $\displaystyle+$ $\displaystyle
i{\bf S}_{i}\cdot({\bf S}_{i-e_{x}}\times{\bf S}_{i-e_{y}})+i{\bf
S}_{i}\cdot({\bf S}_{i-e_{y}}\times{\bf S}_{i+e_{x}})\,.$

After the transformation $i+e_{x}\rightarrow j$, $C_{i}^{\prime}\rightarrow
C_{j}$, we obtain for $j\in B$,

$\displaystyle C_{j}$ $\displaystyle=$
$\displaystyle(S_{j-e_{x}}^{\alpha}-S_{j+e_{x}}^{\alpha})[S_{j}^{\alpha},S_{j}^{\beta}](S_{j-e_{y}}^{\beta}-S_{j+e_{y}}^{\beta})$
(23) $\displaystyle=$ $\displaystyle i{\bf S}_{j}\cdot({\bf
S}_{j+e_{x}}\times{\bf S}_{j+e_{y}})+i{\bf S}_{j}\cdot({\bf
S}_{j+e_{y}}\times{\bf S}_{j-e_{x}})$ $\displaystyle+$ $\displaystyle i{\bf
S}_{j}\cdot({\bf S}_{j-e_{x}}\times{\bf S}_{j-e_{y}})+i{\bf S}_{j}\cdot({\bf
S}_{j-e_{y}}\times{\bf S}_{j+e_{x}})\,.$

Combining the last two equations, we get

$\displaystyle i[h_{x},h_{y}]$ $\displaystyle=$ $\displaystyle-\sum_{i\in
A,B}\,\\{\,{\bf S}_{i}\cdot({\bf S}_{i+e_{x}}\times{\bf S}_{i+e_{y}})+{\bf
S}_{i}\cdot({\bf S}_{i+e_{y}}\times{\bf S}_{i-e_{x}})$ (24)
$\displaystyle+\,{\bf S}_{i}\cdot({\bf S}_{i-e_{x}}\times{\bf
S}_{i-e_{y}})+{\bf S}_{i}\cdot({\bf S}_{i-e_{y}}\times{\bf S}_{i+e_{x}})\,\\}$
$\displaystyle=$ $\displaystyle-\sum_{\square}\,\\{\,{\bf S}_{i}\cdot({\bf
S}_{i+e_{x}}\times{\bf S}_{i+e_{y}})+{\bf S}_{i+e_{x}}\cdot({\bf
S}_{i+e_{x}+e_{y}}\times{\bf S}_{i})$ $\displaystyle+\,{\bf
S}_{i+e_{x}+e_{y}}\cdot({\bf S}_{i+e_{y}}\times{\bf S}_{i+e_{x}})+{\bf
S}_{i+e_{y}}\cdot({\bf S}_{i}\times{\bf S}_{i+e_{x}+e_{y}})\,\\}$
$\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{\square}i(P_{ijkl}-P_{ijkl}^{-1}),$

where the sum over all sites has been rewritten as a sum over all (oriented)
plaquettes $(i,j,k,l)\equiv(i,i+e_{y},i+e_{x}+e_{y},i+e_{x})$ of the square
lattice. Substituting Eq. 24 into Eq. 19 gives the expression of the Floquet
Hamiltonian (7) of the main text.

###  B.2 Finite constant $H_{0}$

When a small constant term $\kappa H_{0}$, $\kappa\ll 1$, is present in the
Hamiltonian a small correction to the Floquet Hamiltonian appears [39], in
addition to the constant term itself,

$\displaystyle\delta H_{F}$ $\displaystyle=$
$\displaystyle\frac{1}{\omega}(-[H_{1},H_{0}]+[H_{-1},H_{0}])\,,$ (25)
$\displaystyle=$ $\displaystyle\kappa\frac{JJ_{1}}{\omega}i[h_{y},h_{0}]\,,$

where $h_{0}=H_{0}/J_{1}$ is the dimensionless uniform (frustrated) Heisenberg
model. By expanding $h_{y}$ and $h_{0}$ one obtains a sum of chiral terms on 2
types of triangles, $(i,i+e_{x},i+e_{y})$ and $(i,i+e_{x},i+2e_{x}+e_{y})$,
with alternating chiralities. Since $\kappa$ is of order $J/\omega$, the
amplitude of $\delta H_{F}$ is of order $\frac{J_{1}J^{2}}{\omega^{2}}$ i.e.
of order $\omega^{-2}$ and, hence, has been neglected in (10).

Note however that, in the case of a ramp, $\delta H_{F}$ has to be multiplied
by $\lambda(t)$, in contrast to $H_{F}^{0}$ multiplied by $\lambda^{2}(t)$ in
(10). Therefore at very short time, $\frac{t}{t_{\rm
ramp}}<\frac{J_{1}}{\omega}$, $\delta H_{F}$ may a priori play some role.
However, we have found numerically, that the Floquet dynamics obtained with
(10) was providing very accurate results compared to those obtained with the
full time-dependent Hamiltonian.

##  Appendix C Derivation of the kick operator

In the stroboscopic picture, the time-evolved state can be written as

$|\Psi(t)\big{>}=\exp{(-iK(t))}\exp{(-iH_{F}t)}|\Psi_{0}\big{>},$

where $K(t)$ is the Hermitian periodic kick operator (of period $T_{\rm
drive}$). The unitary operator $\exp{(-iK(t))}$ corresponds physically to a
change of basis. At all stroboscopic times $t=nT_{\rm drive}$,
$n\in{\mathbb{N}}$, we have $K(t)=0$ and $\exp{(-iK(t))}={\mathbb{I}}_{d}$.
Using Eq. (44) of Ref. [39], we get the leading high-frequency term:

$\displaystyle K(t)$ $\displaystyle=$
$\displaystyle\frac{1}{i\omega}\,[H_{1}(e^{i\omega t}-1)-H_{-1}(e^{-i\omega
t}-1)]\,,$ (26) $\displaystyle=$
$\displaystyle\frac{J}{\omega}\,[(1-\cos{(\omega t)})h_{y}+\sin{(\omega
t)}h_{x}]\,.$

Because of the small prefactor in the expression of $K(t)$, the unitary
remains close to the identity. Most observables are therefore weakly affected
during the micromotion w.r.t. to their behaviors under the effective Floquet
dynamics. This has been confirmed numerically.

_Acknowledgments—_ D.P. acknowledges inspiring discussions with Mari-Carmen
Banũls, Antoine Browaeys, Thierry Giamarchi, Nathan Goldman, David Guery-
Odelin, Vedika Khemani, Michael Kolodrubetz, Thierry Lahaye, Joel Moore, Sen
Niu and Luca Tagliacozzo and support from the TNTOP ANR-18-CE30-0026-01 grant
awarded by the French Research Council. This work was granted access to the
HPC resources of CALMIP center under allocations 2023-P1231 and 2024-P1231.

## References

  * [1] L. Savary and L. Balents,  _Quantum spin liquids: a review_ ,  Reports on Progress in Physics 80(1), 016502 (2016),  10.1088/0034-4885/80/1/016502.
  * [2] V. Kalmeyer and R. B. Laughlin,  _Equivalence of the resonating-valence-bond and fractional quantum Hall states_ ,  Phys. Rev. Lett. 59, 2095 (1987),  10.1103/PhysRevLett.59.2095.
  * [3] B. Bauer, L. Cincio, B. Keller, M. Dolfi, G. Vidal, S. Trebst and A. W. W. Ludwig,  _Chiral spin liquid and emergent anyons in a Kagome lattice Mott insulator_ ,  Nature Communications 5, 5137 (2014),  10.1038/ncomms6137.
  * [4] S.-S. Gong, W. Zhu, J.-X. Zhu, D. N. Sheng and K. Yang,  _Global phase diagram and quantum spin liquids in a spin- $\frac{1}{2}$ triangular antiferromagnet_,  Phys. Rev. B 96, 075116 (2017),  10.1103/PhysRevB.96.075116.
  * [5] A. Wietek and A. M. Läuchli,  _Chiral spin liquid and quantum criticality in extended $s=\frac{1}{2}$ Heisenberg models on the triangular lattice_,  Phys. Rev. B 95, 035141 (2017),  10.1103/PhysRevB.95.035141.
  * [6] J.-Y. Chen, L. Vanderstraeten, S. Capponi and D. Poilblanc,  _Non-abelian chiral spin liquid in a quantum antiferromagnet revealed by an iPEPS study_ ,  Phys. Rev. B 98, 184409 (2018),  10.1103/PhysRevB.98.184409.
  * [7] J.-Y. Chen, S. Capponi, A. Wietek, M. Mambrini, N. Schuch and D. Poilblanc,  _$\mathrm{SU}(3{)}_{1}$ chiral spin liquid on the square lattice: A view from symmetric projected entangled pair states_,  Phys. Rev. Lett. 125, 017201 (2020),  10.1103/PhysRevLett.125.017201.
  * [8] J.-Y. Chen, J.-W. Li, P. Nataf, S. Capponi, M. Mambrini, K. Totsuka, H.-H. Tu, A. Weichselbaum, J. von Delft and D. Poilblanc,  _Abelian ${\mathrm{su}(n)}_{1}$ chiral spin liquids on the square lattice_,  Phys. Rev. B 104, 235104 (2021),  10.1103/PhysRevB.104.235104.
  * [9] S.-S. Gong, W. Zhu and D. N. Sheng,  _Emergent chiral spin liquid: Fractional quantum Hall effect in a Kagome Heisenberg model_ ,  Scientific Reports 4, 6317 (2014),  10.1038/srep06317.
  * [10] Y.-C. He, D. N. Sheng and Y. Chen,  _Chiral spin liquid in a frustrated anisotropic Kagome Heisenberg model_ ,  Phys. Rev. Lett. 112, 137202 (2014),  10.1103/PhysRevLett.112.137202.
  * [11] A. Wietek, A. Sterdyniak and A. M. Läuchli,  _Nature of chiral spin liquids on the Kagome lattice_ ,  Phys. Rev. B 92, 125122 (2015),  10.1103/PhysRevB.92.125122.
  * [12] R. Feynman,  _Simulating physics with computers_ ,  Int. J. Theor. Phys. 21, 467–488 (1982),  10.1007/BF02650179.
  * [13] C. S. Chiu, G. Ji, A. Bohrdt, M. Xu, M. Knap, E. Demler, F. Grusdt, M. Greiner and D. Greif,  _String patterns in the doped Hubbard model_ ,  Science 365(6450), 251 (2019),  10.1126/science.aav3587.
  * [14] M. A. Nichols, L. W. Cheuk, M. Okan, T. R. Hartke, E. Mendez, T. Senthil, E. Khatami, H. Zhang and M. W. Zwierlein,  _Spin transport in a Mott insulator of ultracold fermions_ ,  Science 363(6425), 383 (2019),  10.1126/science.aat4387.
  * [15] S. Smale, P. He, B. A. Olsen, K. G. Jackson, H. Sharum, S. Trotzky, J. Marino, A. M. Rey and J. H. Thywissen,  _Observation of a transition between dynamical phases in a quantum degenerate Fermi gas_ ,  Science Advances 5(8) (2019),  10.1126/sciadv.aax1568.
  * [16] G. Jotzu, M. Messer, R. Desbuquois, M. Lebrat, T. Uehlinger, D. Greif and T. Esslinger,  _Experimental realization of the topological Haldane model with ultracold fermions_ ,  Nature 515(7526), 237 (2014),  10.1038/nature13915.
  * [17] P. Sompet, S. Hirthe, D. Bourgund, T. Chalopin, J. Bibo, J. Koepsell, P. Bojović, R. Verresen, F. Pollmann, G. Salomon, C. Gross, T. A. Hilker _et al._ ,  _Realizing the symmetry-protected haldane phase in Fermi–Hubbard ladders_ ,  Nature 606(7914), 484 (2022),  10.1038/s41586-022-04688-z.
  * [18] N. Goldman, J. Budich and P. Zoller,  _Topological quantum matter with ultracold gases in optical lattices_ ,  Nature Physics 12(1), 639 (2016),  https://doi.org/10.1038/nphys3803.
  * [19] S. Ebadi, T. Wang, H. Levine and et al.,  _Quantum phases of matter on a 256-atom programmable quantum simulator_ ,  Nature 595, 227 (2021),  10.1038/s41586-021-03582-4.
  * [20] G. Semeghini, H. Levine, A. Keesling, S. Ebadi, T. T. Wang, D. Bluvstein, R. Verresen, H. Pichler, M. Kalinowski, R. Samajdar, A. Omran, S. Sachdev _et al._ ,  _Probing topological spin liquids on a programmable quantum simulator_ ,  Science 374(6572), 1242 (2021),  10.1126/science.abi8794,  https://www.science.org/doi/pdf/10.1126/science.abi8794.
  * [21] G. Giudici, M. D. Lukin and H. Pichler,  _Dynamical preparation of quantum spin liquids in Rydberg atom arrays_ ,  Phys. Rev. Lett. 129, 090401 (2022),  10.1103/PhysRevLett.129.090401.
  * [22] B.-Y. Sun, N. Goldman, M. Aidelsburger and M. Bukov,  _Engineering and probing non-abelian chiral spin liquids using periodically driven ultracold atoms_ ,  PRX Quantum 4, 020329 (2023),  10.1103/PRXQuantum.4.020329.
  * [23] M. Kalinowski, N. Maskara and M. D. Lukin,  _Non-abelian Floquet spin liquids in a digital Rydberg simulator_ ,  Phys. Rev. X 13, 031008 (2023),  10.1103/PhysRevX.13.031008.
  * [24] G. Chen, K. R. A. Hazzard, A. M. Rey and M. Hermele,  _Synthetic-gauge-field stabilization of the chiral-spin-liquid phase_ ,  Phys. Rev. A 93, 061601 (2016),  10.1103/PhysRevA.93.061601.
  * [25] A. Mitra,  _Quantum quench dynamics_ ,  Annual Review of Condensed Matter Physics 9(1), 245 (2018),  10.1146/annurev-conmatphys-031016-025451,  https://doi.org/10.1146/annurev-conmatphys-031016-025451.
  * [26] H. Bragança, M. F. Cavalcante, R. G. Pereira and M. C. O. Aguiar,  _Quench dynamics and relaxation of a spin coupled to interacting leads_ ,  Phys. Rev. B 103, 125152 (2021),  10.1103/PhysRevB.103.125152.
  * [27] V. Alba and P. Calabrese,  _Entanglement and thermodynamics after a quantum quench in integrable systems_ ,  Proceedings of the National Academy of Sciences 114(30), 7947 (2017),  10.1073/pnas.1703516114,  https://www.pnas.org/doi/pdf/10.1073/pnas.1703516114.
  * [28] S. Czischek, M. Gärttner and T. Gasenzer,  _Quenches near ising quantum criticality as a challenge for artificial neural networks_ ,  Phys. Rev. B 98, 024311 (2018),  10.1103/PhysRevB.98.024311.
  * [29] S. Trotzky, Y.-A. Chen, A. Flesch, I. P. McCulloch, U. Schollwöck, J. Eisert and I. Bloch,  _Probing the relaxation towards equilibrium in an isolated strongly correlated one-dimensional Bose gas_ ,  Nature Physics 8(4), 325 (2012),  10.1038/nphys2232.
  * [30] P. Titum, J. T. Iosue, J. R. Garrison, A. V. Gorshkov and Z.-X. Gong,  _Probing ground-state phase transitions through quench dynamics_ ,  Phys. Rev. Lett. 123, 115701 (2019),  10.1103/PhysRevLett.123.115701.
  * [31] G. Floquet,  _Sur les équations différentielles linéaires à coefficients périodiques_ ,  Annales scientifiques de l’École Normale Supérieure (Série 2) 12, 47 (1883),  10.24033/asens.220.
  * [32] C. Weitenberg and J. Simonet,  _Tailoring quantum gases by Floquet engineering_ ,  Nat. Phys. 17, 1342–1348 (2021).
  * [33] D. J. Luitz, R. Moessner, S. L. Sondhi and V. Khemani,  _Prethermalization without temperature_ ,  Phys. Rev. X 10, 021046 (2020),  10.1103/PhysRevX.10.021046.
  * [34] N. Goldman and J. Dalibard,  _Periodically driven quantum systems: Effective hamiltonians and engineered gauge fields_ ,  Phys. Rev. X 4, 031027 (2014),  10.1103/PhysRevX.4.031027.
  * [35] P. W. Anderson,  _Resonating valence bonds: A new kind of insulator?_ ,  Materials Research Bulletin 8(2), 153 (1973),  http://dx.doi.org/10.1016/0025-5408(73)90167-0.
  * [36] W. Magnus,  _On the exponential solution of differential equations for a linear operator_ ,  Commun. Pure Appl. Math. 7, 649 (1954),  10.1002/cpa.3160070404.
  * [37] F. Casas, J. A. Oteo and J. Ros,  _Floquet theory: exponential perturbative treatment_ ,  Journal of Physics A: Mathematical and General 34(16), 3379 (2001),  10.1088/0305-4470/34/16/305.
  * [38] S. Rahav, I. Gilary and S. Fishman,  _Effective hamiltonians for periodically driven systems_ ,  Phys. Rev. A 68, 013820 (2003),  10.1103/PhysRevA.68.013820.
  * [39] 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(2), 139 (2015),  10.1080/00018732.2015.1055918,  https://doi.org/10.1080/00018732.2015.1055918.
  * [40] J. Jordan, R. Orús, G. Vidal, F. Verstraete and J. Cirac,  _Classical simulation of infinite-size quantum lattice systems in two spatial dimensions_ ,  Phys. Rev. Letters 101, 250602 (2008),  10.1103/PhysRevLett.101.250602.
  * [41] D. Poilblanc, J. I. Cirac and N. Schuch,  _Chiral topological spin liquids with projected entangled pair states_ ,  Phys. Rev. B 91, 224431 (2015),  10.1103/PhysRevB.91.224431.
  * [42] D. Poilblanc, N. Schuch and I. Affleck,  _$\mathrm{SU}(2{)}_{1}$ chiral edge modes of a critical spin liquid_,  Phys. Rev. B 93, 174414 (2016),  10.1103/PhysRevB.93.174414.
  * [43] M. Suzuki,  _Fractal decomposition of exponential operators with applications to many-body theories and Monte Carlo simulations_ ,  Physics Letters A 146(6), 319 (1990),  https://doi.org/10.1016/0375-9601(90)90962-N.
  * [44] H. F. Trotter,  _On the product of semi-groups of operators_ ,  Proceedings of the American Mathematical Society 10(4), 545 (1959),  10.1090/s0002-9939-1959-0108732-6.
  * [45] D. Poilblanc, N. Schuch, D. Pérez-García and J. I. Cirac,  _Topological and entanglement properties of resonating valence bond wave functions_ ,  Phys. Rev. B 86, 014404 (2012),  10.1103/PhysRevB.86.014404.
  * [46] N. Schuch, D. Poilblanc, J. I. Cirac and D. Pérez-García,  _Resonating valence bond states in the PEPS formalism_ ,  Physical Review B 86(11), 115108 (2012),  10.1103/PhysRevB.86.115108.
  * [47] J.-Y. Chen and D. Poilblanc,  _Topological $\mathbb{Z}_{2}$ resonating-valence-bond spin liquid on the square lattice_,  Phys. Rev. B 97, 161107 (2018),  10.1103/PhysRevB.97.161107.
  * [48] H. Dreyer, L. Vanderstraeten, J.-Y. Chen, R. Verresen and N. Schuch,  _Robustness of critical U(1) spin liquids and emergent symmetries in tensor networks_ ,  arXiv preprint arXiv:2008.04833 (2020),  10.48550/ARXIV.2008.04833.
  * [49] A. Kitaev,  _Fault-tolerant quantum computation by anyons_ ,  Annals of Physics 303(1), 2 (2003),  http://dx.doi.org/10.1016/S0003-4916(02)00018-0.
  * [50] A. E. B. Nielsen, G. Sierra and J. I. Cirac,  _Local models of fractional quantum Hall states in lattices and physical implementation_ ,  Nature Communications 4, 2864 (2013),  10.1038/ncomms3864.
  * [51] D. Poilblanc,  _Investigation of the chiral antiferromagnetic Heisenberg model using projected entangled pair states_ ,  Phys. Rev. B 96, 121118 (2017),  10.1103/PhysRevB.96.121118.
  * [52] J. Yang, Z. Liu and L. Wang,  _Ground state phase diagram and the exotic phases in the spin-1/2 square lattice $J_{1}$-$J_{2}$-$J_{\chi}$ model_ (2024), arXiv:2401.03434.
  * [53] X.-T. Zhang, Y. Huang, H.-Q. Wu, D. N. Sheng and S.-S. Gong,  _Chiral spin liquid and quantum phase diagram of spin- $1/2$ $J_{1}$-$J_{2}$-$J_{\chi}$ model on the square lattice_ (2024), arXiv:2401.07461.
  * [54] J. Hasik, M. Van Damme, D. Poilblanc and L. Vanderstraeten,  _Simulating chiral spin liquids with projected entangled-pair states_ ,  Phys. Rev. Lett. 129, 177201 (2022),  10.1103/PhysRevLett.129.177201.
  * [55] J. I. Cirac, D. Poilblanc, N. Schuch and F. Verstraete,  _Entanglement spectrum and boundary theories with projected entangled-pair states_ ,  Phys. Rev. B 83(24), 245134 (2011),  10.1103/PhysRevB.83.245134.
  * [56] M. S. Rudner, N. H. Lindner, E. Berg and M. Levin,  _Anomalous edge states and the bulk-edge correspondence for periodically driven two-dimensional systems_ ,  Phys. Rev. X 3, 031005 (2013),  10.1103/PhysRevX.3.031005.
  * [57] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Serbyn and Papić,  _Weak ergodicity breaking from quantum many-body scars_ ,  Nat. Phys. 14, 745 (2018).
  * [58] S. Moudgalya, S. Rachel, B. A. Bernevig and N. Regnault,  _Exact excited states of nonintegrable models_ ,  Phys. Rev. B 98, 235155 (2018),  10.1103/PhysRevB.98.235155.
  * [59] C. Gross and I. Bloch,  _Quantum simulations with ultracold atoms in optical lattices_ ,  Science 357(6355), 995 (2017),  10.1126/science.aal3837.
  * [60] W. Hofstetter, J. I. Cirac, P. Zoller, E. Demler and M. D. Lukin,  _High-temperature superfluidity of fermionic atoms in optical lattices_ ,  Phys. Rev. Lett. 89, 220407 (2002),  10.1103/PhysRevLett.89.220407.
  * [61] D. Jaksch and P. Zoller,  _The cold atom hubbard toolbox_ ,  Annals of Physics 315(1), 52 (2005),  https://doi.org/10.1016/j.aop.2004.09.010,  Special Issue.
  * [62] R. Jördens, N. Strohmaier, K. Günter, H. Moritz and T. A. Esslinger,  _Mott insulator of fermionic atoms in an optical lattice_ ,  Nature 455, 204 (2008),  doi: 10.1038/nature07244.
  * [63] K. Yang and Y.-Q. Li,  _Rigorous proof of pseudospin ferromagnetism in two-component bosonic systems with component-independent interactions_ ,  International Journal of Modern Physics B 17(07), 1027 (2003),  10.1142/S0217979203018041,  https://doi.org/10.1142/S0217979203018041.
  * [64] E. Altman, W. Hofstetter, E. Demler and M. D. Lukin,  _Phase diagram of two-component bosons on an optical lattice_ ,  New J. Phys. 5, 113 (2003),  10.1088/1367-2630/5/1/113.
  * [65] A. Browaeys and T. Lahaye,  _Many-body physics with individually-controlled Rydberg atoms_ ,  Nature Physics 16, 132 (2020),  https://doi.org/10.1038/s41567-019-0733-z.
  * [66] T. C. Rutkowski and M. J. Lawler,  _Spin liquid phases of large-spin Mott insulating ultracold bosons_ ,  Phys. Rev. B 93, 094405 (2016),  10.1103/PhysRevB.93.094405.
  * [67] H. Sun, B. Yang, H. Wang and et al.,  _Realization of a bosonic antiferromagnet_ ,  Nat. Phys. 17, 990–994 (2021).
  * [68] A. Mazurenko, C. Chiu, G. Ji, M. Parsons, M. Kanász-Nagy, R. Schmidt, F. Grusdt, E. Demler, D. Greif and M. Greiner,  _A cold-atom Fermi-Hubbard antiferromagnet_ ,  Nature 545, 462 (2017),  10.1038/nature22362.
  * [69] C. S. Chiu, G. Ji, A. Mazurenko, D. Greif and M. Greiner,  _Quantum state engineering of a hubbard system with ultracold fermions_ ,  Phys. Rev. Lett. 120, 243201 (2018),  10.1103/PhysRevLett.120.243201.
  * [70] J. Koepsell, S. Hirthe, D. Bourgund, P. Sompet, J. Vijayan, G. Salomon, C. Gross and I. Bloch,  _Robust bilayer charge pumping for spin- and density-resolved quantum gas microscopy_ ,  Phys. Rev. Lett. 125, 010403 (2020),  10.1103/PhysRevLett.125.010403.
  * [71] Z. Z. Yan, B. M. Spar, M. L. Prichard, S. Chi, H.-T. Wei, E. Ibarra-García-Padilla, K. R. A. Hazzard and W. S. Bakr,  _Two-dimensional programmable tweezer arrays of fermions_ ,  Phys. Rev. Lett. 129, 123201 (2022),  10.1103/PhysRevLett.129.123201.
  * [72] P. Windpassinger and K. Sengstock,  _Engineering novel optical lattices_ ,  Reports on progress in physics. Physical Society (Great Britain) 76, 086401 (2013),  10.1088/0034-4885/76/8/086401.
  * [73] D. Pérez-García, F. Verstraete, J. I. Cirac and M. M. Wolf,  _PEPS as unique ground states of local hamiltonians_ ,  Quant. Inf. Comp. 8, 0650 (2008),  0707.2260.
  * [74] M. Mambrini, R. Orús and D. Poilblanc,  _Systematic construction of spin liquids on the square lattice from tensor networks with SU(2) symmetry_,  Phys. Rev. B 94, 205124 (2016),  10.1103/PhysRevB.94.205124.
  * [75] A. E. B. Nielsen, J. I. Cirac and G. Sierra,  _Laughlin spin-liquid states on lattices obtained from conformal field theory_ ,  Phys. Rev. Lett. 108, 257206 (2012),  10.1103/PhysRevLett.108.257206.