††thanks: Corresponding author<EMAIL_ADDRESS>
# Level spacing distribution of localized phases induced by quasiperiodic
potentials

Chao Yang  Shenzhen Institute for Quantum Science and Engineering, Southern
University of Science and Technology, Shenzhen 518055, China  International
Quantum Academy, Shenzhen 518048, China  Guangdong Provincial Key Laboratory
of Quantum Science and Engineering, Southern University of Science and
Technology, Shenzhen 518055, China  Yucheng Wang  Shenzhen Institute for
Quantum Science and Engineering, Southern University of Science and
Technology, Shenzhen 518055, China  International Quantum Academy, Shenzhen
518048, China  Guangdong Provincial Key Laboratory of Quantum Science and
Engineering, Southern University of Science and Technology, Shenzhen 518055,
China

###### Abstract

Level statistics is a crucial tool in the exploration of localization physics.
The level spacing distribution of the disordered localized phase follows
Poisson statistics, and many studies naturally apply it to the quasiperiodic
localized phase. Here we analytically obtain the level spacing distribution of
the quasiperiodic localized phase, and find that it deviates from Poisson
statistics. Moreover, based on this level statistics, we derive the ratio of
adjacent gaps and find that for a single sample, it is a $\delta$ function,
which is in excellent agreement with numerical studies. Additionally, unlike
disordered systems, in quasiperiodic systems, there are variations in the
level spacing distribution across different regions of the spectrum, and
increasing the size and increasing the sample are non-equivalent. Our findings
carry significant implications for the reevaluation of level statistics in
quasiperiodic systems and a profound understanding of the distinct effects of
quasiperiodic potentials and disorder induced localization.

Introduction.— Quantum localization has consistently been a significant
research area in condensed matter physics. This phenomenon is widely present
in disordered systems, caused by interference from multiply scattered waves
due to system disorder, resulting in the exponential decay of the wave
function and the suppression of transport Anderson1958 ; RMP1 ; RMP2 ;
Kramer1993 . In addition to random disorder, quasiperiodic potentials also
induce localization, and in recent years, they have garnered widespread
interest in both theoretical Soukoulis1981 ; DasSarma1988 ; Biddle2009 ;
XLi2017 ; HYao2019 ; Ganeshan2015 ; Wang1 ; XCZhou ; Wang2022 ; XDeng2019 ;
Ribeiro  and experimental Roati2008 ; Bloch4 ; An2018 ; JiasT ; TXiao2021 ;
Weld ; HengFan  aspects, playing a crucial role in enhancing our understanding
of critical phases TXiao2021 ; Weld ; HengFan ; WangYC2021 ; Kohmoto1990 ,
rich transport behaviors LandiRMP ; Saha ; Dwiputra ; Lacerda ; Jeffrey ,
many-body localization (MBL) IBlochRMP ; BAA ; Schreiber2015 , low dimensional
Anderson transition and mobility edges Soukoulis1981 ; DasSarma1988 ;
Biddle2009 ; XLi2017 ; HYao2019 ; Ganeshan2015 ; Wang1 ; XCZhou ; Wang2022 ;
XDeng2019 ; Ribeiro ; Roati2008 ; Bloch4 ; An2018 ; JiasT .

The energy level distribution of localized phases is completely distinct from
that of extended phases, allowing us to use level spacing statistics to
differentiate between extended and localized phases Shklovskii1993 ;
Mirlin2000 . For disordered systems, the energy levels of the localized phase
are uncorrelated, with no level repulsion, and their distribution follows
Poisson statistics Shklovskii1993 ; Mirlin2000 ; Molcanov1981 . Extending the
statistical patterns of energy levels for disorder-induced localized phases to
quasiperiodic-induced localized phases seems natural. Additionally, the
average of the adjacent gap ratio $\langle r\rangle$ is close to 0.387 XLi2016
; YLiu2020 ; YWang2023 ; Xianlong ; Ray2016 ; JHPixley ; Schiffer2021 ;
Ray2018 ; Khemani2017 ; Roushan2017 ; XiaopengLi ; Modak2015 , which is in
complete agreement with the results predicted by Poisson statistics.
Therefore, the level statistics of quasiperiodic localization systems are
widely accepted to follow Poisson statistics in both single-particle XLi2016 ;
YLiu2020 ; YWang2023 ; Xianlong ; Ray2016 ; JHPixley ; Schiffer2021 ; Ray2018
; Machida1986 ; SNEvangelou ; Roy2019 ; Takada2004 ; YWang2016  and many-body
systems Khemani2017 ; Roushan2017 ; XiaopengLi ; Modak2015 . However, recent
mathematical proof has shown that the distribution of eigenvalues in
quasiperiodic and disordered localized phases exhibits significant differences
Jitomirskaya , implying that the patterns of energy level spacings for the two
cases may also differ. This inspires us to reexamine the level spacing
distribution of the quasiperiodic localized phase and the spectrum’s
differences between this and localized phase induced by disorder.

In this letter, we investigate the level spacing distribution of localized
phases induced by quasiperiodic potentials. We first compare the energy level
distribution $P(E)$, level spacing distribution $P(\delta E)$ and the
distribution of the adjacent gap ratio $P(r)$ for Anderson localization (AL)
phase induced by quasi-periodic potentials with those induced by disorder.
Next, we calculate the number variance of different regions of energy spectrum
in the Aubry-André (AA) model AA . Then, we analytically derive $P(\delta E)$
and $P(r)$ for the AA model’s AL phase. Finally, we calculate $\langle
r\rangle$ and discuss the probability distribution $P(\delta E)$ when
interactions are present.

Figure 1:  (a1-a4) are about disorder-induced localization with $W/J=10$, and
(b1-b4) are about localization induced by quasiperiodic potentials with
$\alpha=F_{22}/F_{23}$, $\theta=0.3\pi$ and $V/J=10$. The size of both systems
is $L=F_{23}=28657$. (a1) and (b1) Eigenenergies in ascending order, and the
index of energy mode $i_{E}$ runs from $1$ to $L$. The green/red/green curves
in (a1) correspond to the lowest/middle/highest $1/3$ energies of the middle
region. (a2) and (b2) The distribution of energies(density of states). The
purple data points in (b2) correspond to $J=0$. (a3) and (b3) Level spacing
statistics with $\delta E=\Delta E/\langle\Delta E\rangle$ and $\langle\Delta
E\rangle$ is the mean level spacing. The red fitting curve in (a3) shows
$P(\delta E)=1.047e^{-1.065\delta E}$. The purple and light blue dashed lines
in (b3) correspond to the fitting results using Eq. (8) (corresponding to
$J=0$) and Eq. (9) (satisfying $p(\delta
E)=-0.07846+\frac{0.5961}{\pi\sqrt{2.592^{2}-\delta E^{2}}}$), respectively.
(a4) and (b4) The distributions of $P(r)$. The red curve in (a4) satisfies
$P(r)=2/(1+r)^{2}$. Inset of (b4): The behavior $P(r)$ with different sizes.

Model and results.—  The AA model is the simplest nontrivial example with a
one-dimensional quasiperiodic potential, described by

$H=J\sum_{j}(c_{j+1}^{\dagger}c_{j}+c_{j}^{\dagger}c_{j+1})+\sum_{j}V_{j}c_{j}^{\dagger}c_{j},$
(1)

where $c_{j}$ ($c^{\dagger}_{j}$) denotes the annihilation (creation) operator
at site $j$, $J$ is the nearest-neighbor hopping coefficient, and
$V_{j}=V\sum_{j}\cos(2\pi\alpha j+\theta)$ with $V$, $\theta$ and $\alpha$
being the quasiperiodic potential amplitude, the phase offset, and an
irrational number, respectively. This model demonstrates self-duality in the
transformation between lattice space and its dual space at $V=2J$, resulting
in the Anderson transition with all eigenstates being extended (localized) for
$V<2J$ ($V>2J$) AA . For simplicity, we fix $J=1$ and set
$\alpha=F_{N-1}/F_{N}$ with $F_{N}$ being the Fibonacci sequence (i.e.,
$F_{1}=1$, $F_{2}=1$, and $F_{N}=F_{N-1}+F_{N-2}$). As $N$ approaches
infinity, $\alpha$ converges to $(\sqrt{5}-1)/2$. Unless otherwise stated, we
take the system size $L=F_{N}$, and use open boundary conditions.

We first compare the level spacing distributions of the localized phase in the
AA model with that induced by random disorder, as shown in Fig. 1. For the
disorder-induced localization, we consider the above Eq. (1), with the onsite
disorder $V_{j}$ being uniformly distributed in the interval $[-W,W]$. We
observe that the energy spectrum of localized phases caused by disorder does
not exhibit significant large gaps [Fig. 1(a1)]. Apart from a decrease in the
density of states (DOS) at the boundaries of the spectrum, the DOS across the
spectrum is uniformly distributed [Fig. 1(a2)]. As a contrast, the energy
spectrum of quasiperiodic localized phases shows two distinct large gaps,
dividing the spectrum into three segments [Fig. 1(b1)]. The numbers of states
in each segment from bottom to top are $F_{N-2}$, $F_{N-3}$, $F_{N-2}$, and at
the boundaries of each segment, the DOS increases [Fig. 1(b2)]. Then we
compare the distribution of energy level spacings, defined as $\Delta
E_{n}=E_{n+1}-E_{n}$, with the eigenvalues $E_{n}$ listed in ascending order.
In the disorder system, the level statistics of localized phases are Poisson:
$P(\delta E)=\frac{1}{\langle\delta E\rangle}exp(-\frac{\delta
E}{\langle\delta E\rangle})$ [Fig. 1(a3)], where $\delta E=\Delta
E/\langle\Delta E\rangle$ and $\langle\delta E\rangle$ is the average of
$\delta E$. Based on the energy level spacing, we can obtain the ratio of
adjacent gaps as $r_{n}=\frac{min(\delta E_{n},\delta E_{n+1})}{max(\delta
E_{n},\delta E_{n+1})}$ DAHuse1 ; DAHuse2 . For Poisson statistics, one can
derive that the distribution of $r$ satisfies $P(r)=2/(1+r)^{2}$ [Fig. 1(a4)],
which gives the average value of $r$ as $\langle
r\rangle=\int_{0}^{1}P(r)rdr=2\ln 2-1\approx 0.387$. However, for the
quasiperiodic localized phase, the energy level spacing noticeably deviates
from Poisson statistics, as indicated by the black data points in Fig. 1(b3).
The distribution $P(r)$ is not $2/(1+r)^{2}$ but rather takes on the form of a
$\delta$ function [Fig. 1(b4)].

Before deriving the distributions $P(\delta E)$ and $P(r)$, we first
investigate the uniformity of level spacings across different regions in the
spectrum. Fig. 2(a) displays three types of spectra, corresponding to the
localized phase in disordered systems and the edge and middle regions in the
middle segment of Fig. 1(b1). The distances between the energy levels of the
disordered system (blue lines) show significant fluctuations and lack
correlation, allowing levels to approach each other arbitrarily closely.
Similar properties are observed in the boundaries of the quasi-periodic
system’s energy spectrum (green lines). However, in the middle region of each
segment of the energy spectrum, level repulsion is observed, expressing the
unlikelihood of levels being degenerate in this system. This suggests that the
energy spectrum of the localized phase in quasi-periodic systems does not lack
correlation between energy levels as observed in disordered systems. To
characterize the uniformity of energy level spacings, we investigate the level
number variance $\Sigma^{2}(M)$, defined as $\Sigma^{2}(M(\epsilon))=\langle
M^{2}(\epsilon)\rangle-\langle M(\epsilon)\rangle^{2}$, where $M(\epsilon)$
quantifies the number of levels within the width $\epsilon$ on the unfolded
scale Dyson1963 ; Guhr1998 ; Bertrand2016 ; WangMBC . For Poisson statistics,
the spectrum exhibits no correlations, resulting in a number variance that is
exactly linear with a slope of one, i.e., $\Sigma^{2}(M)=M$ (blue dashed line
in Fig. 2(b)). Fig. 2(a) shows that the distribution of energy levels in the
middle region of each segment of the energy spectrum is more uniform, leading
to a smaller $\Sigma^{2}$ (red dots in Fig. 2(b)), similar to
$\Sigma^{2}(M)\approx\frac{2}{\pi^{2}}\ln(2\pi M)$ that is obtained from the
Wigner-Dyson distribution. For each segment’s boundary region (green dots in
Fig. 2(b)) and the overall energy spectrum (black dots in Fig. 2(b)), when $M$
is large, the linear slope of their respective $\Sigma^{2}$ is greater than
$1$. This indicates that their energy level distribution is more uneven than
the Poisson distribution. To the best of our knowledge, the spectra that
satisfy the condition of the number variance being linear with a slope greater
than $1$ have not been reported before.

Figure 2:  (a) Examples of spectra. The blue uncorrelated energy levels are
obtained in disordered AL phase. The green and red energy levels respectively
correspond to the edge and middle levels of the middle region in Fig. 1(b1).
Their corresponding number variances are shown in (b). averaged over $30$
samples. The green/red dots are calculated from the lowest/middle one-third
spectrum of the middle region. The black dots and blue curve are calculated
from the full spectrum in disordered and quasiperiodic systems. Here we take
$30$ samples, with each sample specified by choosing an initial phase
$\theta$. The remaining parameters are the same as those in Fig. 1. Figure 3:
Scheme of energy levels for (a) $E_{j}$ and (b) $E_{1,m}$ and $E_{2,m}$. (c)
Level spacings of AA model with $V=10$, $\alpha=F_{N-1}/F_{N}$ and $L=28657$.
The black and purple dots correspond to the hopping amplitudes $J=1$ and
$J=0$, respectively. (d) Scheme of level spacing as a function of $x$.

We now deduce the level spacing distribution $P(\delta E)$ in the AL phase of
the AA model. We set $J=0$, $\alpha=F_{N-1}/F_{N}$, and fix $\theta$, then the
system’s eigenvalues are $E_{j}=V\cos(2\pi jF_{N-1}/F_{N}+\theta)$, with
$j=0,1,2,\cdots,F_{N}-1$ [Fig. 3(a)]. We introduce $n$, setting it equal to
$jF_{N-1}$ mod($F_{N}$), then $E_{n}=-V\cos[2\pi n/F_{N}+\theta+\pi]$, and it
is easy to verify that the range of $n$ is $n=0,1,2,...,F_{N}-1$. Shifting the
labels of energies $n=m-n_{0}$ with
$n_{0}=(\theta+\pi-\theta_{1})/(\frac{2\pi}{F_{N}})$, one can obtain

$E_{m}=-V\cos(2\pi\frac{m}{F_{N}}+\theta_{1}),\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ m=1,2,\cdots,F_{N}.$ (2)

By selecting the appropriate value for $n_{0}$, one can make the range of
$\theta_{1}$ satisfy $\theta_{1}\in[-2\pi/F_{N},0)$ explain1 . We then
separate the energy levels into two parts, as shown in Fig. 3(b). For
$m=1,\cdots,F_{N}/2$ we denote the energies by $E_{1,m}$; for
$m=F_{N}/2+1,F_{N}/2+2,\cdots,F_{N}$ we relabel $m\rightarrow F_{N}+1-m$ and
denote the energies by $E_{2,m}$, hence $E_{1,m}=-V\cos(2\pi
m/F_{N}+\theta_{1})$ and $E_{2,m}=-V\cos[2\pi(m-1)/F_{N}-\theta_{1}]$, with
$m=1,2,\cdots,F_{N}/2$. It is convenient to introduce variables
$x_{m}=2\pi(m-1/2)/F_{N}$ and
$y=\pi/F_{N}+\theta_{1}\in[-\pi/F_{N},\pi/F_{N})$, and then the energies
become

$\begin{cases}E_{1}(x_{m})=-V\cos(x_{m}+y),\\\
E_{2}(x_{m})=-V\cos(x_{m}-y),\end{cases}$ (3)

The energies are naturally ordered:

$E_{1}(x_{m})<E_{1}(x_{m+1}),\leavevmode\nobreak\ \leavevmode\nobreak\
E_{2}(x_{m})<E_{2}(x_{m+1}).$ (4)

If $0<y\leq\frac{\pi}{F_{N}}$, the total energies are ordered by
$E_{2}(x_{1})<E_{1}(x_{1})<E_{2}(x_{2})<E_{1}(x_{2})<\cdots$, thus $\Delta
E_{1}(x_{m})=E_{1}(x_{m})-E_{2}(x_{m})$ and $\Delta
E_{2}(x_{m})=E_{2}(x_{m+1})-E_{1}(x_{m})$. Combining Eq. (3), we can obtain
that

$\begin{cases}\Delta E_{1}(x_{m})=2Vb_{1}\sin(x_{m}),\\\ \Delta
E_{2}(x_{m})=2Vb_{2}\sin(x_{m}+\frac{\pi}{F_{N}}).\end{cases}\\\ $ (5)

where $b_{1}=\sin y$ and $b_{2}=\sin(\pi/F_{N}-y)$. When $-\pi/F_{N}<y<0$,
$E_{1}(x_{1})<E_{2}(x_{1})<E_{1}(x_{2})<E_{2}(x_{2})<\cdots$, one can obtain
Eq. (5), and it still holds true, with the only difference being that
$b_{1}=-\sin y$ and $b_{2}=\sin(\pi/F_{N}+y)$.

We consider the limit $L\rightarrow\infty$ and hence $x_{m}\rightarrow x$. For
the sake of convenience and comparison, we consider $\delta
E_{1}(x_{m})=\frac{\Delta E_{1}(x_{m})}{\langle\Delta
E_{1}(x_{m})\rangle}=\frac{\Delta E_{1}(x_{m})}{2V/L}$, and combining Eq. (5),
we obtain

$\delta E_{1}(x_{m})=b_{1}L\sin(x_{m}).$ (6)

Similarly,

$\delta E_{2}(x_{m})=b_{2}L\sin(x_{m}+\frac{\pi}{F_{N}}).$ (7)

We note that $b_{1}$ and $b_{2}$ are of the order of $1/F_{N}$, so $b_{1}L$
and $b_{2}L$ are of the order of $1$. Therefore, the distribution of $\delta
E$ consists of two branches, as shown in Fig. 3(c), and these two branches
satisfy Eq. (6) and Eq. (7), respectively. Then one can calculate that the
total number of states for the energy smaller than $\delta E$ is $P_{1}(\delta
E_{1}\leq\delta E)=2P(0\leq x\leq\arcsin(\frac{\delta
E}{b_{1}L}))=\frac{2\arcsin(\frac{\delta E}{b_{1}L})}{\pi}$ [see Fig. 3(d)].
Hence the probability distribution is

$p_{1}(\delta E)=\frac{dP_{1}(\delta E_{1}\leq\delta E)}{d\delta
E}=\frac{2}{\pi\sqrt{(b_{1}L)^{2}-\delta E^{2}}}.$ (8)

Similarly, one can obtain $p_{2}(\delta
E)=\frac{2}{\pi\sqrt{(b_{2}L)^{2}-\delta E^{2}}}$. The total probability
distribution should be considered as the sum of the two.

Figure 4:  (a) level spacing’s distributions and (b) $P(r)$ of AA model with
$\alpha=(\sqrt{5}-1)/2$. The remaining parameters are the same as those in
Fig. 3(c) and Fig. 1(b4). (c) For $\alpha=(\sqrt{5}-1)/2$ (black data points)
and $\alpha=F_{N-1}/F_{N}$ (red and green data points), $r$ varies with
different initial phases $\theta$ and the strength of the quasiperiodic
potential $V$. The green/red dots are calculated from the lowest/middle one-
third of the spectrum in the central region [see Fig. 1(b1)], and the black
dots are obtained from the entire energy spectrum. (d) For
$\alpha=(\sqrt{5}-1)/2$ (red data) and $\alpha=F_{N-1}/F_{N}$ (blue data),
$\langle r\rangle$ averaged over $30$ samples and all energy levels varies
with changes in size.

We previously discussed the case of $V\gg J$. For the general case, it is
challenging to derive its expression. However, considering the statistically
unified form $P(\delta E)=Aexp(-B\frac{\delta E}{\langle\delta E\rangle})$ for
the level spacing distribution induced by disorder, with the fitting
parameters $A$ and $B$ changing with increasing disorder strength, we believe
that for AL induced by quasi-periodic potential, the level spacing
distribution should satisfy a unified form

$p(\delta E)=\sum_{\kappa}\frac{A_{\kappa}\Theta((b_{\kappa}L)^{2}-\delta
E^{2})}{\pi\sqrt{(b_{\kappa}L)^{2}-\delta E^{2}}}+C,$ (9)

where $\Theta$ is the step function explaintheta , $\kappa$ represents the
number of branches, $A_{\kappa}$ and $C$ are the undetermined parameters that
describe the scaling and translation of $P(\delta E)$ with $J=0$. In Fig.
1(b3), the purple and light blue dashed lines correspond to the fitting
results using Eq. (8) (corresponding to $J=0$) and Eq. (9), respectively. We
observe that Eq. (9) can effectively describe the level spacing distribution
$P(\delta E)$ of quasiperiodic AL.

We further derive the distribution of the adjacent gap ratio $r$ through the
use of the defining equation $P(r)=\int d(\delta E_{n},\delta
E_{n+1})\delta(r-\frac{\min\\{\delta E_{n},\delta E_{n+1}\\}}{\max\\{\delta
E_{n},\delta E_{n+1}\\}})p(\delta E_{n},\delta E_{n+1})$. For the case of
$J=0$, $\delta E_{n}$ and $\delta E_{n+1}$ respectively correspond to the two
purple lines in Fig. 3(c), which are described by Eq. (6) and Eq. (7).
Considering $\pi/F_{N}\rightarrow 0$, $\delta(r-\frac{\min\\{\delta
E_{n},\delta E_{n+1}\\}}{\max\\{\delta E_{n},\delta
E_{n+1}\\}})=\delta(r-\frac{\min\\{b_{1},b_{2}\\}}{\max\\{b_{1},b_{2}\\}})$
can be brought outside the integral, so

$P(r)=\delta(r-\frac{\min\\{b_{1},b_{2}\\}}{\max\\{b_{1},b_{2}\\}}).$ (10)

Thus, the distribution $P(r)$ is $\delta$ function, as shown in Fig. 1(b4),
which is clearly different from the $P(r)$ given by Poisson statistics. When
$\alpha=(\sqrt{5}-1)/2$, the number of branches of $\delta E$ (i.e., $\kappa$
in Eq. (9)) will exceed $2$ (see Supplementary Material SM ), as shown in Fig.
4(a). Consequently, multiple peaks appear in $P(r)$, as depicted in Fig. 4(b).
For the case of $J\neq 0$, using perturbation theory, we can demonstrate that
$b_{1}$ and $b_{2}$ in Eq. (5) need to be multiplied by the same factor SM .
Therefore, from Eq. (10), $P(r)$ is independent of both $J$ and $V$ in the AL
phase. Additionally, the expressions for $b_{1}$ and $b_{2}$ include the
initial phase $\theta$, hence the peak positions of $P(r)$ depend on $\theta$.
From Fig. 4(c), we observe that the values of $r$ are independent of $V$ and
the position in the energy spectrum. However, they depend on $\theta$ and on
whether $\alpha$ takes the value $(\sqrt{5}-1)/2$ or $F_{N-1}/F_{N}$.

Then we consider the sample average of $r$, which is equivalent to average
over $y$, i.e., $\langle
r\rangle=\frac{F_{N}}{\pi}\int_{y_{1}}^{y_{2}}dy\frac{\min\\{b_{1},b_{2}\\}}{\max\\{b_{1},b_{2}\\}}$
with $y_{1}=0$ and $y_{2}=\pi/F_{N}$ or $y_{1}=-\pi/F_{N}$ and $y_{2}=0$.
Substituting the expressions for $b_{1}$ and $b_{2}$, one easily obtains
$\langle r\rangle=2\ln 2-1\approx 0.387$ SM , as shown in Fig. 4(d). We note
that although the result of $\langle r\rangle$ is the same as that given by
Poisson distribution, it is not caused by Poisson statistics. When
$\alpha=(\sqrt{5}-1)/2$, $\langle r\rangle$ is a value close to $0.4$,
distinct from $0.387$ [Fig. 4(d)]. When the interaction is added, even with a
fixed $\theta$, the level spacing distribution of the quasiperiodic MBL phase
still follows Poisson statistics SM .

Conclusion and Discussion.— We have derived the level-spacing distribution of
quasiperiodic AL phase, which satisfies Eq. (9) and found that it does not
follow Poisson statistics, whereas quasiperiodic MBL’s spectrum conforms to
Poisson statistics. In addition, we found more differences in the spectrum
between the quasiperiodic and disordered AL phases. Specifically: (1) the
former exhibits different degrees of uniformity in level spacing across
different spectral regions, and the overall distribution is even more uneven
than a Poisson distribution, whereas the latter shows a relatively uniform
level distribution across different spectral regions; (2) the distribution of
$P(r)$ for the former is a $\delta$ function dependent only on the initial
phase, while the distribution of $P(r)$ for the latter follows
$P(r)=2/(1+r)^{2}$; (3) the sample-averaged value $\langle r\rangle$ for the
former depends on whether $\alpha=F_{N-1}/F_{N}$ or $\alpha=(\sqrt{5}-1)/2$.
Although, for the case of $\alpha=F_{N-1}/F_{N}$, the obtained $\langle
r\rangle$ is the same as that obtained from Poisson statistics, it does not
originate from Poisson statistics. Further, there are spatial correlations in
quasiperiodic systems, indicating that increasing the number of samples is not
equivalent to increasing the size, which is in contrast to Poisson statistics.
The energy spectrum of quasiperiodic systems can be experimentally determined
in various systems, such as semiconductor quantum dots Kiczynski2022  or
superconducting qubits Roushan2017 .

###### Acknowledgements.

We thank Qi Zhou for making us aware of Ref. Jitomirskaya . This work is
supported by National Key R&D Program of China under Grant No.2022YFA1405800,
the National Natural Science Foundation of China (Grant No.12104205), the Key-
Area Research and Development Program of Guangdong Province (Grant No.
2018B030326001), Guangdong Provincial Key Laboratory (Grant
No.2019B121203002).

## References

  * (1) P. W. Anderson, Phys. Rev. 109, 1492 (1958).
  * (2) P. A. Lee and T. V. Ramakrishnan, Disordered electronic systems, Rev. Mod. Phys. 57, 287 (1985).
  * (3) F. Evers and A. D. Mirlin, Anderson transitions, Rev. Mod. Phys. 80, 1355 (2008).
  * (4) B. Kramer and A. MacKinnon, Localization: theory and experiment, Rep. Prog. Phys. 56, 1469 (1993).
  * (5) C. M. Soukoulis and E. N. Economou, Localization in One-Dimensional Lattices in the Presence of Incommensurate Potentials, Phys. Rev. Lett. 48, 1043 (1981).
  * (6) S. Das Sarma, S. He, and X. C. Xie, Mobility edge in a model one-dimensional potential, Phys. Rev. Lett. 61, 2144 (1988).
  * (7) J. Biddle, B. Wang, D. J. Priour Jr, and S. Das Sarma, Localization in one-dimensional incommensurate lattices beyond the Aubry-André model, Phys. Rev. A 80, 021603 (2009); J. Biddle and S. Das Sarma, Predicted mobility edges in one-dimensional incommensurate optical lattices: an exactly solvable model of Anderson localization, Phys. Rev. Lett. 104, 070601 (2010).
  * (8) X. Li, X. Li, and S. Das Sarma, Mobility edges in one dimensional bichromatic incommensurate potentials, Phys. Rev. B 96, 085119 (2017); D. Vu and S. Das Sarma, Generic mobility edges in several classes of duality-breaking one-dimensional quasiperiodic potentials, Phys. Rev. B 107, 224206 (2023).
  * (9) H. Yao, H. Khoudli, L. Bresque, and L. Sanchez-Palencia, Critical behavior and fractality in shallow one-dimensional quasiperiodic potentials, Phys. Rev. Lett. 123, 070405 (2019).
  * (10) S. Ganeshan, J. H. Pixley, and S. Das Sarma, Nearest neighbor tight binding models with an exact mobility edge in one dimension, Phys. Rev. Lett. 114, 146601 (2015).
  * (11) Y. Wang, X. Xia, L. Zhang, H. Yao, S. Chen, J. You, Q. Zhou, and X.-J. Liu, One dimensional quasiperiodic mosaic lattice with exact mobility edges, Phys. Rev. Lett. 125, 196604 (2020).
  * (12) X.-C. Zhou, Y. Wang, T.-F. J. Poon, Q. Zhou, and X.-J. Liu, Exact new mobility edges between critical and localized states, Phys. Rev. Lett. 131, 176401 (2023).
  * (13) Y. Wang, L. Zhang, W. Sun, T.-F. J. Poon, and X.-J. Liu, Quantum phase with coexisting localized, extended, and critical zones, Phys. Rev. B 106, L140203 (2022).
  * (14) X. Deng, S. Ray, S. Sinha, G. Shlyapnikov, and L. Santos, One-dimensional quasicrystals with power-law hopping, Phys. Rev. Lett. 123, 025301 (2019).
  * (15) M. Gonçalves, B. Amorim, E. V. Castro, and P. Ribeiro, Hidden dualities in 1D quasiperiodic lattice models, SciPost Phys. 13, 046 (2022); M. Gonçalves, B. Amorim, E. V. Castro, and P. Ribeiro, Renormalization-Group Theory of 1D quasiperiodic lattice models with commensurate approximants, Phys. Rev. B 108, L100201 (2023); M. Gonçalves, B. Amorim, E. V. Castro, and P. Ribeiro, Critical phase dualities in 1D exactly-solvable quasiperiodic models, Phys. Rev. Lett. 131, 186303 (2023).
  * (16) G. Roati, C. D’Errico, L. Fallani, M. Fattori, C. Fort, M. Zaccanti, G. Modugno, M. Modugno, and M. Inguscio, Anderson localization of a non-interacting Bose-Einstein condensate, Nature (London) 453, 895 (2008).
  * (17) H. P. Lüschen, S. Scherg, T. Kohlert, M. Schreiber, P. Bordia, X. Li, S. D. Sarma, and I. Bloch, Single-particle mobility edge in a one-dimensional quasiperiodic optical lattice, Phys. Rev. Lett. 120, 160404 (2018);T. Kohlert, S. Scherg, X. Li, H. P. Lüschen, S. D. Sarma, I. Bloch, and M. Aidelsburger, Observation of many-body localization in a one-dimensional system with single-particle mobility edge, Phys. Rev. Lett. 122, 170403 (2019).
  * (18) F. A. An, E. J. Meier, and B. Gadway, Engineering a flux-dependent mobility edge in disordered zigzag chains, Phys. Rev. X 8, 031045 (2018); F. A. An, K. Padavić, E. J. Meier, S. Hegde, S. Ganeshan, J. H. Pixley, S. Vishveshwara, and B. Gadway, Observation of tunable mobility edges in generalized Aubry-André lattices, Phys. Rev. Lett. 126, 040603 (2021).
  * (19) Y. Wang, J.-H. Zhang, Y. Li, J. Wu, W. Liu, F. Mei, Y. Hu, L. Xiao, J. Ma, C. Chin, and S. Jia, Observation of Interaction-Induced Mobility Edge in an Atomic Aubry-André Wire, Phys. Rev. Lett. 129, 103401 (2022).
  * (20) T. Xiao, D. Xie, Z. Dong, T. Chen, W. Yi, and B. Yan, Observation of topological phase with critical localization in a quasi-periodic lattice, Science Bulletin 66, 2175 (2021).
  * (21) T. Shimasaki, M. Prichard, H. E. Kondakci, J. Pagett, Y. Bai, P. Dotti, A. Cao, T.-C. Lu, T. Grover, and D. M. Weld, Anomalous localization and multifractality in a kicked quasicrystal, arXiv:2203.09442.
  * (22) H. Li, Y.-Y. Wang, Y.-H. Shi, K. Huang, X. Song, G.-H. Liang, Z.-Y. Mei, B. Zhou, H. Zhang, J.-C. Zhang, et al., Observation of critical phase transition in a generalized Aubry-André-Harper model with superconducting circuits, npj Quantum Information 9, 40 (2023).
  * (23) Y. Wang, L. Zhang, S. Niu, D. Yu, X.-J. Liu, Realization and detection of non-ergodic critical phases in optical Raman lattice, Phys. Rev. Lett. 125, 073204 (2020).
  * (24) Y. Hatsugai and M. Kohmoto, Energy spectrum and the quantum Hall effect on the square lattice with next-nearest-neighbor hopping, Phys. Rev. B 42, 8282 (1990); J. H. Han, D. J. Thouless, H. Hiramoto, and M. Kohmoto, Critical and bicritical properties of Harper’s equation with next-nearest-neighbor coupling, Phys. Rev. B 50, 11365 (1994).
  * (25) G. T. Landi, D. Poletti, and G. Schaller, Nonequilibrium boundary-driven quantum systems: Models, methods, and properties, Rev. Mod. Phys. 94, 045006 (2022).
  * (26) M. Saha, S. K. Maiti, and A. Purkayastha, Anomalous transport through algebraically localized states in one dimension, Phys. Rev. B 100, 174201 (2019); A. Purkayastha, A. Dhar, and M. Kulkarni, Nonequilibrium phase diagram of a one-dimensional quasiperiodic system with a single-particle mobility edge, Nonequilibrium phase diagram of a one-dimensional quasiperiodic system with a single-particle mobility edge, Phys. Rev. B 96, 180204(R) (2017); M. Saha, B. P. Venkatesh, and B. K. Agarwalla, Quantum transport in quasiperiodic lattice systems in the presence of Büttiker probes, Phys. Rev. B 105, 224204 (2022).
  * (27) D. Dwiputra and F. P. Zen, Environment-assisted quantum transport and mobility edges, Phys. Rev. A 104, 022205 (2021).
  * (28) A. M. Lacerda, J. Goold, and G. T. Landi, Dephasing enhanced transport in boundary-driven quasiperiodic chains, Phys. Rev. B 104, 174203 (2021); V. Balachandran, S. R. Clark, J. Goold, and D. Poletti, Energy Current Rectification and Mobility Edges, Phys. Rev. Lett. 123, 020603 (2019); C. Chiaracane, M. T. Mitchison, A. Purkayastha, G. Haack, and J. Goold, Quasiperiodic quantum heat engines with a mobility edge, Phys. Rev. Research 2, 013093 (2020).
  * (29) T.-F. J. Poon, Y. Wan, Y. Wang, and X.-J. Liu, Anomalous quantum transport in 2D asymptotic quasiperiodic system, arXiv:2312.04349.
  * (30) D. A. Abanin, E. Altman, I. Bloch, and M. Serbyn, Colloquium: Many-body localization, thermalization, and entanglement, Rev. Mod. Phys. 91, 021001 (2019).
  * (31) D. M. Basko, I. L. Aleiner, and B. L. Altshuler, Metal-insulator transition in a weakly interacting many-electron system with localized single-particle states, Ann. Phys. 321, 1126 (2006).
  * (32) M. Schreiber, S. S. Hodgman, P. Bordia, H. P. Lüschen, M. H. Fischer, R. Vosk, E. Altman, U. Schneider, and I. Bloch, Observation of many-body localization of interacting fermions in a quasirandom optical lattice, Science 349, 842 (2015).
  * (33) B. I. Shklovskii, B. Shapiro, B. R. Sears, P. Lambrianides, and H. B. Shore, Statistics of spectra of disordered systems near the metal-insulator transition, Phys. Rev. B 47, 11487 (1993).
  * (34) A. D. Mirlin, Statistics of energy levels and eigenfunctions in disordered systems, Phys. Rep. 326, 259 (2000).
  * (35) S. A. Molčanov, The Local Structure of the Spectrum of the One-Dimensional Schrödinger Operator, Commun. Math. Phys. 78, 429 (1981).
  * (36) X. Li, J. H. Pixley, D.-L. Deng, S. Ganeshan, and S. Das Sarma, Quantum nonergodicity and fermion localization in a system with a single-particle mobility edge, Phys. Rev. B 93, 184204 (2016).
  * (37) Y. Liu, X.-P. Jiang, J. Cao, and S. Chen, Non-Hermitian mobility edges in one-dimensional quasicrystals with parity-time symmetry, Phys. Rev. B 101, 174205 (2020).
  * (38) Y. Wang, L. Zhang, Y. Wan, Y. He, and Y. Wang, Two-dimensional vertex-decorated Lieb lattice with exact mobility edges and robust flat bands, Phys. Rev. B 107, L140201 (2023).
  * (39) S. Cheng, R. Asgari, and G. Xianlong, From topological phase to transverse Anderson localization in a two-dimensional quasiperiodic system, Phys. Rev. B 108, 024204 (2023).
  * (40) S. Ray, M. Pandey, A. Ghosh and S. Sinha, Localization of weakly interacting Bose gas in quasiperiodic potential, New J. Phys. 18, 013013 (2016).
  * (41) J. H. Pixley, J. H. Wilson, D. A. Huse, and S. Gopalakrishnan, Weyl Semimetal to Metal Phase Transitions Driven by Quasiperiodic Potentials, Phys. Rev. Lett. 120, 207604 (2018).
  * (42) S. Schiffer, X.-J. Liu, H. Hu, and J. Wang, Anderson localization transition in a robust PT-symmetric phase of a generalized Aubry-André model, Phys. Rev. A 103, L011302 (2021).
  * (43) S. Ray, A. Ghosh, and S. Sinha, Drive-induced delocalization in the Aubry-André model, Phys. Rev. E 97, 010101(R) (2018).
  * (44) V. Khemani, D. N. Sheng, and D. A. Huse, Two Universality Classes for the Many-Body Localization Transition, Phys. Rev. Lett. 119, 075702 (2017).
  * (45) P. Roushan, C. Neill, J. Tangpanitanon, V.M. Bastidas, A. Megrant, R. Barends, Y. Chen, Z. Chen, B. Chiaro, A. Dunsworth, A. Fowler, B. Foxen, M. Giustina, E. Jeffrey, J. Kelly, E. Lucero, J. Mutus, M. Neeley, C. Quintana, D. Sank, A. Vainsencher, J. Wenner, T. White, H. Neven, D. G. Angelakis, J. Martinis, Spectroscopic signatures of localization with interacting photons in superconducting qubits, Science 358, 1175 (2017).
  * (46) X. Li, S. Ganeshan, J. H. Pixley, and S. D. Sarma, Many-Body Localization and Quantum Nonergodicity in a Model with a Single-Particle Mobility Edge, Phys. Rev. Lett. 115, 186601 (2015).
  * (47) R. Modak and S. Mukerjee, Many-Body Localization in the Presence of a Single-Particle Mobility Edge, Phys. Rev. Lett. 115, 230401 (2015).
  * (48) K. Machida and M. Fujita, Quantum energy spectra and one-dimensional quasiperiodic systems, Phys. Rev. B 34, 7367 (1986); M. Fujita and K. Machida, Spectral Properties of One-Dimensional Quasi-Crystalline and Incommensurate Systems, J. Phys. Soc. Jpn. 56, 1470 (1987).
  * (49) S. N. Evangelou and E. N. Economou, Spectral density correlations and eigenfunction fluctuations in one-dimensional quasi-periodic systems, J. Phys.: Condens. Matter 3 5499 (1991).
  * (50) N. Roy and A. Sharma, Study of counterintuitive transport properties in the Aubry-André-Harper model via entanglement entropy and persistent current, Phys. Rev. B 100, 195143 (2019).
  * (51) Y. Takada, K. Ino, and M. Yamanaka, Statistics of spectra for critical quantum chaos in one-dimensional quasiperiodic systems, Phys. Rev. E 70, 066203 (2004).
  * (52) Y. Wang, Y. Wang, and S. Chen, Spectral statistics, finite-size scaling and multifractalanalysis of quasiperiodic chain with p-wave pairing, Eur. Phys. J. B 89, 254 (2016).
  * (53) S. Jitomirskaya, Critical phenomena, arithmetic phase transitions, and universality: some recent results on the almost Mathieu operator, Current Developments in Mathematics (2019); A. Avila and S. Jitomirskaya, In preparation.
  * (54) S. Aubry and G. André, Analyticity breaking and Anderson localization in incommensurate lattices, Ann. Israel Phys. Soc. 3, 133 (1980).
  * (55) V. Oganesyan and D. A. Huse, Localization of interacting fermions at high temperature, Phys. Rev. B 75, 155111 (2007).
  * (56) A. Pal and D. A. Huse, Many-body localization phase transition, Phys. Rev. B 82, 174411 (2010).
  * (57) F. J. Dyson and M. L. Mehta, Statistical Theory of the Energy Levels of Complex Systems, J. Math. Phys. 4, 701 (1963).
  * (58) T. Guhr, A. Müller-Groeling, and H. A. Weidenmüller, Random-matrix theories in quantum physics: common concepts, Phys. Rep. 299, 189 (1998).
  * (59) C. L. Bertrand and A. M. García-García, Anomalous Thouless energy and critical statistics on the metallic side of the many-body localization transition, Phys. Rev. B 94, 144201 (2016).
  * (60) Y. Wang, C. Cheng, X.-J. Liu, and D. Yu, Many-body critical phase: extended and nonthermal, Phys. Rev. Lett. 126, 080602 (2021).
  * (61) When $\theta_{1}$ satisfy $\theta_{1}\in[-2\pi/F_{N},0)$, one can easily find that there always exists an integer $n_{0}$ such that $n_{0}=(\theta+\pi-\theta_{1})/(\frac{2\pi}{F_{N}})$ holds true.
  * (62) $\Theta(x)=0$ if $x\leq 0$ and $\Theta(x)=1$ if $x>0$. As shown in Fig. 3(d), when the line of $\delta E$ intersects both the red and blue curves, $P(\delta E)$ is the sum calculated from both branches. However, when $\delta E$ only intersects the blue curve, for the red curve, $(b_{\kappa}L)^{2}-(\delta E)^{2}<0$. Therefore, $P(\delta E)$ only includes the blue branch in this case.
  * (63) See Supplemental Material for details on (I) deriving the level-spacing distribution of the AA model with $\alpha=(\sqrt{5}-1)/2$; (II) Calculating level spacings using perturbation theory; (III) Calculating the sample-averaged value of the adjacent gap ratio; (IV) Level statistics of many-body quasiperiodic systems.
  * (64) M. Kiczynski, S. K. Gorman, H. Geng, M. B. Donnelly, Y. Chung, Y. He, J. G. Keizer, and M. Y. Simmons, Engineering topological states in atom-based semiconductor quantum dots, Nature (London) 606, 694 (2022).

## Supplementary Material:
Level spacing distribution of localized phases caused by quasiperiodic
potentials

In the Supplementary Materials, we first derive the level-spacing distribution
of the Aubry-André (AA) model with $\alpha=(\sqrt{5}-1)/2$. Then, we calculate
level spacings for the case of $J\neq 0$ using perturbation theory and
determine the sample-averaged value of the adjacent gap ratio. Finally, we
study the level statistics of many-body quasiperiodic systems.

For convenience, we here rewrite the AA model:

$H=J\sum_{j}(c_{j+1}^{\dagger}c_{j}+c_{j}^{\dagger}c_{j+1})+\sum_{j}V\sum_{j}\cos(2\pi\alpha
j+\theta)c_{j}^{\dagger}c_{j}.$ (S1)

###  .1 Level-spacing distribution of AA model with $\alpha=(\sqrt{5}-1)/2$

In this section we show the energies levels of the AA model for fixed $\theta$
in the limit $J=0$. The energies are

$E_{j}^{(0)}=V\cos(2\pi\alpha j+\theta),\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ j=0,1,2,\cdots,F_{N}-1.$ (S2)

with $\alpha=\frac{\sqrt{5}-1}{2}$ and $\theta\in[-\pi,\pi)$. Separating
$\alpha$ into two parts:

$\alpha=\frac{F_{N-1}}{F_{N}}+\beta_{N}\equiv\alpha_{N}+\beta_{N},$ (S3)

with $\alpha_{N}=\frac{F_{N-1}}{F_{N}}$ is the ratio of Fibonacci sequence,
and $\beta_{N}$ can be evaluated by

$\beta_{N}=\alpha-\frac{F_{N-1}}{F_{N}}=\alpha-\frac{\alpha^{1-N}-(-\alpha)^{N-1}}{\alpha^{-N}-(-\alpha)^{N}}\sim\frac{1+\alpha^{2}}{\alpha}(-1)^{N+1}\alpha^{2N},$
(S4)

so $\beta_{N}\rightarrow 0$ exponentially as $N\rightarrow\infty$. Hence the
energies in Eq. S2 can be represented as

$E_{j}^{(0)}=-V\cos(2\pi\frac{\gamma_{j}}{F_{N}}+2\pi\beta_{N}j+\theta+\pi).$
(S5)

with $\gamma_{j}=jF_{N-1}(mod\leavevmode\nobreak\ F_{N})$. Define
$n=jF_{N-1}(mod\leavevmode\nobreak\ F_{N})$, then we have
$j=nF_{N-1}(mod\leavevmode\nobreak\ F_{N})$ for even $N$. Hence the energies
can be relabeled by

$E_{n}^{(0)}=-V\cos[2\pi\frac{n}{F_{N}}+2\pi\beta_{N}\gamma_{n}+\theta+\pi],\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ n=0,1,\cdots,F_{N}-1.$ (S6)

Shifting the labels of energies $n=m-n_{0}$ with

$n_{0}=\frac{\theta+\pi-\theta_{1}}{\frac{2\pi}{F_{N}}},\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ n_{0}=1,\cdots,F_{N},$ (S7)

where $\theta_{1}$ is

$\theta_{1}=[\theta+\pi](mod\leavevmode\nobreak\
\frac{2\pi}{F_{N}})-\frac{2\pi}{F_{N}},\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \theta_{1}\in[-\frac{2\pi}{F_{N}},0).$ (S8)

Then we see

$E_{m}^{(0)}=-V\cos(2\pi\frac{m}{F_{N}}+2\pi\beta_{N}\gamma_{m-n_{0}}+\theta_{1}),\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ m=1,2,\cdots,F_{N}.$ (S9)

Now separating the energies levels into two parts. For
$m=1,\cdots,\frac{F_{N}}{2}$ we denote the energies by $E_{1,m}^{(0)}$; for
$m=\frac{F_{N}}{2}+1,\frac{F_{N}}{2}+1,\cdots,F_{N}$ we relabel $m\rightarrow
F_{N}+1-m$ and denote the energies by $E_{2,m}^{(0)}$, hence

$\begin{cases}E_{1,m}^{(0)}=-V\cos(2\pi\frac{m}{F_{N}}+2\pi\beta_{N}\gamma_{m-n_{0}}+\theta_{1}),\\\
E_{2,m}^{(0)}=-V\cos(2\pi\frac{m-1}{F_{N}}-2\pi\beta_{N}\gamma_{1-m-n_{0}}-\theta_{1}),\end{cases}$
(S10)

with $m=1,2,\cdots,\frac{F_{N}}{2}$. It is convenient to define the variable
$x_{m}=2\pi\frac{m-\frac{1}{2}}{F_{N}}$,
$y=\frac{\pi}{F_{N}}+\theta_{1}\in[-\frac{\pi}{F_{N}},\frac{\pi}{F_{N}})$ and
$z_{m}=2\pi\beta_{N}\gamma_{m}$, the energies are

$\begin{cases}E_{1}^{(0)}(x_{m})=-V\cos(x_{m}+y+z_{m-n_{0}}),\\\
E_{2}^{(0)}(x_{m})=-V\cos(x_{m}-y-z_{1-m-n_{0}}).\end{cases}$ (S11)

Now we first consider the limit $\beta_{N}=0(z=0)$ and suppose $F_{N}$ is even
for simplification. The energies are naturally ordered:

$E_{1}^{(0)}(x_{m})<E_{1}^{(0)}(x_{m+1}),\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ E_{2}^{(0)}(x_{m})<E_{2,m}^{(0)}(x_{m+1}).$ (S12)

Hence the full energies is ordered by

$\begin{cases}if&\leavevmode\nobreak\ \leavevmode\nobreak\ 0\leq
y\leq\frac{\pi}{F_{N}},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
E_{2}^{(0)}(x_{1})<E_{1}^{(0)}(x_{1})<E_{2}^{(0)}(x_{2})<E_{1}^{(0)}(x_{2})<\cdots.\\\
if&\leavevmode\nobreak\ \leavevmode\nobreak\
-\frac{\pi}{F_{N}}<y<0,\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\
E_{1}^{(0)}(x_{1})<E_{2}^{(0)}(x_{1})<E_{1}^{(0)}(x_{2})<E_{2}^{(0)}(x_{2})<\cdots.\end{cases}$
(S13)

The level spacing for $0\leq y\leq\frac{\pi}{F_{N}}$ satisfy

$\begin{split}if&\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ 0\leq
y\leq\frac{\pi}{F_{N}},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\begin{cases}\Delta
E_{1}^{(0)}(x_{m})=E_{1}^{(0)}(x_{m})-E_{2}^{(0)}(x_{m}),\\\ \Delta
E_{2}^{(0)}(x_{m})=E_{2}^{(0)}(x_{m+1})-E_{1}^{(0)}(x_{m}).\end{cases}\\\
if&\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ -\frac{\pi}{F_{N}}<y<0,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\begin{cases}\Delta
E_{1}^{(0)}(x_{m})=E_{2}^{(0)}(x_{m})-E_{1}^{(0)}(x_{m}),\\\ \Delta
E_{2}^{(0)}(x_{m})=E_{1}^{(0)}(x_{m+1})-E_{2}^{(0)}(x_{m}).\end{cases}\end{split}$
(S14)

After a straightforward calculation, we can get

$\begin{cases}\delta E_{1}^{(0)}(x_{m})=b_{1}F_{N}\sin(x_{m}),\\\ \delta
E_{2}^{(0)}(x_{m})=b_{2}F_{N}\sin(x_{m}+\frac{\pi}{F_{N}}).\end{cases}$ (S15)

with $\delta E=F_{N}\Delta E/2V$. We see that the level spacings lie on two
sine functions with amplitude $2Vb_{1/2}$ and

$\begin{split}if&\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ 0\leq
y\leq\frac{\pi}{F_{N}},\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ b_{1}=\sin y,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ b_{2}=\sin(\frac{\pi}{F_{N}}-y),\\\
if&\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ -\frac{\pi}{F_{N}}<y<0,\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ b_{1}=-\sin
y,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
b_{2}=\sin(\frac{\pi}{F_{N}}+y).\end{split}$ (S16)

We see the amplitudes satisfy

$b_{1}+b_{2}=\frac{\pi}{F_{N}}+O(F_{N}^{-3}).$ (S17)

Figure S1: (a)Scheme of energies levels for non vanishing $\beta_{N}$.
(b)Level spacing for $\alpha=\frac{\sqrt{5}-1}{2}$. The amplitude $b_{j}$s are
calculated from Eq. S19.

For the case $\beta_{N}\neq 0$, there is another phase shift
$z_{m}=2\pi\beta_{N}\gamma_{m}$ in Eq. S10 except $\theta_{1}$. In the large
$N$ limit, the Fibonacci number $F_{N}=\frac{\alpha^{-N}}{\sqrt{5}}$, then the
phase shift can be evaluated from Eq. S4 that
$|2\pi\beta_{N}F_{N}|\sim\frac{1+\alpha^{2}}{5\alpha}(-1)^{N+1}\frac{2\pi}{F_{N}}\sim
0.8944\frac{\pi}{F_{N}}$, which is in the same order of $\frac{\pi}{F_{N}}$.
Hence the order of energies Eq. S12 remains while the order of full energies
Eq. S13 can break as shown in Fig. S1(a). Here we show an example for
$\theta=0.4\pi$ and $N=18$. Then we can derive $n_{0}=1809$ and
$y=0.6\frac{\pi}{F_{N}}$. All the energy level spacings lie on

$\begin{split}|\Delta
E_{1}|=&|E_{1}(x_{m})-E_{2}(x_{m})|=2Vb_{1}\sin(x_{m}+\frac{z_{m-n_{0}}+z_{1-m-n_{0}}}{2}),\\\
|\Delta
E_{2}|=&|E_{2}(x_{m+1})-E_{1}(x_{m})|=2Vb_{2}\sin(x_{m}+\frac{\pi}{F_{N}}+\frac{z_{m-n_{0}}+z_{-m-n_{0}}}{2}),\\\
|\Delta
E_{3}|=&|E_{1}(x_{m+1})-E_{2}(x_{m})|=2Vb_{3}\sin(x_{m}+\frac{\pi}{F_{N}}+\frac{z_{1+m-n_{0}}-z_{1-m-n_{0}}}{2}),\\\
|\Delta
E_{4}|=&|E_{1}(x_{m+1})-E_{1}(x_{m})|=2Vb_{4}\sin(x_{m}+y+\frac{z_{1+m-n_{0}}+z_{m-n_{0}}}{2}),\\\
|\Delta
E_{5}|=&|E_{2}(x_{m+1})-E_{2}(x_{m})|=2Vb_{5}\sin(x_{m}+\frac{\pi}{F_{N}}-y-\frac{z_{-m-n_{0}}+z_{1-m-n_{0}}}{2}),\end{split}$
(S18)

where the amplitude

$\begin{split}b_{1}=&|\sin(y+\frac{z_{m-n_{0}}+z_{1-m-n_{0}}}{2})|,\\\
b_{2}=&|\sin(\frac{\pi}{F_{N}}-y-\frac{z_{-m-n_{0}}+z_{m-n_{0}}}{2})|,\\\
b_{3}=&|\sin(\frac{\pi}{F_{N}}+y+\frac{z_{1+m-n_{0}}+z_{1-m-n_{0}}}{2})|,\\\
b_{4}=&|\sin(\frac{\pi}{F_{N}}+\frac{z_{1+m-n_{0}}-z_{m-n_{0}}}{2})|,\\\
b_{5}=&|\sin(\frac{\pi}{F_{N}}+\frac{z_{1-m-n_{0}}-z_{-m-n_{0}}}{2})|.\end{split}$
(S19)

It can be seen that $b_{4}=b_{5}$ is independent of $y$ and $\theta$. From the
identity

$z_{m+a}+z_{-m+b}=2\pi\beta_{N}kF_{N}+z_{a+b},\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ k=1,2,$ (S20)

where $z_{a+b}$ is a constant independent of $m$, so each $b_{j}$ has two
values, and the amplitude can be several of these $b_{j}$s. In Fig. S1(b) we
show the level spacing as a function of $x$, and see that there are five
relative $b_{j}$s. As the adjacent gap $P(r)$ is the ratio of nearest level
spacing, more delta peaks emerge as shown in the main text. After averaging
over $y$, the mean adjacent gap is around $0.4$ from numerical simulation.

###  .2 Energies levels of AA model from perturbative perspective

In this section, we discuss the energies levels for finite hopping amplitude
from perturbation theory. Separating the Hamiltonian into two parts

$H=H_{J}+H_{V},$ (S21)

with $H_{J/V}$ denoting the kinetic/potential energies with
$|H_{J}|\ll|H_{V}|$. Up to the second order of $J/V$, the energies are

$E_{j}=E_{j}^{(0)}+\sum_{l\neq j}\frac{|\langle
j|H_{J}|l\rangle|^{2}}{E_{j}^{(0)}-E_{l}^{(0)}}=E_{j}^{(0)}+\frac{J^{2}E_{j}^{(0)}}{[E_{j}^{(0)}]^{2}-V^{2}\cos^{2}(\pi\alpha)}.$
(S22)

where $E_{j}^{(0)}$ is the unperturbated energy in Eq. S2. We see two gaps are
opened at $E_{gap}=\pm V\cos(\pi\alpha)$ and the perturbation breaks down near
the band edge. Near the center the band $E_{j}\sim 0$, the energies are
renormalized by $1-\frac{J^{2}}{V^{2}\cos^{2}(\pi\alpha)}$. In Fig.
S2(a1)-(a3) we show energies levels for different $V$. The perturbation works
well near the center of the band even for small qausi-periodic potential.

Figure S2: (a1)-(a3)The perturbative energy levels. (b1)-(b3)The perturbative
level spacings. The hopping amplitude $J=1$, the phase factor $\theta=0.2\pi$,
and the system size $L=1597$. The strength of qausi-periodic potential is
$V=10$ for (a1)(b1), $V=7$ for (a2)(b2), and $V=4$ for (a3)(b3) respectively.

Now we move to the behavior of level spacing. Substituting Eq. S22 into Eq.
S14 we see the leading order of level spacing are

$\begin{split}\delta E_{1}(x_{m})=&b_{1}F_{N}\sin
x_{m}\bigg{[}1-\frac{J^{2}}{V^{2}}\frac{\cos^{2}(\pi\alpha)+\cos^{2}x_{m}}{[\cos^{2}(\pi\alpha)-\cos^{2}x_{m}]^{2}}\bigg{]},\\\
\delta E_{2}(x_{m})=&b_{2}F_{N}\sin
x_{m}\bigg{[}1-\frac{J^{2}}{V^{2}}\frac{\cos^{2}(\pi\alpha)+\cos^{2}x_{m}}{[\cos^{2}(\pi\alpha)-\cos^{2}x_{m}]^{2}}\bigg{]}.\end{split}$
(S23)

We see that the level spacing is modified by a function dependent on $x_{m}$,
and could divergence near $\cos x_{m}=\pm\cos(\pi\alpha)$, but the ratio
$\frac{\delta E_{1}}{\delta E_{2}}=\frac{b_{1}}{b_{2}}$ is a constant for
fixed phase. Hence the adjacent gap $P(r)$ is almost $J$ independent. In Fig.
S2(b1)-(b3) we show the level spacing as a function of $x$, we see that the
perturbation works well even near the band edge.

###  .3 The derivation of the average value of $r$

As mentioned in the main text, the sample average of $r$ is equivalent to
averaging over $y$. We set $\alpha=\frac{F_{N-1}}{F_{N}}$ and suppose $0\leq
y\leq\frac{\pi}{F_{N}}$. When $L\rightarrow\infty$, $b_{1}=\sin y\sim y$ and
$b_{2}=\sin(\pi/F_{N}-y)\sim\pi/F_{N}-y$, so

$\begin{split}\langle
r\rangle=&\frac{F_{N}}{\pi}\int_{0}^{\frac{\pi}{F_{N}}}dy\frac{\min\\{b_{1},b_{2}\\}}{\max\\{b_{1},b_{2}\\}}\\\
=&\frac{F_{N}}{\pi}(\int_{0}^{\frac{\pi}{2F_{N}}}dy\frac{y}{\frac{\pi}{F_{N}}-y}+\int_{\frac{\pi}{2F_{N}}}^{\frac{\pi}{F_{N}}}dy\frac{\frac{\pi}{F_{N}}-y}{y})\\\
=&2\ln 2-1.\end{split}$ (S24)

We observe that the result is consistent with the Poisson distribution, but it
is not caused by Poisson statistics.

###  .4 level statistics of many-body quasiperiodic systems

On the AA model, we introduce nearest-neighbor interactions
$U\sum_{j}n_{j}n_{j+1}$ with a fixed initial phase $\theta$. From Fig. S3, it
can be observed that its level-spacing distribution follows Poisson
statistics. Therefore, the localized phases in many-body quasiperiodic systems
exhibit Poisson statistics.

Figure S3:  The distribution of adjacent gap for $J=1$, $V=5$, $\theta=0.3\pi$
and $U=2$. The system size $L=16$ and the number of particles $N=8$. The
fitting function $P(\delta E)=1.742e^{-1.903\delta E}$.


# SuperCaustics: Real-time, open-source simulation of transparent objects for
deep learning applications

1st Mehdi Mousavi  Department of Computer Science
Georgia State University
Atlanta, United States
<EMAIL_ADDRESS>2nd Rolando Estrada
Department of Computer Science
Georgia State University
Atlanta, United States
<EMAIL_ADDRESS>
###### Abstract

Transparent objects are a very challenging problem in computer vision. They
are hard to segment or classify due to their lack of precise boundaries, and
there is limited data available for training deep neural networks. As such,
current solutions for this problem employ rigid synthetic datasets, which lack
flexibility and lead to severe performance degradation when deployed on real-
world scenarios. In particular, these synthetic datasets omit features such as
refraction, dispersion and caustics due to limitations in the rendering
pipeline. To address this issue, we present SuperCaustics, a real-time, open-
source simulation of transparent objects designed for deep learning
applications. SuperCaustics features extensive modules for stochastic
environment creation; uses hardware ray-tracing to support caustics,
dispersion, and refraction; and enables generating massive datasets with
multi-modal, pixel-perfect ground truth annotations. To validate our proposed
system, we trained a deep neural network from scratch to segment transparent
objects in difficult lighting scenarios. Our neural network achieved
performance comparable to the state-of-the-art on a real-world dataset using
only 10% of the training data and in a fraction of the training time. Further
experiments show that a model trained with SuperCaustics can segment different
types of caustics, even in images with multiple overlapping transparent
objects. To the best of our knowledge, this is the first such result for a
model trained on synthetic data. Both our open-source code and experimental
data are freely available online.

###### Index Terms:

Deep Learning, Computer Vision, Synthetic Data, Dataset, AI Simulation,
Virtual Environment, Caustics

##  I Introduction

Detecting transparent objects is one of the most challenging problems in
computer vision, because these objects do not form disjoint boundaries with
their environment. In addition, different types of transparent objects have
unique characteristics based on the materials they’re made of, and they come
in varied thicknesses and densities that affects their appearance. Light
passing through a transparent object will bounce multiple times before it
reaches our eyes (or a camera), causing effects such as refraction, specular
reflections, caustics and dispersion; and the dynamic nature of light makes
them appear in radically different ways when exposed to light from different
angles.

Figure 1: Images rendered in SuperCaustics. Changing image features in real
time is trivial in SuperCaustics. Note the changes in lighting, sharpness of
caustics, realistic backdrop and specular reflections on the transparent
objects.

Deep learning has flourished into a rapidly thriving field over the last
decade. In particular, it has achieved state-of-the-art results in problems
for which large amounts of (usually labelled) data is available. However,
there has been relatively little work on applying deep learning to transparent
object detection in part because gathering and labeling a sufficiently large
dataset for this problem is cumbersome and difficult. Most recent publications
on transparent object detection have utilized small, custom-made datasets of
real data, gathered and labelled by the authors themselves to validate or test
their algorithm or DCNN model. The actual training data consists primarily of
large volumes of synthetic images created in 3D modelling software (e.g.,
Blender [1]).

An important limitation of these rigid image datasets is their static nature.
For instance, once an image has been generated, one cannot change its
illumination, or replace the reflection maps in the environment. For real
images, labeling the ground truth is also a significant challenge, especially
for cameras with high pixel counts.  Depth, in particular, is almost
impossible to label by hand, so one is limited to active depth sensing (e.g.,
RGB-D sensors) or passive methods like disparity matching. However, depth
cameras often register errors when encountering transparent objects due to
their reflective and refractive properties [2]. In computer vision
applications where transparent objects are involved —such as robotic grasping—
precise labeling is an integral part of the entire pipeline, and having fault-
free data is crucial to the function of a robotics grasp algorithm. For
example, the Cleargrasp pipeline [2], the current state-of-the-art for deep-
learning based transparent object detection, uses a segmentation map to
complement and remove faulty depth signals from a depth image.

Moreover, no matter how extensive and large a pre-made training dataset is, it
can only cover so many lighting conditions. The positions of objects, type of
surfaces, camera angles, and placement of light sources will only cover a
small subset of all possible combinations. This is a challenge for transparent
objects because their appearance is highly dependent on how they are lit. For
example, they might cast shadows, caustics, or a combination of the two simply
by varying the angle or intensity of the light sources. As such, models
trained on synthetic datasets with only a small range of lightning scenarios
suffer performance degradation when deployed on real-world data.

More specifically, as we review in Section II, clever deep learning solutions
are often stifled by rigid datasets that cannot represent the full possibility
space of the target domain, are difficult to change after acquisition, require
significant human effort to produce or use proprietary technology for which
the source-code and project files used to generate them are rarely released
for free use by the public.

To address this gap, we propose SuperCaustics, an easy to use, highly
customizable, easily extendable, open-source real-time simulation tool
designed for generating dynamic, massive datasets in highly complex scenarios.
As we detail further in Section III, SuperCaustics is a modular set of
customizable classes developed within NVIDIA’s RTX Branch of Unreal Engine 4
(Epic Games, USA)[3, 4], allowing even those researchers who are unfamiliar
with Unreal Engine to create, maintain and iterate on their own photorealistic
customized datasets for vision tasks involving transparent objects.
SuperCaustics features four modules:

  1. 1.

A stochastic scene generation system that creates the geometry of a scene
based on parameters like object type, number of objects, camera type,
lighting, and background.

  2. 2.

A prop manager that adds or removes customizable background items to increase
scene complexity and variety.

  3. 3.

A ground truth core that automatically overlays a wide variety of accurate,
pixelwise ground-truth annotations.

  4. 4.

A data ablation core [5] that enables changing characteristics of the objects,
background, environment, lighting, or rendering without affecting scene
composition.

Figure 2: Showcase of some of the available ground truth images in
SuperCaustics. First row (from left to right): caustics-enabled image, camera-
space normals, world-space normals. Second row (from left to right): caustics
segmentation, depth image, segmentation mask. (Image best viewed on screen.)

Our proposed system allows for fine granular control over small details of the
scene before and after generation. Each module has a set of customizable and
extendable controls, including tools for data ablation. As detailed in [5],
data ablation is the study of the effects of isolated changes in the features
of data on the performance of a deep learning network (e.g. having the exact
same image under a different lighting condition). Using the Data Ablation
Core, we are able to change the characteristics of a scene in real-time, e.g.,
different lighting color, presence of sharp or soft caustics, light source
angle, texture or shape of the backdrop, or changes in roughness/color of the
glass surfaces (see Fig. 4). The tools in the Data Ablation Core can also be
used to increase variety in the captured images.

A virtual dataset generated from SuperCaustics contains the features of the
target domain, free from acquisition errors or labeling bias. More
importantly, as our experiments confirm, a small, targeted, carefully crafted
data-set can efficiently match the performance achieved by state-of-the-art
methods that train their systems on massive datasets.

Figure 3: As illustrated, our system has four components: (1) Scene Generator
Module; (2) Prop Manager Module; (3) Ground-truth Core; (4) Data Ablation
Core. Each of these modules comes with its own sub-modules that enable
generation of a customizable scene in real-time. We use Blueprint, Unreal’s
visual scripting language, for the ground-truth annotations and data ablation
controls. We use a separate Python interaction module—Controller—for data
collection.

##  II Related Work

The old and the new - A problem that persists: Recognizing transparent objects
has been a difficult problem to solve in vision. Early solutions used general
features of transparent objects, e.g., specular reflections, for object
recognition [6]. More sophisticated approaches include estimating bounding
boxes based on the refraction distortions that transparent objects’ edges
create on the background [7]. Newer research trains deep convolutional neural
networks on massive datasets to reach a pixel-wise labeling of transparent
objects [2, 8]. Also, some researchers have devised clever methods to exploit
sensor failures in depth images to localize transparent objects [9].
Regardless, often times researchers are forced to use low-quality real data
for testing their algorithm or training a neural network [2, 7], and other
times researchers resort to creating rigid virtual datasets to train large
models, such as a DCNN. [8, 9, 2].

Synthetic Data:  Synthetic data has proven useful in various computer vision
tasks, such as depth estimation [10, 2], surface normal estimation [5, 2],
robotic grasping [11], and object segmentation [12] by multiple independent
researchers [13, 12, 14, 5]. Modern computer graphics can achieve near-
photorealism, so synthetic data has become a viable alternative in situations
where acquiring or labeling real data is difficult. However, there are still
few datasets that meet the image quality required to be considered photo-
realistic [8], and those that do, rarely include transparent objects because
of the complexities they introduce in rendering and generating correct ground-
truth labels. More importantly, to the best of our knowledge, there is no
photorealistic, synthetic dataset that includes the tools and project files
needed to generate the data itself. More often, these synthetic datasets are
static and use proprietary code and technologies that does not allow for
reproduction or change in the dataset. We believe flexibility in creating the
image data is crucial since it allows for greater degrees of freedom when
evaluating a specific method, and allows researchers to make corrections in
the data and address weak spots in the model’s performance.

The need for reliable data: One common factor in transparent object detection
is the lack of reliable data. In older research, we see the use of personally-
gathered real data that suffers from lack of variety, noise, or lack of
labels. To get past this obstacle, we see various works generate their own
synthetic datasets using tools like Blender [2, 1] or NVIDIA Omniverse [8,
15]. As a result, these datasets consist of huge collections of static images.
In addition, the developers of these datasets rarely release the project files
or source code used to generate the data. Nevertheless, in most cases the
target domain of the final algorithm is the real world. As such, to be able to
generalize well into the target domain these datasets have to have tens of
thousands of images, at least, making the training very expensive and time
consuming. Despite all this data, we still see significant performance
degradation when shifting from crisp, high-resolution synthetic data to
blurry, grainy real images that are lit entirely differently. These synthetic
datasets are non-changeable (static) and more often than not, they do not
include graphically intense phenomena that happen in real-world scenarios
(e.g., dispersion or caustics), which makes them very narrow in domain and
makes the trained algorithms sensitive to these lighting conditions. [8, 2].

SuperCaustics vs. Cleargrasp: Arguably, the most significant effort in
modeling transparent objects is Cleargrasp [2]. In this study, researchers
trained deep learning models to estimate depth, surface normals, and occlusion
boundaries using a mix of synthetic and real data. Their results indicate that
in addition to a significant performance loss when shifting from a synthetic
to a real domain, the models consistently misclassify caustics. In particular,
they wrongly classify caustics as separate transparent objects. The authors of
Cleargrasp speculated that this was because the same amount of sharp caustics
were not present in the synthetic dataset due to limitations in the rendering
software. A more flexible simulation tool, such as our proposed system, allows
for addition or removal of such features. In particular, SuperCaustics can
produce images with sharp caustics. More generally, our system can incorporate
features of the real domain in the training data to make trained models more
robust to such features in the target domain.

Overall, prior work on transparent object detection has been hampered by a
lack of reliable, flexible data. To address this gap, we have developed
SuperCaustics, an extendable, open-source, user-friendly simulation tool that
allows for fine controls over every variable in the simulation. Our system can
generate millions of unique images out of the box. In addition, it allows
researchers to import 3D models of their choice, set their desired parameters,
and capture synthetic data from fully customizable virtual environments. We
describe our system in more detail, below.

##  III SuperCaustics

SuperCaustics is a simulation tool made in Unreal Engine for generating
massive computer vision datasets that include transparent objects. Unreal
Engine is the engine of choice for projects with high-resolution, real-time 3D
graphics. It is free for both commercial and non-commercial use and its source
code is publicly available and extended by the community. We use one such
extended version created for hardware raytracing by NVIDIA RTX graphics cards
[4].

Figure 4: SuperCaustics allows fine controls over fundamental characteristics
of every variable in the simulation. e.g., controlling the softness of the
light-source and its caustics. Left: Very soft caustics, Middle: Moderately
soft, Right: Very Sharp.

The key contributions of SuperCaustics are: 1. Open-source, free set of user-
friendly, extendable tools for creating custom datasets. 2. Customizable
photorealistic scenes using real-time hardware ray-tracing, 3. Automatic,
multimodal, accurate ground-truth generation for visual tasks, 4. Python and
UE4 interplay for data acquisition and processing.

Below, we detail each major module of SuperCaustics, and provide explanations
of the different ground-truth annotations available in SuperCaustics.

###  III-A Generator Module

The Generator Module creates a scene based on given parameters like item type,
shape, number of items and/or range of items. It leverages the physics engine
inside UE4 to drop the items on the scene within a customizable distance
offset. A random physics impulse can be added to the objects upon spawning to
increase randomness and add variety to how they come to rest. The intensity of
the impulse is also customizable.

###  III-B Props Manager

Similar to the Generator Module, a separate module manages the background
items in the simulation. This module accepts different 3D objects as prop
inputs, and manages their position and rotation within the first few frames of
the simulation until they come to rest by the physics engine.

###  III-C Data Ablation Core

Data Ablation is a powerful method for determining weak points and improving
robustness in a neural network’s performance [5]. The Data Ablation Core in
SuperCaustics allows for such experiments within the simulation. Additionally,
the features of the Data Ablation Core can be used to further increase the
variety in the final dataset. We can adjust the cameras (camera matrix),
lighting, backdrops, reflection profiles (HDRI maps), and properties of the
glass material (color, specular, thickness, roughness, texture). Additionally,
all properties of the ray-tracing engine are exposed for experimentation
within the Data Ablation module.

###  III-D Ground Truth Core

SuperCaustics also includes scripts for automatically estimating depth,
surface normals (world space and camera space), object masks, outlines,
occlusion boundaries and caustics segmentation. The Ground Truth Core can be
readily extended by adding additional scripts (or modifying the available
scripts). We use post-processing shaders (dubbed materials in UE4), to overlay
these properties over the image, enabling pixel-perfect alignment between the
data and the ground-truth labels. (See Figure 2) We detail each ground-truth
mode below:

####  III-D1 Depth Ground Truth

We calculate the normalized distance between each pixel that belongs to a
specific object and the camera. By default, we set the real-life range of
depth to 10 meters, which covers the entire environment. This range is
customizable in the Ground Truth Core.

####  III-D2 Surface Normals Ground Truth

We estimate the normal vector w.r.t to each 3D surface, then color each pixel
to indicate the vector’s direction. We use 6 main colors to show 6 axis of
direction (positive and negative xyz). Surface normals can be calculated in 2
different spaces: camera-space and world-space. Normals of a pixel in world-
space are defined with respect to the world Cartesian coordinate system. In
world-space, each normal vector points towards a fixed axis in space,
regardless of where it is being observed from. Camera-space is the space in
which points are defined with respect to the camera coordinate system. Each
normal vector has the transformation of the camera itself applied to it,
making the camera the origin point of the coordinate system. This means the
transformations of the camera in world-space is implicitly applied to every
normal vector in the image. This is a simplification compared to world-space
surface normals. (See Figure 2.)

TABLE I: Models trained with SuperCaustics Data generalize well to novel
scenarios

Experiment | Scores Obtained
---|---
Training | Testing | Accuracy | F1 Score | Precision | Recall | IOU
SuperCaustics | SuperCaustics (Known Objects) | 0.9924 | 0.9361 | 0.9975 | 0.8818 | 0.8801
SuperCaustics-R | Cleargrasp-Real (Novel) | 0.9560 | 0.6941 | 0.5627 | 0.9056 | 0.5316
ClearGrasp | Cleargrasp-Real (Novel) | NR | NR | NR | NR | 0.5800

####  III-D3 Object Mask Ground Truth

In Unreal Engine, it is easy to map visible pixels to their corresponding 3D
objects. Our Blueprint script uses this mapping to automatically detect and
overlay pixel-perfect semantic labels on the objects in the scene.

####  III-D4 Outlines Ground Truth

The edges of transparent objects are sometimes hard to detect visually. In
this Ground Truth pass, we can highlight the edges of objects with a
customizable thickness (e.g., 1px or 20px).

####  III-D5 Occlusion Boundaries Ground Truth

This Ground Truth pass shows surfaces where transparent objects are touching
another surface (like the backdrop). Specifically, this pass shows where there
is a depth discontinuity. The logic of this ground truth pass is adapted from
Cleargrasp [2].

####  III-D6 Local & Non-local Caustics Mask Ground Truth

In SuperCaustics, we have fine control over the characteristics of the image,
and we can create the exact same image under different settings. We use the
difference in the luminosity between two rendered images (one with caustics,
and another without) to reach a pixel-wise segmentation of local caustics
(inside transparent object) and non-local caustics (projected unto another
surface) in the image.

##  IV Evaluation

To validate the usefulness of SuperCaustics, we trained a deep convolutional
neural network (DCNN) to do a pixel-wise labeling of transparent objects with
both soft and sharp caustics. We chose pixel-wise labeling because it was more
suitable given the limited hardware we had available. To demonstrate the
flexibility of SuperCaustics, we aimed to test our model on the real test set
released by Cleargrasp. To do that, we set the angles of the simulated Intel
real-sense cameras to resemble similar perspectives and based the parameters
of our simulation on the characteristics of the Cleargrasp’s real-world
dataset (i.e., occasional sharp caustics from hard lights, general soft
lighting, presence of a tote box and wooden and cloth backdrops among
synthetic backgrounds). Below, we compare our network’s performance to the
network proposed in Cleargrasp [2]. Our results indicate that our network
achieved comparable results with much fewer training data.

###  IV-A Experimental Setup

####  IV-A1 Hardware

We conducted all our experiments in a Dell Precision 7920R server with two
Intel Xeon Silver 4110 CPUs, two GeForce GTX 1080 Ti graphics cards, and 128
GBs of RAM.

####  IV-A2 Image acquisition

Using the features described in Section III, we generated 9000
1920$\times$1080 synthetic images from 650 randomized scenarios. In each
scenario, a reflection map and a backdrop was randomly chosen from a bank of
33 HDRI mappings and 33 backdrops. Then, a number of transparent objects were
dropped into the scene from off-camera, so that they came to rest naturally
using UE4’s physics engine. Then, the Prop Manager module added the input
props in random locations inside the scene. For gathering the images, 12 Intel
Realsense cameras were simulated in various locations and angles. In every
camera angle, the Data Ablation Core rotated the main light of the scene to
cast shadows and caustics at various angles before capturing an image and its
ground-truth values. For training, the images were split in the following
ratio: 8000 for training and 1000 for validation.

####  IV-A3 Deep Neural Networks

We used an implementation of U-net [16] as our segmentation neural network. In
particular, we used the tools provided by the EasyTorch Library to prototype
and run the experiments. EasyTorch is a Pytorch-based deep learning library
used in [17, 18, 19]. We evaluated our performance on mean intersection-over-
union (IOU) for the test set. For the loss function, we used the negative log
likelihood loss (NLL). We trained our network for 30 epochs on purely
synthetic data from SuperCaustics. Afterwards, we used the SuperCaustics model
as a pre-trained starting point (labeled SuperCaustics-R in our results),
which was trained on Cleargrasp Known Objects for an additional four epochs.

TABLE II: Caustics Segmentation

Experiment | Scores Obtained
---|---
Training | Testing | Accuracy | F1 Score | Precision | Recall | IOU
Caustic Segmentation | Caustic Segmentation (Known Objects) | 0.9981 | 0.9386 | 0.9541 | 0.9237 | 0.8844
Caustic Segmentation | Caustic Segmentation (Novel Objects) | 0.9958 | 0.8102 | 0.8682 | 0.7594 | 0.6809

###  IV-B Results & Discussion

Table I shows the object segmentation performance of our network trained on
SuperCaustics data, and Table II shows the equivalent results when segmenting
caustics directly. In particular, the first table shows that our
SuperCaustics-R model (which was also trained on a subset of real data from
Cleargrasp) generalizes well to novel objects. Our network trained only with
synthetic data achieves 88% IOU on the same domain, and it is robust to
challenging image features like sharp caustics and segmenting multiple
overlapping objects. Moreover, our SuperCaustics-R model achieves a
performance of 53% IOU (comparable to Cleargrasp’s 58%) on the Novel Objects
test set. Below, we analyze these results in more detail.

Figure 5: Image (left), Prediction (middle), Ground-truth (right) First two
rows show SuperCaustics Model on SuperCaustics data. Second two rows show
SuperCaustics-R on Cleargrasp-Novel.

Performance, Efficiency and Scalability: Compared to [2], our training
procedure is much more efficient. By producing higher quality, wider-domain
data, we can achieve roughly the same performance on the Novel Objects test
set using only 10% of the amount of data and in a fraction of the training
time used by Cleargrasp (34 epochs on about 8000 images, compared to 100
epochs on 80000 images). Moreover, our training setup did not require 8 Tesla
v100 cards as used by Cleargrasp. Achieving comparable results more
efficiently supports our assertion that for difficult, complex tasks like
labeling shiny transparent objects with caustics, using a rigid, massive
dataset is not ideal, and one can always benefit from the availability of a
better curated, higher quality dataset. That is why open-source, user-friendly
simulations like SuperCaustics are important, because for every use-case
scenario, a well-made custom dataset is going to achieve comparable results
quicker and cheaper than traditional, massive datasets.

Robustness to Caustics:  Fig. 5 shows examples of images with sharp caustics
that were produced by our SuperCaustics system. Our trained model does not
confuse caustics with separate transparent objects and is quite robust to
caustics from various angles and light sources. This is due to the presence of
natural caustics in the dataset.

Performance on the Real Domain:  As Fig. 5 shows, our trained model seems to
have learned the general characteristics of transparent objects from the
SuperCaustics synthetic data. However, as opposed to Cleagrasp where known
objects were present in the synthetic data, our model was not exposed to the
shape of these novel objects. Performance in the real domain was also affected
by the grainy, noisy nature of the real-camera. Prior work has shown that
segmentation models are very sensitive to fidelity [5], thus, a change in the
general fidelity of the domain had a significant, but expected impact on the
performance of the model.

A Solution to Cleargrasp’s main problem — Segmenting Non-local Caustics:  In
Cleargrasp’s original experiments, their trained models falsely classify
caustics as separate transparent objects [2]. One possible solution for this
problem is having a separate system that only classifies caustics, and because
predictions in these regions tend to be invalid, masking the caustics out of
the predictions. This setup is straightforward with SuperCaustics. To test
this, we generated 2300 1280$\times$720 images using the caustics segmentation
mask ground-truth. We then trained a segmentation model to identify sharp,
refractive, and reflective non-local caustics for 50 epochs. We then generated
60 images using novel objects (not seen before by the model) to test the
generalization of the model. As Table II shows, our model was able to segment
caustics cast by novel objects at 68% IOU rate. As such, we believe that
training a model with our SuperCaustics system will prove a feasible solution
for masking out the erroneous pixels in Cleargrasp.

Figure 6: Image (left), Ground-truth (middle), Prediction (right),
SuperCaustics Model on Segmenting Caustics cast by Novel transparent objects.

Sensitivity to Brightness: Our trained model seems to be very sensitive to
brightness in the features of the objects. In our tests with synthetic data,
where the camera’s exposure underexposes shadows, our model seems to make
mistakes in segmenting transparent objects.

Side-Effect of the Patch Based Approach: In some of the images, there is some
noise visible in the labeled pixels. We believe this is due to the patch-based
approach used by the U-Net neural network; specifically, the noise is a result
of the artifacts from stitching these patches together in the final image. In
most cases, patches are surrounded by uniform neighboring pixels, which could
be used as information in a context-based post processing algorithm to get rid
of the noise.

##  V Conclusion & Future Work

Conclusion: In this paper, we present SuperCaustics, an extendable, open-
source, and user-friendly simulation tool for transparent objects that allows
for fine controls over every variable in the simulation. Our system allows
users to create customized datasets using its flexible simulation tools. With
SuperCaustics, we hope to address the dearth of high-quality data for
transparent objects, particularly data that can be modified post-acquisition.
We demonstrated that by using a better curated data, a network can achieve
performance comparable to when training with much larger rigid datasets. In
particular, we were able to achieve comparable performance to the state-of-
the-art networks trained on Cleargrasp’s data using a dataset of only 9000 RGB
images. Specifically, we found that high-quality data is linked to better
performance in segmentation, as our SuperCaustics-R model was able to approach
the performance of state-of-the-art (cleargrasp) with only 34 epochs of
training and 10% of the data, on far more reasonable hardware.

We believe that SuperCaustics will open exciting avenues for using synthetic
data in machine learning, particularly for vision research related to
transparent objects.

Future Work:  Potential research possibilities enabled by the availability of
readily-expandable, open-source simulations are endless. However, we are
working on improving and developing our own deep learning pipeline for other
tasks such as surface normal estimation and occlusion-boundary detection. As
we have shown, an exciting avenue that can make full use of SuperCaustics is
classifying projected local and non-local caustics in images of transparent
objects. This can be useful in complex pipelines that are sensitive to
caustics and need a masking pass to make sure only transparent objects are
processed when estimating surface normals, depth, or segmentation boundaries.
Other avenues for future research are performing data ablation studies, where
isolated features in the data are changed to determine sensitive points in the
performance of a specific neural network, and using the class based features
developed in SuperCaustics to help solve problems such as active learning,
lifelong learning and class incremental learning [20, 21, 22, 23]. Finally, we
will continue developing tools to expand the capabilities of our SuperCaustics
system, as well as improving its stability and compatibility with different
target domains.

## References

  * [1] B. O. Community, _Blender - a 3D modelling and rendering package_ , Blender Foundation, Stichting Blender Foundation, Amsterdam, 2018. [Online]. Available: http://www.blender.org
  * [2] S. S. Sajjan, M. Moore, M. Pan, G. Nagaraja, J. Lee, A. Zeng, and S. Song, “Cleargrasp: 3d shape estimation of transparent objects for manipulation,” 2019\.
  * [3] Epic Games, “Unreal engine.” [Online]. Available: https://www.unrealengine.com
  * [4] NVIDIA, “NvRTX branch of Unreal Engine,” 2021. [Online]. Available: https://github.com/NvRTX/UnrealEngine/tree/NvRTX\\_Caustics-4.26
  * [5] M. Mousavi, A. Khanal, and R. Estrada, “Ai playground: Unreal engine-based data ablation tool for deep learning,” in _International Symposium on Visual Computing_. Springer, 2020, pp. 518–532.
  * [6] Osadchy, Jacobs, and Ramamoorthi, “Using specularities for recognition,” in _Proceedings Ninth IEEE International Conference on Computer Vision_ , 2003, pp. 1512–1519 vol.2.
  * [7] M. Fritz, G. Bradski, S. Karayev, T. Darrell, and M. Black, “An additive latent feature model for transparent object recognition,” _Advances in Neural Information Processing Systems_ , vol. 22, pp. 558–566, 2009.
  * [8] L. Zhu, A. Mousavian, Y. Xiang, H. Mazhar, J. van Eenbergen, S. Debnath, and D. Fox, “RGB-D local implicit function for depth completion of transparent objects,” _CoRR_ , vol. abs/2104.00622, 2021. [Online]. Available: https://arxiv.org/abs/2104.00622
  * [9] V. Seib, A. Barthen, P. Marohn, and D. Paulus, “Friend or foe: exploiting sensor failures for transparent object localization and classification,” in _2016 International Conference on Robotics and Machine Vision_ , A. V. Bernstein, A. Olaru, and J. Zhou, Eds., vol. 10253, International Society for Optics and Photonics. SPIE, 2017, pp. 94 – 98. [Online]. Available: https://doi.org/10.1117/12.2266255
  * [10] K. Rematas, I. Kemelmacher-Shlizerman, B. Curless, and S. M. Seitz, “Soccer on your tabletop,” _CoRR_ , vol. abs/1806.00890, 2018. [Online]. Available: http://arxiv.org/abs/1806.00890
  * [11] J. Mahler, J. Liang, S. Niyaz, M. Laskey, R. Doan, X. Liu, J. A. Ojea, and K. Goldberg, “Dex-net 2.0: Deep learning to plan robust grasps with synthetic point clouds and analytic grasp metrics,” 2017.
  * [12] V. Veeravasarapu, C. Rothkopf, and R. Visvanathan, “Model-driven simulations for computer vision,” in _2017 IEEE Winter Conference on Applications of Computer Vision (WACV)_ , 2017, pp. 1063–1071.
  * [13] V. Haltakov, C. Unger, and S. Ilic, “Framework for generation of synthetic ground truth data for driver assistance applications,” in _GCPR_ , 2013.
  * [14] A. Gaidon, Q. Wang, Y. Cabon, and E. Vig, “Virtual worlds as proxy for multi-object tracking analysis,” 2016.
  * [15] NVIDIA. Nvidia omniverse platform. [Online]. Available: https://developer.nvidia.com/nvidia-omniverse-platform
  * [16] O. Ronneberger, P. Fischer, and T. Brox, “U-net: Convolutional networks for biomedical image segmentation,” in _International Conference on Medical image computing and computer-assisted intervention_. Springer, 2015, pp. 234–241.
  * [17] A. Khanal, “Easy torch,” https://github.com/sraashis/easytorch, 2020.
  * [18] A. Khanal and R. Estrada, “Dynamic deep networks for retinal vessel segmentation,” _Frontiers in Computer Science_ , vol. 2, p. 35, 2020. [Online]. Available: https://www.frontiersin.org/article/10.3389/fcomp.2020.00035%7D,DOI={10.3389/fcomp.2020.00035},ISSN={2624-9898}
  * [19] S. Motevali, A. Khanal, and R. Estrada, “Optic disc segmentation using disk-centered patch augmentation,” 2021.
  * [20] J. K. Mandivarapu, B. Camp, and R. Estrada, “Self-net: Lifelong learning via continual self-modeling,” _Frontiers in Artificial Intelligence_ , vol. 3, p. 19, 2020. [Online]. Available: https://www.frontiersin.org/article/10.3389/frai.2020.00019
  * [21] J. Mandivarapu, B. Camp, and R. Estrada, “Deep active learning via open set recognition,” _CoRR_ , vol. abs/2007.02196, 2020. [Online]. Available: https://arxiv.org/abs/2007.02196
  * [22] B. Camp, J. K. Mandivarapu, and R. Estrada, “Continual learning with deep artificial neurons,” 2020.
  * [23] J. K. Mandivarapu, E. Bunch, Q. You, and G. Fung, “Efficient document image classification using region-based graph neural network,” 2021.