# Exploring quantum criticality and ergodicity-breaking dynamics in spin-1
Kitaev chains via single-ion anisotropies

Wen-Yi Zhang  College of Physics, Nanjing University of Aeronautics and
Astronautics, Nanjing, 211106, China  Key Laboratory of Aerospace Information
Materials and Physics (NUAA), MIIT, Nanjing 211106, China  Qing-Min Hu
College of Physics, Nanjing University of Aeronautics and Astronautics,
Nanjing, 211106, China  Key Laboratory of Aerospace Information Materials and
Physics (NUAA), MIIT, Nanjing 211106, China  Jie Ren  Department of Physics,
Changshu Institute of Technology, Changshu 215500, China  Liangsheng Li
National Key Laboratory of Scattering and Radiation, Beijing 100854, China
Wen-Long You<EMAIL_ADDRESS>College of Physics, Nanjing University of
Aeronautics and Astronautics, Nanjing, 211106, China  Key Laboratory of
Aerospace Information Materials and Physics (NUAA), MIIT, Nanjing 211106,
China

###### Abstract

We investigate the topological gauge-theory terms and quantum criticality in
spin-1 Kitaev chain with generic single-ion anisotropies (SIAs). The ground-
state phase diagram, including Kitaev spin liquid (KSL) phase and gapless
dimer phase, is determined by the infinite time evolving block decimation
(iTEBD) method. By varying the strength of uniaxial SIA, there is quantum
phase transition between the KSL phase and the dimer phase, which is
equivalent to the confinement-deconfinement transition in the lattice
Schwinger model with a topological $\theta$-angle of $\pi$. Here, we
demonstrate an added rhombic SIA to the $\mathbb{Z}_{2}$ symmetric model
shifts topological angle $\theta$ away from $\pi$, leading to the emergence of
$y$-ferroquadrupole and $x$-ferroquadrupole phases for negative and positive
values of rhombic SIAs, respectively. By adjusting these rhombic SIAs, the
system from the $y$-ferroquadrupole to the $x$-ferroquadrupole phases
undergoes either a crossover through the KSL phase or a genuine phase
transition through the deconfined line, thus providing an example for
unnecessary phase transitions. We find the spin-1 Hamiltonian can be exactly
mapped to the spin-1/2 effective extended PXP Hamiltonian coupled to all
positive $\mathbb{Z}_{2}$ gauge charges, where the uniaxial SIA is equivalent
to the uniform detuning while the rhombic SIA plays a role of staggered
detuning. We then examine the hierarchical fragmentation of the Hilbert space
and its associated dynamics in specific limiting cases. The quantum many-body
scar (QMBS) arises from the Fibonacci constraint under a weak field of SIAs.
When the staggered detuning is sufficiently large, the resultant slow dynamics
can be nominally understood using a second-order effective Hamiltonian derived
by the Schrieffer-Wolff transformation. Our work introduces a theoretical
framework for simulating topological $\theta$-angles and ergodicity-breaking
dynamics and explores higher-spin generalizations of the scarred model,
establishing a compelling link between solid-state-based quantum systems and
lattice gauge theories that can be implemented using state-of-the-art cold-
atom quantum simulators.

##  I Introduction

Recent advancements in synthetic quantum matter have propelled the exploration
of non-equilibrium dynamics in isolated quantum systems. The eigenstate
thermalization hypothesis (ETH) that sets the foundation of statistical
physics, still emerges as a central framework in quantum physics [1, 2, 3, 4].
Rydberg atoms, with their high principal quantum numbers, exhibit strong
repulsive interactions when close together, preventing two such atoms from
occupying nearby space [5]. The quantum simulator utilizing Rydberg atoms is
capable of emulating quantum many-body dynamics and probing unprecedented
ground-state quantum phases [6, 7, 8]. A seminal experiment involving a
Rydberg atom chain, initiated from a density wave state, showed persistent
quantum revivals [9]. This observation elucidated that embedded within the
predominantly thermal spectrum, a minority of eigenstates significantly
contribute the atypical dynamics. These zero-measure non-thermal eigenstates
delineate as the quantum many-body scar (QMBS). Recognized as a paradigm of
weak ergodicity breaking in disorder-free nonintegrable models, QMBSs have
garnered widespread attention and have now become an entire area of active
research. A kinetically-constrained PXP model was derived to effectively
elucidate the mechanisms responsible for the observed phenomena in Rydberg
atom experiments [10]. Subsequently, the oscillatory dynamics for certain
initial states was also observed across various experimental platforms [11,
12]. Meanwhile, multiple approaches have been proposed for constructing
Hamiltonians with small invariant subspaces that weakly violate the
ergodicity. Analytic expressions for QMBS states have also been identified in
the spin-1 Affleck-Kennedy-Lieb-Tasaki (AKLT) model [13]. Other systems that
display QMBSs include the transverse field Ising ladder [14], the two-
dimensional PXP model [15, 16], the spin-1 XY model [17, 18, 19], and the non-
Hermitian many-body systems [20]. The approaches for constructing Hamiltonians
are predominantly divided into three categories: the projector embedding
method [21, 22, 23], the spectrum generating algebra technique [24, 25, 26],
and the Krylov restricted thermalization strategy [27, 28, 29].

While researchers have explored QMBSs across diverse models, the seminal PXP
model, where these states were first identified, has never faded from view.
Notably, an explicit formulation of certain QMBS states within the PXP model
has been elucidated [30]. By introducing perturbations to the PXP model [31],
one can either engineer perfect revivals [32] or achieve integrability [33].
Recent studies have demonstrated that the PXP model is intricately embedded in
the spin-1 Kitaev model [34, 35], in which the uniaxial single-ion anisotropy
(SIA) plays an equivalent role of the static uniform detuning [36].
Specifically, in the context of achieving targeted Rydberg excitations, the
detuning can be meticulously adjusted at any designated site within an atom
array using external Raman lasers [37, 38]. Most prior research has assumed
uniform detuning in Rydberg atoms [39, 40, 41, 42, 43, 44], but latest studies
have increasingly focused on the effects of staggered detuning [45]. Recent
advancements in experiments with Rydberg atom arrays have demonstrated that
constrained Hamiltonian dynamics can be connected to U(1) lattice gauge theory
(LGT) [46, 47], creating a mathematical equivalence between the PXP model and
the lattice Schwinger model, which incorporates a topological angle
$\theta\\!=\\!\pi$ [48, 49], and the deviation $\theta\\!-\\!\pi$ of the
lattice Schwinger model is equivalent to an additional staggered field applied
to the Rydberg system. Notably, Pan et al. have successfully implemented a
tunable $\theta$-angle in U(1) LGT with an atom-number-resolved optical
lattice quantum simulator, allowing for the direct observation of the
confinement-deconfinement transition [50]. In the paper, we propose a method
for realizing topological terms in spin-1 Kitaev chain with generic SIAs. We
illustrate that the staggered detuning present in the extended PXP model can
be effectively mimicked by the incorporation of rhombic SIAs in spin-1 Kitaev
model. Notably, both the staggered terms in the spin-1/2 PXP model and the
rhombic SIA terms in spin-1 Kitaev model can be equivalently mapped to the
topological $\theta$-angle term in gauge theory. Through the demonstrations,
the intriguing correspondence between these different models is established.

The structure of this paper is as follows: In Sec. II, we introduce a one-
dimensional spin-1 Kitaev model with tunable three-component SIAs and employ
both the exact diagonalization (ED) method and the infinite time-evolving
block decimation (iTEBD) algorithm [51] to determine the ground-state phase
diagram. In Sec. III, we scrutinize the quantum criticality in the spin-1
Kitaev model by varying both uniaxial and rhombic SIAs. Furthermore, we
utilize a Schrieffer-Wolff transformation in Sec. IV to develop a perturbative
description of the spin-1 Kitaev chain in the large limit of rhombic SIA. In
Sec. V, we demonstrate the spin-1 Hamiltonians can be exactly mapped to the
spin-1/2 extended PXP Hamiltonians in the ground-state manifold characterize
by all positive $\mathbb{Z}_{2}$ gauge charges. In addition to the gauge
constraint, the Hilbert space exhibits a further fragmentation in the strong-
confinement limit. We compare the quench dynamics for various values of
rhombic SIAs. The glassy dynamics can be well understood by the effective
Schrieffer-Wolff Hamiltonian from the perspective of Hilbert space
fragmentation. Finally, a summary is provided in Sec. VI.

##  II Spin-1 Kitave model with three-component single-ion anisotropies

Figure 1:  Phase diagram of Hamiltonian (1) with
$A_{x}=A_{r}\sin\theta\cos\phi$, $A_{y}=A_{r}\sin\theta\sin\phi$, and
$A_{z}=A_{r}\cos\theta$, calculated by the iTEBD method with bond dimension
$\chi=120$. (a) Representation of the phase diagram in three-dimensional
spherical coordinates. The orange (blue) line denotes the cross-section for
$\phi=\pi/4$ ($\phi=3\pi/4$), which will be detailed in (b) [(c)]. (b) Two-
dimensional cross-section for $\phi=\pi/4$, i.e., $A_{x}=A_{y}$ and
$A_{r}\in[0,5]$. The grey semi-infinite shaded region indicates the gapless
dimer phase, with the grey dashed line corresponding to the phase boundaries
as $A_{r}\to\infty$, where $\theta_{0}=0.3041\pi$ denotes the angle with the
$z$-axis. (c) Two-dimensional cross-section for $\phi=3\pi/4$, i.e.,
$A_{x}=-A_{y}$, and $A_{r}\in[0,1.5]$. In this case, Eq. (1) reduces to
Hamiltonian (5). The dimer phase degrades into a ray line, represented by the
purple line. Blue arrows indicate where the strength of the rhombic SIA
$A_{x}\to\infty$ and the uniaxial SIA $A_{z}=0$. Dashed lines represent
$A_{r}=1$, and light blue crosses indicate that the system undergoes a
crossover rather than a phase transition as the lowest excitation gap never
closes. See more details in the main text.

We start with the one-dimensional spin-1 Kitaev models with tunable three-
component SIAs. The Hamiltonian is given by

$\displaystyle\hat{H}_{\rm
KA}\\!=\\!\sum_{j=1}^{N/2}K_{j}\\!(S_{2j-1}^{x}S_{2j}^{x}\\!+\\!S_{2j}^{y}S_{2j+1}^{y})\\!+\\!\sum_{\alpha,j=1}^{N}\\!A_{\alpha}\\!(S_{j}^{\alpha})^{2}.$
(1)

Here, $A_{\alpha}~{}(\alpha=x,y,z)$ denotes the amplitude of the
$\alpha$-component SIA, and $K_{j}$ parameterizes the strength of the bond-
dependent Kitaev exchange coupling. $S_{j}^{\alpha}$ ($\alpha=x,y,z$)
represents the $\alpha$-component of the spin-1 operator at the $j$-th site
among total $N$ sites, obeying the SU(2) algebra as
$[S_{i}^{\alpha},S_{j}^{\beta}]=i\delta_{ij}\epsilon_{\alpha\beta\gamma}S_{j}^{\gamma}$,
where $\epsilon_{\alpha\beta\gamma}$ is the antisymmetric tensor. Hereafter,
we assume the uniform coupling with $K_{j}=1$ unless stated otherwise. We
adopt a specific spin-1 representation defined as follows:

$\displaystyle|x\rangle\\!$ $\displaystyle=$
$\displaystyle\\!\frac{1}{\sqrt{2}}(|-1\rangle\\!-\\!|1\rangle),|y\rangle\\!=\\!\frac{i}{\sqrt{2}}(|-1\rangle\\!+\\!|1\rangle),|z\rangle\\!=\\!|0\rangle,\quad$
(2)

where $|m\rangle$ represents the eigenstate of the spin operator $S^{z}$ with
eigenvalues $m=-1$, $0$, $1$. In this representation, we have
$S_{\beta\gamma}^{\alpha}=-i\epsilon_{\alpha\beta\gamma}$. Notably, the spin
operators $S^{\alpha}$ ($\alpha=x,y,z$) exhibit interdependence due to the
constraint $(\bm{S}_{j})^{2}=S(S+1)=2$. To parameterize the amplitudes
${A_{\alpha}}$ for the $\alpha$ component of SIA, we employ spherical
coordinates:

$\displaystyle
A_{x}=A_{r}\sin\theta\cos\phi,A_{y}=A_{r}\sin\theta\sin\phi,A_{z}=A_{r}\cos\theta.\quad$
(3)

For convenience, we here define the parameterization of the Hamiltonian in Eq.
(1) using the polar angle $\theta\in[0,2\pi)$, the azimuthal angle
$\phi\in[0,\pi)$, and the radius $A_{r}$. Considering the phase transition
from the gapless dimer phase to the Kitaev spin liquid (KSL) phase occurs at
$A_{zc}=-0.655$ for $\theta=\pi$ [36], a more general expression for the
critical points is established as follows:

$\displaystyle
A_{rc}\\!=\\!-\\!\frac{A_{zc}\sin\theta_{0}}{\sin(\theta-\theta_{0})}\textrm{
and }\theta_{0}\\!<\\!\theta\\!<\\!\pi\\!+\\!\theta_{0},\textrm{ for
}\phi=\pi/4,$ (4)

where $\theta_{0}\equiv\arctan(\sin^{-1}\phi)=0.3041\pi$. To effectively
delineate the phase structure, we present the comprehensive phase diagram of
$\hat{H}_{\rm KA}$ in Eq. (1) through the lens of three-dimensional
$A_{\alpha}$ coordinates, as depicted in Fig. 1(a). Here, the orange and blue
lines represent the phase planes for $\phi=\pi/4$ and $\phi=3\pi/4$,
respectively, with their detailed behaviors explored via cross-sectional
views. For $\phi=\pi/4$ in Fig. 1(b), the emergence of the dimer phase is
distinctly marked by the grey shaded region, where the grey dashed line that
indicates $\theta_{0}$ as the pivotal angle in relation to the $z$-axis
delineates the phase boundary in the limiting case, i.e., $A_{r}\to\infty$.
Figure 1(c) depicts the phase diagram of the Hamiltonian given in Eq. (5) for
$\phi=3\pi/4$. At $\theta=\pi$, a phase transition is triggered by increasing
$A_{r}$ from the KSL phase to the gapless dimer phase. Interestingly, the KSL-
dimer transition is equivalent to the confinement-deconfinement transition.
[54]. Moreover, the yellow dot highlights the 2D $\mathbb{Z}_{2}$ tricritical
point, known as Coleman’s quantum phase transition point [55]. Notably, the
yellow dot is also dubbed as unnecessary quantum critical point [52, 53].
Since the phases on either side of this critical point are identical, the
transition can be circumvented by following certain symmetric paths within the
entire parameter space. In fact, the phase diagram is completely consistent
with that of the $U(1)$ QLM with topological $\theta$-angle [56]. The dimer
phase is represented along the purple line, while the KSL phase encompasses
the region beyond this boundary. Notably, within the KSL phase, the area for
$\theta$ ranging from $0$ to $\pi$ is distinguished as the $x$-ferroquadrupole
phase. Conversely, for $\theta$ spanning $(\pi,2\pi)$, the phase identified as
the $y$-ferroquadrupole phase. Importantly, the light blue crosses mark the
crossover phenomena, rather than phase transitions. This interpretation is
supported by the lack of gap closure in the lowest excitation spectrum and the
continuity observed in the quadrupole correlation functions.

##  III Effect of tunable single-ion anisotropy interactions for $\phi=3\pi/4$

Next, we focus on the scenario where $\phi=3\pi/4$ (i.e., $A_{x}=-A_{y}$), the
model simplifies to a spin-1 Kitaev chain with both tunable rhombic and
uniaxial single-ion anisotropies. The Hamiltonian in such case is given by

$\displaystyle\hat{H}$ $\displaystyle=$
$\displaystyle\sum_{j=1}^{N/2}(S_{2j-1}^{x}S_{2j}^{x}+S_{2j}^{y}S_{2j+1}^{y})$
(5) $\displaystyle-$
$\displaystyle\frac{A_{r}\sin\theta}{\sqrt{2}}\sum_{j=1}^{N}\left[(S_{j}^{x})^{2}-(S_{j}^{y})^{2}\right]+A_{r}\cos\theta\sum_{j=1}^{N}(S_{j}^{z})^{2},\quad\;$

where the second term represents the strength of rhombic SIA and the third
term denotes the amplitude of uniaxial SIA.

In our investigation of quantum phases in one-dimensional spin-1 Kitaev models
endowed with tunable three-component SIAs, we initially examine the ground
state local magnetizations, $\langle S_{j}^{\alpha}\rangle=0$
$(\alpha=x,y,z)$, across the model’s parameter space. This absence of local
magnetization, indicating the conservation of time-reversal symmetry.
Subsequently, we focus on calculating the two-point spin-spin correlations
$C^{\alpha\beta}_{i,j}=\langle S^{\alpha}_{i}S^{\beta}_{j}\rangle$, where
$\alpha,\beta=x,y,z$. For $\phi=\pi/4$, the Hamiltonian, as defined in Eq.
(1), exhibits invariance under an operation that combines a $\pi/2$-rotation
about the $z$-axis and a single-site translation. This symmetry leads to
equalities in the finite-size system expressed as

$\displaystyle
C^{xx}_{1,2}=C^{yy}_{2,3},C^{yy}_{1,2}=C^{xx}_{2,3},C^{zz}_{1,2}=C^{zz}_{2,3}.$
(6)

However, our previous research has shown that the dimer phase will emerge upon
turning on three-component SIAs $A_{r}>A_{rc}$ with
$\phi=\pi/4,\theta_{0}<\theta<\pi+\theta_{0}$. The dimer phase is distinctly
marked by a pattern of alternating nearest-neighbor spin-spin correlations and
thus breaks the translational symmetry, which is characterized by the dimer
order parameter defined by
$O_{D}=|\langle\bm{S}_{1}\cdot\bm{S}_{2}\rangle-\langle\bm{S}_{2}\cdot\bm{S}_{3}\rangle|$.
For a more detailed analysis, it is insightful to examine the dimer order
parameter into its component-wise contributions along the $x$, $y$, and $z$
axes, such as

$\displaystyle O_{D}^{x}$ $\displaystyle=$
$\displaystyle|\langle{S}_{1}^{x}{S}_{2}^{x}\rangle-\langle{S}_{2}^{y}{S}_{3}^{y}\rangle|,$
$\displaystyle O_{D}^{y}$ $\displaystyle=$
$\displaystyle|\langle{S}_{1}^{y}{S}_{2}^{y}\rangle-\langle{S}_{2}^{x}{S}_{3}^{x}\rangle|,$
$\displaystyle O_{D}^{z}$ $\displaystyle=$
$\displaystyle|\langle{S}_{1}^{z}{S}_{2}^{z}\rangle-\langle{S}_{2}^{z}{S}_{3}^{z}\rangle|.$
(7)

Meanwhile, we consider the magnetic order in the KSL phase. The KSL is
characterized by the quadrupolar order in which the spin-rotational symmetry
is broken whereas both translational and time-reversal symmetries are
retained. A symmetric and traceless rank-2 quadrupole tensor operator can be
given by

$\displaystyle Q_{j}^{\alpha\beta}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\left(S_{j}^{\alpha}S_{j}^{\beta}+S_{j}^{\beta}S_{j}^{\alpha}\right)-\frac{1}{3}\bm{S}_{j}^{2}\delta_{\alpha\beta},$
(8)

with $\alpha,\beta\in\\{x,y,z\\}$ and the Kronecker delta
$\delta_{\alpha\beta}$ at site $j$. We find that all off-diagonal components
of quadrupole moments are zero for the whole parameter space for the ground
states, i.e., $\langle Q_{j}^{xy}\rangle$ =$\langle
Q_{j}^{yz}\rangle$=$\langle Q_{j}^{zx}\rangle$=0. We calculate the two
remaining types of quadrupole orders as follows:

$\displaystyle Q_{j}^{x^{2}-y^{2}}$ $\displaystyle=$ $\displaystyle
Q_{j}^{xx}-Q_{j}^{yy},\quad Q_{j}^{3z^{2}-r^{2}}=\sqrt{3}Q_{j}^{zz},$ (9)

and their respective quadrupole correlation functions:

$\displaystyle\mathcal{Q}_{i,j}^{\gamma}$ $\displaystyle=$
$\displaystyle\langle Q_{i}^{\gamma}Q_{j}^{\gamma}\rangle-\langle
Q_{i}^{\gamma}\rangle\langle
Q_{j}^{\gamma}\rangle,\gamma=x^{2}-y^{2},3z^{2}-r^{2}.\quad$ (10)

Nonzero quadrupole moments characterize ground states that break spin-
rotational symmetry through anisotropic spin fluctuations, while maintaining
both translational and time-reversal symmetries. This phenomenon mirrors the
magnetic counterpart of liquid crystals, indicating that our spin-1 systems
are in spin nematic quadrupole phases.

In this study, we employ the infinite time-evolving block decimation (iTEBD)
technique with the bond dimension set to $\chi=120$ [51], to compute the order
parameters for the Hamiltonian Eq. (5). We investigate the quadrupole order
parameters $\langle Q_{j}^{\gamma}\rangle$ and their correlation functions
$\mathcal{Q}_{j,j+1}^{\gamma}$, particularly at $A_{r}=1$—as delineated by the
dashed line in Fig. 1(c)—as a function of the angle $\theta$. These findings
are illustrated in Fig. 2(a). Notably, the analysis of quadrupole correlations
on both odd and even bonds reveals a pronounced discontinuity at $\theta=\pi$,
indicating the occurrence of first-order quantum phase transitions. Adjacent
to $\theta=\pi/2$, the anisotropy terms in Eq. (5) are predominantly
characterized by a finite rhombic SIA with a negative coefficient, while the
uniaxial SIA becomes negligible, which leads to $\langle
Q_{j}^{x^{2}-y^{2}}\rangle$ approaches unity, heralding the ground state is
characterized by $|S_{j}^{y}=0\rangle$ and dominant spin fluctuations are
within the $xz$ plane. This uniaxial spin nematic state is defined by
$y$-ferroquadrupole ordering and persists over the range
$\theta\in(\pi/2,\pi)$. For the interval $\theta\in(\pi,3\pi/2)$, the system’s
evolution can be understood by a dual mapping
$U^{\prime}=\prod_{j}\exp\left(i\frac{\pi}{2}S_{j}^{z}\right)$, leading to a
reorientation of the $x$ and $y$ axes. This transformation coincides with the
presence of a finite rhombic SIA bearing a positive sign, under which
conditions $\langle Q_{j}^{x^{2}-y^{2}}\rangle$ tends toward -1. Consequently,
the system’s ground state transitions to a configuration characterized by
$|S_{j}^{x}=0\rangle$, denoting $x$-ferroquadrupole ordering. Although the
quadrupole order parameter $\langle Q_{j}^{x^{2}-y^{2}}\rangle$ is continuous
across the $\theta=\pi$, the corresponding correlation functions exhibit a
jump. One can find that as $\theta$ varies from $\pi/2$ to $3\pi/2$,
$\mathcal{Q}_{1,2}^{3z^{2}-r^{2}}$ changes from a finite value to zero while
$\mathcal{Q}_{1,2}^{3z^{2}-r^{2}}$ shows an opposite trend.
$\mathcal{Q}_{i,j}^{x^{2}-y^{2}}$ also displays discontinuity at the critical
point, yet remains nonzero in both phases. This demonstrates that the two-
point $\mathcal{Q}_{i,j}^{3z^{2}-r^{2}}$ correlation serves as an effective
order parameter for identifying uniaxial spin nematic states. Figure. 2(b)
illustrates that when $\theta=0$, only the dominant term of the correlation
functions undergo a continuous change, indicating that the transition from
$x$-ferroquadrupole to $y$-ferroquadrupole is an intra-phase shift, signifying
a crossover rather than a phase transition.

##  IV Effective Hamiltonian for $\phi=3\pi/4$

Figure 2:  Order parameters of the Hamiltonian (5). (a) Quadrupole correlation
functions $\mathcal{Q}_{1,2}^{x^{2}-y^{2}}$,
$\mathcal{Q}_{2,3}^{x^{2}-y^{2}}$, $\mathcal{Q}_{1,2}^{3z^{2}-r^{2}}$ and
$\mathcal{Q}_{2,3}^{3z^{2}-r^{2}}$ at $A_{r}=1$ [see Fig. 1(c) dark red dashed
line] as a function of $\theta$. Inset shows quadrupole operators $\langle
Q_{1}^{x^{2}-y^{2}}\rangle$ and $\langle Q_{1}^{3z^{2}-r^{2}}\rangle$ at
$A_{r}=1$ [see Fig. 1(c) dark red dashed line] as a function of $\theta$. (b)
Quadrupole correlation functions at $A_{r}=1$ [see Fig. 1(c) green dashed
line] as a function of $\theta$. Inset shows quadrupole operators at $A_{r}=1$
[see Fig. 1(c) green dashed line] as a function of $\theta$. (c) The three
components of dimer order parameter $O_{D}$ as a function of $A_{r}$. From the
KSL phase to the dimer phase is at $\theta=\pi,A_{r}=0.655$. Here we use the
iTEBD method and the bond dimension is set as $\chi=120$.

To enhance the understanding of this model, let us consider a unitary
transformation on the even sites:

$\displaystyle U$ $\displaystyle=$ $\displaystyle\prod_{j}\exp(i\pi
S_{2j}^{x})\exp(i\frac{\pi}{2}S_{2j}^{z}),$ (11)

which will rewrite the biperiodic Kitave interactions in Eq. (5) into a
translationally invariant form. After the unitary transformation, the
Hamiltonian Eq. (5) can be rewritten:

$\displaystyle\tilde{H}$ $\displaystyle=$ $\displaystyle
K\sum_{j=1}^{N}S_{j}^{x}S_{j+1}^{y}+\frac{A_{r}\sin\theta}{\sqrt{2}}\sum_{j=1}^{N}(-1)^{j}\left[(S_{j}^{x})^{2}-(S_{j}^{y})^{2}\right]$
(12) $\displaystyle+$ $\displaystyle
A_{r}\cos\theta\sum_{j=1}^{N}(S_{j}^{z})^{2}.$

It is easy to see that the uniaxial SIA term remains in its original form, but
the rhombic SIA term could get a staggered phase factor. Interestingly, under
the rotation (11) the local bond parity operators are defined by

$\displaystyle\hat{W}_{j}=\Sigma_{j}^{y}\Sigma_{j+1}^{x},$ (13)

where the site parity operaotr is defined as $\Sigma_{j}^{a}\equiv e^{i\pi
S_{j}^{a}}$. Notably, these operators are diagonal within the bases defined in
Eq. (2). One can readily find that $\hat{W}_{j}$ is invariant by inspecting
$[\hat{W}_{j},\tilde{H}]=0$. As the eigenvalues of $\Sigma_{j}^{a}$ are $\pm
1$, the eigenvalues of $\hat{W}_{j}$ in Eq. (13) are related to
$\mathbb{Z}_{2}$-valued invariants, i.e., $w_{j}=\pm 1$. It is straightforward
to deduce that for a pair of nearest neighbor sites $\langle j,j{+}1\rangle$,
total $3{\times}3=9$ allowed states can be distinguished into the $w_{j}=1$
sector spanned by $|xy\rangle$, $|xz\rangle$, $|yx\rangle$, $|zy\rangle$,
$|zz\rangle$ and the $w_{j}=-1$ sector spanned by $|xx\rangle$, $|yy\rangle$,
$|yz\rangle$, $|zx\rangle$. The Hilbert space $\mathcal{H}$ decomposes into
$2^{N}$ dynamically disconnected Krylov subspaces of varying sizes,
$\mathcal{H}=\bigoplus_{n}\mathcal{K}_{n}$, each characterized by vectors
$\vec{w}=\\{w_{1},w_{n},\ldots,w_{N}\\}$.

###  IV.1 The north-south direction

Let us further consider the various limiting cases for the Hamiltonian in Eq.
(12). In the simplest limiting case of $\theta=0,A_{r}\to\infty$ ($A_{z}\gg
K$), all spins are confined to the state $|0\rangle$ and the ground state is
obviously the KSL phase. For $\theta=\pi,A_{r}\to\infty$ ($A_{z}\ll-K$), all
spins are allowed to be in the states $|\pm 1\rangle$. Therefore the ground
state is expected to undergo a phase transition. Three components of dimer
order parameter $O_{D}$ as a function of $A_{r}$ at $\theta=\pi$ in Fig. 2(c).
It is observed that the joint symmetry defined in Eq. (6) is broken after
$A_{r}>|A_{zc}|$. The $x$ component of the dimer order parameter increases
smoothly from zero to a finite value as the parameter $A_{r}$ increase across
the critical value, indicating a second-order transition occurs at
$\theta=\pi,A_{r}=0.655$. The presence of nonvanishing dimer correlations
before the critical point in dimer phase can be attributed to the limitations
imposed by the finite bond dimension. The dimer orders are associated with the
spontaneous breaking of translational symmetry of $O_{D}$ in an infinite
system. The emergence of the dimer ordering is distinct from the general
mechanism for the formation of dimerized phases, which is typically induced by
inherent bond alternation and the resulting breaking of translational
symmetry. The ground state is two-fold degenerate for dimer phase in the
thermodynamic limit, which is in contrast to the gapped ground state for the
KSL phase.

###  IV.2 The east-west direction

Next, we consider another limit $A_{r}\to\infty$ (equivalently $K\to 0$). We
analyze the scenario where $\theta=3\pi/2$, where the Hamiltonian (12)
simplifies to a form where the staggered rhombic SIA ($A_{x}=-A_{y}\gg K$ and
$A_{z}=0$).

$\displaystyle\tilde{H}$ $\displaystyle=$
$\displaystyle\sum_{j=1}^{N}S_{j}^{x}S_{j+1}^{y}+\frac{A_{r}}{\sqrt{2}}\sum_{j=1}^{N}(-1)^{j+1}\left[(S_{j}^{x})^{2}-(S_{j}^{y})^{2}\right].\quad$
(14)

Note the condition $\theta=\pi/2$ is obtained through a direct spatial
translation related to $\theta=3\pi/2$. In the limit of $A_{r}\to\infty$, the
system can be studied perturbatively using an effective Hamiltonian approach
based on the Schrieffer-Wolff transformation [57]. Considering
$\Sigma_{j}^{a}=1-2(S_{j}^{a})^{2}$ $(a=x,y,z)$, we rewrite $\tilde{H}$ in Eq.
(14) as $\tilde{H}=\tilde{H}_{0}+\tilde{V}$ , where

$\displaystyle\tilde{H}_{0}$ $\displaystyle=$
$\displaystyle\frac{A_{r}}{2\sqrt{2}}\sum_{j=1}^{N}(-1)^{j+1}(\Sigma_{j}^{y}-\Sigma_{j}^{x}),$
(15)

and

$\displaystyle\tilde{V}=\sum_{j=1}^{N}S_{j}^{x}S_{j+1}^{y}.$ (16)

To obtain the effective Hamiltonian, a unitary transformation is applied,
which is generated by an anti-Hermitian operator $S$:

$\displaystyle\tilde{H}_{\rm{eff}}=e^{S}\tilde{H}e^{-S}\equiv\tilde{H}_{0}+\sum_{n=1}^{\infty}\tilde{H}_{\rm{eff}}^{(n)},$
(17)

where each term $\tilde{H}_{\rm{eff}}^{(n)}$ is of order $(1/A_{r})^{n}$. The
Schrieffer-Wolff generator $S$ is formulated to ensure
$[\tilde{H}_{\rm{eff}}^{(n)},\tilde{H}_{0}]=0$ for any $n$. This is achieved
by expanding $S=\sum_{m=1}^{{\infty}}S^{(m)}$, where each $S^{(m)}$ takes the
order $(1/A_{r})^{m}$. Typically, selecting $S\equiv S^{(1)}$ simplifies the
first order term $\tilde{H}_{\rm{eff}}^{(1)}=0$, allowing for the computation
of higher-order terms based on $S^{(1)}$ and $\tilde{H}_{0}$. In this case, we
can obtain

$\displaystyle\tilde{H}_{\rm{eff}}\\!$ $\displaystyle=$
$\displaystyle\tilde{H}_{0}+\mathcal{P}([S^{(1)},\tilde{H}_{0}]+\tilde{V})\mathcal{P}$
(18)
$\displaystyle+\\!\mathcal{P}([S^{(1)},\\!\tilde{V}]\\!+\\!\frac{1}{2!}[S^{(1)},\\![S^{(1)},\\!\tilde{H}_{0}]])\mathcal{P}\\!+\\!\cdots\\!,$

where $\mathcal{P}$ is the projection operator that maps each term in the
expression above onto the eigenbasis of $\tilde{H}_{0}$. The requirement for
$\tilde{H}_{\rm{eff}}^{(1)}=0$ implies that $S^{(1)}$ must satisfy

$\displaystyle[S^{(1)},\tilde{H}_{0}]+\tilde{V}=0.$ (19)

Direct verification shows that Eq. (19) is satisfied by:

$\displaystyle
S^{(1)}=\frac{1}{\sqrt{2}A_{r}}\sum_{j=1}^{N}(-1)^{j+1}S_{j}^{x}\Sigma_{j}^{z}S_{j+1}^{y}.$
(20)

It is readily find that $[S^{(1)},\hat{W}_{j}]=0$. The leading-order effective
Hamiltonian is then determined as up to the second order:

$\displaystyle\tilde{H}_{\rm{eff}}^{(2)}$ $\displaystyle=$
$\displaystyle\frac{1}{4\sqrt{2}A_{r}}\sum_{j=1}^{N}(-1)^{j}(\Sigma_{j}^{z}-\Sigma_{j}^{y})(1-\Sigma_{j+1}^{y}).$
(21)

See Appendix A for details. In this case, when expanded up to the second
order, the effective Hamiltonian exhibits a diagonal structure in the
specified bases (2), representing an emergent conserved quantity enforced by
the Schrieffer-Wolff transformation. The rationale behind this will become
clearer in the spin-1/2 representation discussed in the following section.

Figure 3:  (a) Illustration of the topological mapping between the spin-1
Kitaev chain and the spin-1/2 PXP model. The prime lattice defined by the
spin-1 sites $\\{j\\}$ is mapped to the dual lattice with spin-1/2 degrees of
freedom located at the bond centers $\\{j+1/2\\}$. The upper figure shows how
the spin-1 state $|\cdots yxyxzzy\cdots\rangle$ is transformed into the
spin-1/2 state
$|\cdots\downarrow\uparrow\downarrow\uparrow\downarrow\downarrow\downarrow\uparrow\cdots\rangle$,
and the lower figure shows the mapping of the spin-1 state $|\cdots
xyxyxzz\cdots\rangle$ to the spin-1/2 state
$|\cdots\uparrow\downarrow\uparrow\downarrow\uparrow\downarrow\downarrow\downarrow\cdots\rangle$.
(b) The Hilbert space graph of $\hat{H}_{\rm{sdPXP}}$ in Eq. (25) with a
sufficiently large value of $A_{r}$ for $N=6$ sites with periodic boundary
conditions. Each cluster is the system’s eigenfunction, and the number of wave
functions represents the degeneracy of the system at this energy level. Nodes
of the graph are grouped according to the staggered number $n_{\rm s}$, from
the $|\mathbb{Z}_{2}\rangle$ or $|\mathbb{Z}^{\prime}_{2}\rangle$ states.

##  V Mapping to spin-1/2 extended PXP model

In the context of a spin-1 Kitaev chain, it has been established that the
ground state lies within the flux-free sector, characterized by
$\vec{w}=\\{1,1,\cdots,1\\}$. It is interesting to observe that such fractal
dimension-based subspace is the largest among all the sectors. Within this
constrained Hilbert space, there exists a one-to-one mapping between
permissible two-site configurations, $\langle j,j+1\rangle$, on a prime
lattice and the corresponding dual lattice, which features a spin-$1/2$ degree
of freedom located at the bond center, denoted as $i\equiv j+1/2$. This
correspondence is described by the following rule:

$\displaystyle|\cdots
zz\cdots\rangle_{j,j+1}\leftrightarrow|\cdots\downarrow\downarrow\downarrow\cdots\rangle_{j-\frac{1}{2},j+\frac{1}{2},j+\frac{3}{2}},$
$\displaystyle|\cdots
yx\cdots\rangle_{j,j+1}\leftrightarrow|\cdots\downarrow\uparrow\downarrow\cdots\rangle_{j-\frac{1}{2},j+\frac{1}{2},j+\frac{3}{2}},$
$\displaystyle|\cdots
zy\cdots\rangle_{j,j+1}\leftrightarrow|\cdots\downarrow\downarrow\uparrow\cdots\rangle_{j-\frac{1}{2},j+\frac{1}{2},j+\frac{3}{2}},$
$\displaystyle|\cdots
xz\cdots\rangle_{j,j+1}\leftrightarrow|\cdots\uparrow\downarrow\downarrow\cdots\rangle_{j-\frac{1}{2},j+\frac{1}{2},j+\frac{3}{2}},$
$\displaystyle|\cdots
xy\cdots\rangle_{j,j+1}\leftrightarrow|\cdots\uparrow\downarrow\uparrow\cdots\rangle_{j-\frac{1}{2},j+\frac{1}{2},j+\frac{3}{2}}.$
(22)

A schematic illustration of the mapping is presented in Fig. 3(a). Within the
flux-free sector, the Hamiltonian of Eq. (12) can be exactly mapped to a
compact form in the spin-1/2 representation:

$\displaystyle\hat{H}_{\rm ePXP}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}P_{i-1}X_{i}P_{i+1}+2A_{r}\cos\theta\sum_{i=1}^{N}P_{i-1}n_{i}P_{i+1}$
(23) $\displaystyle-$
$\displaystyle\sqrt{2}A_{r}\sin\theta\sum_{i=1}^{N}(-1)^{i}P_{i-1}n_{i}P_{i+1},$

where
$X=|\downarrow\rangle\langle\uparrow|+|\uparrow\rangle\langle\downarrow|$ and
$P=|\downarrow\rangle\langle\downarrow|$ acts as a projector that ensures
neighboring sites are not both in the excited state, and
$n=1-P=|\uparrow\rangle\langle\uparrow|$. It has been established the one-
dimensional lattice Schwinger model that contains confinement and
deconfinement phases, has been realized by using ultracold bosons in optical
lattices and the Rydberg-atom array by virtue of the PXP model.

In the extended PXP model in Eq. (23), particularly along the axis where
$\theta=\pi$, a large tends to favor atoms in the excited state. However, the
Rydberg blockade effect prohibits neighboring atoms from simultaneously
occupying the Rydberg state. Consequently, the energetically favorable
configuration under these constraints is for half of the atoms to be in the
excited state, with the Rydberg atoms occupying every other site in an
alternating sequence. This scenario results in an antiferromagnetic spin
configuration in the ground state, which breaks the lattice translational
symmetry and exhibits double degeneracy. Such degeneracy is disrupted when
$\theta$ deviates from the value of $\pi$. In the extended PXP model, the
staggered potential that differentiates on-site energies between even and odd
sites selectively stabilizes one of the degenerate states. This asymmetry
inhibits the occurrence of a quantum phase transition. In this framework,it is
pertinent to consider our Kitaev model as an ideal candidate for simulating
the lattice Schwinger model, where the $\theta$ parameter analogously serves
as the topological angle. Remarkably, both the mass of the matter field and
the topological angle can be experimentally modulated in the Kitaev model,
allowing exploration of the phase diagram from both directions. Furthermore,
there exists an exact mapping between the Kitaev model and the PXP model
within the ground-state manifold, where the $\mathbb{Z}_{2}$ gauge charge
$w_{j}=1$ for all $j$-th bonds. The automatic satisfaction of the
$\mathbb{Z}_{2}$ gauge condition in the ground-state manifold eliminates the
need for an explicit gauge constraint for achieving equivalence between the
LGT and the PXP model. Therefore, the KSL-dimer transition is reminiscence of
celebrated confinement-deconfinement transition, highlighting a unique
interplay between quantum simulation and topological phenomena.

The Rydberg blockade introduces kinetic constraints that are crucial for the
atypical dynamics observed in QMBS states. At the initialization time $t=0$,
the system is prepared in product states
$|\psi(0)\rangle=|\mathbb{Z}_{k}\rangle$, which are defined as follows:

$\displaystyle|\mathbb{Z}_{k}\rangle=|\cdots\underbrace{\uparrow\cdots\uparrow\downarrow}_{k}\underbrace{\uparrow\cdots\uparrow\downarrow}_{k}\cdots\underbrace{\uparrow\cdots\uparrow\downarrow}_{k}\cdots\rangle.$
(24)

Subsequently, the evolution of the system is governed by the PXP Hamiltonian
denoted as $\hat{H}_{\rm{PXP}}=\sum_{i=1}^{N}P_{i-1}X_{i}P_{i+1}$, with the
state at any time $t$ given by $|\psi(t)\rangle=\exp(-i\hat{H}_{\rm
PXP}t)|\psi(0)\rangle$. It was noted that the quantum quench from either
$|\mathbb{Z}_{2}\rangle$ or $|\mathbb{Z}_{3}\rangle$ exhibits periodic
revivals in the quantum fidelity $F(t)=|\langle\psi(0)|\psi(t)\rangle|^{2}$,
while $|\mathbb{Z}_{1}\rangle$ or $|\mathbb{Z}_{4}\rangle$ thermalize under
time evolution. The observed oscillations and apparent nonergodic dynamics are
due to the existence of equal spacing of the QMBS eigenstates. The von Neumann
entanglement entropy for a subsystem consisting of $L$ sites is defined as
$S_{L}=-{\rm tr}(\rho_{L}^{(n)}\log_{2}\rho_{L}^{(n)})$, where
$\rho_{L}^{(n)}$ is the reduced density matrix obtained by tracing out the
$N-L$ sites of the complement system from the density matrix
$\rho^{(n)}=|\chi_{n}\rangle\langle\chi_{n}|$ of each $n$th eigenstate
$|\chi_{n}\rangle$.

Figure 4: Half-chain entanglement entropies $S_{N/2}$ as a function of energy
$E$ for the staggered detuned PXP model in Eq. (25). (a) $A_{r}=1$, (b)
$A_{r}=5$, and (c) $A_{r}=10$. All data presented are for systems with length
$N=20$.  Figure 5: Energy level distribution for $\hat{H}_{\rm{sdPXP}}$ in Eq.
(25) with $N=20$ sites. The upper inset shows two sets of energy levels spaced
by $\sqrt{2}A_{r}$, corresponding to the primary clusters. The lower inset
shows the energy level spacing $1/\sqrt{2}A_{r}$ of nearly degenerate states
within each cluster, corresponding to the secondary clusters.

As discussed in Eq. (23), the spin-1 Kitaev model with rhombic and uniaxial
single-ion anisotropies can be perfectly mapped to the spin-1/2 PXP model with
staggered detuning and uniform detuning within the flux-free subspace [36].
The detuning term is ubiquitous in actual experiments. For example, the
uniform detuning of the driving laser from the excited state can be finely
tuned in cold-atom platforms. The static uniform detuning has been discussed
in detail in our previous work, and some interesting results have been
obtained. For example, the quench dynamics starting from the product states is
symmetric between positive and negative values of the static uniform detuning.
We find that the coherent oscillations of quantum fidelity and certain local
observables are sustained against small static uniform detuning perturbations
in a quantum quench from special initial states. While the oscillation
amplitudes of these observables decay with time as the static uniform detuning
strength is increased, the system completely thermalizes upon approaching the
critical point. Next we will focuses on the effect of the staggered detuning,
specifically at $\theta=3\pi/2$. In this case, the Hamiltonian for the
staggered detuned PXP model is given by

$\displaystyle\hat{H}_{\rm{sdPXP}}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\\!P_{i-1}X_{i}P_{i+1}\\!+\\!\sqrt{2}A_{r}\\!\sum_{i=1}^{N}\\!(-1)^{i}\\!P_{i-1}n_{i}P_{i+1}.\quad\quad$
(25)

Figure 4 shows the half-chain entanglement entropy $S_{N/2}$ as a function of
energy $E$ for different values of the staggered detuning $A_{r}$. Figure 4(a)
demonstrates that with a low amplitude of $A_{r}$, the spectrum is populated
by a large number of non-degenerate states. As $A_{r}$ is increased, the
energy levels restructure towards a uniformly spaced distribution [58]. As
shown in Fig. 4(b), there exists a mix of states including anomalous states
characterized by low $S_{N/2}$ and typical thermal state with high $S_{N/2}$
within each primary cluster. As $A_{r}$ reaches sufficiently high values, the
eigenstates of the effective Hamiltonian modeled in Eq. (25) organize into
nearly equidistant clusters, as depicted in Fig. 4(c). Figure 5 depicts the
energy level distribution for $\hat{H}_{\rm{sdPXP}}$ in Eq. (25) at a moderate
value of $A_{r}$, which clearly identifies several primary sets of energy
levels, each separated by $\sqrt{2}A_{r}$. For systems with an even number of
sites $N$, the energy levels are organized into $N+1$ distinct sectors. This
classification hinges on the staggered occupancy number defined by:

$\displaystyle n_{\rm s}=\sum_{j=1}^{N}(-1)^{j}n_{j}.$ (26)

The staggered occupancy number plays a critical role, analogous to the Hamming
distance, in differentiating states across these sectors. This concept is
illustrated in Fig. 3(b), where $n_{s}$ effectively acts as a metric for
spatial coordinates, defining physical boundaries between
$|\mathbb{Z}_{2}\rangle$ and $|\mathbb{Z}^{\prime}_{2}\rangle$ states, where
$|\mathbb{Z}_{2}^{\prime}\rangle=|\downarrow\uparrow\downarrow\uparrow\downarrow\uparrow\cdots\rangle$
is obtained by translating one lattice spacing on $|\mathbb{Z}_{2}\rangle$.
Note that the concept of a non-Hermitian scar has been utilized, where the
Hamming distance serves as a proxy for spatial coordinates, and these two
states function as physical boundaries [59]. Within these primary clusters,
there is secondary clustering of states, spaced at intervals of
$1/\sqrt{2}A_{r}$. As $A_{r}$ increases towards infinity, the degeneracy of
the system markedly increases. Despite partial resolution of these
degeneracies by the Schrieffer-Wolff Hamiltonian, the spectrum continues to
exhibit a considerable degree of degeneracy at finite $A_{r}$ values.

The presence of clustering of eigenstates in Eq. (25) for strong staggered
detuning $A_{r}$, which leads to violation of ETH at intermediate frequencies,
can be comprehended via the effective spin-1/2 Hamiltonian. According to the
mapping rule (V), Eq. (15) can be rewritten as

$\displaystyle\hat{H}_{\rm{eff}}^{(0)}\\!$ $\displaystyle=$
$\displaystyle\\!\sqrt{2}A_{r}\sum_{i=1}^{N}(-1)^{i}P_{i-1}n_{i}P_{i+1}\equiv\sqrt{2}A_{r}\hat{H}_{\rm{sPnP}},\quad$
(27)

and Eq. (21) is recast into

$\displaystyle\hat{H}_{\rm{eff}}^{(2)}\\!=\\!\frac{1}{\sqrt{2}A_{r}}\\!\sum_{i=1}^{N}(-1)^{i}\\!P_{i-1}\\!Z_{i}P_{i+1}\\!\equiv\\!\frac{1}{2\sqrt{2}A_{r}}\\!\hat{H}_{\rm{sPZP}}.$
(28)

Meanwhile, the Schrieffer-Wolff generator $S^{(1)}$ to the spin-1/2 system can
also be mapped as

$\displaystyle
S^{(1)}\\!=\\!\frac{i}{\sqrt{2}A_{r}}\\!\sum_{i=1}^{N}(-1)^{i}P_{i-1}Y_{i}P_{i+1}\\!\equiv\\!\frac{i}{\sqrt{2}A_{r}}\hat{H}_{\rm{sPYP}}.$
(29)

One finds the forming of primary clusters is a result of near integrable
nature of the effective Hamiltonian in Eq. (27) for large staggered detuning.
We show that these secondary clusters at a much smaller quasienergy scale
result from the second-order effective Hamiltonian $\hat{H}_{\rm{eff}}^{(2)}$
in Eq. (28) at a large $A_{r}$. In fact, the diagonal nature can be attributed
to the roles of Eqs. (27), (28), and (29) as the generators ”X”, ”Z”, and ”Y”
of a nearly-SU(2) algebra. Specifically, the commutation relations are given
by:

$\displaystyle\left[\hat{H}_{\rm{sPnP}},\hat{H}_{\rm{PXP}}\right]$
$\displaystyle=$ $\displaystyle i\hat{H}_{\rm{sPYP}},$ (30a)
$\displaystyle\left[\hat{H}_{\rm{sPnP}},\hat{H}_{\rm{sPYP}}\right]$
$\displaystyle=$ $\displaystyle-i\hat{H}_{\rm{PXP}},$ (30b)
$\displaystyle\left[\hat{H}_{\rm{PXP}},\hat{H}_{\rm{sPYP}}\right]$
$\displaystyle=$ $\displaystyle i\hat{H}_{\rm{sPZP}}.$ (30c)

It is interesting to note that the right-hand side of Eq. (30c) does not
adhere strictly to the SU(2) algebra. This deviation is encapsulated by the
expression:

$\displaystyle\hat{H}_{\rm{sPZP}}=(\hat{H}_{\rm{sPnP}}+O_{zzz}),$ (31)

where $O_{zzz}=\frac{1}{2}\sum_{i=1}^{N}(-1)^{i}Z_{i-1}Z_{i}Z_{i+1}$
represents a residual contribution. Such nearly-SU(2) algebra was studied in
the context of QMBS in the PXP model in Refs. [32, 60, 61, 62]. The algebraic
counterpart for Eq. (14) is detailed in Appendix A.

Figure 6: Dynamics of quantum fidelity for the staggered detuned PXP model in
Eq. (25) and the effective Hamiltonian
$\hat{H}_{\rm{eff}}^{(0)}+\hat{H}_{\rm{eff}}^{(2)}$ [cf. Eqs. (27-28)] with
$N=28$ sites starting from initial state
$(|\mathbb{Z}_{2}\rangle+|\mathbb{Z}_{2}^{\prime}\rangle)/\sqrt{2}$. (a)
$A_{r}=1$, (b) $A_{r}=5$, and (c) $A_{r}=10$. Here, the first three periods
are chosen for comparison.

The observed clustering behaviors leading to violations of ETH in the
staggered Rydberg chain are not attributable to conventional Hilbert space
fragmentation, as Eq. (25) does not exhibit dynamically disconnected sectors
for any finite $A_{r}$. We explore the dynamics at short timescales, focusing
on an effective Hamiltonian up to second order, induced by deformations as
detailed in Eq. (25). Figure 6 shows the quantum fidelity between the original
Hamiltonian $\hat{H}_{\rm{sdPXP}}$ in Eq. (25) and the effective Hamiltonian,
represented by $\hat{H}_{\rm{eff}}^{(0)}+\hat{H}_{\rm{eff}}^{(2)}$ [cf. Eqs.
(27-28)], for various values of $A_{r}$ at $N=28$. The initial state used for
these measurements is the cat state
$(|\mathbb{Z}_{2}\rangle+|\mathbb{Z}_{2}^{\prime}\rangle)/\sqrt{2}$. In Fig.
6(a) the effective Hamiltonian reveals no decay in fidelity over the short-
time evolution, indicating a significant display of slow dynamics. At very
small values of $A_{r}$, the fidelity evolution of the effective Hamiltonian
diverges clearly from that of the original model in Eq. (25), due to the
absence of off-diagonal terms of the effective Hamiltonian. As $A_{r}$
increases, the trajectories of the two models begin to converge, yet a
noticeable deviation emerges after the second oscillation period as
illustrated in Fig. 6(b). Notably, Fig. 6(c) shows that for $A_{r}=10$, the
periodic revivals in the quantum fidelity for the effective Hamiltonian,
starting from the superposed initial state
$(|\mathbb{Z}_{2}\rangle+|\mathbb{Z}_{2}^{\prime}\rangle)/\sqrt{2}$, nearly
perfectly align with those of the original Hamiltonian $\hat{H}_{\rm{sdPXP}}$
in Eq. (25). Consequently, it is concluded that when $A_{r}$ reaches an order
of magnitude around 10, the effective Hamiltonian provides an accurate
description for the dynamics of the original model.

The glassy dynamics with memory effects which imply that the system can retain
memory of its initial state for longer times due to the restricted dynamics
within isolated subspaces, observed in Fig. 6 highlights the role of Hilbert
space fragmentation (HSF) in shaping the system’s behavior. We have previously
identified a hierarchical fragmentation within the spin-1 Kitaev model’s
Hilbert space. The primary level of fragmentation divides the entire Hilbert
space into $2^{N}$ dynamically isolated Krylov subspaces. These subspaces vary
in size and are characterized by symmetry sector
$\vec{w}=\\{w_{1},w_{n},\cdots,w_{N}\\}$ of gauge charge. The largest of these
subspace, imposed by all positive $\mathbb{Z}_{2}$ gauge charges, is termed
the Fibonacci Hilbert space [35]. Connectivity graphs of the spin-1/2
representation further delineate the hierarchical fragmentation structure in
Fig. 3(b), where the subsectors are further differentiated by their degeneracy
levels and clustering, as shown in Table 1.

$N$ | Degeneracy | $n_{\rm clustering}$ | Dimension
---|---|---|---
$4$ |  |  |  |  | $1$ | $2$ | $1$ | $2$ | $1$ |  |  |  |  | $5$ | $7$
$6$ |  |  |  | $1$ | $3$ | $3$ | $4$ | $3$ | $3$ | $1$ |  |  |  | $7$ | $18$
$8$ |  |  | $1$ | $4$ | $6$ | $8$ | $9$ | $8$ | $6$ | $4$ | $1$ |  |  | $9$ | $47$
$10$ |  | $1$ | $5$ | $10$ | $15$ | $20$ | $21$ | $20$ | $15$ | $10$ | $5$ | $1$ |  | $11$ | $123$
$12$ | $1$ | $6$ | $15$ | $26$ | $39$ | $48$ | $52$ | $48$ | $39$ | $26$ | $15$ | $6$ | $1$ | $13$ | $322$
Table 1: The degenerate manifold encoded in Eq. (25) for a sufficiently large
value of $A_{r}$, clustering number and the dimension of Hilbert space for
different even system sizes $N$.

One can find the exponential growth in both the number of clusters and the
dimensionality of the Hilbert space with increasing system size $N$. In
particular, the largest symmetry sector within the Fibonacci space is
identified by $n_{s}=0$. As this hierarchy is distinguished by how the number
of disconnected subspaces within a symmetry sector scales with the system size
$N$, it is tempting to understand the form of fragmentation. In systems
exhibiting strong fragmentation, the largest connected subspace within any
given symmetry sector is exponentially smaller than the entire sector. The
ratio of the dimension of the largest connected subspace $\mathcal{D}_{\rm
max\ subsector}$ to the dimension of the full symmetry sector
${\mathcal{D}_{\rm full\ sector}}$ decreases exponentially with the system
size $N$ [63, 64]

$\displaystyle\frac{\mathcal{D}_{\rm max\ subsector}}{\mathcal{D}_{\rm full\
sector}}\sim{\rm exp}(-\kappa N),$ (32)

where $\kappa>0$. Otherwise, the system is said to be weakly fragmented.
Figure 7 shows both the spin-1 Kitaev chain in Eq. (14) and the staggered
detuned PXP model in Eq. (25) with a sufficiently large value of $A_{r}$
display strong HSF. Remarkably, the decay rate in the former case is more
pronounced, suggesting that the Hilbert space of the spin-1 Kitaev chain is
indeed “more fragmented” than the staggered detuned PXP model.

Figure 7: $\mathcal{D}_{\rm max\ subsector}/\mathcal{D}_{\rm full\ sector}$ as
the function of $N$ for different spin. Strong fragmentation for (a) the
spin-1 Kitaev chain in Eq. (14) and (b) the spin-1/2 Hamiltonian
$\hat{H}_{\rm{sdPXP}}$ in Eq. (25) with a sufficiently large value of $A_{r}$.

##  VI Summary

In summary, we study a one-dimensional spin-1 Kitaev models with three-
component single-ion anisotropies (SIAs). To effectively describe the effects
of SIAs on ground-state properties and far-from-equilibrium dynamics, we use a
three-dimensional coordinate, i.e., $A_{x}\\!=\\!A_{r}\sin\theta\cos\phi$,
$A_{y}\\!=\\!A_{r}\sin\theta\sin\phi$ and $A_{z}\\!=\\!A_{r}\cos\theta$. At
zero SIA field ($A_{r}=0$), the local bond operators of the pure Kitaev chain
are conserved with all eigenvalues $w_{j}=\pm 1$. Even in the presence of
SIAs, these local bond operators retain their conservation, playing a role of
the gauge charge. In terms of the infinite time evolving block decimation
(iTEBD) method, the comprehensive phase diagram is given. In the zero SIA
field limit, the Kitaev chain hosts the Kitaev spin liquid (KSL) ground state,
characterized by extremely short-ranged spin-spin correlations. At
$\theta=\pi$, where the $\phi$ is irrelevant, a phase transition is triggered
by increasing $A_{r}$ from the KSL phase to the dimer phase. Due to the
spontaneous breaking of translational symmetry, the dimer phase, which is two-
fold degenerate, is analogous to the deconfined phase in the lattice Schwinger
model, which incorporates a celebrated topological $\theta$-angle of $\pi$. To
this end, the KSL-dimer transition is equivalent to the confinement-
deconfinement transition. When a translation-symmetry-breaking term is added,
the topological angle $\theta$ starts to deviate from $\pi$. In this case,
such Coleman’s transition no longer exists and the deconfined line disappears
as a consequence of charge-parity symmetry breaking. We here find that the
polar angle $\theta$ can act as the topological angle in the U(1) lattice
gauge theory. When $\phi=\pi/4$, the deconfined line expands into a semi-
infinite region, while for other $\phi$ values, it remains a simple deconfined
line.

We focus on the scenario where $\phi=3\pi/4$, where the SIAs can be recast
into the mixture of uniaxial and rhombic SIAs. Under these circumstances, such
two types of SIAs can be tuned independently. Notably, the varying uniaxial
SIA induces the KSL-dimer transition at $\theta=\pi$ and a finite rhombic SIA
can modify $\theta$ away from $\pi$. Considering the local quadrupole order
parameter (8) coupled with the rhombic SIA in Hamiltonian (5), the spin
nematic order, which breaks the spin-rotational symmetry, will be induced when
the rhombic SIA becomes predominant. Consequently, a positive rhombic SIA
leads to the emergence of the $x$-ferroquadrupole phase, whereas a negative
sign of the rhombic SIA stabilizes the $y$-ferroquadrupole phase.
Interestingly, as the rhombic SIAs are tuned from negative to positive values,
the system from the $y$-ferroquadrupole to the $x$-ferroquadrupole phases can
either pass through the deconfined line via a first-order phase transition or
undergo a crossover through the KSL phase, which gives us an example of
unnecessary criticality.

It has been recognized the spin-1/2 PXP model is embedded in the ground-state
manifold of the spin-1 Kitaev model, where the topological charges are $+1$.
The quantum quench from either $|\mathbb{Z}_{2}\rangle$ or
$|\mathbb{Z}_{2}^{\prime}\rangle$ states exhibits periodic revivals in the
quantum fidelity at the PXP model, referred to as quantum many-body scars
(QMBSs). To investigate the dynamical properties far from the pure Kitaev
limit, we utilize a Schrieffer-Wolff transformation to develop a perturbative
description of the spin-1 Kitaev chain in the large limit of rhombic SIA. We
find that for a sufficiently large value of $A_{r}$, the effective
Hamiltonian, expanded up to the second order, remains diagonal. In order to
provide a clearer explication, we map the spin-1 Hamiltonian (12) to the
spin-1/2 extended PXP Hamiltonian (25) exactly in the ground-state manifold
characterized by all positive $\mathbb{Z}_{2}$ gauge charges, where the
uniaxial SIA is equivalent to the uniform detuning while the rhombic SIA plays
a role of staggered detuning. For the PXP model with spatially uniform
detuning, the mapping to U(1) lattice gauge theory is that the system hosts a
regime exhibiting quark-antiquark confinement, which provides a route to
investigate confinement in continuous gauge theories with quantum simulators.
Intriguingly, for the case of the PXP model with staggered detuning, a
Schrieffer-Wolff transformation adopted in the large staggered detuning
enforces an additional conservation law for the staggered occupation number
$n_{s}$ (26). We show that when the staggered detuning is sufficiently large,
the eigenspectrum is strongly clustered. In the large staggered detuning
limit, an additional conservation law for the staggered occupation number
$n_{s}$ emerges and the energy level exhibit a equally spaced distributions.

Returning to the spin-1 Kitaev chain, due to the existence of the conserved
gauge charges $\hat{W}_{j}$ (13), a strong fragmentation of the Hilbert space
into exponentially many disconnected sectors. The largest subspace with all
positive gauge charges corresponds to the Fibonacci Hilbert space, reminiscent
of a spin-1/2 system with Rydberg blockade. We trace the origin of the
emergent conserved quantity (28) to the nearly-SU(2) algebra, which is also
believed to underlie QMBS in the PXP model at zero staggered field.
Furthermore, the emergent symmetry (28) leads, at the second order in the
effective Hamiltonian description, to a substantially more fragmented Hilbert
space. The interplay of these two conservation laws results in hierarchical
Hilbert space fragmentation (HSF), leading to glassy dynamics in specific
limiting cases, which reveals no decay in fidelity over the short-time
evolution. Despite the surprising algebraic connection, such slow dynamics
belongs to a scenario distinct from the QMBS arising from the Fibonacci
constraint. Our work not only introduces a theoretical perspective in
simulating the topological $\theta$-angle and ergodicity-breaking dynamics,
but also paves the way for exploring higher-spin generalizations of the
scarred model. Such solid-state-based quantum systems could be experimentally
realized by state-of-the-art cold-atom quantum simulators.

###### Acknowledgements.

We thank Gaoyong Sun, Chengjie Zhang, and Xuejia Yu for valuable discussions
on related work. This work is supported by the National Natural Science
Foundation of China (NSFC) under Grant No. 12174194, Postgraduate Research &
Practice Innovation Program of Jiangsu Province, under Grant No. KYCX23_0347,
Opening Fund of the Key Laboratory of Aerospace Information Materials and
Physics (Nanjing University of Aeronautics and Astronautics), MIIT, Top-notch
Academic Programs Project of Jiangsu Higher Education Institutions (TAPP), and
stable supports for basic institute research under Grant No. 190101.

##  Appendix A Schrieffer-Wolff transformation and nearly-SU(2) algebra in
spin-1 systems

The derivation process is as follows. First, we define the projection operator
that projects the Hamiltonian to the flux-free sector
$\vec{w}{=}\\{1,1,\cdots,1\\}$

$\displaystyle\mathcal{P}=\prod_{j}(\frac{\Sigma_{j}^{y}\Sigma_{j+1}^{x}+\mathbbm{1}}{2}).$
(33)

The requirement for $\tilde{H}_{\rm{eff}}^{(1)}=0$ implies that $S^{(1)}$ must
satisfy $[S^{(1)},\tilde{H}_{0}]+\tilde{V}=0.$ Next, we calculate the
commutation relations in the constrained Hilbert space

$\displaystyle\mathcal{P}[S^{(1)},H_{0}]\mathcal{P}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{j=1}^{N}\left(\mathcal{P}\left[(-1)^{j+1}S_{j}^{x}\Sigma_{j}^{z}S_{j+1}^{y},(-1)^{j+1}(S_{j}^{x})^{2}-(S_{j}^{y})^{2}\right]\right.\mathcal{P}$
$\displaystyle+\left.\mathcal{P}\left[(-1)^{j+1}S_{j}^{x}\Sigma_{j}^{z}S_{j+1}^{y},(-1)^{j+2}(S_{j+1}^{x})^{2}-(S_{j+1}^{y})^{2}\right]\mathcal{P}\right)$
$\displaystyle=$
$\displaystyle\frac{1}{2}\left(-S_{j}^{x}S_{j+1}^{y}-S_{j}^{x}S_{j+1}^{y}\right)$
$\displaystyle=$ $\displaystyle-\tilde{V}.$ (34)

Therefore, the selected $S^{(1)}$ satisfies $[S^{(1)},H_{0}]+V=0$. The second-
order effective Hamiltonian can be calculated as

$\displaystyle\tilde{H}_{\rm{eff}}^{(2)}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\mathcal{P}[S^{(1)},V]\mathcal{P}$ (35)
$\displaystyle=$
$\displaystyle\frac{1}{2}\mathcal{P}\left[\frac{1}{\sqrt{2}A_{r}}\sum_{j=1}^{N}(-1)^{j+1}S_{j}^{x}\Sigma_{j}^{z}S_{j+1}^{y},\sum_{j=1}^{N}S_{j}^{x}S_{j+1}^{y}\right]\mathcal{P}$
$\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}A_{r}}\sum_{j=1}^{N}(-1)^{j}(S_{j}^{x})^{2}\Sigma_{j}^{z}(S_{j+1}^{y})^{2}$
$\displaystyle=$
$\displaystyle\frac{1}{4\sqrt{2}A_{r}}\sum_{j=1}^{N}(-1)^{j}(\Sigma_{j}^{z}-\Sigma_{j}^{y})(1-\Sigma_{j+1}^{y}).$

The diagonality of $\tilde{H}_{\rm{eff}}^{(2)}$ in Eq. (21) is a consequence
of the fact that $\tilde{H}_{0}$ in Eq. (15) and $\tilde{V}$ in Eq. (16) can
be viewed as generators of a “nearly-SU(2)” algebra. The nearly-SU(2) algebra
is generated by the operators $\tilde{H}_{\rm n}$, $\tilde{H}_{\rm X}$, and
$\tilde{H}_{\rm Y}$, defined by

$\displaystyle\tilde{H}_{\rm n}$ $\displaystyle\equiv$
$\displaystyle\frac{\tilde{H}_{0}}{\sqrt{2}A_{r}}=\frac{1}{4}\sum_{j=1}^{N}(-1)^{j+1}(\Sigma_{j}^{y}-\Sigma_{j}^{x}),$
(36) $\displaystyle\tilde{H}_{\rm X}$ $\displaystyle\equiv$
$\displaystyle\tilde{V}=\sum_{j=1}^{N}S_{j}^{x}S_{j+1}^{y},$ (37)
$\displaystyle\tilde{H}_{\rm Y}$ $\displaystyle\equiv$
$\displaystyle-i\sqrt{2}A_{r}S^{(1)}=i\sum_{j=1}^{N}(-1)^{j}S_{j}^{x}\Sigma_{j}^{z}S_{j+1}^{y},$
(38)

and

$\displaystyle\tilde{H}_{\rm Z}\\!$ $\displaystyle\equiv$
$\displaystyle\\!2\\!\sqrt{2}\\!A_{r}\\!\tilde{H}_{\rm{eff}}^{(2)}\\!=\\!\frac{1}{2}\\!\sum_{j=1}^{N}\\!(-1)^{j}\\!(\\!\Sigma_{j}^{z}\\!-\\!\Sigma_{j}^{y}\\!)\\!(\\!1\\!-\\!\Sigma_{j+1}^{y}\\!).$
(39)

For periodic boundary condition, the commutation relations among these
generators are given by

$\displaystyle\left[\tilde{H}_{\rm n},\tilde{H}_{\rm X}\right]$
$\displaystyle=$ $\displaystyle i\tilde{H}_{\rm Y},$ (40)
$\displaystyle\left[\tilde{H}_{\rm n},\tilde{H}_{\rm Y}\right]$
$\displaystyle=$ $\displaystyle-i\tilde{H}_{\rm X},$ (41)
$\displaystyle\left[\tilde{H}_{\rm X},\tilde{H}_{\rm Y}\right]$
$\displaystyle=$ $\displaystyle i\tilde{H}_{\rm Z}.$ (42)

From the above, we see that $\tilde{H}_{0}$ and $\tilde{V}$ are the ”$Z$” and
”$X$” generators of the nearly-SU(2) algebra. The Schrieffer-Wolff generator
$S^{(1)}$ in Eq. (20) that eliminates off-diagonal matrix elements of Eq. (14)
to leading order is related to the ”$Y$” generator of this algebra.

## References

  * Srednicki [1999] M. Srednicki, The approach to thermal equilibrium in quantized chaotic systems, Journal of Physics A: Mathematical and General 32, 1163 (1999).
  * D’Alessio and Rigol [2014] L. D’Alessio and M. Rigol, Long-time behavior of isolated periodically driven interacting lattice systems, Phys. Rev. X 4, 041048 (2014).
  * Luca D’Alessio and Rigol [2016] A. P. Luca D’Alessio, Yariv Kafri and M. Rigol, From quantum chaos and eigenstate thermalization to statistical mechanics and thermodynamics, Advances in Physics 65, 239 (2016).
  * Deutsch [2018] J. M. Deutsch, Eigenstate thermalization hypothesis, Reports on Progress in Physics 81, 082001 (2018).
  * Browaeys and Lahaye [2020] A. Browaeys and T. Lahaye, Many-body physics with individually controlled Rydberg atoms, Nature Physics 16, 132 (2020).
  * Ho  _et al._ [2019] W. W. Ho, S. Choi, H. Pichler, and M. D. Lukin, Periodic orbits, entanglement, and quantum many-body scars in constrained models: Matrix product state approach, Phys. Rev. Lett. 122, 040603 (2019).
  * Giudici  _et al._ [2022] G. Giudici, M. D. Lukin, and H. Pichler, Dynamical preparation of quantum spin liquids in Rydberg atom arrays, Phys. Rev. Lett. 129, 090401 (2022).
  * Yao  _et al._ [2022] Z. Yao, L. Pan, S. Liu, and P. Zhang, Bounding entanglement entropy using zeros of local correlation matrices, Phys. Rev. Res. 4, L042037 (2022).
  * Bernien  _et al._ [2017] H. Bernien, S. Schwartz, A. Keesling, H. Levine, A. Omran, H. Pichler, S. Choi, A. S. Zibrov, M. Endres, M. Greiner, V. Vuletić, and M. D. Lukin, Probing many-body dynamics on a 51-atom quantum simulator, Nature 551, 579 (2017).
  * Turner  _et al._ [2018] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Serbyn, and Z. Papić, Weak ergodicity breaking from quantum many-body scars, Nature Physics 14, 745 (2018).
  * Su  _et al._ [2023] G.-X. Su, H. Sun, A. Hudomal, J.-Y. Desaules, Z.-Y. Zhou, B. Yang, J. C. Halimeh, Z.-S. Yuan, Z. Papić, and J.-W. Pan, Observation of many-body scarring in a Bose-Hubbard quantum simulator, Phys. Rev. Res. 5, 023010 (2023).
  * Zhang  _et al._ [2023a] P. Zhang, H. Dong, Y. Gao, L. Zhao, J. Hao, J.-Y. Desaules, Q. Guo, J. Chen, J. Deng, B. Liu, W. Ren, Y. Yao, X. Zhang, S. Xu, K. Wang, F. Jin, X. Zhu, B. Zhang, H. Li, C. Song, Z. Wang, F. Liu, Z. Papić, L. Ying, H. Wang, and Y.-C. Lai, Many-body Hilbert space scarring on a superconducting processor, Nature Physics 19, 120 (2023a).
  * Moudgalya  _et al._ [2018] S. Moudgalya, S. Rachel, B. A. Bernevig, and N. Regnault, Exact excited states of nonintegrable models, Phys. Rev. B 98, 235155 (2018).
  * van Voorden  _et al._ [2020] B. van Voorden, J. c. v. Minář, and K. Schoutens, Quantum many-body scars in transverse field Ising ladders and beyond, Phys. Rev. B 101, 220305 (2020).
  * Michailidis  _et al._ [2020] A. A. Michailidis, C. J. Turner, Z. Papić, D. A. Abanin, and M. Serbyn, Stabilizing two-dimensional quantum scars by deformation and synchronization, Phys. Rev. Res. 2, 022065 (2020).
  * Lin  _et al._ [2020] C.-J. Lin, V. Calvera, and T. H. Hsieh, Quantum many-body scar states in two-dimensional Rydberg atom arrays, Phys. Rev. B 101, 220304 (2020).
  * Schecter and Iadecola [2019] M. Schecter and T. Iadecola, Weak ergodicity breaking and quantum many-body scars in spin-1 XY magnets, Phys. Rev. Lett. 123, 147201 (2019).
  * Chattopadhyay  _et al._ [2020] S. Chattopadhyay, H. Pichler, M. D. Lukin, and W. W. Ho, Quantum many-body scars from virtual entangled pairs, Phys. Rev. B 101, 174308 (2020).
  * Gotta  _et al._ [2023] L. Gotta, S. Moudgalya, and L. Mazza, Asymptotic quantum many-body scars, Phys. Rev. Lett. 131, 190401 (2023).
  * Chen  _et al._ [2023] Q. Chen, S. A. Chen, and Z. Zhu, Weak ergodicity breaking in non-Hermitian many-body systems, SciPost Phys. 15, 052 (2023).
  * Omiya and Müller [2023a] K. Omiya and M. Müller, Fractionalization paves the way to local projector embeddings of quantum many-body scars, Phys. Rev. B 108, 054412 (2023a).
  * Omiya and Müller [2023b] K. Omiya and M. Müller, Quantum many-body scars in bipartite Rydberg arrays originating from hidden projector embedding, Phys. Rev. A 107, 023318 (2023b).
  * Logarić  _et al._ [2024] L. Logarić, S. Dooley, S. Pappalardi, and J. Goold, Quantum many-body scars in dual-unitary circuits, Phys. Rev. Lett. 132, 010401 (2024).
  * Moudgalya  _et al._ [2020] S. Moudgalya, N. Regnault, and B. A. Bernevig, $\eta$-pairing in hubbard models: From spectrum generating algebras to quantum many-body scars, Phys. Rev. B 102, 085140 (2020).
  * Sanada  _et al._ [2023] K. Sanada, Y. Miao, and H. Katsura, Quantum many-body scars in spin models with multibody interactions, Phys. Rev. B 108, 155102 (2023).
  * Deng and Yang [2023] W. Deng and Z.-C. Yang, Using models with static quantum many-body scars to generate time-crystalline behavior under periodic driving, Phys. Rev. B 108, 205129 (2023).
  * Moudgalya  _et al._ [2022] S. Moudgalya, B. A. Bernevig, and N. Regnault, Quantum many-body scars and Hilbert space fragmentation: a review of exact results, Reports on Progress in Physics 85, 086501 (2022).
  * Francica and Dell’Anna [2023] G. Francica and L. Dell’Anna, Hilbert space fragmentation in a long-range system, Phys. Rev. B 108, 045127 (2023).
  * Andreadakis and Zanardi [2023] F. Andreadakis and P. Zanardi, Coherence generation, symmetry algebras, and Hilbert space fragmentation, Phys. Rev. A 107, 062402 (2023).
  * Lin and Motrunich [2019] C.-J. Lin and O. I. Motrunich, Exact quantum many-body scar states in the Rydberg-blockaded atom chain, Phys. Rev. Lett. 122, 173401 (2019).
  * Mark  _et al._ [2020] D. K. Mark, C.-J. Lin, and O. I. Motrunich, Exact eigenstates in the lesanovsky model, proximity to integrability and the PXP model, and approximate scar states, Phys. Rev. B 101, 094308 (2020).
  * Choi  _et al._ [2019] S. Choi, C. J. Turner, H. Pichler, W. W. Ho, A. A. Michailidis, Z. Papić, M. Serbyn, M. D. Lukin, and D. A. Abanin, Emergent SU(2) dynamics and perfect quantum many-body scars, Phys. Rev. Lett. 122, 220603 (2019).
  * Khemani  _et al._ [2019a] V. Khemani, C. R. Laumann, and A. Chandran, Signatures of integrability in the dynamics of Rydberg-blockaded chains, Phys. Rev. B 99, 161101 (2019a).
  * You  _et al._ [2020] W.-L. You, G. Sun, J. Ren, W. C. Yu, and A. M. Oleś, Quantum phase transitions in the spin-1 Kitaev-Heisenberg chain, Phys. Rev. B 102, 144437 (2020).
  * You  _et al._ [2022] W.-L. You, Z. Zhao, J. Ren, G. Sun, L. Li, and A. M. Oleś, Quantum many-body scars in spin-1 Kitaev chains, Phys. Rev. Res. 4, 013103 (2022).
  * Zhang  _et al._ [2023b] W.-Y. Zhang, Y.-N. Wang, D. Liu, J. Ren, J. Li, N. Wu, A. M. Oleś, and W.-L. You, Quantum many-body scars in spin-1 Kitaev chain with uniaxial single-ion anisotropy, Phys. Rev. B 108, 104411 (2023b).
  * Ebadi  _et al._ [2021] S. Ebadi, T. T. Wang, H. Levine, A. Keesling, G. Semeghini, A. Omran, D. Bluvstein, R. Samajdar, H. Pichler, W. W. Ho, S. Choi, S. Sachdev, M. Greiner, V. Vuletić, and M. D. Lukin, Quantum phases of matter on a 256-atom programmable quantum simulator, Nature 595, 227 (2021).
  * Bluvstein  _et al._ [2021] D. Bluvstein, A. Omran, H. Levine, A. Keesling, G. Semeghini, S. Ebadi, T. T. Wang, A. A. Michailidis, N. Maskara, W. W. Ho, S. Choi, M. Serbyn, M. Greiner, V. Vuletić, and M. D. Lukin, Controlling quantum many-body dynamics in driven Rydberg atom arrays, Science 371, 1355 (2021).
  * Ghosh  _et al._ [2018] R. Ghosh, A. Sen, and K. Sengupta, Ramp and periodic dynamics across non-Ising critical points, Phys. Rev. B 97, 014309 (2018).
  * Whitsitt  _et al._ [2018] S. Whitsitt, R. Samajdar, and S. Sachdev, Quantum field theory for the chiral clock transition in one spatial dimension, Phys. Rev. B 98, 205118 (2018).
  * Keesling  _et al._ [2019] A. Keesling, A. Omran, H. Levine, H. Bernien, H. Pichler, S. Choi, R. Samajdar, S. Schwartz, P. Silvi, S. Sachdev, P. Zoller, M. Endres, M. Greiner, V. Vuletić, and M. D. Lukin, Quantum Kibble–Zurek mechanism and critical dynamics on a programmable Rydberg simulator, Nature 568, 207 (2019).
  * Bull  _et al._ [2019] K. Bull, I. Martin, and Z. Papić, Systematic construction of scarred many-body dynamics in 1D lattice models, Phys. Rev. Lett. 123, 030601 (2019).
  * Khemani  _et al._ [2019b] V. Khemani, C. R. Laumann, and A. Chandran, Signatures of integrability in the dynamics of Rydberg-blockaded chains, Phys. Rev. B 99, 161101 (2019b).
  * Mukherjee  _et al._ [2020] B. Mukherjee, A. Sen, D. Sen, and K. Sengupta, Dynamics of the vacuum state in a periodically driven Rydberg chain, Phys. Rev. B 102, 075123 (2020).
  * Surace  _et al._ [2020] F. M. Surace, P. P. Mazza, G. Giudici, A. Lerose, A. Gambassi, and M. Dalmonte, Lattice gauge theories and string dynamics in Rydberg atom quantum simulators, Phys. Rev. X 10, 021041 (2020).
  * Zache  _et al._ [2019] T. V. Zache, N. Mueller, J. T. Schneider, F. Jendrzejewski, J. Berges, and P. Hauke, Dynamical topological transitions in the massive schwinger model with a $\theta$ term, Phys. Rev. Lett. 122, 050403 (2019).
  * Zache  _et al._ [2022] T. V. Zache, M. Van Damme, J. C. Halimeh, P. Hauke, and D. Banerjee, Toward the continuum limit of a $(1+1)\mathrm{D}$ quantum link schwinger model, Phys. Rev. D 106, L091502 (2022).
  * Desaules  _et al._ [2023] J.-Y. Desaules, A. Hudomal, D. Banerjee, A. Sen, Z. Papić, and J. C. Halimeh, Prominent quantum many-body scars in a truncated Schwinger model, Phys. Rev. B 107, 205112 (2023).
  * Wang  _et al._ [2023] H.-Y. Wang, W.-Y. Zhang, Z. Yao, Y. Liu, Z.-H. Zhu, Y.-G. Zheng, X.-K. Wang, H. Zhai, Z.-S. Yuan, and J.-W. Pan, Interrelated thermalization and quantum criticality in a lattice gauge simulator, Phys. Rev. Lett. 131, 050401 (2023).
  * [50] W.-Y. Zhang, Y. Liu, Y. Cheng, M.-G. He, H.-Y. Wang, T.-Y. Wang, Z.-H. Zhu, G.-X. Su, Z.-Y. Zhou, Y.-G. Zheng, H. Sun, B. Yang, P. Hauke, W. Zheng, J. C. Halimeh, Z.-S. Yuan, and J.-W. Pan, Observation of microscopic confinement dynamics by a tunable topological $\theta$-angle, arXiv:2306.11794 (2023) .
  * Vidal [2007] G. Vidal, Classical simulation of infinite-size quantum lattice systems in one spatial dimension, Phys. Rev. Lett. 98, 070201 (2007).
  * Bi and Senthil [2019] Z. Bi and T. Senthil, Adventure in topological phase transitions in $3+1$-D: Non-Abelian deconfined quantum criticalities and a possible duality, Phys. Rev. X 9, 021034 (2019).
  * Prakash  _et al._ [2023] A. Prakash, M. Fava, and S. A. Parameswaran, Multiversality and unnecessary criticality in one dimension, Phys. Rev. Lett. 130, 256401 (2023).
  * Buyens  _et al._ [2016] B. Buyens, J. Haegeman, H. Verschelde, F. Verstraete, and K. Van Acoleyen, Confinement and string breaking for ${\mathrm{QED}}_{2}$ in the hamiltonian picture, Phys. Rev. X 6, 041040 (2016).
  * Coleman [1976] S. Coleman, More about the massive schwinger model, Annals of Physics 101, 239 (1976).
  * Halimeh  _et al._ [2022] J. C. Halimeh, I. P. McCulloch, B. Yang, and P. Hauke, Tuning the topological $\theta$-angle in cold-atom quantum simulators of gauge theories, PRX Quantum 3, 040316 (2022).
  * Schrieffer and Wolff [1966] J. R. Schrieffer and P. A. Wolff, Relation between the Anderson and Kondo Hamiltonians, Phys. Rev. 149, 491 (1966).
  * Mukherjee  _et al._ [2022] B. Mukherjee, A. Sen, and K. Sengupta, Periodically driven Rydberg chains with staggered detuning, Phys. Rev. B 106, 064305 (2022).
  * [59] R. Shen, F. Qin, J.-Y. Desaules, Z. Papić, and C. H. Lee, Enhanced many-body quantum scars from the non-Hermitian Fock skin effect, arXiv:2403.02395 (2024) .
  * Iadecola  _et al._ [2019] T. Iadecola, M. Schecter, and S. Xu, Quantum many-body scars from magnon condensation, Phys. Rev. B 100, 184312 (2019).
  * Bull  _et al._ [2020] K. Bull, J.-Y. Desaules, and Z. Papić, Quantum scars as embeddings of weakly broken Lie algebra representations, Phys. Rev. B 101, 165139 (2020).
  * Chen and Iadecola [2021] I.-C. Chen and T. Iadecola, Emergent symmetries and slow quantum dynamics in a Rydberg-atom chain with confinement, Phys. Rev. B 103, 214304 (2021).
  * Sala  _et al._ [2020] P. Sala, T. Rakovszky, R. Verresen, M. Knap, and F. Pollmann, Ergodicity breaking arising from hilbert space fragmentation in dipole-conserving hamiltonians, Phys. Rev. X 10, 011047 (2020).
  * Yang  _et al._ [2020] Z.-C. Yang, F. Liu, A. V. Gorshkov, and T. Iadecola, Hilbert-space fragmentation from strict confinement, Phys. Rev. Lett. 124, 207602 (2020).