Probing the edge between integrability and quantum chaos
in interacting few-atom systems

arXiv:2104.12934v2 [quant-ph] 24 Jun 2021

Thomás Fogarty1 , Miguel Ángel Garcı́a-March2 , Lea F. Santos3 , and N.L. Harshman4
1

Quantum Systems Unit, Okinawa Institute of Science and Technology Graduate University, Onna, Okinawa 904-0495, Japan

2

Instituto Universitario de Matemática Pura y Aplicada, Universitat Politècnica de València, E-46022 València, Spain

3

Department of Physics, Yeshiva University, New York, New York 10016, USA

4

Department of Physics, American University, 4400 Massachusetts Ave. NW, Washington, DC 20016, USA

Interacting quantum systems in the
chaotic domain are at the core of various ongoing studies of many-body physics,
ranging from the scrambling of quantum
information to the onset of thermalization.
We propose a minimum model
for chaos that can be experimentally realized with cold atoms trapped in onedimensional multi-well potentials. We explore the emergence of chaos as the number of particles is increased, starting withas few as two, and as the number of wells isincreased, ranging from a double well to a
multi-well Kronig-Penney-like system. In
this way, we illuminate the narrow boundary between integrability and chaos in ahighly tunable few-body system. We showthat the competition between the particle
interactions and the periodic structure of
the confining potential reveals subtle indications of quantum chaos for 3 particles, while for 4 particles stronger signatures are seen. The analysis is performed
for bosonic particles and could also be extended to distinguishable fermions.

1 Introduction
The interest in quantum chaos, especially whencaused by the interactions between particles,
has grown significantly in the last few years
due to its relationship with several questionsof current experimental and theoretical researchthat arise in atomic, molecular, optical, condensed matter, and high energy physics, as wellas in quantum information science. In interacting many-body quantum systems, quantum
chaos ensures thermalization and the validity ofthe eigenstate thermalization hypothesis (ETH)[17, 28, 114], hinders localization [74, 79, 98],
facilitates the phenomenon of many-body quantum scarring [8, 107], leads to diffusive transport [12, 55] and causes the fast spread of quantum information [75, 91].
Quantum chaos refers to properties of the spectrum and eigenstates that appear in the quantum domain when the classical counterpart ofthe system is chaotic in the sense of mixing andpositive Lyapunov exponent. The features are
similar to what we find in random matrix theory [76], namely the eigenvalues are strongly correlated [44] and the eigenstates in the meanfield basis are close to random vectors [114].
This quantum-classical correspondence is well established for systems with few degrees of freedom [103], such as billiards, the kicked rotor, andthe Dicke model, where the source of chaos is respectively the shape of the billiard, the strengthof the kicks, and the collective interaction between light and matter. In the case of systems
with many interacting particles, the semiclassical analysis is challenging and sometimes not welldefined, so the common approach has been to refer to many-body quantum systems that present
the above mentioned properties of spectrum andeigenstates as chaotic, even when the classicallimit is not analyzed.
The purpose of this work is to identify a minimum model of interacting particles that is chaotic
and that can be experimentally studied with cold
atoms. There are theoretical examples in theliterature of quantum systems with only 3 or 4interacting particles that already exhibit chaotic
properties. They include the cesium atom, whichhas 4 valence electrons [33]; systems composedof 4 particles of unequal masses in a harmonic

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

1

trap [49] and 3 particles with unequal masses
on a ring [54]; 3 or 4 excitations in spin-1/2
chains with short-range [94] or long-range couplings [115]; and even spin-1/2 chains with only
3 sites [77]. In the context of thermalization due
to chaos, we also find works that obtained the
Fermi-Dirac distribution in systems with only 4particles [31, 32, 34, 57, 97].
We consider a one-dimensional system with N
identical particles that is split into wells separated by delta-function barriers. A single barrier defines a double-well system and many barriers results in the finite Kronig–Penney model[65, 89]. We focus on the sector of states symmetric under particle exchange and spatial parity,and that the particles interact via contact interactions which are modeled with a delta function
as in the Lieb-Liniger Hamiltonian [72, 81]. The
system is integrable when there are finite barriers or finite interactions (not both). It is alsosolvable in the limiting cases of infinite barrier
strength and infinite interaction strength. However, we provide numerical evidence that when
the interaction strength and the barrier strength
are simultaneously finite, integrability is broken. In this case, strong signatures of quantum
chaos emerge for N = 4 particles in the presence of just one barrier (double-well system). We
also demonstrate that the signatures of quantum
chaos get enhanced as we increase the number ofbarriers, in which case strong level repulsion isverified for as few as N = 3 particles. Our analysis is done for bosons, but can also be extended tosystems with a small number of distinguishable
fermions.

2 Experimental realization
The experimental realization of our model can bedone with a controllable number of interacting
atoms trapped in one or several one-dimensional
traps. Like experiments based on atoms in optical lattices [14, 70], such systems allow precisioncontrol and find potential applications in quantum engineering and quantum technologies [102].
They are also suitable test beds for addressing
the question of the transition from few- to manybody systems [15, 16, 102]. Our choice of the
Kronig-Penney potential is also motivated by recent experiments that have achieved coherent optical lattices with a sub-wavelength spatial struc-ture that consists of ultra-narrow barriers [109].
This technique can create sharp-box potentials
and the heights of individual barriers can be varied, which adds a further tunable parameter tothe system.
Our results apply to few-body atomic systems with identical particles that possess spatial wave functions symmetric under particle exchange. Wave functions with this spatial symmetry are clearly relevant for the description ofbosons, but our results are equally valid for distinguishable fermions, i.e. identical fermions withinternal degrees of freedom like spin. For a system of distinguishable fermions, the antisymmetry required by particle statistics can be carriedby the spin or internal wave function. All permitted symmetries for three fermions or bosons (distinguishable or indistinguishable) are discussed
in [36, 46, 50]; for more particles in [47]; for few
particles in double or few wells in [37, 48].
Due to losses via three-body re-combination,
experiments with a few bosons trapped in the ultracold regime have an additional difficulty whencompared with those with a few fermions. However, for the order of tens of bosons, it was shownin [18] that one can successfully load dipole traps
by means of evaporative cooling. Smaller numberof atoms can be loaded in optical lattices [112], inarrays of double wells [25], or in a two-site opticalring, which can be appropriately reshaped into a
Gaussian trap [52]. An alternative experimental route leading to cooled atoms trapped in several wells is that of few atoms in optical tweezers,
which for a single atom in the ground state wasaccomplished in [60]. In subsequent papers it was
experimentally demonstrated the trapping of two
87 Rb bosonic atoms in two wells [61] or in uniformly filled arrays of traps [68]. Recent advances
allow the laser cooling of atoms in optical tweezers [88]; the trapping of individual atoms in optical tweezer arrays [93]; and even the loading ofatoms one by one [4, 30, 99] in a one-dimensional
array [30].
Another experimental context for the resultsin this paper is the ground breaking experimentsthat showed the accessibility and versatility ofsystems with a very small number of interacting
fermions. Reference [100] demonstrated that adeterministic number of ultracold fermions could
be extracted from a larger ensemble by applying a tightly confined one-dimensional dimple po-

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

2

tential. In this experiment, a few 6 Li atoms
in the two lowest-energy Zeeman substates weretrapped. The strength of the interactions between atoms in different spin states could becontrolled via Feshbach resonance. In a subsequent experiment by the same group, they considered one atom of one species (an impurity)interacting with an increasing number of identical fermions, being able to build a small Fermi
sea adding fermions one by one in a controllable
manner [110]. One can have more than two components in a few fermion system, as in the experiment reported in [84], where a one dimensional
system with a tunable number of spin components was realized.
In addition to the controllable number ofatoms, another ingredient required for the realization of our model is the ability to change
the trapping potential, creating double, triple orgenerally multiwell potentials. The same groupthat realized few trapped fermions in [100] wasalso able to trap few fermions in one dimensional
double-well [7, 78] and multi-well systems [6].

3 ModelIn the simplest case of equally-spaced barriers,
the Hamiltonian describing our system takes theform
1
H(N, W, τ, γ) = T N + τ V N,W + γU N
1

(1a)

where
TN = −
V N,W =

N
X
∂2

∂x2i
i=1

N W
−1
X
X

δ(xi − πk/W )

(1b)
(1c)

i=1 k=1

UN =

X
hiji

δ(xi − xj ).

(1d)interactions, as it is the most common scenario
for ultracold bosons.
The goal of our analysis of the Hamiltonian
H(N, W, τ, γ) is to understand how the signatures of quantum chaos scale with the numberof particles N , the number of wells W , the barrier strength τ and the interaction strength γ.
The particularly simple form (1) of the Hamiltonian means that it is amenable to analytic andnumeric calculations.

3.1 Solvable and integrable limiting cases
For a fixed N and W , consider the parameter
space window τ ∈ [0, ∞) and γ ∈ [0, ∞), depicted in Fig. 1. An interesting feature of the
Hamiltonian H(N, W, τ, γ) is that the model has
exact solutions for four special limiting cases:
1. H(N, W, 0, 0) = T N : N non-interacting
identical particles in a one-dimensional infinite square well with width L.
2. H(N, W, 0, ∞): N hard-core identical particles in a one-dimensional infinite square wellwith width L.
3. H(N, W, ∞, 0): N non-interacting identical particles distributed in W identical onedimensional infinite square wells with width
L/W .
4. H(N, W, ∞, ∞):N hard-core identicalparticles distributed in W identical onedimensional infinite square wells with width
L/W .
In all of these four cases, the configuration space
is sectioned into one or more N -dimensional polytopes with high symmetry. Solving for the spectrum is equivalent to solving the Schrödinger
equation with Dirichlet boundary conditions.
Exact solutions for the Schrödinger equation inthese polytopes can be constructed from symmetrized combinations of one-particle states using methods of Refs. [59, 64, 82, 108].
Beyond these four exactly solvable special
cases, the Hamiltonian H(N, W, τ, γ) is also integrable along the four edges of the (τ, γ) parameter space window (cf. Fig. 1). Referring to
the four ‘corner’ models denoted above, the integrable limits are:This model realizes a system with N identicalinteracting particles of mass m trapped in a onedimensional box of length L that is disrupted
by W − 1 delta-barriers. In the equation above,xi ∈ [0, π] are the positions of the particles scaledby the length L/π, the energy scale is providedby 1 = ~2 π 2 /(2mL2 ) (henceforth set to unity),and τ ≥ 0 and γ ≥ 0 are unitless parametersdescribing the barrier strength, and interaction 1 ↔ 2 H(N, W, 0, γ): Without barriers, the Hamiltonian (1) for N identical particles in anstrength, respectively. Here we consider repulsive
Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

3

integrable

(0,∞)

1

(0,0)

integrable
CHAOS?
(t,g)

integrable

4

(∞,∞)

integrable2

be different partitions of interacting particles and the size of the well is L/W .
Note that these four integrable models allow us
to establish an adiabatic map between the energy levels for the solvable (and superintegrable)models at the corners of parameter space.

3

(∞,0)

Figure 1: Depiction of the (compactified) (τ, γ) parameter space window τ ∈ [0, ∞) and γ ∈ [0, ∞). At the four
corners of the window, the Hamiltonian H(N, W, τ, γ)
has exact solutions and is superintegrable. The fouredges connecting these corners are integrable models
where exact solutions can be found by the solution of
coupled transcendental equations. Except for possibly
the case of N = 2, the H(N, W, τ, γ) does not appear
integrable or solvable for arbitrary (τ, γ) in this parameter space

infinite square well with delta-function interactions is solvable by coordinate Bethe
ansatz [39, 40, 80]. The energies are the solution of coupled transcendental equations
that depend on N and γ.
1 ↔ 3 H(N, W, τ, 0): N non-interacting particles
in an infinite square well with W −1 barriers.
The Hamiltonian separates into N identical one-dimensional sub-Hamiltonians. The
spectrum for each sub-Hamiltonian is obtained via solution to transcendental equations that depend on τ [89].
2 ↔ 4 H(N, W, τ, ∞): In this limit of finite wellswith hard-core interactions, the solutions
are Tonks-Girardeau constructions derived
from the symmetrized Slater determinantsof the solutions of the non-interacting case
H(N, W, τ, 0) [41].
3 ↔ 4 H(N, W, ∞, γ): With infinite barriers, theconfiguration space fractures into W N configurations, i.e. each particle i ∈ {1, . . . , N }is in a well j ∈ {1, . . . , W } with width L/W .
Each of these particle-well configurations is
an independent system, since with τ → ∞there is no tunneling among wells. Withineach configuration, the solutions are similar
to the coordinate Bethe ansatz solutions ofH(N, W, 0, γ), except that now there may

3.2 Symmetries and degeneracies
To analyze signatures of chaos, the Hilbert space
H must first be decomposed into subspaces with
fixed symmetry. The Hamiltonian (1) has two
symmetries for any interaction strength and barrier strength in the (τ, γ) parameter window.
First, it is symmetric under particle permutations SN , i.e. any permutation p ∈ SN is represented as a linear transformation in configuration
space (x1 , x2 , . . . , xN ) → (xp1 , xp2 , . . . , xpN ) thatleaves the Hamiltonian (1) invariant. The second symmetry is total parity inversion Π about
the center of the well, implemented as the affine
linear transformation xi → π − xi for all i.
These two symmetries allow eigenstates to beclassified by irreducible representations of SN
and of parity Π and reduce the total Hilbert
space H into subspaces with a given symmetry(cf. [47]). In this work, we focus on the sector of
Hilbert space H[N ]+ ⊂ H containing states withbosonic symmetry under particle exchange andpositive parity. Note that a special property ofdelta-function operators V N,W and U N is thatthey vanish on certain subspaces of H, because
the support of the delta-functions coincides withnodal lines of symmetrized eigenstates of T N . In
particular, there is a subspace of H[N ]+ uponwhich V N,W vanishes when N is even or both N
and W are odd. Also, in the limits τ → ∞ andγ → ∞, both V N,W and U N must vanish on any
states with finite energy; in this limit the wave
functions must have nodal surfaces that coincidewith the support of these operators.
Degeneracies in the spectrum either originate
with the symmetries of the Hamiltonian or they
are designated ‘accidental’. Because the only
symmetries for generic (τ, γ) are particle exchanges SN and parity Π, there should be nodegeneracies originating in symmetry in the interior of parameter space. This is because theirreducible representations [N ]+ for totally symmetric bosonic states are one-dimensional. However, additional symmetries arise in the limiting

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

4cases of zero and infinite strengths τ or γ, including the symmetries of separability and systemdecoupling due to infinite barriers [48]. Thesesymmetries explain the integrability of the edge
models and the superintegrability of the cornermodels. When τ and γ are both finite, then these
additional symmetries are broken and we expect
these energy levels that cross at certain parameters along the edges of the model space to repel.
Level repulsion, which is a signature of chaos, isthus expected in the interior of this model space(see Fig. 1), which is indeed what we verify forN ≥ 3.
Besides degeneracies originating in symmetries, several other kinds of accidental degeneracies should be possible in our system. The
variation of two parameters is sufficient to allowfor ‘diabolic points,’ topologically-stable conic
degeneracies between two energy levels [9, 11].
These degeneracies can be distinguished from‘near misses’ by looking at how the wave function transforms when small loops in control space
are taken and in subsequent work we will consider such loops. There are also degeneracies
for the solvable corner models arising from thenumber theory of decomposing integers into sumsof squares of integers. Also called Pythagorean
degeneracies [48, 101], the density of such degeneracies grows slowly with energy but eventually comes to dominate the spectrum [9]. Finally, there are other degeneracies characteristicof Bethe-ansatz solvable systems [53] that should
be relevant for some of the integrable edge models. It is not clear how these additional ‘accidental’ degeneracies at the corners and edges affectthe spectrum of the interior of parameter space,but they inform our interpretation of the levelstatistics presented below for the N = 2 case.

3.3 Density of states
The density of states of the Hamiltonian
H(N, W, τ, γ) is used to interpret some of the numerical results below. Perhaps surprisingly, it isindependent of τ and γ to leading order in the energy E. The density of states in the sector H[N ]+isρ[N ]+ (E) =1
π N/2
E N/2−1
2N +2 (N + 1)! Γ(N/2 + 1)
+O [N, W, τ, γ] (E N/2−3/2 ),(2)

Figure 2: Density of states for N = 2 (blue points),
N = 3 (red points) and N = 4 (yellow points). Solidlines show the leading term in Eq. (2) for the respective
particle numbers. Symbols denote different system parameters: γ = 10 and τ = 0 (circles), γ = τ = 10 and
W = 2 (triangles) and γ = τ = 10 and W = 10 (stars).

where we have used the notation O [N, W, τ, γ]
to indicate that generally, the coefficient in frontof the subleading term proportional to E N/2−3/2
will depend on all the parameters of the Hamiltonian.
We derive Eq. (2) in the Appendix for the foursolvable ‘corner’ models, and infer that it holds
for the entire parameter window. More generally, adding one-dimensional delta-functions toan otherwise free problem should not change thedensity of states. By Weyl’s Law [56], the volume of phase space can be related to the density of states, and boundary conditions (whether
Dirichlet, Neumann or hybrid) do not change thephase space volume to leading order.
In Fig. 2 we show that the leading term in
Eq. (2) agrees with the scaling of the density ofstates regardless of number of wells W , barrier
strength τ or interaction strength γ. Because
the leading term in Eq. (2) is the same for theentire parameter window for a given N , at least
heuristically we expect that the variation in levelstatistics with the parameters (τ, γ) is governedby the properties of the subleading term (or evenlower order terms). In the Appendix, we support
this hypothesis by looking at how energy levels
change as the model parameters are adiabaticallytuned along integrable edge models to connect
the spectra of the solvable corner models.
Note that for N = 2, the density of states is

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

5

a constant and the subleading correction is actually decreasing with energy as E −1/2 . In Fig. 2,
all data points for N = 2 lie on top of the line
given by the leading term in Eq. (2). For N = 3,
the density of states increases as E 1/2 and thesubleading term is constant. Only for N = 4 areboth the leading and the subleading term of thedensity of states growing with energy. If the subleading term is important for understanding thedensity of level crossings as we propose, then thishelps to explain why N = 4 is the threshold when
the level statistics conform to the expectations ofrandom matrix theory across such a wide rangeof parameter space, whereas the N = 3 only in alimited range.

4 Indicators of chaos
In the parameter window of interest, finite interactions and finite barrier heights, the Hamiltonian (1) is diagonalized numerically using a basis consisting of N -particle non-interacting eigen+states of T N in the Hilbert space H[N ] [46]. In
the following we discuss the different indicatorsof chaos that are used in this work along with thenumerical results.

4.1 Energy spacing distribution
To study the degree of short-range correlations
between the eigenvalues, we use the distribution
P (s) of the spacings between neighboring levels obtained after unfolding the spectrum. In
generic integrable models, the level spacing distribution is Poissonian, PP (s) = e−s , when theenergy levels are uncorrelated and not prohibitedfrom crossing, although different shapes emerge
for “picket-fence”-kind of spectra [10, 38, 85]
and systems with an excessive number of degeneracies [113]. For chaotic systems with realand symmetric Hamiltonian matrices, as the one
considered here, the level spacing distribution
follows the Wigner-Dyson distribution[44, 76],
2
PWD (s) = (πs/2) exp −πs /4 , which indicatesthat the eigenvalues are correlated and repel eachother.
First we investigate the minimal case for chaos
in our system, where only one barrier is inserted
centrally in the square well, thereby creating adouble well potential (W = 2). In Fig. 3 weshow the level spacing distributions for 2, 3 andFigure 3: Energy level spacing distributions for (a,d,g)
N = 2 particles, (b,e,h) N = 3 particles and (c,f,i)
N = 4 particles with different interaction strengths γ.
The barrier height is fixed to τ = 10 and the number
of wells is W = 2. The red solid line is the Poissonian
distribution and the blue solid line is the Wigner-Dyson
distribution.

4 particles for different interactions and at a fixedbarrier height of τ = 10.
For weak interactions, γ = 1, a picket
fence pattern is noticeable for N = 2 particles[Fig. 3 (a)], which indicates non-generic correlations in the energy spectrum [10]. Additionally,
the peak at s = 0 signals excessive degeneracies
typical for solvable systems [53, 113]. For N = 3,
the picket fence structure vanishes and the distribution is closer to Poissonian with some remaining evidence of additional degeneracies at s = 0[Fig. 3 (b)]. For N = 4, the distribution is also
close to Poissonian, except for a slight decreaseat s = 0 that provides evidence for level repulsion
already at weak interactions [Fig. 3 (c)].
At stronger interactions (γ ≥ 10), the distributions for N = 2 particles in Fig. 3 (d) and
Fig. 3 (g) become close to Poissonian and level
repulsion is not evident. This Poissonian distribution suggests that for the N = 2 and W = 2
the Hamiltonian is integrable (or effectively so)for all values of (τ, γ); see the conclusion for a further discussion on this point. In contrast, somedegree of level repulsion is already noticeable for
N = 3 [Fig. 3 (e) and Fig. 3 (h)] and the WignerDyson distribution is visible for N = 4 particles[Fig. 3 (f) and Fig. 3 (i)].
To explore the crossover of the energy spectrum from Poissonian to Wigner-Dyson and therole of interactions and the potential barrier, we

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

6

Figure 4: Brody distribution parameter β as a functionof τ and γ for (a,d) N = 2, (b,e) N = 3, and (c,f)
N = 4 particles. The number of wells are (a-c) W = 2
and (d-f) W = 10. For W = 2 wells the maximum
βmax = 0.8269 is found for N = 4 particles at γ = 20
and τ = 12.5, while for W = 10, we find βmax = 0.9512
for N = 4 particles at γ = 45 and τ = 12.5. Black areas
show β < 0 indicating degeneracies and picket-fence
spectra.

fit our numerical results to the Brody distribution [22] (see an alternative in [58]),
Pβ (s) = (β + 1)bsβ exp(−bsβ+1 ),
 

b= Γ

β+2
β+1

(3)

β+1

.

For the Wigner-Dyson distribution, β ∼ 1, while
the Poissonian distribution leads to β ∼ 0.
Close to the integrability limits, vanishing interactions (γ ≈ 0) with finite barriers (τ 6= 0) or
finite interactions (γ 6= 0) with vanishing barrierheight (τ ≈ 0), the energy spacing distributions
can display a picket-fence or quasi-Poissonian distribution as discussed in Fig. 3 for N = 2, which
results in β < 0. We highlight these non-genericregions as black in Figs. 4 (a)-(f), showing thatthey are more prevalent at smaller particle numbers, while their footprint in the parameter space
almost vanishes for N = 4. As we increase
the interactions (γ  0) and the barrier heights(τ  0), the degeneracies of the N = 2 system
are destroyed and the level spacing distributions
become closer to Poissonian for a large regionof the parameter space with β remaining below
0.35.
For more particles, the distributions are different, with indications of energy level repulsion

emerging for N = 3 in the region 15 . γ . 45
and 5 . τ . 35 where β ∼ 0.5 [Fig. 4 (b)]. For
N = 4 the fitting parameter β in almost theentire parameter space is larger than 0.5, withthe areas of integrability confined to the edges ofthe (τ, γ) parameter space [Fig. 4 (c)]. In fact,
we find the maximum to be around β ∼ 0.83 at
γ = 20 and τ = 12.5, indicating that the energy
spacing distribution gets close to Wigner-Dyson.
In Figs. 4 (d)-(f) we also consider multiple
wells, with W = 10. For more wells, the densityof states does not change to leading order in τ ,but a greater proportion of the eigenstates of T N
feel the effect of the barriers more acutely. For
example, for N = 2 and W = 2, the barrier potential V N,W vanishes on half of the eigenstatesof T N in H[N ]+ , whereas for N = 2 and W = 10,
the barrier potential V N,W only vanishes on one
tenth of the eigenstates of T N in H[N ]+ .
Comparing the case of W = 10 to W = 2
in Fig. 4, we see that larger values of β are indeed found for W = 10 and the region of theparameter space where β is large has increased
significantly for both N = 3 and N = 4. Thissuggests that in systems with a large number ofwells, the energy level repulsion is enhanced andchaos could be observed for systems with as few
as N = 3 particles. However, for N = 2, we finda maximum of β ≈ 0.5 and this value does not
approach the chaotic limit by increasing the number of wells. As before, the case of only two particles resists the transition to chaos. In Sect. 4.4
below, we discuss in more detail how the signatures of chaos change as the number of barriers
is increased.

4.2 Off-diagonal ETH
Quantum chaos ensures the validity of the ETH,
so we can also use the indicators of ETH to detect
the transition to chaos. Two conditions need tobe satisfied for a few-body observable O, evolving
according toO(t) =X∗Cm
Cn e−i(En −Em )t Omn

m6=n

+

X

|Cm |2 Omm ,(4)mto reach thermal equilibrium. In the equation
above, Cm = hm|Ψ(0)i is the overlap between
the eigenstate |mi of the Hamiltonian that describes the system and the initial state |Ψ(0)i,

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

7

viation of the distribution. For Gaussian distributions the kurtosis is given by KÔ = 3.

Figure 5: Inverse of the kurtosis KT−1 as a function of τand γ for (a,d) N = 2 (b,e) N = 3 and (c,f) N = 4
particles. The number of wells are (a-c) W = 2 and (df) W = 10. For W = 2 the minimum kurtosis is 4.3592
at γ = 30 and τ = 15, while for W = 10 the minimum
kurtosis is 3.3274 at γ = 20 and τ = 7.5.

and Omn = hm|Ô|ni. The first condition isthat of equilibration, which depends on the firstterm on the right-hand side of Eq. (4). At longtimes, due to the lack of degeneracies of chaotic
Hamiltonians, the small values of the coefficients
Cm ’s obtained when the systems is quenched far
from equilibrium, and the small values of the offdiagonal elements Omn ’s caused by the chaotic
eigenstates, the first term in Eq. (4) leads tosmall fluctuations that decrease with system sizeand cancel out on average. So apart from small
fluctuations, the observable reaches its infiniteP
time average m |Cm |2 Omm . The second condition is that this infinite-time average approaches
the thermodynamic average as the system size increases, confirming that the equilibrium is indeed
thermal. These two steps are usually referred toas off-diagonal- and diagonal-ETH, respectively.
Here, we consider the off-diagonal-ETH to detect the transition to chaos. The distribution of
Omn in chaotic (thermalizing) systems is Gaussian [13, 19–21, 66, 92], reflecting the chaotic
structure of the eigenstates, while other forms
emerge for integrable models [66]. To quantify
how close the distribution of Omn is to a Gaussian, we use the kurtosis,

In our calculations of KÔ , we consider the kinetic energy operator T N (in the following wedrop the superscript to simplify the notation).
To ensure no effects from the bottom edge of thespectrum, we choose a energy window far fromthe ground state and include only states within
Em ∈ [Emid − ∆E, Emid + ∆E], where the center of the energy window is high in the spectrum
at Emid = 700 and its width is ∆E = 100. In
Figs. 5 (a)-(f) we show the inverse of the kurtosis, 1/KT , with the maximal value 1/3 indicating
a Gaussian distribution and therefore the presence of chaos.
For N = 4 particles and W = 2 [Fig. 5 (c)],the distribution is close to Gaussian for barrierheights 1 . τ . 30 and interactions 20 . γ . 60,
which coincides with the region of β > 0.8 in
Fig. 4 (c). We choose the minimum value ofthe kurtosis in Fig. 5 (c) and show the Gaussian probability distribution of Tmn in Fig. 6 (c).
In contrast, for lower particle numbers, N = 2, 3,
the distribution of the off-diagonal elements Tmn
is sharply peaked, as seen in Fig. 6 (a) and
Fig. 6 (b), which is consistent with the kurtosis
having values KT  3 in Fig. 5 (a) and Fig. 5 (b).

(5)

The kurtosis shows a similar enhancement due
to the presence of more barriers, attaining a minimum of KT ≈ 3.3 for N = 4 particles in W = 10
wells [Fig. 5 (f)]. Indeed, for 1 . τ < 20 and abroad range of interactions, the kurtosis is close
to 3. Increasing the number of barriers also
moves the band of minimal kurtosis to lower values of τ , as lower barrier heights are necessary
to retain the competition with the inter-particle
interactions. In a similar region of the parameter space, there is also a visible minimum of the
kurtosis for N = 3 particles [Fig. 5 (e)]. In fact,
when taking the interaction and barrier height
which give the minimum value of kurtosis [γ = 10
and τ = 7.5 for N = 4], there are distinct Gaussian probability distributions for both N = 4 andN = 3 particles [see Figs. 6 (f,g)]. However, for
N = 2 the off-diagonal elements of the kinetic energy operator Tmn do not indicate the presenceof chaos in either the kurtosis [Fig. 5 (d)] or the
probability distribution [Fig. 6 (e)].

where h.i indicates the average over all pairs of
eigenstates |mi =6 |ni and σ is the standard de-To quantify how the off-diagonal elements Tmn
behave as a function of the energy difference ω =KÔ =h(Omn − hOmn i)4 i,
σ4Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

8

Figure 6: Off-diagonal elements of the kinetic energy operator Tmn for (a-d) W = 2 wells and (e-h) W = 10 wells.
The probability distributions in the respective potentials are shown for (a,e) N = 2 (b,f) N = 3 and (c,g) N = 4.
(d,h) The ratio ΓT as a function of ω = |Em − En | in an energy window at Emid = 700 and width ∆E = 100. Theinteraction strengths and barrier heights are chosen at the point of minimum kurtosis for 4 particles in Fig.5(c,f),namely γ = 30 and τ = 15 for W = 2, and γ = 20 and τ = 7.5 for W = 10.

|Em − En | we also show in Fig. 6 the ratioΓT =|Tmn |2

2,(6)|Tmn |

which is equal to π/2 for a Gaussian distribution.
In Fig. 6 (d) and Fig. 6 (h) this ratio is shown
for W = 2 and W = 10 wells, respectively, forthe same Hamiltonian parameters that gave thelowest values of the kurtosis discussed in the preceding paragraph. For W = 2 wells, the valuesof ΓT for N = 4 sit very close to π/2 over a large
range of ω, confirming that the Gaussianity of the
distribution is preserved at different energy spacings. The ratio for N = 2, 3 has a large variancewith the majority of the points being far from
π/2, which is indicative of the peaked distributions in Fig. 6 (a) and Fig. 6 (b). For W = 10
wells, we find that not only ΓT ≈ π/2 for N = 4
particles, but also the N = 3 system approaches
this result. This suggests that the N = 3 system
can be tuned between integrability and chaos bychanging the trapping potential.

“ramp” [27]. The correlation hole correspondsP
to a dip below hS∞ i = h n |Cn |4 i, which is
the infinite-time average (saturation value) ofhSP (t)i. This dip emerges also in the spectralPform factor h n,m e−i(En −Em )t i, but contrary tothis one, the survival probability is a true dynamical quantity. In the equation above, h.i indicatesaverages. The survival probability is non-selfaveraging [87, 96], so the correlation hole is not
visible unless averages are performed. They canbe done over initial states, disorder realizations,or, as in our case, they correspond to moving time
averages. The correlation hole detects short- andlong-range correlations in the spectrum, and in
addition, it does not require unfolding the spectrum or separating it by symmetries [29, 92].
In cold atom systems the survival probability iscommonly used to probe the non-equilibrium dynamics of few- [23, 62] and many-body systems[26, 35, 42, 63], and can be experimentally measured using interferometric techniques [24].
We take an initial state that has a homogeneous probability distribution centered at Emid
inside an energy window of width ∆E

4.3 Survival Probability(

Spectral correlations get manifested also in theevolution of the survival probability,2hSP (t)i = h|hΨ(0)|Ψ(t)i| i=hX|Cn |2 |Cm |2 e−i(En −Em )t i, (7)n,min the form of what is known as the correlation hole [2, 29, 43, 51, 67, 69, 71, 73, 86, 95,
104–106, 111], recently referred to also as theρ(E) =

1
2∆E

0

for E ∈ [Emid − ∆E, Emid + ∆E]otherwise.

In [67, 95, 106] a general analytic solution wasderived for the survival probability for chaotic
systems. For the square distribution used here,
this solution is given by"

1 − S∞ sin2 (∆E t)∆E t
SP (t) =η− b2
η−1(∆E t)2πη

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.



#

+S∞ .
(8)
9

Figure 7: Survival probability of (a,d) 2 particles (b,e), 3 particles, and (c,f) 4 particles. The number of wells is (a-c)
W = 2 and (d-f) W = 10. The initial state is a square distribution centered at Emid = 700 and of width ∆E = 100.
Black lines are the moving time averages on the logarithmic timescale with window [log10 t − ∆t, log10 t + ∆t] and∆t = 0.02. Green horizontal dashed line indicates the infinite time average of the survival probability S∞ . For
N = 3 and 4, red lines represent the analytic solution. For N = 4, the grey points represent the unaveraged value
for comparison. The parameters used for each row correspond to those that give the minimum kurtosis for 4 particles
in Fig.5, namely γ = 30 and τ = 15 for W = 2, and γ = 20 and τ = 7.5 for W = 10.

Here, η is the number of energy eigenvectors in
the energy window ∆E. The initial decay of thesurvival probability is captured by the first termin the brackets in Eq. (8), after which the dynamics is described by the two-level form factor

cay and subsequent ramp of SP (t) matching precisely the analytics. For N = 3 and W = 10
[Fig. 7 (e)], partial revivals still obscure the minimum of the correlation hole, but the ramp towards saturation can be seen to follow the analytic results.

b2 (t̄) = 1 − 2t̄ + t̄ ln(2t̄ + 1) Θ(1 − t̄)"+ t̄ ln2t̄ + 12t̄ − 1!

#

− 1 Θ(t̄ − 1),with Θ(t̄) the Heaviside step function.
In Fig. 7 the survival probability is shown fordifferent number of particles and number of wells(grey points). The moving time average (blacklines) smooths the data and allows us to identify
the correlation hole with its ramp towards S∞ .
For W = 2 and N = 2, 3 [Figs. 7 (a,b)] there
is no obvious indication of a dip below S∞ that
could be described by the two-level form factor.
For N = 4 [Fig. 7 (c)], a noticeable dip manifests
below the saturation point and follows closely the
analytic solution in Eq. (8) (red lines). This is indicative of chaotic behavior. In the W = 10 wellsystem the correlation hole for N = 4 is evenmore pronounced [Fig. 7 (f)] with the initial de-

4.4 Dependence on the number of wells
For a more systematic analysis of the onset ofchaos as the number of wells in the system isincreased, we show in Fig. 8 (a) the minimumof the kurtosis Kmin = min[KT ] in the rangeγ, τ ∈ [0, 100]. Here we fix the size of the box
L and focus on N = 3 and N = 4 particles. For
both cases the number of wells dictates the emergence of chaos, but to different degrees. For thetrivial case of W = 1 (no barriers) the kurtosis islarge and both N = 3 and N = 4 are integrable.
For W > 1 the kurtosis of N = 4 takes low values, Kmin ∼ 3, indicating chaos, essentially irrespective of the number of wells when W . 13.
However, for N = 3 the inclusion of more barriers causes a more subtle change to the kurtosis,
which decreases slowly as more barriers are intro-

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

10

duced, attaining a minimum of Kmin ∼ 4 in theregion of 6 . W . 13. It is in this region that
the N = 3 system displays chaotic signatures as
discussed in the previous sections. Interestingly,
increasing the number of wells further (W > 13)
results in an increase of the kurtosis and indications of chaos are diminished. A similar tendency
is seen for N = 4, albeit in a less drastic manner,as the kurtosis increases at a lower rate. For both
N = 3 and N = 4 we expect that as the limit of
N/W → 0 is approached the contest between theordering of the particles in the wells and their interactions is reduced and that the system slowly
returns to what we would see in the continuum:
particles in a box [90].
In Fig. 8 (b) and Fig. 8 (c) we show the optimal interactions and barrier heights for achieving
the minimum kurtosis shown in Fig. 8 (a). Figure 8 (c) shows that the optimal barrier strengths
τ for N = 3 and N = 4 have close agreement
for all W , and that these values decrease withincreasing number of wells (for W ≤ 20). This
reduction in τ is necessary to preserve the competition between barriers and interactions, as when
the number of barriers is increased the impactof τ is magnified. This effect can be seen in
the shift of the chaotic region in Fig. 5 (c) and(f). Fig. 8 (b) shows that the optimal interaction
strengths γ found for both N = 3 and N = 4 converge to a similar value in the region 7 . W . 20,
which encompasses the areas of low kurtosis in
Fig. 8 (a). This suggests that both N = 3 andN = 4 particles are chaotic in the same regions of
the parameter space [τ, γ]. Increasing the number of wells beyond W = 20 breaks this trend andthe parameters for N = 3 and N = 4 diverge as
the indications of chaos are lost.

5 Discussion and Conclusion
To summarize, between N = 2 and N = 4 ourmodel makes a clear transition. For N = 4 particles and just W = 2 wells, there are clear signatures of chaos when both the barrier and theinteractions have finite strength. The density ofstates grows rapidly with energy, and numericalanalysis gives evidence of the highly-correlated
spectrum typical of random matrices, of the validity of the off-diagonal ETH, and of a clear
ramp in the survival probability. As the number of wells is increased, evidence for the onset

Figure 8: (a) Minimum of the kurtosis as a functionof the number of wells for N = 3 (red triangles) andN = 4 (yellow circles). Kurtosis is calculated in anenergy window of width ∆E = 100 around Emid = 700.
(b) Interaction and (c) barrier height corresponding to
the points of minimum kurtosis in (a).

of chaos becomes more robust, even for weakerbarriers and interactions.
In contrast, for N = 2, the spectrum is numerically indistinguishable from an integrable systemthroughout the parameter window for W = 2,
and deviates only slightly from this as the number of wells is increased. Several possible origins
for the non-chaotic behavior of N = 2 can be
hypothesized, including some undiscovered integral of motion, some sort of partial integrabilityof a subset of eigenstates of the Hamiltonian, aproliferation of accidental degeneracies from several sources that create the appearance of integrability, or a combination of all these. Thesepossibilities will be considered in more detail ina subsequent work. For now, we note with someirony that integrability is much less generic andmore complicated than chaos.
The intermediate case of N = 3 stands at the
ragged edge where symmetry and integrability
dissolve into chaos and randomness. By tuning
the interaction and barrier parameters and increasing the number of barriers, the full gamutof possibilities can be realized on the same small
atomic system. The N = 3 case has some common aspects with the N = 2 case and other features similar to N = 4. For example, the densityof states grows with the energy like the N = 4
and unlike the N = 2 case, although it grows
sublinearly while the N = 4 case grows linearly.
On the other hand, the probability distributions
for the off-diagonal elements of the kinetic energy
are closer to the N = 2 than to the N = 4 caseAccepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

11

for only two wells, but when increasing the number of wells, it gets closer to the Gaussian distribution characteristic of the N = 4 case. The
variety of possible scenarios for N = 3 atoms
is most clear in Fig. 8 (a), where the degree ofchaoticity measured with the minimum kurtosis
Kmin makes sweeping changes as the number ofwells is increased.
As we discuss above, current experiments withsubwavelength lattice potentials and deterministically prepared few-body states provide the ideal
platform to probe this boundary between integrability and chaos. Furthermore, these small
atomic systems have been proposed as the working units of larger quantum information processing devices and protocols. Therefore, understanding their information, control, and entanglement properties in these different regimes becomes important.
For example, this model lends itself naturally
to “digitization”. The number of particles in eachwell becomes a useful observable in the infinite
barrier limit; they are integrals of motion in fact.
Similarly, coherent superpositions of eigenstatesof these or other integrals of motion in the limiting edge models could be used for storing quantum information (c.f. [83]). A quench from an
integrable limit to the chaotic parameter regime
would break these integrals of motion and effectively scramble the information held in the initialstate.
An important extension of this work that is
experimentally relevant is determining how sensitive our results are to the idealizations of deltabarriers, precisely symmetric positions and uniform barrier height. We expect that small deviations of the periodic Kronig-Penney lattice,such as finite width barriers, would not significantly alter the chaotic regions of our system.
The versatility of our model also allows us to explore the possibility of emergent integrability [1]when more disorder is introduced into the system via non-regularly spaced barriers or barriersof different heights. Both aspects are inspiring
and we leave them for future research. Anotherinteresting scenario which we have not exploredis the case in which the interactions are attractive, whereby the system is now furnished witha bound state and whose limit at infinite interactions is the so-called super Tonks-Girardeaugas [3, 5, 45].

Acknowledgments
The authors thank M. Olshanii, T. Busch, A.
Fabra and M. Boubakour for insights on integrability and conversations about chaos. TF
acknowledges support from JSPS KAKENHI21K13856 and the Okinawa Institute of Scienceand Technology Graduate University. We are
grateful for the help and support provided by
the Scientific Computing and Data Analysis section of Research Support Division at OIST. LFS
was supported by the NSF grant No. DMR1936006. M.A.G.M. acknowledges funding fromthe Spanish Ministry of Education and Vocational Training (MEFP) through the Beatriz
Galindo program 2018 (BEAGAL18/00203) andSpanish Ministry MINECO (FIDEUA PID2019106901GBI00/10.13039/501100011033).

A Density of states derivation
To derive the density of states, we first calculate
the total number of states (with any symmetryor parity) and energy less than E:
N (E) =

1
π N/2
E N/2
2N Γ(N/2 + 1)

(9)+O [N, W, τ, γ] (E (N −1)/2 ).
The density of states in Eq. (2) is the derivativeof this with respect to energy.
To establish this result (9), first consider thesimplest limiting case H(N, W, 0, 0). There is asolution of H(N, W, 0, 0) for every set of nonnegative integers n = {n1 , . . . , nN } with energy
P2
En = Ni=1 ni . The space of solutions therefore
is a (hyper)cubic lattice in the all-positive ‘quadrant’ (really 2N -rant) of RN . Since each statetakes up a unit volume in this quantum number
space, to find the number of states NN (E) withenergy less than E, one √takes the volume of an
N -ball with radius r = E and divides by 2N
to account for the all-positive condition, giving
the leading term in Eq. (9). The first correction
term comes from the N -sphere boundary of the
N -ball, which has one dimension lower.
The spectrum of H(2, 2, 0, 0) is depicted asmodel 1 in Fig. 9. In this simple case, stateswith energy less than E√ lie within the quartercircle with radius r = E. That quarter disk
therefore has area πE/4, agreeing with (9).

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

12

2

4

(0, )( ,  )1

3(0, 0)(, 0)

Figure 9: These four diagrams depict the spectrum forthe four solvable cases of N = 2 particle and W = 2
wells; c.f. Fig. 1 and Sec. 3.1. Model 1: no interactions
or barriers; model 2: infinite interactions no barriers;model 3 no interactions, infinite barriers; model 4 infinite interactions and barriers. Each dot represents anenergy eigenstate with energy E = n21 + n22 given by thesum of the squares of the integer coordinates (n1 , n2 )of the point. States with n1 ≥ n2 represent the symmetrized states with positive parity in black and with
with negative parity in cyan. States with n1 < n2 represented by magenta and yellow dots are antisymmetricand have negative and positive parity, respectively. Inmodels 1 and 2, all states are two-fold degenerate except
in model 1 (no barriers, no interactions) when n1 = n2 .
Additionally, models 3 and 4 have additional two-folddegeneracies (in model 4, states along the diagonal represented by split dots with twice the area), four-fold (inmodel 4, represented by quartered dots along the diagonal) and eight-fold (in models 3 and 4, represented bypairs of quartered dots with exchanged integer quantum
numbers). The dashed√ quarter-circle is the boundary
E = 150 with radius 150 and is included to aid visualization of the spectral flow.

Parallel arguments give the same leading terms
for the other three solvable corner models. For
example, in the limit of no barriers and infinite
interactions H(N, W, 0, ∞), we follow the construction of Girardeau and find that there are N !
states for every strictly increasing set of integers{n|n1 < n2 < · · · nN }; c.f. [47]. Therefore, infinite interactions exclude cases when two or morequantum numbers are the same. This is depicted
for H(2, 2, 0, ∞) in Fig. 9, where the spectrum is
missing the diagonal states with n1 = n2 . For
N = 2, the number of states missing from thearea estimation in Eq. (9) is therefore propor√
tional to the length of this boundary r = E.
More generally for N particles, the geometrical structures in the quantum number ‘quadrant’corresponding to all numbers different remains
N -dimensional even in the presence of infinite
interactions. However, the structures with twoquantum numbers equal correspond to interior
boundaries of the quadrant and have dimension
N − 1. Therefore corrections to account for the‘missing states’ appear at the subleading orderrN −1 = E (N −1)/2 . Further corrections to thenumber of states appear at next-to-subleading order E (N −2)/2 when either two pairs of quantum
numbers are the same or three quantum numbers
are the same.
Note that in the example with N = 2 in Fig. 9,
the Tonks-Girardeau map shifts the symmetricstates (n1 , n2 ) in model 1 to (n1 + 1, n2 ) in model
2. Since an integrable model connects these twolimiting cases, this establishes a one-to-one adiabatic mapping between the two spectra. From
this mapping the number of level crossings that
occur as γ is tuned from 0 to ∞ can be explicitly
calculated without actually solving for the spectrum on the integrable model that links these twocases.
Similarly, for the case H(N, W, ∞, 0) of infinite barriers and no interactions, there are W N
solutions for every set of non-negative integers nwhere all integers ni are multiples of the number
of wells W . This increases the volume associated with each set of quantum numbers from 1
to W N and that factor of 1/W N exactly cancels
the degeneracy factor W N giving the same leading term. This is depicted for the simplest caseof N = 2 and W = 2 in model 3 of Fig. 9. As
before, we do expect the coefficient on the subleading term to depend on the intricate combina-

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

13

torics of putting N identical particles in W wells.
A parallel argument holds for the fourth solvablemodel H(N, W, ∞, ∞).
Note again that these four models with exact
solutions are connected by integral models andexact spectral maps can be constructed for allof these cases. Because these four extreme cases
all have the same leading term in the density ofstates that depends only on N , we assume that
all the models that lie in this region of parameter space have the same property, and this isfurther supported by the phase space argument
presented in Sec. 4.1.
When we elect to consider only one symmetry
sector, e.g. bosons with positive parity, that reduces the number of states (9) by a factor of 1/2
for parity and 1/N ! for symmetrization at leading
order. At subleading order, the correction coefficient depends on the combinatorics of N particles in W wells and on the parameters (τ, γ).
For example, for the simplest case of N = 2 andW = 2 depicted, we see the importance of statesalong the diagonal n1 √= n2 , which would result
in corrections of order E for the length of that
diagonal.
As a final comment, if the leading term of thedensity of states is independent of the parameters (τ, γ), then the subleading term contains information about the density of level crossings forthe integrable models which presumably becomeslevel repulsions in the non-integrable parameterregion. Future work will investigate this connection between spectral flow, level crossings, andenergy level statistics.

References[1] Dmitry A. Abanin, Ehud Altman, Immanuel Bloch, and Maksym Serbyn.
Colloquium:
Many-body localization,
thermalization, and entanglement. Rev.
Mod. Phys., 91:021001, May 2019. DOI:
10.1103/RevModPhys.91.021001.
URL
https://link.aps.org/doi/10.1103/
RevModPhys.91.021001.
[2] Y. Alhassid and R. D. Levine. Spectral
autocorrelation function in the statisticaltheory of energy levels. Phys. Rev. A,
46:4650–4653, 1992. DOI: 10.1103/PhysRevA.46.4650. URL https://link.aps.
org/doi/10.1103/PhysRevA.46.4650.

[3] G. E. Astrakharchik,J. Boronat,
J. Casulleras, and S. Giorgini.
BeyondtheTonks-Girardeau
gas:
Strongly correlated regime in quasione-dimensional bose gases.
Phys.
Rev. Lett., 95:190407, Nov 2005. DOI:
10.1103/PhysRevLett.95.190407.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.95.190407.
[4] Daniel Barredo, Sylvain de Léséleuc,
Vincent Lienhard, Thierry Lahaye, andAntoine Browaeys.
An atom-by-atom
assembler of defect-free arbitrary twodimensional atomic arrays. Science, 354(6315):1021–1023, 2016. DOI: 10.1126/science.aah3778.
[5] M T Batchelor, M Bortz, X W Guan,
and N Oelkers. Evidence for the superTonks-Girardeau gas.
J. Stat. Mech.,
2005(10):L10001–L10001, oct 2005. DOI:
10.1088/1742-5468/2005/10/l10001. URL
https://doi.org/10.1088/1742-5468/
2005/10/l10001.
[6] J. H. Becher, E. Sindici, R. Klemt,
S. Jochim, A. J. Daley, and P. M.
Preiss.
Measurement of identical particle entanglement and the influenceof antisymmetrization.
Phys. Rev.
Lett., 125:180402, Oct 2020.
DOI:
10.1103/PhysRevLett.125.180402.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.125.180402.
[7] Andrea Bergschneider,
Vincent M
Klinkhamer, Jan Hendrik Becher, Ralf
Klemt, Lukas Palm, Gerhard Zürn, Selim
Jochim, and Philipp M Preiss. Experimental characterization of two-particle
entanglement through position and momentum correlations. Nat. Phys., 15(7):640–644, 2019. DOI: 10.1038/s41567-0190508-6.
[8] Hannes Bernien, Sylvain Schwartz, Alexander Keesling, Harry Levine, Ahmed Omran, Hannes Pichler, Soonwon Choi,
Alexander S. Zibrov, Manuel Endres,Markus Greiner, Vladan Vuletić, andMikhail D. Lukin.
Probing manybody dynamics on a 51-atom quantum
simulator.
Nature, 551(7682):579–584,
2017. ISSN 1476-4687. DOI: 10.1038/nature24622. URL https://doi.org/10.

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

141038/nature24622.
[9] M. V. Berry.
Quantizing a classically ergodic system:
Sinai’s billiard and the KKR method.
Ann.
Phys., 131(1):163–216, January 1981.
ISSN 0003-4916.
DOI: 10.1016/00034916(81)90189-5.
URL
https://www.sciencedirect.com/science/
article/pii/0003491681901895.
[10] M. V. Berry and M. Tabor. Level clustering in the regular spectrum. Proc. R.
Soc. Lond. A, 356:375 – 394, 1977. DOI:
10.1098/rspa.1977.0140.
[11] M. V. Berry and M. Wilkinson. Diabolical Points in the Spectra of Triangles. Proc. R. Soc. Lond. A, 392(1802):
15–43, 1984.
ISSN 0080-4630.
DOI:
10.1098/rspa.1984.0022. URL https://www.jstor.org/stable/2397740.
[12] B. Bertini, F. Heidrich-Meisner, C. Karrasch, T. Prosen, R. Steinigeweg, andM. Žnidarič. Finite-temperature transport
in one-dimensional quantum lattice models. Rev. Mod. Phys., 93:025003, May 2021.
DOI:
10.1103/RevModPhys.93.025003.
URL
https://link.aps.org/doi/10.
1103/RevModPhys.93.025003.
[13] Wouter Beugeling, Roderich Moessner, and Masudul Haque. Off-diagonal
matrix elements of local operators in
many-body quantum systems.
Phys.
Rev. E, 91:012144, Jan 2015.
DOI:
10.1103/PhysRevE.91.012144.
URL
https://link.aps.org/doi/10.1103/
PhysRevE.91.012144.
[14] Immanuel Bloch, Jean Dalibard, and Wilhelm Zwerger. Many-body physics with ultracold gases. Rev. Mod. Phys., 80:885–
964, Jul 2008. DOI: 10.1103/RevModPhys.80.885. URL https://link.aps.
org/doi/10.1103/RevModPhys.80.885.
[15] D Blume.
Few-body physics with ultracold atomic and molecular systems
in traps.
Rep. Prog. Phys., 75(4):046401, mar 2012. DOI: 10.1088/00344885/75/4/046401. URL https://doi.
org/10.1088/0034-4885/75/4/046401.
[16] Doerte Blume. Jumping from two andthree particles to infinitely many. Physics,
3:74, 2010. DOI: 10.1103/Physics.3.74.
[17] F. Borgonovi, F. M. Izrailev, L. F. Santos,

and V. G. Zelevinsky. Quantum chaos andthermalization in isolated systems of interacting particles. Phys. Rep., 626:1, 2016.
DOI: 10.1016/j.physrep.2016.02.005.
[18] R. Bourgain, J. Pellegrino, A. Fuhrmanek,
Y. R. P. Sortais, and A. Browaeys.
Evaporative cooling of a small number ofatoms in a single-beam microscopic dipole
trap. Phys. Rev. A, 88:023428, Aug 2013.
DOI: 10.1103/PhysRevA.88.023428. URL
https://link.aps.org/doi/10.1103/
PhysRevA.88.023428.
[19] Marlon Brenes, John Goold, and Marcos Rigol.
Low-frequency behavior ofoff-diagonal matrix elements in the integrable XXZ chain and in a locally
perturbed quantum-chaotic XXZ chain.
Phys. Rev. B, 102:075127, Aug 2020.
DOI: 10.1103/PhysRevB.102.075127. URL
https://link.aps.org/doi/10.1103/
PhysRevB.102.075127.
[20] Marlon Brenes, Tyler LeBlond, John
Goold, and Marcos Rigol.
Eigenstate thermalization in a locally perturbed integrable system.
Phys. Rev.
Lett., 125:070605, Aug 2020.
DOI:
https://doi.org/10.1103/PhysRevLett.125.070605.
URL
https://link.aps.org/doi/10.
1103/PhysRevLett.125.070605.
[21] Marlon Brenes, Silvia Pappalardi, Mark T.
Mitchison, John Goold, and Alessandro
Silva. Out-of-time-order correlations andthe fine structure of eigenstate thermalisation. arXiv:2103.01161, 2021. URL https://arxiv.org/abs/2103.01161.
[22] T. A. Brody, J. Flores, J. B. French,
P. A. Mello, A. Pandey, and S. S. M.
Wong. Random-matrix physics: spectrum
and strength fluctuations. Rev. Mod. Phys.,
53:385, 1981.
DOI: 10.1103/RevModPhys.53.385. URL https://link.aps.
org/doi/10.1103/RevModPhys.53.385.
[23] Steve Campbell, Miguel Ángel Garcı́aMarch, Thomás Fogarty, and Thomas
Busch. Quenching small quantum gases:
Genesis of the orthogonality catastrophe. Phys. Rev. A, 90:013617, Jul 2014.
DOI: 10.1103/PhysRevA.90.013617. URL
https://link.aps.org/doi/10.1103/
PhysRevA.90.013617.
[24] Marko Cetina, Michael Jag, Rianne S.

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

15

Lous, Isabella Fritsche, Jook T. M.
Walraven, Rudolf Grimm, Jesper Levinsen, Meera M. Parish, Richard Schmidt,
Michael Knap, and Eugene Demler. Ultrafast many-body interferometry of impurities coupled to a Fermi sea. Science, 354(6308):96–99, 2016. ISSN 0036-8075. DOI:
10.1126/science.aaf5134.
[25] P. Cheinet, S. Trotzky, M. Feld, U. Schnorrberger, M. Moreno-Cardoner, S. Fölling,and I. Bloch. Counting atoms using interaction blockade in an optical superlattice.
Phys. Rev. Lett., 101:090404, Aug 2008.
DOI:
10.1103/PhysRevLett.101.090404.
URL
https://link.aps.org/doi/10.
1103/PhysRevLett.101.090404.
[26] Aurélia Chenu, Javier Molina-Vilaplana,and Adolfo del Campo.
Work statistics, Loschmidt echo and information
scrambling in chaotic quantum systems.
Quantum, 3:127, March 2019.
ISSN
2521-327X. DOI: 10.22331/q-2019-03-04127. URL https://doi.org/10.22331/
q-2019-03-04-127.
[27] J. Cotler, N. Hunter-Jones, J. Liu,and B. Yoshida.
Chaos, complexity,
and random matrices. J. High Energy
Phys., 2017(11):48, Nov 2017.
DOI:
10.1007/JHEP11(2017)048. URL https://doi.org/10.1007/JHEP11(2017)048.
[28] L. D’Alessio, Y. Kafri, A. Polkovnikov,
and M. Rigol.
From quantum chaosand eigenstate thermalization to statistical mechanics and thermodynamics.
Adv. Phys., 65(3):239–362, 2016. DOI:
10.1080/00018732.2016.1198134.
URL
https://doi.org/10.1080/00018732.
2016.1198134.
[29] Javier de la Cruz, Sergio LermaHernández, and Jorge G. Hirsch. Quantum
chaos in a system with high degree of symmetries. Phys. Rev. E, 102:032208, Sep
2020. DOI: 10.1103/PhysRevE.102.032208.
URL
https://link.aps.org/doi/10.
1103/PhysRevE.102.032208.
[30] Manuel Endres, Hannes Bernien, Alexander Keesling, Harry Levine, Eric R. Anschuetz, Alexandre Krajenbrink, Crystal
Senko, Vladan Vuletic, Markus Greiner,and Mikhail D. Lukin. Atom-by-atom assembly of defect-free one-dimensional cold

atom arrays. Science, 354(6315):1024–
1027, 2016. DOI: 10.1126/science.aah3752.
[31] Hans Feldmeier and Jürgen Schnack.
Molecular dynamics for fermions. Rev.
Mod. Phys., 72:655–688, Jul 2000. DOI:
10.1103/RevModPhys.72.655.
URL
https://link.aps.org/doi/10.1103/
RevModPhys.72.655.
[32] V. V. Flambaum and F. M. Izrailev. Statistical theory of finite Fermi systems
based on the structure of chaotic eigenstates. Phys. Rev. E, 56:5144, 1997. DOI:
10.1103/PhysRevE.56.5144.
[33] V. V. Flambaum, A. A. Gribakina, G. F.
Gribakin, and M. G. Kozlov. Structure ofcompound states in the chaotic spectrumof the Ce atom: Localization properties,
matrix elements, and enhancement of weakperturbations. Phys. Rev. A, 50:267–296,
1994. DOI: 10.1103/PhysRevA.50.267.
[34] V. V. Flambaum, F. M. Izrailev, and
G. Casati.
Towards a statistical theory of finite Fermi systems and compound states: Random two-body interaction approach. Phys. Rev. E, 54:2136–
2139, Aug 1996.
DOI: 10.1103/PhysRevE.54.2136. URL https://link.aps.
org/doi/10.1103/PhysRevE.54.2136.
[35] Thomás Fogarty, Sebastian Deffner,Thomas Busch, and Steve Campbell.
Orthogonality catastrophe as a consequence of the quantum speed limit. Phys.
Rev. Lett., 124:110601, Mar 2020. DOI:
10.1103/PhysRevLett.124.110601.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.124.110601.
[36] M. A. Garcı́a-March, B. Juliá-Dı́az,
G. E. Astrakharchik, J. Boronat, and
A. Polls. Distinguishability, degeneracy,
and correlations in three harmonicallytrapped bosons in one dimension. Phys.
Rev. A, 90:063605, Dec 2014.
DOI:
10.1103/PhysRevA.90.063605.
URL
https://link.aps.org/doi/10.1103/
PhysRevA.90.063605.
[37] M. A. Garcı́a-March, A. Yuste, B. JuliáDı́az, and A. Polls.
Mesoscopic superpositions of Tonks-Girardeau states
and the Bose-Fermi mapping.
Phys.
Rev. A, 92:033621, Sep 2015.
DOI:
10.1103/PhysRevA.92.033621.
URL

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

16

https://link.aps.org/doi/10.1103/
PhysRevA.92.033621.
[38] M A Garcia-March, S van Frank, M Bonneau, J Schmiedmayer, M Lewenstein,and Lea F Santos. Relaxation, chaos,and thermalization in a three-mode modelof a Bose–Einstein condensate. New J.
Phys., 20(11):113039, nov 2018. DOI:
10.1088/1367-2630/aaed68. URL https://doi.org/10.1088/1367-2630/aaed68.
[39] M. Gaudin. Boundary energy of a Bose gas
in one dimension. Physical Review A, 4(1):386–394, July 1971. DOI: 10.1103/PhysRevA.4.386. URL http://link.aps.org/
doi/10.1103/PhysRevA.4.386.
[40] Michel Gaudin.
The Bethe Wavefunction.
Cambridge University Press,
Trans. Caux, Jean-Sébastien edition, April
2014. ISBN 978-1-107-04585-9. DOI:
10.1017/CBO9781107053885.
[41] Marvin Girardeau. Relationship between
Systems of Impenetrable Bosons andFermions in One Dimension.
Journalof Mathematical Physics,
1(6):516–
523, November 1960.
ISSN 00222488.
DOI: 10.1063/1.1703687.
URL http://jmp.aip.org/resource/1/jmapaq/v1/
i6/p516_s1?isAuthorized=no.
[42] J. Goold, T. Fogarty, N. Lo Gullo, M. Paternostro, and Th. Busch. Orthogonality
catastrophe as a consequence of qubit
embedding in an ultracold Fermi gas.
Phys. Rev. A, 84:063632, Dec 2011.
DOI: 10.1103/PhysRevA.84.063632. URL
https://link.aps.org/doi/10.1103/
PhysRevA.84.063632.
[43] T. Guhr and H.A. Weidenmüller. Correlations in anticrossing spectra and scattering
theory. Analytical aspects. Chem. Phys.,
146:21 – 38, 1990. DOI: 10.1016/03010104(90)90003-R.
URL
http:
//www.sciencedirect.com/science/
article/pii/030101049090003R.
[44] T. Guhr, A. Mueller-Gröeling, and H. A.
Weidenmüller.
Random matrix theories in quantum physics: Common concepts. Phys. Rep., 299:189, 1998. DOI:
10.1016/S0370-1573(97)00088-4.
[45] Elmar Haller, Mattias Gustavsson, Manfred J. Mark, Johann G. Danzl, Russell
Hart, Guido Pupillo, and Hanns-ChristophNägerl. Realization of an excited, strongly
correlated quantum gas phase. Science, 325(5945):1224–1227, 2009. DOI: 10.1126/science.1175850.
[46] N. L. Harshman. Symmetries of threeharmonically trapped particles in one
dimension. Phys. Rev. A, 86:052122, Nov
2012. DOI: 10.1103/PhysRevA.86.052122.
URL
https://link.aps.org/doi/10.
1103/PhysRevA.86.052122.
[47] N. L. Harshman. One-dimensional traps,
two-body interactions, few-body symmetries: II. N particles. Few-Body Systems,
57(1):45–69, 2016.
ISSN 0177-7963,
1432-5411. DOI: 10.1007/s00601-015-10255.
URL http://link.springer.com.
proxyau.wrlc.org/article/10.1007/
s00601-015-1025-5.
[48] N. L. Harshman. Infinite barriers andsymmetries for a few trapped particles in one dimension.
Phys. Rev.
A, 95(5):053616, May 2017.
DOI:
10.1103/PhysRevA.95.053616.
URL
https://link.aps.org/doi/10.1103/
PhysRevA.95.053616.
[49] N. L. Harshman, Maxim Olshanii,
A. S. Dehkharghani, A. G. Volosniev,
Steven Glenn Jackson, and N. T. Zinner.
Integrable families of hard-core
particles with unequal masses in a onedimensional harmonic trap. Phys. Rev. X,
7:041001, Oct 2017. DOI: 10.1103/PhysRevX.7.041001. URL https://link.aps.
org/doi/10.1103/PhysRevX.7.041001.
[50] N.L. Harshman. One-dimensional traps,
two-body interactions, few-body symmetries: I. One, two, and three particles. FewBody Systems, 57(1):11–43, 2016. DOI:
10.1007/s00601-015-1024-6.
[51] U. Hartmann, H.A. Weidenmüller, andT. Guhr.
Correlations in anticrossing
spectra and scattering theory: Numerical simulations.
Chem. Phys., 150(3):311 – 320, 1991. ISSN 0301-0104.
DOI:
http://dx.doi.org/10.1016/03010104(91)87105-5.
URL
http:
//www.sciencedirect.com/science/
article/pii/0301010491871055.
[52] Xiaodong He, Peng Xu, Jin Wang, andMingsheng Zhan. High efficient loadingof two atoms into a microscopic optical

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

17

trap by dynamically reshaping the trapwith a spatial light modulator. Opt. Express, 18(13):13586–13592, 2010. DOI:
10.1364/OE.18.013586.
[53] Ole J. Heilmann and Elliott H. Lieb.
Violation of the noncrossing rule: The
Hubbard Hamiltonian for benzene.
Ann. N.Y. Acad. Sci., 172(15):584–
617, 1971.
ISSN 1749-6632.
DOI:
10.1111/j.1749-6632.1971.tb34956.x. URL
https://nyaspubs.onlinelibrary.
wiley.com/doi/abs/10.1111/j.
1749-6632.1971.tb34956.x.
[54] David Huber, Oleksandr V. Marchukov,
Hans-Werner Hammer, and Artem G.
Volosniev.
Morphology of three-bodyquantum states from machine learning.
New Journal of Physics, 2021.
DOI: 10.1088/1367-2630/ac0576.
URL
http://iopscience.iop.org/article/
10.1088/1367-2630/ac0576.
[55] Marko Žnidarič.
Weak integrability
breaking: Chaos with integrability signature in coherent diffusion.
Phys.
Rev. Lett., 125:180605, Oct 2020. DOI:
10.1103/PhysRevLett.125.180605.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.125.180605.
[56] Victor Ivrii.
100 years of Weyl’s
law.
Bull. Math. Sci., 6(3):379–452,
October 2016. ISSN 1664-3607, 16643615.
DOI: 10.1007/s13373-016-0089y. URL https://link.springer.com/
article/10.1007/s13373-016-0089-y.
[57] F. M. Izrailev. Quantum-classical correspondence for isolated systems of interacting particles: Localization and ergodicity in energy space.
Phys. Scr.,
T90:95–104, 2001. DOI: 10.1238/Physica.Topical.090a00095.
[58] Felix M. Izrailev.
Simple modelsof quantum chaos:Spectrum andeigenfunctions.
Phys. Rep., 196(5):299 – 392, 1990.
ISSN 0370-1573.
DOI:
https://doi.org/10.1016/03701573(90)90067-C.
URL
http:
//www.sciencedirect.com/science/
article/pii/037015739090067C.
[59] Sudhir R. Jain. Exact solution of the
Schrödinger equation for a particle in a
regular N-simplex. Phys. Lett. A, 372(12):1978–1981, March 2008. ISSN 0375-9601.
DOI:
10.1016/j.physleta.2007.11.016.
URL
https://www.sciencedirect.
com/science/article/pii/
S0375960107016027.
[60] A. M. Kaufman, B. J. Lester, and C. A.
Regal. Cooling a single atom in an optical tweezer to its quantum ground state.
Phys. Rev. X, 2:041014, Nov 2012. DOI:
10.1103/PhysRevX.2.041014.
[61] A. M. Kaufman, B. J. Lester, C. M.
Reynolds, M. L. Wall, M. Foss-Feig,
K. R. A. Hazzard, A. M. Rey, and C. A.
Regal.
Two-particle quantum interference in tunnel-coupled optical tweezers.
Science, 345(6194):306–309, 2014. DOI:
10.1126/science.1250057.
[62] Tim Keller and Thomás Fogarty. Probing the out-of-equilibrium dynamics
of two interacting atoms.
Phys.
Rev. A, 94:063620, Dec 2016.
DOI:
10.1103/PhysRevA.94.063620.
URL
https://link.aps.org/doi/10.1103/
PhysRevA.94.063620.
[63] Michael Knap, Aditya Shashi, Yusuke
Nishida, Adilet Imambekov, Dmitry A.
Abanin, and Eugene Demler.
Timedependent impurity in ultracold fermions:
Orthogonality catastrophe and beyond.
Phys. Rev. X, 2:041020, Dec 2012.
DOI: 10.1103/PhysRevX.2.041020. URL
https://link.aps.org/doi/10.1103/
PhysRevX.2.041020.
[64] H. R. Krishnamurthy, H. S. Mani, and
H. C. Verma.
Exact solution of theschrodinger equation for a particle ina tetrahedral box. J. Phys. A, 15(7):2131–2137, July 1982. ISSN 0305-4470.
DOI: 10.1088/0305-4470/15/7/024. URL
https://doi.org/10.1088/0305-4470/
15/7/024.
[65] M. Lacki, M. A. Baranov, H. Pichler,
and P. Zoller. Nanoscale “dark state”
optical potentials for cold atoms. Phys.
Rev. Lett., 117:233001, Nov 2016. DOI:
10.1103/PhysRevLett.117.233001.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.117.233001.
[66] Tyler LeBlond, Krishnanand Mallayya,
Lev Vidmar, and Marcos Rigol. Entanglement and matrix elements of observables

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

18

in interacting integrable systems. Phys.
Rev. E, 100:062134, Dec 2019.
DOI:
10.1103/PhysRevE.100.062134.
URL
https://link.aps.org/doi/10.1103/
PhysRevE.100.062134.
[67] S. Lerma-Hernández, D. Villaseñor, M. A.
Bastarrachea-Magnani, E. J. TorresHerrera, L. F. Santos, and J. G. Hirsch.
Dynamical signatures of quantum chaos
and relaxation time scales in a spin-boson
system. Phys. Rev. E, 100:012218, Jul
2019. DOI: 10.1103/PhysRevE.100.012218.
URL
https://link.aps.org/doi/10.
1103/PhysRevE.100.012218.
[68] Brian J. Lester, Niclas Luick, Adam M.
Kaufman, Collin M. Reynolds, andCindy A. Regal. Rapid production of
uniformly filled arrays of neutral atoms.
Phys. Rev. Lett., 115:073003, Aug 2015.
DOI:
10.1103/PhysRevLett.115.073003.
URL
https://link.aps.org/doi/10.
1103/PhysRevLett.115.073003.
[69] Luc Leviandier, Maurice Lombardi, Rémi
Jost, and Jean Paul Pique. Fourier transform: A tool to measure statistical levelproperties in very complex spectra. Phys.
Rev. Lett., 56:2449–2452, Jun 1986. DOI:
10.1103/PhysRevLett.56.2449.
[70] Maciej Lewenstein,
Anna Sanpera,
and Veronica Ahufinger.
Ultracold
Atoms in Optical Lattices: Simulating
quantum many-body systems.
Oxford University Press, 2012.
DOI:
10.1093/acprof:oso/9780199573127.001.0001.
[71] F Leyvraz, A Garcı́a, H Kohler, andT H Seligman. Fidelity under isospectral
perturbations: a random matrix study.
J. Phys. A, 46(27):275303, 2013. DOI:
10.1088/1751-8113/46/27/275303.
URL
http://stacks.iop.org/1751-8121/46/
i=27/a=275303.
[72] Elliott H. Lieb and Werner Liniger.
Exact analysis of an interacting Bosegas. I. The general solution and theground state.
Phys. Rev., 130:1605–
1616, May 1963.
DOI: 10.1103/PhysRev.130.1605. URL https://link.aps.
org/doi/10.1103/PhysRev.130.1605.
[73] M. Lombardi and T. H. Seligman. Universal and nonuniversal statistical properties of levels and intensities for chaotic

Rydberg molecules. Phys. Rev. A, 47:
3571–3586, May 1993. DOI: 10.1103/PhysRevA.47.3571. URL http://link.aps.
org/doi/10.1103/PhysRevA.47.3571.
[74] D. Luitz and Y. Bar Lev. The ergodic sideof the many-body localization transition.
Ann. Phys.(Berlin), 529(7):1600350, 2017.
DOI: 10.1002/andp.201600350.
[75] Juan Maldacena, Stephen H. Shenker,
and Douglas Stanford.
A bound on
chaos. J. High Energy Phys., 2016(8):106,
2016.
DOI: 10.1007/JHEP08(2016)106.
URL
http://dx.doi.org/10.1007/
JHEP08(2016)106.
[76] M. L. Mehta. Random Matrices. Academic
Press, Boston, USA, 1991.
[77] Nicolás Mirkin and Diego Wisniacki.
Quantum chaos, equilibration, and control
in extremely short spin chains. Phys.
Rev. E, 103:L020201, Feb 2021. DOI:
10.1103/PhysRevE.103.L020201.
URL
https://link.aps.org/doi/10.1103/
PhysRevE.103.L020201.
[78] Simon Murmann, Andrea Bergschneider,
Vincent M Klinkhamer, Gerhard Zürn,
Thomas Lompe, and Selim Jochim. Two
fermions in a double well: Exploring afundamental building block of the Hubbard model.
Phys. Rev. Lett., 114(8):080402, 2015.
DOI: 10.1103/PhysRevLett.114.080402.
[79] R. Nandkishore and D.A. Huse. Manybody localization and thermalization inquantum statistical mechanics. Annu. Rev.
Condens. Matter Phys., 6:15, 2015. DOI:
10.1146/annurev-conmatphys-031214014726.
[80] N. Oelkers, M. T. Batchelor, M. Bortz,
and X.-W. Guan. Bethe ansatz studyof one-dimensional Bose and Fermi gaseswith periodic and hard wall boundary conditions. J. Phys. A, 39(5):1073, 2006.
ISSN 0305-4470.
DOI: 10.1088/03054470/39/5/005.
URL http://stacks.
iop.org/0305-4470/39/i=5/a=005.
[81] M. Olshanii. Atomic scattering in the presence of an external confinement and a gasof impenetrable bosons. Phys. Rev. Lett.,
81:938–941, Aug 1998. DOI: 10.1103/PhysRevLett.81.938. URL https://link.aps.
org/doi/10.1103/PhysRevLett.81.938.

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

19

[82] Maxim Olshanii and Steven G. Jackson.
An exactly solvable quantum four-bodyproblem associated with the symmetriesof an octacube. New J. Phys., 17(10):105005, 2015. ISSN 1367-2630. DOI:
10.1088/1367-2630/17/10/105005.
URL
http://stacks.iop.org/1367-2630/17/
i=10/a=105005.
[83] Maxim Olshanii, Thibault Scoquart,
Dmitry Yampolsky, Vanja Dunjko, andSteven Glenn Jackson. Creating entanglement using integrals of motion. Phys.
Rev. A, 97(1):013630, January 2018.
DOI: 10.1103/PhysRevA.97.013630. URL
https://link.aps.org/doi/10.1103/
PhysRevA.97.013630.
[84] Guido Pagano, Marco Mancini, Giacomo Cappellini, Pietro Lombardi, Florian Schäfer, Hui Hu, Xia-Ji Liu, Jacopo Catani, Carlo Sias, Massimo Inguscio, and Leonardo Fallani. A onedimensional liquid of fermions with tunable spin. Nat. Phys., 10(3):198–201, 2014.
DOI: 10.1038/nphys2878. URL https://doi.org/10.1038/nphys2878.
[85] Akhilesh Pandey and Ramakrishna Ramaswamy. Level spacings for harmonicoscillator systems.
Phys. Rev. A, 43:4237–4243, Apr 1991. DOI: 10.1103/PhysRevA.43.4237. URL http://link.aps.
org/doi/10.1103/PhysRevA.43.4237.
[86] J. P. Pique, Y. Chen, R. W. Field,
and J. L. Kinsey.
Chaos and dynamics on 0.5–300 ps time scales invibrationally excited acetylene: Fouriertransform of stimulated-emission pumping spectrum.
Phys. Rev. Lett., 58:
475–478, Feb 1987. DOI: 10.1103/PhysRevLett.58.475. URL http://link.aps.
org/doi/10.1103/PhysRevLett.58.475.
[87] R. E. Prange.
The spectral form
factor is not self-averaging.
Phys.
Rev. Lett., 78:2280, 1997.
DOI:
10.1103/PhysRevLett.78.2280.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.78.2280.
[88] Cindy Regal. Bringing order to neutral
atom arrays. Science, 354(6315):972–973,
2016. DOI: 10.1126/science.aaj2145.
[89] Irina Reshodko, Albert Benseny, Judit
Romhányi, and Thomas Busch. Topolog-

ical states in the Kronig–Penney modelwith arbitrary scattering potentials. New
J. Phys., 21(1):013010, jan 2019. DOI:
10.1088/1367-2630/aaf9bf. URL https://doi.org/10.1088/1367-2630/aaf9bf.
[90] Marcos Rigol and Alejandro Muramatsu.
Ground-state properties of hard-corebosons confined on one-dimensional optical lattices. Phys. Rev. A, 72:013604, Jul
2005. DOI: 10.1103/PhysRevA.72.013604.
URL
https://link.aps.org/doi/10.
1103/PhysRevA.72.013604.
[91] Daniel A. Roberts, Douglas Stanford,
and Leonard Susskind. Localized shocks.
JHEP, 2015(3):51, 2015.
ISSN 10298479.
DOI: 10.1007/JHEP03(2015)051.
URL
http://dx.doi.org/10.1007/
JHEP03(2015)051.
[92] Lea F. Santos, Francisco Pérez-Bernal,
andE.
JonathanTorres-Herrera.
Speck of chaos.
Phys. Rev. Research, 2:043034, Oct 2020.
DOI:
10.1103/PhysRevResearch.2.043034.
URL
https://link.aps.org/doi/10.
1103/PhysRevResearch.2.043034.
[93] S. Saskin, J. T. Wilson, B. Grinkemeyer,
and J. D. Thompson. Narrow-line cooling and imaging of ytterbium atoms inan optical tweezer array.
Phys. Rev.
Lett., 122:143002, Apr 2019.
DOI:
10.1103/PhysRevLett.122.143002.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.122.143002.
[94] M. Schiulaz, M. Távora, and L. F. Santos. From few- to many-body quantum
systems. Quantum Sci. Technol., 3:044006,
2018. DOI: 10.1088/2058-9565/aad913.
[95] Mauro Schiulaz, E. Jonathan TorresHerrera, and Lea F. Santos.
Thouless and relaxation time scales inmany-body quantum systems.
Phys.
Rev. B, 99:174313, May 2019.
DOI:
10.1103/PhysRevB.99.174313.
URL
https://link.aps.org/doi/10.1103/
PhysRevB.99.174313.
[96] Mauro Schiulaz, E. Jonathan TorresHerrera, Francisco Pérez-Bernal, and
Lea F. Santos.
Self-averaging in
many-body quantum systems out of
equilibrium: Chaotic systems.
Phys.
Rev. B, 101:174312, May 2020. DOI:

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

20

10.1103/PhysRevB.101.174312.
URL
https://link.aps.org/doi/10.1103/
PhysRevB.101.174312.
[97] J. Schnack and H. Feldmeier. Statisticalproperties of fermionic molecular dynamics. Nuc. Phys. A, 601(2):181 – 194, 1996.
ISSN 0375-9474.
DOI: 10.1016/03759474(95)00505-6.
URL
http:
//www.sciencedirect.com/science/
article/pii/0375947495005056.
[98] M. Schreiber, S. S. Hodgman, Pr. Bordia, H. P. Lüschen, M. H. Fischer,
R. Vosk, E. Altman, U. Schneider,and I. Bloch.
Observation of manybody localization of interacting fermionsin a quasirandom optical lattice. Science, 349(6250):842–845, 2015.
ISSN
0036-8075. DOI: 10.1126/science.aaa7432.
URL http://science.sciencemag.org/
content/349/6250/842.
[99] Kai-Niklas Schymik, Vincent Lienhard,
Daniel Barredo, Pascal Scholl, Hannah
Williams, Antoine Browaeys, and Thierry
Lahaye. Enhanced atom-by-atom assemblyof arbitrary tweezer arrays. Phys. Rev. A,
102:063107, Dec 2020. DOI: 10.1103/PhysRevA.102.063107.
[100] Friedhelm Serwane,
Gerhard Zürn,
Thomas Lompe, TB Ottenstein, AN Wenz,
and S Jochim. Deterministic preparationof a tunable few-fermion system. Science, 332(6027):336–338, 2011.
DOI:
10.1126/science.1201351.
[101] G. B. Shaw. Degeneracy in the particlein-a-box problem. J. Phys. A, 7(13):1537–
1546, September 1974. ISSN 0301-0015.
DOI: 10.1088/0305-4470/7/13/008. URL
https://doi.org/10.1088/0305-4470/
7/13/008.
[102] Tomasz Sowiński and Miguel Ángel Garcı́aMarch. One-dimensional mixtures of several ultracold atoms: a review. Rep.
Progr. Phys., 82(10):104401, 2019. DOI:
10.1088/1361-6633/ab3a80.
[103] H-J. Stöckmann.
Quantum Chaos:an introduction. Cambridge University
Press, Cambridge, UK, 2006.
DOI:
10.1017/CBO9780511524622.
[104] E. J. Torres-Herrera and Lea F. Santos. Extended nonergodic states in disordered many-body quantum systems. Ann.

Phys. (Berlin), 529:1600284, 2017. ISSN
1521-3889. DOI: 10.1002/andp.201600284.
URL http://dx.doi.org/10.1002/andp.
201600284.
[105] E. J. Torres-Herrera and Lea F. Santos. Dynamical manifestations of quantum
chaos: correlation hole and bulge. Philos.
Trans. Royal Soc. A, 375(2108):20160434,
2017. DOI: 10.1098/rsta.2016.0434.
[106] E. J. Torres-Herrera, Antonio M. Garcı́aGarcı́a, and Lea F. Santos.
Generic
dynamical features of quenched interacting
quantum systems: Survival probability,
density imbalance, and out-of-time-orderedcorrelator. Phys. Rev. B, 97:060303, Feb
2018. DOI: 10.1103/PhysRevB.97.060303.
URL
https://link.aps.org/doi/10.
1103/PhysRevB.97.060303.
[107] C. J. Turner, A. A. Michailidis, D. A.
Abanin, M. Serbyn, and Z. Papić.
Weak ergodicity breaking from quantum many-body scars.
Nat. Phys.,
14(7):745–749, 2018.
ISSN 17452481.
DOI: 10.1038/s41567-018-01375.
URL https://doi.org/10.1038/
s41567-018-0137-5.
[108] J. W. Turner.
On the quantum
particle in a polyhedral box.
J.
Phys. A, 17(14):2791–2797, October 1984.
ISSN 0305-4470.
DOI: 10.1088/03054470/17/14/022. URL https://doi.org/
10.1088/0305-4470/17/14/022.
[109] Y. Wang, S. Subhankar, P. Bienias,
M. Lacki, T-C. Tsui, M. A. Baranov, A. V.
Gorshkov, P. Zoller, J. V. Porto, and S. L.
Rolston. Dark state optical lattice witha subwavelength spatial structure. Phys.
Rev. Lett., 120:083601, Feb 2018. DOI:
10.1103/PhysRevLett.120.083601.
URL
https://link.aps.org/doi/10.1103/
PhysRevLett.120.083601.
[110] AN Wenz, G Zürn, Simon Murmann,
I Brouzos, T Lompe, and S Jochim. From
few to many: Observing the formation of a
Fermi sea one atom at a time. Science, 342(6157):457–460, 2013. DOI: 10.1126/science.1240516.
[111] Joshua Wilkie and Paul Brumer. Timedependent manifestations of quantum
chaos. Phys. Rev. Lett., 67:1185–1188, Sep
1991. DOI: 10.1103/PhysRevLett.67.1185.

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

21

URL
https://link.aps.org/doi/10.
1103/PhysRevLett.67.1185.
[112] Sebastian Will, Thorsten Best, Ulrich
Schneider, Lucia Hackermüller, Dirk-SörenLühmann, and Immanuel Bloch. Timeresolved observation of coherent multibody interactions in quantum phase revivals. Nature, 465(7295):197–201, May
2010. ISSN 1476-4687. DOI: 10.1038/nature09036. URL https://doi.org/10.
1038/nature09036.
[113] Pablo R. Zangara, Axel D. Dente, E. J.
Torres-Herrera, Horacio M. Pastawski,
A. Iucci, and Lea F. Santos. Time fluctuations in isolated quantum systems of interacting particles. Phys. Rev. E, 88:032913,
2013. DOI: 10.1103/PhysRevE.88.032913.
[114] V. Zelevinsky, B. A. Brown, N. Frazier,and M. Horoi. The nuclear shell model asa testing ground for many-body quantum
chaos. Phys. Rep., 276:85–176, 1996. DOI:
10.1016/S0370-1573(96)00007-5.
[115] Guy Zisling, Lea F. Santos, andYevgeny Bar Lev. How many particles
make up a chaotic many-body quantum
system?
SciPost Phys., 10:88, 2021.
DOI:
10.21468/SciPostPhys.10.4.088.
URL https://scipost.org/10.21468/
SciPostPhys.10.4.088.

Accepted in Quantum 2021-06-21, click title to verify. Published under CC-BY 4.0.

22