#  Reconstruction of A Scale–Dependent Gravitational Phase Transition

Marzieh Farhang<EMAIL_ADDRESS>Department of Physics, Shahid Beheshti
University, 1983969411, Tehran, Iran
Nima Khosravi  Department of Physics, Shahid Beheshti University, 1983969411,
Tehran, Iran
<EMAIL_ADDRESS>
###### Abstract

In this work we extend our earlier phenomenological model for a gravitational
phase transition (GPT) (Farhang & Khosravi, 2021) and its generalization to
early times (Khosravi & Farhang, 2022) by letting the modifications in the
linearly–perturbed Einstein equations be scale-dependent. These modifications
are characterized as deviations of the parameters $\mu(z,k)$ and $\gamma(z,k)$
from their values in general relativity (GR). The scale-dependent amplitudes
of modified $\mu(z,k)$ and $\gamma(z,k)$ and the parameters defining the phase
transition, along with the standard cosmological parameters, are measured by
various data combinations. Out of the perturbation parameters, we construct
gravity eigenmodes which represent patterns of perturbations best detectable
by data. We detect no significant deviation from GR in these parameters.
However, the larger parameter space produced due to the new degrees of freedom
allows for the reconciliation of various datasets which are in tension in
$\Lambda$CDM. In particular, we find $H_{0}=71.9\pm 9.2$ from anisotropies of
the Cosmic Microwave Background as measured by Planck (Planck Collaboration et
al., 2020) and various measurements of the Baryonic Acoustic Oscillations, in
agreement with local Hubble measurements (Riess et al., 2022). We also find
that the $\sigma_{8}$ tension between the measurements of Dark Energy Survey
(Kazin et al., 2014) and Planck is reduced to less than $1\sigma$.

##  1 Introduction

The standard model of cosmology, aka the $\Lambda$CDM scenario, in spite of
its incredible achievements in explaining the observed sky, suffers from
serious limitations. Certain fundamental questions are left unanswered in this
scenario such as the nature of dark matter and dark energy and the
cosmological constant problem. Moreover, the model has been challenged by
observations in recent years. Among the most severe is the inconsistency of
the Hubble constant as directly measured by local data (Riess et al., 2022;
Murakami et al., 2023) and $\Lambda$CDM-based inferences from the Cosmic
Microwave Background (CMB) anisotropies (Planck Collaboration et al., 2020).
The Hubble tension is accompanied by milder ones including the $\sigma_{8}$
tension (Abbott et al., 2018; Heymans et al., 2020), the deviation of the
measured CMB lensing amplitude from unity, the parameter inconsistencies in
the CMB low and high $\ell$ multipoles, and the CMB spatial anomalies (Planck
Collaboration et al., 2020). All these anomalies can be due to systematics in
observational and analysis pipeline. In this paper, however, we follow the
viewpoint that the tensions might be due to theoretical insufficiencies. We
mainly focus on $H_{0}$ and $\sigma_{8}$ tensions here.

In the attempts to address the $H_{0}$ tension through modifications to the
$\Lambda$CDM model, a very natural place to look for is dark energy. See,
e.g., Di Valentino et al. (2017); Zhao et al. (2017); Khosravi et al. (2019);
Yang et al. (2018b, a); Banihashemi et al. (2019); Li & Shafieloo (2019);
Keeley et al. (2019); Raveri (2020); Yan et al. (2020); Di Valentino et al.
(2020c); Braglia et al. (2020); Gómez-Valent et al. (2020); Lucca & Hooper
(2020); Di Valentino et al. (2020a, b); Banihashemi et al. (2020, 2021) and
for a comprehensive review see Di Valentino et al. (2021). Modifications in
the early universe are also able to lessen the $H_{0}$ tension by reducing the
sound horizon $r_{s}$ (Poulin et al., 2019; Knox & Millea, 2020; Keeley et
al., 2020; Vagnozzi, 2021). These two quantities are degenerate and an early
phase of dark energy can help relax the $H_{0}$ tension, even in the presence
of the baryonic acoustic oscillations, or BAO data (see, however, the
discussions in, e.g., Clark et al. (2023); Allali et al. (2021) on swelling
the $S_{8}$ tension). The $H_{0}$ tension can also be reduced by going beyond
the slow roll inflation (Keeley et al., 2020). However, these models often
fail to explain the various datasets simultaneously (Jedamzik et al., 2021).
In an attempt for simultaneously addressing various tensions and explaining
the different datasets in one model, the authors studied a phenomenological
model of modified gravity with some physical roots in the physics of critical
phenomena (Banihashemi et al., 2020, 2019, 2021). In these works we modelled a
gravitational phase transition (GPT) at both background and linear
perturbation levels by a smooth step function in redshift (Farhang & Khosravi,
2021; Khosravi & Farhang, 2022, referred to as the GPTs throughout this work).
The degrees of freedom in the GPT scenario could simultaneously explain the
Hubble and $\sigma_{8}$ tension, as well as the internal CMB inconsistencies
of low-high $\ell$ and lensing amplitudes. However, the BAO data were found to
heavily push the parameters toward their $\Lambda$CDM values and thus be in
tension with direct Hubble measurements.

In this paper, we extend the GPT model by letting the transition be
scale–dependent. More specifically, the modifications to the gravity sector in
perturbation equations are allowed to happen with different amplitudes at
different scales. This scale-dependent gravity can be motivated by, e.g.,
physics of critical phenomena where there is a possibility to have
scale–dependent forces due to nontrivial physics of microstructures. The
scale-dependency is included in the analysis in a data–driven fashion which
allows for the construction of gravity eigenmodes, or GEMs. GEMs describe
patterns of deviations from general relativity, or GR, that data are most
sensitive to, and in general depend on the datasets used. We then investigate
how the newly introduced degrees of freedom impact the tension and whether
they point to deviations from GR through reconstructing the scale-dependent
perturbation parameters.

The paper is organized as follows. We first introduce the model for gravity
eigenmodes in Section 2 and the data used to constrain it in Section 3. The
results of the confrontation of the model with data is presented in Section 4.
We conclude the paper with discussions of the results and remarks on future
work in Section 5.

##  2 Gravity Eigenmodes

In the gravitational eignemode scenario, we assume that gravity undergoes a
phase transition, modelled phenomenologically by a "tanh", as introduced and
expanded in in the GPTs. The transition affects both background expansion and
perturbation equations, and being only a function of redshift, it is assumed
to be scale–independent in the GPTs. In this paper we extend the model further
to allow for scale-dependent modifications to perturbation equations.

We assume a flat, homogenous and isotropic universe whose background evolution
is described by

$\displaystyle
H^{2}(z)=H_{0}^{2}\bigg{[}\Omega_{r}(1+z)^{4}+\Omega_{m}(1+z)^{3}+\Omega_{\Lambda}(z)\bigg{]}.$

Scalar perturbations in the linear regime are introduced to the background
metric through the functions $\Psi$ and $\Phi$ in the Newtonian gauge,

$\displaystyle ds^{2}=a^{2}(\tau)[$ $\displaystyle-$
$\displaystyle(1+2\Psi(\tau,\vec{x}))d\tau^{2}+(1-2\Phi(\tau,\vec{x}))d\vec{x}^{2}].$

The perturbations $\Psi$ and $\Phi$ must satisfy the linearly–perturbed
Einstein equations. These equations in Fourier space can be modified
phenomenologically as

$\displaystyle k^{2}\Psi=-\mu(z,k)\,4\pi
Ga^{2}\,\left[\rho\Delta+3(\rho+P)\sigma\right],$ $\displaystyle
k^{2}\left[\Phi-\gamma(z,k)\Psi\right]=\mu(z,k)\,12\pi
Ga^{2}\,(\rho+P)\sigma.$

Deviations in $\mu(z,k)$ and $\gamma(z,k)$ from unity would then embrace any
general redshift– and scale–dependent modifications to GR.

Following the GPTs, we model the gravitational transition in the background by

$\displaystyle\Omega_{\Lambda}(z)=\Omega_{\Lambda\rm{CDM}}+\Delta_{\Lambda}\,\frac{1+\tanh\big{[}\alpha\,(z_{\rm
t}-z)\big{]}}{2},$ (1)

where $z_{\rm t}$ and $\alpha$ represent the redshift and inverse of the width
of the transition. $\Omega_{\Lambda\rm{CDM}}$ corresponds to the cosmological
constant in $\Lambda$CDM, so $\Delta_{\Lambda}=0$ would be equivalent to
$\Lambda$CDM. For the flat universe assumption we impose $\Omega_{\rm
r}+\Omega_{\rm m}+\Omega_{\Lambda}(z=0)=1$, where
$\Omega_{\Lambda}(z=0)=\Omega_{\Lambda\rm{CDM}}+\Delta_{\Lambda}[1+\tanh(\alpha
z_{\rm t})]/2$. Similarly, for the modifications to perturbation equations, we
assume a transition in $\mu$ and $\gamma$,

$\displaystyle\mu(z,k)$ $\displaystyle=$
$\displaystyle\mu_{\rm{GR}}+\Delta_{\mu}(k)\frac{1+\tanh\big{[}\alpha\,(z_{\rm
t}-z)\big{]}}{2},$ $\displaystyle\gamma(z,k)$ $\displaystyle=$
$\displaystyle\gamma_{\rm{GR}}+\Delta_{\gamma}(k)\frac{1+\tanh\big{[}\alpha\,(z_{\rm
t}-z)\big{]}}{2},$ (2)

where $\mu_{\rm{GR}}=\gamma_{\rm{GR}}=1$ and
$\Delta_{\mu}(k)=\Delta_{\gamma}(k)=0$ correspond to GR. It should be noted
that in this model the redshift dependences are assumed to be the same for
modifications in the background (Eq. 1), and perturbation equations (Eq. 2).
Scale-dependence of $\mu(z,k)$ and $\gamma(z,k)$ happens through
$\Delta_{\mu}$ and $\Delta_{\gamma}$ in a model-independent way. We divide the
$k$-range of interest $[k_{\rm{min}},k_{\rm{max}}]$ into $n_{k}$ bins.
Deviations from GR are parameterized as amplitudes of top–hats in the $k$-bins
of $\mu$ and $\gamma$, i.e., through non-zero $\Delta_{\mu_{i}}$ and
$\Delta_{\gamma_{i}}$, $i=1,..,n_{k}$. Therefore, in addition to the standard
$\Lambda$CDM parameters, the parameter space of this phenomenological model of
cosmology would consist of $\\{z_{\rm
t},\alpha,\Delta_{\Lambda},\Delta_{\mu_{1}},...,\Delta_{\mu_{n_{k}}},\Delta_{\gamma_{1}},...,\Delta_{\gamma_{n_{k}}}\\}$
(see Table 1). We modify the code MGCAMB111https://github.com/sfu-
cosmo/MGCAMB.
MGCAMB is an implementation of CAMB, https://camb.info/, for modified gravity
models.(Hojjati et al., 2011), to take into account our desired changes. The
new parameters are then fed into the Cosmological MonteCarlo code,
CosmoMC222https://cosmologist.info/cosmomc/, to be measured with various
dataset combinations.

In general, the correlation between these perturbation parameters could be
large, in particular for adjacent bins, due to similar traces they leave on
data. On the other hand, adding $\sim 2n_{k}$ new parameters is a quite
expensive approach to describe possible deviations around the astonishingly
simple $\Lambda$CDM. One could therefore consider constructing $2n_{k}$
uncorrelated linear combinations of $\Delta_{\mu_{i}}$’s and
$\Delta_{\gamma_{i}}$’s and keeping only the most tightly constrained ones.
The new parameters are generated through the eigen-decomposition of the
perturbation blocks of the covariance matrix, $\bf{C}_{\mu}$ and ${\bf
C}_{\gamma}$. They are often referred to as the principal components or
eigenmodes of the covariance matrix, hence the name gravity eigenmodes.
Estimates of their errors are given by the square roots of the eigenvalues of
$\bf{C}_{\mu}$ and ${\bf C}_{\gamma}$. The eigenmodes and their measured means
and estimated uncertainties can then be used to reconstruct the overall
patterns of observed deviation from unity for $\mu$ and $\gamma$, at scales of
interest.

parameter | description
---|---
$\Delta_{\Lambda}$ | transition amplitude for $\Omega_{\Lambda}$
$z_{\rm t}$ | transition redshift
$\alpha$ | inverse of the width of transition
$\Delta_{\mu_{i}}$ | transition amplitude for $\mu$, $i\in\\{1,...,n_{k}\\}$
$\Delta_{\gamma_{i}}$ | transition amplitude for $\gamma$, $i\in\\{1,...,n_{k}\\}$
Table 1: Model parameters.

##  3 datasets

In this work we use Planck CMB measurements of temperature and polarization
anisotropies and its lensing, referred to as P18 (Planck Collaboration et al.,
2020), the local measurements of the Hubble constant, labeled as R21 (Riess et
al., 2022), BAO measurements consisting of 6DFGS (Beutler et al., 2011), SDSS
MGS (Ross et al., 2015) and BOSS DR12 (Alam et al., 2017), and the
measurements of cosmic shear and galaxy clustering by Dark Energy Survey,
referred to as DES (Kazin et al., 2014). We also consider a mock futuristic
and extremely high–precision measurement of the Hubble constant with
$H_{0}=73.30\pm 0.01$, to investigate the implications of such a measurement
for the construction of $\Delta_{\mu}(k)$ and $\Delta_{\gamma}(k)$. We refer
to this Hubble measurement as $H_{0}73$ for short.

##  4 results

| $\Omega_{\rm b}h^{2}$ | $\Omega_{\rm c}h^{2}$ | $H_{0}$ | $\tau$ | $\ln(10^{10}A_{\rm s})$ | $n_{\rm s}$
---|---|---|---|---|---|---
P18+ |  | GEM |  |  |  |
BAO+R21 | $0.0224\pm 0.0003$ | $0.119\pm 0.002$ | $70.3\pm 6.8$ | $0.053\pm 0.007$ | $3.040\pm 0.016$ | $0.967\pm 0.007$
BAO | $0.0225\pm 0.0002$ | $0.119\pm 0.002$ | $71.9\pm 9.2$ | $0.053\pm 0.008$ | $3.037\pm 0.017$ | $0.969\pm 0.006$
R21 | $0.0225\pm 0.0002$ | $0.119\pm 0.002$ | $72.3\pm 2.6$ | $0.054\pm 0.009$ | $3.039\pm 0.017$ | $0.967\pm 0.005$
|  |  $\Lambda$CDM |  |  |  |
BAO | $0.0224\pm 0.0001$ | $0.119\pm 0.001$ | $67.7\pm 0.4$ | $0.0561\pm 0.0071$ | $3.047\pm 0.014$ | $0.967\pm 0.004$
Table 2: Standard cosmological parameters in GEM and $\Lambda$CDM measured by
various datasets.

The goal is to find how small deviations from GR, parametrized in a semi-blind
way, would alleviate some of the cosmological tensions. We take $n_{k}=20$
logarithmically–spaced bins in the range
$[k_{\rm{min}},k_{\rm{max}}]=[10^{-5},1]$. For the the redshift of transition
$z_{\rm t}$ we first search in $[10^{-4},10^{4}]$ to cover a large redshift
range of cosmological interest. We explore this range logarithmically for a
relatively fair coverage of different redshift scales. The transition is also
allowed to occur both ways, from $\Lambda$CDM at late times to the perturbed
model at early times and vice versa by allowing $\alpha$ to be positive as
well as negative. However, we find no support from data for such a transition
at high $z$. We therefore limit the search to a linear search in the range
$z\in[10^{-4},20]$ for higher resolution. We are interested to see whether the
Hubble tension between R21 and P18+BAO datasets and the $\sigma_{8}$ tension
between P18 and DES are reduced in this scenario.

Figure 1: The posterior probability distribution of $H_{0}$ for three data
combinations. The shaded green regions illustrate the R21 measurement of
$H_{0}$. Top: The background evolution is assumed to be the same as in
$\Lambda$CDM. Bottom: The background evolves according to Eq. 1.

###  4.1 Hubble tension

We first investigate how various combinations of P18, BAO and R21 constrain
standard and GEM parameters. We compare the results for two cases, with and
without modifications in the background evolution. This would help better
distinguish the impact of $\mu$ and $\gamma$ perturbations from background
modifications on observables.

Figure 2: The first three GEMs, or $k$-dependent gravity eigenmodes, for $\mu$
and $\gamma$, with various datasets. The background is allowed to experience a
phase transition as in Eq. 1.

Figure 1 illustrates the posterior probability of $H_{0}$ for these two cases.
In both plots the likelihoods for P18 and P18+BAO agree, confirming the
consistency of the two datasets as far as $H_{0}$ is concerned. It is,
however, evident from the top panel that modifications in $\mu$ and $\gamma$
alone, i.e., with the background assumed to be the same as $\Lambda$CDM,
cannot resolve the Hubble tension between P18 (or P18+BAO) and R21, and there
still exists a $\sim 2-3\sigma$ tension between the two measurements. The
joint P18+R21 likelihood, despite the inconsistency of P18 and R21 in this
scenario, is plotted to see how far P18 would push the bestfit $H_{0}$ of R21
toward its preferred lower values. Allowing the background to also deviate
from $\Lambda$CDM would largely increase the uncertainty in the $H_{0}$
measurement, and let $H_{0}$ adopt higher values even with P18, as shown in
the bottom panel. With this enhanced uncertainty of $H_{0}$ measurement, the
Hubble tension is fully removed. It is interesting to note that the joint
P18+R21 measurement of $H_{0}$, which are now justifiably combined, yields a
quite large bestfit for $H_{0}$, unlike the P18+R21 measurement in the fixed-
background scenario. Also note that background–only modifications, and with no
new degrees of freedom in perturbation equations, would not relax the Hubble
tension, as was explored and fully discussed in Khosravi & Farhang (2022).

Table 2 compares the mean values and errors for the standard cosmological
parameters in GEM and $\Lambda$CDM scenarios, for various dataset
combinations. Note that the datasets R21 and P18+BAO are incompatible in
$\Lambda$CDM and cannot be combined. Moreover, we find that $\Delta_{\Lambda}$
is basically unconstrained and hits the prior boundaries $[-1,1]$ for all data
combinations used here. The transition is preferred to occur at $z_{t}<1$ and
in both directions (determined by the sign of $\alpha$). However, a transition
from $\Lambda$CDM at high redshifts to GEM at lower redshifts is slightly
favoured by data. We note that although the degrees of freedom in $\mu$ and
$\gamma$ are important to relax the tension, these parameters are found to be
consistent with GR.

The obscuration of possible deviation in the new parameters from their GR and
$\Lambda$CDM values could be explained by the large number of perturbation
parameters, as discussed in Section 2. We therefore use the $n_{k}\times
n_{k}$ correlation matrices of perturbation amplitudes $\bf{C}_{\mu}$ and
$\bf{C}_{\gamma}$ to generate the gravity eigenmodes. The amplitudes of the
modes and their errors can be calculated respectively from the measurements of
the top–hat amplitudes and their covariance matrices. This method of eigenmode
construction is data–driven in the sense that the modes can be rank–ordered so
that the first few ones represent perturbation patterns with most detectable
imprints on data. Therefore one expects smaller errors on the amplitudes of
these modes, compared to the top–hats we started with, and a higher chance of
detecting semi–blind, scale–dependent deviations from GR.

The first three perturbation eigenmodes of $\Delta_{\mu}(k)$ and
$\Delta_{\gamma}(k)$ are presented in Figure 2. The $\Delta_{\mu}$ eigenmodes
have their main features in $k\in\\{0.001,0.1\\}$, extending to higher $k$’s
as the mode number increases. The modes for the full data combination
P18+BAO+R21 are the smoothest with fewer wiggles, and constructed with fewer,
more dominant bumps. For $\gamma$ the case is different and the modes seem to
be extended in a large $k$-range. It should be stated that for some
$\Delta_{\gamma_{i}}$’s, unlike $\Delta_{\mu_{i}}$’s, the errors are dominated
by the prior ranges. On the other hand, the priors cannot be increased farther
due to numerical limitations. This would make the patterns of $\gamma$
eigenmodes affected by noise and numerical insufficiencies.

The estimated errors for these modes are listed in Table 3. All the mode
amplitudes are found to be consistent with zero, and therefore no hint $\geq
1\sigma$ is found for deviation from GR. However, it is important to stress
that the Hubble tension is fully resolved in our phenomenological extension of
GR. It is interesting to see the $\Delta_{\mu}(k)$ and $\Delta_{\gamma}(k)$
patterns leading to this resolution. We reconstruct these perturbation
patterns using the first five eigenmodes. We use the linear dependance of the
modes on the $\mu$ and $\gamma$ top-hats, and their measured amplitudes to
infer the means and errors of the GEMs. The resulting patterns of deviations
are shown in Figure 3 as solid blue curves. The shaded blue areas around the
mean patterns illustrate the region of $1\sigma$ uncertainty and are found
based on the estimated errors of the mode amplitudes and their pairwise
independence. Therefore the price to resolve the Hubble tension in this model
has been the (insignificant) deviations from GR in $\mu$ and $\gamma$, mainly
for $k\in\\{0.001,0.1\\}$. The fact that GEM can resolve the Hubble tension,
yet no deviations from $\Lambda$CDM are detected, persuades one to explore how
a very tight $H_{0}$ measurement, in substantial tension with early
$\Lambda$CDM-inferred measurements of the Hubble constant, would derive GEM
parameters away from their GR values. The green solid lines and shaded areas
in Figure 3 illustrate the means and $1\sigma$ regions of uncertainty for the
reconstructed $\Delta_{\mu}(k)$ and $\Delta_{\gamma}(k)$ patterns with this
$H_{0}73$ mock dataset (see Section 3). To compensate for the tightly measured
high value of $H_{0}$, $\mu$ and $\gamma$ need to differ from $1$ at more than
$1\sigma$ at various scales. With this mock dataset we also find a higher than
$5\sigma$ detection of $\Delta_{\Lambda}$ with $\Delta_{\Lambda}=0.44\pm
0.08$, and with a transition at $z_{\rm t}=0.14\pm 0.11$ with $\alpha=2.2\pm
0.8$.

###  4.2  Tension in $\sigma_{8}$

We use the data from DES to constrain parameters in the GEM scenario and
investigate how they compare to P18 predictions. The measurements of
$\sigma_{8}$ are of prime interest in this case as they happen to show a mild,
$\sim 2\sigma$ disagreement in the $\Lambda$CDM framework for P18 and DES. For
the DES case we have imposed a relatively large prior on $H_{0}$, i.e.,
$H_{0}\in[55,95]$. In the GEM scenario we find

$\displaystyle\sigma_{8}=0.78\pm 0.11~{}~{}~{}{\rm DES}$ (3)
$\displaystyle\sigma_{8}=0.85\pm 0.07~{}~{}~{}~{}{\rm P18}$ (4)

showing that the $\sigma_{8}$ tension is reduced to less than one sigma.

| $\sigma_{\mu_{1}}$ | $\sigma_{\mu_{2}}$ | $\sigma_{\mu_{3}}$ | $\sigma_{\gamma_{1}}$ | $\sigma_{\gamma_{2}}$ | $\sigma_{\gamma_{3}}$
---|---|---|---|---|---|---
P18+BAO | $0.38$ | $0.40$ | $0.43$ | $0.45$ | $0.47$ | $0.49$
P18+R21 | $0.35$ | $0.37$ | $0.38$ | $0.46$ | $0.46$ | $0.47$
P18+R21+BAO | $0.34$ | $0.35$ | $0.37$ | $0.44$ | $0.45$ | $0.46$
P18+$H_{0}$73 | 0.04 | 0.04 | 0.07 | 0.04 | 0.06 | 0.07
Table 3: The estimated uncertainty of the eigenmodes plotted in Figure 2.

Figure 3: The reconstructed perturbations in $\mu(k)$ and $\gamma(k)$ in the
GEM scenario, with the first five eigenmodes, for the two data combinations
P18+BAO+R21 and P18+$H_{0}73$.

##  5 Discussion

In this work we explored a phenomenological extension to general relativity,
with perturbation equations allowed to deviate from GR in a scale–dependent
way. Based on the two earlier GPT works, this deviation was assumed to
accompany a phase transition in redshift, both at the background and
perturbation levels. The scale-dependency was investigated in a semi–blind way
by generating gravity eigen-modes, as patterns of modifications to GR with the
highest impact on data. Our results showed that the enhanced parameter space
could consistently accommodate CMB, BAO, local Hubble measurements and data
from large scale structure. However, it should be noted that no significant
departure from $\Lambda$CDM and GR was observed.

To test the model against a more restrictive scenario, we assumed a mock
measurement of the Hubble constant $H_{0}=73.30\pm 0.01$. This led to tighter
constraints on our free parameters. In particular, our forecasts showed the
inconsistency of $\Delta_{\mu}(k)$, which represents modification to Poisson
equation, from GR predictions if the Hubble constant insisted on the large
value $H_{0}=73.30$ (see Figure 3). On the other hand, $\Delta_{\gamma}(k)$
which modifies the lensing potential, did not show any significant departure
from GR. We therefore expect future data with their higher resolution to be
able to distinguish the GEM model from GR, with tighter constraints on GEM
parameters. The $k$–dependent freedom in $\mu$ and $\gamma$ amplitudes may
also have the potential to relax the low-high $\ell$ and lensing tensions in
CMB. We leave this to future work.

The $k$-reconstruction is a plausible approach to look for modified gravity
signatures. However, introducing many degrees of freedom to the extremely
simple and neat $\Lambda$CDM is a quite high price to relax the tensions. Yet,
it should be noted that this could be an issue for a theoretical model, while
the GEM scenario is phenomenological in redshift space and data–driven in
$k$-space. For such a theoretically–blind approach, the ultimate goal could be
adding as many degrees of freedom as numerically possible to allow the data to
freely choose the ones best describing it. This is simply not feasible due to
limited available information in the data and computational power. Hence one
can cut the mode–hierarchy by keeping the few modes with the tightest
constraints and leave the rest out of the analysis.

It should also be noted that this search for non-trivial $k$-dependency could
be motivated by the physics of critical phenomena where it is possible for the
system to behave differently in different scales, associated with a transition
in time. That is why we restricted the gravitational modification in this work
to be accompanied with a phase transition. This assumption can however be
relaxed to explore a more general reconstruction scenario of gravity in
redshift and scale, without assuming a restrictive theoretical prior (see
Raveri, 2020, for reconstructed modifications of gravity in time). This
approach would suffer from the many degrees of freedom introduced to cover
possible modification in large spans of redshifts and scales. It is therefore
limited by numerical insufficiencies, and calls for a clever way to look at.
We leave this to future work.

##  6 Acknowledgement

Part of the numerical computation of this work was carried out on the
computing cluster of the Canadian Institute for Theoretical Astrophysics
(CITA), University of Toronto.

## References

  * Abbott et al. (2018) Abbott, T., et al. 2018, Phys. Rev. D, 98, 043526, doi: 10.1103/PhysRevD.98.043526
  * Alam et al. (2017) Alam, S., et al. 2017, Mon. Not. Roy. Astron. Soc., 470, 2617, doi: 10.1093/mnras/stx721
  * Allali et al. (2021) Allali, I. J., Hertzberg, M. P., & Rompineve, F. 2021, Phys. Rev. D, 104, L081303, doi: 10.1103/PhysRevD.104.L081303
  * Banihashemi et al. (2021) Banihashemi, A., Khosravi, N., & Shafieloo, A. 2021, JCAP, 06, 003, doi: 10.1088/1475-7516/2021/06/003
  * Banihashemi et al. (2019) Banihashemi, A., Khosravi, N., & Shirazi, A. H. 2019, Phys. Rev. D, 99, 083509, doi: 10.1103/PhysRevD.99.083509
  * Banihashemi et al. (2020) —. 2020, Phys. Rev. D, 101, 123521, doi: 10.1103/PhysRevD.101.123521
  * Beutler et al. (2011) Beutler, F., Blake, C., Colless, M., et al. 2011, Mon. Not. Roy. Astron. Soc., 416, 3017, doi: 10.1111/j.1365-2966.2011.19250.x
  * Braglia et al. (2020) Braglia, M., Ballardini, M., Emond, W. T., et al. 2020, Phys. Rev. D, 102, 023529, doi: 10.1103/PhysRevD.102.023529
  * Clark et al. (2023) Clark, S. J., Vattis, K., Fan, J., & Koushiappas, S. M. 2023, Phys. Rev. D, 107, 083527, doi: 10.1103/PhysRevD.107.083527
  * Di Valentino et al. (2017) Di Valentino, E., Melchiorri, A., & Mena, O. 2017, Phys. Rev. D, 96, 043503, doi: 10.1103/PhysRevD.96.043503
  * Di Valentino et al. (2020a) Di Valentino, E., Melchiorri, A., Mena, O., & Vagnozzi, S. 2020a, Phys. Dark Univ., 30, 100666, doi: 10.1016/j.dark.2020.100666
  * Di Valentino et al. (2020b) —. 2020b, Phys. Rev. D, 101, 063502, doi: 10.1103/PhysRevD.101.063502
  * Di Valentino et al. (2020c) Di Valentino, E., Mukherjee, A., & Sen, A. A. 2020c.  https://arxiv.org/abs/2005.12587
  * Di Valentino et al. (2021) Di Valentino, E., Mena, O., Pan, S., et al. 2021, Class. Quant. Grav., 38, 153001, doi: 10.1088/1361-6382/ac086d
  * Efstathiou (2021) Efstathiou, G. 2021, Mon. Not. Roy. Astron. Soc., 505, 3866, doi: 10.1093/mnras/stab1588
  * Farhang & Khosravi (2021) Farhang, M., & Khosravi, N. 2021, Phys. Rev. D, 103, 083523, doi: 10.1103/PhysRevD.103.083523
  * Gómez-Valent et al. (2020) Gómez-Valent, A., Pettorino, V., & Amendola, L. 2020, Phys. Rev. D, 101, 123513, doi: 10.1103/PhysRevD.101.123513
  * Heymans et al. (2020) Heymans, C., Tröster, T., Asgari, M., et al. 2020, arXiv e-prints, arXiv:2007.15632.  https://arxiv.org/abs/2007.15632
  * Hojjati et al. (2011) Hojjati, A., Pogosian, L., & Zhao, G.-B. 2011, JCAP, 08, 005, doi: 10.1088/1475-7516/2011/08/005
  * Jedamzik et al. (2021) Jedamzik, K., Pogosian, L., & Zhao, G.-B. 2021, Commun. in Phys., 4, 123, doi: 10.1038/s42005-021-00628-x
  * Kazin et al. (2014) Kazin, E. A., et al. 2014, Mon. Not. Roy. Astron. Soc., 441, 3524, doi: 10.1093/mnras/stu778
  * Keeley et al. (2019) Keeley, R. E., Joudaki, S., Kaplinghat, M., & Kirkby, D. 2019, JCAP, 12, 035, doi: 10.1088/1475-7516/2019/12/035
  * Keeley et al. (2020) Keeley, R. E., Shafieloo, A., Hazra, D. K., & Souradeep, T. 2020, JCAP, 09, 055, doi: 10.1088/1475-7516/2020/09/055
  * Khosravi et al. (2019) Khosravi, N., Baghram, S., Afshordi, N., & Altamirano, N. 2019, Phys. Rev. D, 99, 103526, doi: 10.1103/PhysRevD.99.103526
  * Khosravi & Farhang (2022) Khosravi, N., & Farhang, M. 2022, Phys. Rev. D, 105, 063505, doi: 10.1103/PhysRevD.105.063505
  * Knox & Millea (2020) Knox, L., & Millea, M. 2020, Phys. Rev. D, 101, 043533, doi: 10.1103/PhysRevD.101.043533
  * Li & Shafieloo (2019) Li, X., & Shafieloo, A. 2019, Astrophys. J. Lett., 883, L3, doi: 10.3847/2041-8213/ab3e09
  * Lucca & Hooper (2020) Lucca, M., & Hooper, D. C. 2020.  https://arxiv.org/abs/2002.06127
  * Murakami et al. (2023) Murakami, Y. S., Riess, A. G., Stahl, B. E., et al. 2023.  https://arxiv.org/abs/2306.00070
  * Planck Collaboration et al. (2020) Planck Collaboration, Aghanim, N., Akrami, Y., et al. 2020, A&A, 641, A6, doi: 10.1051/0004-6361/201833910
  * Poulin et al. (2019) Poulin, V., Smith, T. L., Karwal, T., & Kamionkowski, M. 2019, Phys. Rev. Lett., 122, 221301, doi: 10.1103/PhysRevLett.122.221301
  * Raveri (2020) Raveri, M. 2020, Phys. Rev. D, 101, 083524, doi: 10.1103/PhysRevD.101.083524
  * Riess et al. (2022) Riess, A. G., et al. 2022, Astrophys. J. Lett., 934, L7, doi: 10.3847/2041-8213/ac5c5b
  * Ross et al. (2015) Ross, A. J., Samushia, L., Howlett, C., et al. 2015, Mon. Not. Roy. Astron. Soc., 449, 835, doi: 10.1093/mnras/stv154
  * Vagnozzi (2021) Vagnozzi, S. 2021.  https://arxiv.org/abs/2105.10425
  * Yan et al. (2020) Yan, S.-F., Zhang, P., Chen, J.-W., et al. 2020, Phys. Rev. D, 101, 121301, doi: 10.1103/PhysRevD.101.121301
  * Yang et al. (2018a) Yang, W., Mukherjee, A., Di Valentino, E., & Pan, S. 2018a, Phys. Rev. D, 98, 123527, doi: 10.1103/PhysRevD.98.123527
  * Yang et al. (2018b) Yang, W., Pan, S., Di Valentino, E., et al. 2018b, JCAP, 09, 019, doi: 10.1088/1475-7516/2018/09/019
  * Zhao et al. (2017) Zhao, G.-B., et al. 2017, Nature Astron., 1, 627, doi: 10.1038/s41550-017-0216-z


supp.pdf

††thanks<EMAIL_ADDRESS>
# Universal Aspect of Relaxation Dynamics in Random Spin Models

Tian-Gang Zhou  Institute for Advanced Study, Tsinghua University, Beijing,
100084, China  Wei Zheng  Hefei National Research Center for Physical Sciences
at the Microscale and School of Physical Sciences, University of Science and
Technology of China, Hefei 230026, China  CAS Center for Excellence in Quantum
Information and Quantum Physics, University of Science and Technology of
China, Hefei 230026, China  Hefei National Laboratory, University of Science
and Technology of China, Hefei 230088, China  Pengfei Zhang  Department of
Physics, Fudan University, Shanghai, 200438, China

###### Abstract

The concept of universality is a powerful tool in modern physics, allowing us
to capture the essential features of a system’s behavior using a small set of
parameters. In this letter, we unveil universal spin relaxation dynamics in
anisotropic random Heisenberg models with infinite-range interactions at high
temperatures. Starting from a polarized state, the total magnetization can
relax monotonically or decay with long-lived oscillations, determined by the
sign of a universal single function
$A=-\xi_{1}^{2}+\xi_{2}^{2}-4\xi_{2}\xi_{3}+\xi_{3}^{2}$. Here
$(\xi_{1},\xi_{3},\xi_{3})$ characterizes the anisotropy of the Heisenberg
interaction. Furthermore, the oscillation shows up only for $A>0$, with
frequency $\Omega\propto\sqrt{A}$. To validate our theory, we compare it to
numerical simulations by solving the Kadanoff-Baym (KB) equation with a melon
diagram approximation and the exact diagonalization (ED). The results show our
theoretical prediction works in both cases, regardless of a small system size
$N=8$ in ED simulations. Our study sheds light on the universal aspect of
quantum many-body dynamics.

_Introduction.–_ A complete description of realistic many-body systems always
contains a large number of parameters. For example, typical solid-state
material contains complicated interactions between electrons, phonons, nuclei,
and impurities. However, properties that are of physical interest can usually
be captured by simple toy models with few parameters. This is a remarkable
consequence of universality. The universality states that microscopically
different systems can share the same physics at large scales. It usually
emerges in the low-energy limit. For example, the phase transition of many-
body systems can be classified into universality classes determined only by
the symmetry and dimension of systems [1, 2]. Low-energy scattering between
atoms can be well described by a single parameter, the scattering length
$a_{s}$, despite details of underlying microscopic interaction potentials [3].
Aiming at deepening our understanding of realistic systems, discovering new
universalities becomes an important subject in modern quantum many-body
physics.

Recent years have witnessed a great breakthrough in understanding novel
quantum dynamics in many-body systems both theoretically [4, 5, 6, 7, 8, 9,
10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28,
29, 30] and experimentally [31, 32, 33, 34, 35, 36, 37]. In the previous
studies on relaxation, most universal dynamics emerges in the low temperature
or long time regime. That reflects the microscopic details of models are
smoothed out in the low energy scale. However at high temperature and short
time, it is common believed that most microscopic details are involved in the
evolution. Such that the evolution is highly model dependent and hard to
observe a universal dynamics. In this work, we unveil that a universal aspect
of relaxation dynamics which shows up in an anisotropic Heisenberg model with
all-to-all interactions even at high temperatures and short time. The
Hamiltonian reads:

$\hat{H}=\sum\limits_{1\leq i<j\leq
N}J_{ij}(\xi_{1}\hat{S}^{x}_{i}\hat{S}^{x}_{j}+\xi_{2}\hat{S}^{y}_{i}\hat{S}^{y}_{j}+\xi_{3}\hat{S}^{z}_{i}\hat{S}^{z}_{j})-h(t)\sum\limits_{1\leq
i\leq N}\hat{S}^{x}_{i}.$ (1)

This model with different anisotropy parameters $(\xi_{1},\xi_{2},\xi_{3})$
has been realized in cold molecules [38, 39], NV centers [40, 41], trapped
fermions [42], Rydberg atoms [43, 44], high spin atoms [45], and solid-state
NMR systems [46, 47]. A schematic figure is presented in FIG. 1 (a). Because
of random locations or complicated spatial wavefunctions of spin carriers,
$J_{ij}$ is usually modeled as independent random Gaussian variables with
expectation $\overline{J_{ij}}=\bar{J}/N$ and variance $\overline{\delta
J_{ij}^{2}}=4J^{2}/N$.

Figure 1: (a). Schematics of the random spin model with random (anisotropic)
Heisenberg interactions $J_{ij}$ in the magnetic field $h$. (b). Different
dynamical behaviors of the system for different anisotropy parameters
$(1,\xi_{2},\xi_{3})$.The boundary line is determined by
$A=-\xi_{1}^{2}+\xi_{2}^{2}-4\xi_{2}\xi_{3}+\xi_{3}^{2}=0$, which is symmetric
under the reflection along $\xi_{2}=\pm\xi_{3}$. OR and NOR denote the
oscillating regime and non-oscillating relaxation regime respectively,
distinguished by features of the magnetization relaxation process.

We focus on the following protocol: The system is prepared at high
temperatures with a polarization field $h(t<0)=h$, which induces a
magnetization in the $x$ direction. We then monitor the relaxation of the
total magnetization after turning off $h$ suddenly at $t=0$. We find the total
magnetization decays either monotonically or with long-lived oscillations,
depending on $A=-\xi_{1}^{2}+\xi_{2}^{2}-4\xi_{2}\xi_{3}+\xi_{3}^{2}$. The
oscillation only appears for $A>0$, in which case the frequency satisfies
$\Omega\propto J\sqrt{A}$. Importantly, this phenomenon should be understood
as a universal property of the relaxation dynamics since the criterion only
contains a specific combination of anisotropic parameters, instead of full
details of the microscopic model (1). To validate our theoretical prediction,
we further perform numerical simulations based on the Kadanoff-Baym (KB)
equation with melon-diagram approximations, and the exact diagonalization
(ED). Numerical results show the theoretical prediction works in both cases,
although we are limited to a small system size $N=8$ in the ED. Our work also
provides a novel theoretical framework to analyze the dynamics of randomly
interacting quantum spin models.

_Theoretical Analysis.–_ We are interested in the relaxation dynamics of total
magnetization. Our theoretical analysis is based on the path-integral approach
on the Keldysh contour, as elaborated in [48, 49]. To begin with, we observe
that the random spin model can be written in terms of Abrikosov fermion
operators $\hat{c}_{i,s}$ with spins $s=\uparrow,\downarrow$ in the single
occupation subspace. Explicitly, we have
$\hat{S}^{\alpha}_{i}=\sum_{ss^{\prime}}\frac{1}{2}\hat{c}^{\dagger}_{i,s}(\sigma^{\alpha})_{ss^{\prime}}\hat{c}_{i,s^{\prime}}$,
where $\alpha=x,y,z$ and $\sigma^{\alpha}$ denote the corresponding Pauli
matrices. Since the Hamiltonian (1) exhibits $\pi$ rotation symmetries along
the $x$ axis, the total magnetization can only be along the $x$ axis. We thus
introduce $m(t)\equiv\langle\hat{S}^{x}(t)\rangle$. After imposing the
symmetry constraints [50], the magnetization can be computed by real-time
Green’s functions of fermion operators:

$m(t)=-iG^{>}_{\uparrow\downarrow}(t,t)=-iG^{<}_{\uparrow\downarrow}(t,t),$
(2)

where we have defined
$G^{>}_{ss^{\prime}}(t_{1},t_{2})\equiv-i\sum_{l}\left<c_{l,s}(t_{1})c_{l,s^{\prime}}^{\dagger}(t_{2})\right>/N$
and $G^{<}_{ss^{\prime}}(t_{1},t_{2})\equiv
i\sum_{l}\left<c_{l,s^{\prime}}^{\dagger}(t_{2})c_{l,s}(t_{1})\right>/N$.

The relaxation dynamics of $m(t)$ can then be computed once we obtain the
Green’s functions $G^{\gtrless}(t_{1},t_{2})$. It is known that the evolution
of $G^{\gtrless}(t_{1},t_{2})$ is governed by the Kadanoff-Baym equation,
which can be derived by the Schwinger-Dyson equation on the Schwinger-Keldysh
contour.

$\displaystyle i\partial_{t_{1}}$ $\displaystyle
G^{\gtrless}+\frac{1}{2}h_{\text{eff}}(t_{1})\sigma^{x}G^{\gtrless}=\Sigma^{R}\circ
G^{\gtrless}+\Sigma^{\gtrless}\circ G^{A},$ (3)
$\displaystyle-i\partial_{t_{2}}$ $\displaystyle
G^{\gtrless}+\frac{1}{2}h_{\text{eff}}(t_{2})G^{\gtrless}\sigma^{x}=G^{R}\circ\Sigma^{\gtrless}+G^{\gtrless}\circ\Sigma^{A}.$

Here we have introduced self-energies $\Sigma^{\gtrless}$ and $\Sigma^{R/A}$.
We define the operation $\circ$ for functions with two time variables as
$f\circ g\equiv\int dt_{3}~{}f(t_{1},t_{3})g(t_{3},t_{2})$. The retarded and
advanced Green’s functions $G^{R/A}$ are related to $G^{\gtrless}$ by
$G^{R/A}=\pm\Theta\left(\pm t_{12}\right)\left(G^{>}-G^{<}\right)$, where
$\Theta(t)$ is the Heaviside step function. Similar relations work for self-
energies $\Sigma^{R/A}$. $h_{\text{eff}}(t)=h(t)+\bar{J}m(t)$ is the effective
magnetic field, which includes the mean-field contribution from $\bar{J}$. For
$t<0$, the system is prepared in thermal equilibrium. Consequently, we have
$G^{\gtrless}(t_{1},t_{2})=G^{\gtrless}_{\beta}(t_{12})$ for $t_{1},t_{2}<0$.
For either $t_{1}>0$ or $t_{2}>0$, the Green’s functions evolve due to the
quantum quench and should be obtained by solving Eq. (3) after the self-energy
is specified.

The approximation comes in when we try to relate the self-energies to Green’s
functions. After transforming into the Abrikosov fermion representation, the
random Heisenberg interaction takes the form of random fermion scatterings.
Interestingly, such random interaction terms is a close analog of the
celebrated complex Sachdev-Ye-Kitaev (SYK) model [51, 52, 53, 54, 55, 56, 57].
Motivated by this observation, here we make the melon diagram approximation
for the fermion self-energy. A formal argument to control errors is to
generalize the Hamiltonian (1) into large-$M$ spins, as in the seminal work by
Sachdev and Ye [58]. Consequently, we have

$\displaystyle\Sigma^{\gtrless}$
$\displaystyle(t_{1},t_{2})=\frac{J^{2}}{4}\sum_{\alpha,\alpha^{\prime}}\xi_{\alpha}\xi_{\alpha^{\prime}}$
(4)
$\displaystyle\sigma^{\alpha^{\prime}}G^{\gtrless}(t_{1},t_{2})\sigma^{\alpha}\text{Tr}\left[\sigma^{\alpha^{\prime}}G^{\gtrless}(t_{1},t_{2})\sigma^{\alpha}G^{\lessgtr}(t_{2},t_{1})\right],$

Here $\alpha,\alpha^{\prime}\in\\{1,2,3\\}$ and we omit spin indices for
conciseness. We have introduced the anisotropy vector
$\bm{\xi}=(\xi_{1},\xi_{2},\xi_{3})$. The melon diagram approximation may fail
in the low-temperature limit if the system exhibits spin glass orders [59]. In
this work, we avoid this problem by focusing on the high-temperature regime
with $\beta J\ll 1$. Combining Eq. (2), (3), and (4) leads to a set of closed
equations which determines the relaxation of the magnetization.

Figure 2: The numerical result for the evolution of the magnetization $m(t)$
by numerically solving: (a) the Kadanoff-Baym equation (3) and (b) exact
diagonalization with system size $N=8$. Initially, the system is in thermal
equilibrium with $\beta J=0.04$, $\bar{J}=0$ and $h/J=10$. We take $\xi_{1}=1$
and consider four different anisotropy parameters
$(\xi_{2},\xi_{3})=(0.8,-1.5)$, $(0.8,1.5)$, $(1,-2)$, and $(1,2)$, which
corresponds to $A=6.69$, $-2.91$, $12$, and $-4$. The results show that the
relaxation of $m(t)$ is monotonic/oscillating if $A<0$/$A>0$. These two
numerical results match each other to good precision despite a small $N$.

Typical numerical results for $m(t)$ obtained by two methods are shown in FIG.
2. Here we consider examples with $\xi_{1}=1$ and
$(\xi_{2},\xi_{3})=(0.8,-1.5)$, $(0.8,1.5)$, $(1,-2)$, and $(1,2)$. We set the
initial temperature $\beta J=0.04$, the polarization field $h/J=10$ and
$\bar{J}=0$. In the long-time limit, the system exhibits the quantum
thermalization to the thermal ensemble with $h=0$. In this case, $\pi$
rotations along $y$ or $z$ also become the symmetry of the Hamiltonian, which
makes $m(\infty)=0$. According to the relaxation process, different anisotropy
parameters can be divided into two groups, under which $m(t)$ relaxes
monotonically (for $(\xi_{2},\xi_{3})=(0.8,1.5)$ and $(1,2)$) or with long-
lived oscillations (for $(\xi_{2},\xi_{3})=(0.8,-1.5)$ and $(1,-2)$).
Furthermore, we numerically checked that the presence of the oscillation is
stable against deformations of parameters. As a result, we propose the
Hamiltonian (1) with different $(\xi_{1},\xi_{2},\xi_{3})$ can be separated in
parameter regimes with oscillating relaxation (OR) versus non-oscillating
relaxation (NOR), as shown in FIG. 1. In supplementary material [50], we
verify that the difference in the dynamical behavior can not be detected in
equilibrium, for example, via spin susceptibility.

_Oscillation mode.–_ After evolving for a long time, the total magnetization,
as well as off-diagonal components of Green’s functions, becomes very small.
Consequently, we can perform a linearized analysis of the KB equation to
reveal the mechanism for the oscillation and determine the criterion for
different dynamical behaviors. A differential equation that governs the long-
time evolution of the magnetization can be derived following a few steps:

Step 1.– The linearized analysis can be largely simplified after the Keldysh
rotation. We introduce the standard Keldysh Green’s function of fermions as
$G^{K}=G^{>}+G^{<}$. The total magnetization can be expressed as its off-
diagonal component:

$m(t)=-iG^{K}_{\uparrow\downarrow}(t,t)/2.$ (5)

We can further combine equations in (3) to derive the equation for $G^{K}$.
The result reads

$G^{K}=G^{R}\circ\Sigma^{K}\circ G^{A},\ \ \ \ \ \text{with}\
\Sigma^{K}=\Sigma^{>}+\Sigma^{<}.$ (6)

Step 2.– We linearize Eq. (6) around the equilibrium solution in the long-time
limit after the quantum thermalization. We expand
$G^{a}(t_{1},t_{2})=G^{a,\beta_{f}}(t_{12})+\delta G^{a}(t_{1},t_{2})$, where
$G^{a,\beta_{f}}(t)$ is the equilibrium Green’s function on the final state.
Leaving details into the supplementary material [50], the off-diagonal element
of (6) reads

$\displaystyle\delta G^{K}_{\uparrow\downarrow}$
$\displaystyle=G^{R,\beta_{f}}_{\uparrow\uparrow}\circ\delta\Sigma^{K}_{\uparrow\downarrow}\circ
G^{A,\beta_{f}}_{\uparrow\uparrow},$ (7)
$\displaystyle\delta\Sigma^{K}_{\uparrow\downarrow}$
$\displaystyle=\frac{1}{4}J^{2}A\left((G^{>,\beta_{f}}_{\uparrow\uparrow})^{2}+(G^{<,\beta_{f}}_{\uparrow\uparrow})^{2}\right)\delta
G^{K}_{\uparrow\downarrow}.$

where we have introduced
$A=-\xi_{1}^{2}+\xi_{2}^{2}-4\xi_{2}\xi_{3}+\xi_{3}^{2}$ as in the
introduction. Since $G^{K,\beta_{f}}_{\uparrow\downarrow}=0$, Eq. (5) is
equivalent to $m(t)=-i\delta G^{K}_{\uparrow\downarrow}(t,t)/2$.

Step 3.– To proceed, we need to obtain approximations for
$G^{a,\beta_{f}}_{\uparrow\uparrow}$. In thermal equilibrium with $h=0$, the
self-energies (4) can be simplified as

$\Sigma^{\gtrless,\beta_{f}}_{ss^{\prime}}(t)=-\frac{J^{2}\bm{\xi}^{2}}{2}G^{\gtrless,\beta_{f}}_{ss}(t)^{3}\delta_{ss^{\prime}},$
(8)

where we have used $G^{>,\beta_{f}}_{ss}(t)=-G^{<,\beta_{f}}_{ss}(-t)$ due to
the particle-hole symmetry. Eq. (8) then matches the self-energy of the
Majorana SYK4 model with effective coupling constant $J|\bm{\xi}|/\sqrt{2}$.
It is known that at high temperatures $\beta J\ll 1$, the SYK model can be
described by weakly interacting quasi-particles [60]. Taking the Lorentzian
approximation, we have

$\displaystyle G_{\uparrow\uparrow}^{R/A,\beta_{f}}(t)\approx\mp i\Theta(\pm
t)e^{-\Gamma|t|/2},\ \ \ \ \ \
G_{\uparrow\uparrow}^{\gtrless,\beta_{f}}(t)\approx\mp ie^{-\Gamma|t|/2}/2.$
(9)

with quasi-particle decay rate $\Gamma\propto J$.

Step 4.– Finally, combining (7) and (9), $\delta G^{K}_{\uparrow\downarrow}$
satisfies the differential equation

$\left(\partial_{t_{1}}+\frac{\Gamma}{2}\right)\left(\partial_{t_{2}}+\frac{\Gamma}{2}\right)\delta
G^{K}_{\uparrow\downarrow}=-\frac{A}{8}J^{2}e^{-\Gamma|t_{12}|}\delta
G^{K}_{\uparrow\downarrow}.$ (10)

We have inversed the retarded/advanced Green’s functions (9) using
$\left(\partial_{t}+\frac{\Gamma}{2}\right)G^{R,\beta_{f}}_{\uparrow\uparrow}(t)=\left(\partial_{t}-\frac{\Gamma}{2}\right)G^{A,\beta_{f}}_{\uparrow\uparrow}(t)=-i\delta(t)$.

Eq. (10) is the starting point for analyzing the relaxation dynamics. Since it
is invariant under time translations, we separate out the center-of-mass time
dependence by introducing $\delta
G^{K}_{\uparrow\downarrow}(t_{1},t_{2})=\text{Re}~{}e^{-\lambda\frac{t_{1}+t_{2}}{2}}\varphi(t_{12})$.
The relaxation is oscillatory only if $\lambda$ is complex. Interestingly,
$\varphi(t_{12})$ then satisfies the 1D Schr$\ddot{\text{o}}$dinger equation

$-\frac{(\Gamma-\lambda)^{2}}{4}\varphi(t_{12})=-\partial_{t_{12}}^{2}\varphi(t)+\frac{A}{8}J^{2}e^{-\Gamma|t_{12}|}\varphi(t_{12}),$
(11)

where $-\frac{(\Gamma-\lambda)^{2}}{4}$ plays the role of the energy $E$ and
$\frac{A}{8}J^{2}e^{-\Gamma|t_{12}|}$ plays the role of potential $V$. Eq.
(11) suggests the boundary line between the oscillating regime and the non-
oscillating regime is at $A=0$: For $A<0$, the potential energy is negative.
It is known that in 1D any attractive potential exhibits at least one bound
state. Denoting the energy of the ground state as $-|E_{0}|$, we can solve
$\lambda=\Gamma-2\sqrt{|E_{0}|}$, which is real. Consequently, we expect the
magnetization relaxes monotonically. For $A>0$, the potential is repulsive.
The eigenstates of the (11) are scattering modes with continuous positive
energy $E$. We find $\lambda=\Gamma\pm 2i\sqrt{E}$, which is complex. This
leads to oscillations in the relaxation process.

To further determine the typical oscillation frequency $\Omega$, we need to
determine the typical energy $E$ that contributes to the quench dynamics.
According to Eq. (5), the magnetization probes the decay of the wave function
at $t_{12}=0$, where the potential energy is $\sim AJ^{2}$. For $E\ll AJ^{2}$,
the eigenstate has exponentially small weight near $t_{12}=0$. As a result,
the corresponding contribution to $m(t)$ can be neglected. We can approximate

$m(t)\sim\int_{AJ^{2}}dE~{}c(E)e^{-\Gamma t-2i\sqrt{E}t}.$ (12)

Here $c(E)$ is some smooth function determined by the initial condition. We
then expect $\Omega\approx c_{0}\sqrt{A}J$, with some $O(1)$ constant $c_{0}$
which does not depend on parameters in the Hamiltonian (1) and should be
extracted using numerics. Interestingly, the result predicts the oscillation
period $T=2\pi/\Omega$ diverges as we approach $A=0$, which can be viewed as
an analog of the divergence of the correlation length in traditional phase
transition described by order parameters.

We comment that our results unveil the universality of relaxation dynamics in
random spin models. Although the microscopic model in (1) contains several
parameters, the criterion for the different relaxation behaviors, as well as
the oscillation frequency, only depends on a specific combination $A$. This is
a direct analog of universality in the scattering theory, where for a
complicated potential, the low-energy scattering problem can only depend on a
specific combination of microscopic parameters, which is the scattering
length.

Figure 3: A comparison between the theoretical prediction
$\Omega_{Th}/J=c_{0}\sqrt{A}$ and numerical simulations. Here we choose
$c_{0}=1$. In each case, we randomly choose the anisotropic parameter
$(\xi_{1},\xi_{2},\xi_{3})$. Initially, the system is also in thermal
equilibrium with $\beta J=0.04$, $\bar{J}=0$, and $h/J=10$. The numerical data
are obtained by fitting numerics based on the KB equation and ED method. The
error bars correspond to standard deviations when concerning the different
fitting regions.

_Numeric tests.–_ Now we compare the theoretical predictions in the last
section to numerics. We can firstly test our criteria using numerics by
solving the KB equation. After fixing $\xi_{1}=1$, we find $A=6.69$, $-2.91$,
$12$, and $-4$ respectively for $(\xi_{2},\xi_{3})=(0.8,-1.5)$, $(0.8,1.5)$,
$(1,-2)$, and $(1,2)$. Our theory then predict oscillating relaxation for
$(0.8,-1.5),(1,-2)$ and non-oscillating relaxation for $(0.8,1.5),(1,2)$. This
is consistent with the evolution of $m(t)$ as shown in FIG. 2(a). To test our
theoretical prediction beyond the melon diagram approximation, we also perform
numerical simulations using the ED for a system size $N=8$ and 1000 random
realizations. The system was initially prepared in thermal equilibrium with
$\beta J=0.04$ and $h/J=10$. We choose anisotropy $\xi_{1}=1$,
$(\xi_{2},\xi_{3})=(0.8,-1.5)$, $(0.8,1.5)$, $(1,-2)$, and $(1,2)$, as in the
previous section. Despite using a small $N$, FIG. 2(b) shows the
oscillating/non-oscillating relaxation in the ED matches the prediction of the
KB equation.

We further compare our prediction of the oscillation frequency $\Omega\approx
c_{0}\sqrt{A}J$ to numerical results. We obtain $\Omega$ in numerics by
fitting $m(t)=m_{0}\cos(\Omega t+\theta)e^{-\Gamma t}+m_{\text{offset}}$. Here
$m_{0}$ is the amplitude, $\theta$ is the phase, $\Gamma$ is the quasi-
particle decay rate, and $m_{\text{offset}}$ is the offset which is
significant in the finite $N$ ED numerics. The fitting particularly focuses on
the matching in the small $m(t)$ region. Hence the detailed fitting region and
the error bars caused by such ambiguity are left to the supplementary material
[50]. The results are shown in FIG. 3. We randomly choose the anisotropic
parameters, and the first 15 cases correspond to $A>0$, and the last 4 cases
to $A<0$[50]. Among the $A>0$ cases, the mean ratio between the numerical data
and the polynomial $A$ reads
$\overline{\Omega_{\text{KB}}/(J\sqrt{A})}=0.995\pm 0.018$ and
$\overline{\Omega_{\text{ED}}/(J\sqrt{A})}=0.94\pm 0.04$. Therefore, we set
$c_{0}=1$ for theoretical predictions in FIG 3. Although the error bars for ED
numerics are significantly larger than KB numerics since the calculation is
based on the finite $N=8$ system, we find the theoretical prediction of the
oscillation frequency almost matches the KB results and the ED results, up to
the error bars. From FIG. 3, most notably, the OR and NOR relaxation are
sharply distinguished by the $A>0$ or $A<0$ criterion, which is perfectly
aligned with our theoretical analysis.

_Discussions.–_ In this work, we show that the random Heisenberg model with
all-to-all interactions exhibits universal relaxation dynamics governed by a
single parameter $A=-\xi_{1}^{2}+\xi_{2}^{2}-4\xi_{2}\xi_{3}+\xi_{3}^{2}$.
Unlike traditional examples where the universality emerges in the low-energy
limit, here the universal physics appears at high temperatures. For $A<0$, the
magnetization decays monotonically after we turn off the polarization field.
For $A>0$, long-lived oscillation appears during the relaxation process, with
a frequency $\Omega\propto J\sqrt{A}$. Our theoretical analysis is based on
the path-integral approach on the Keldysh contour, which is verified by
comparing our theory to numerical simulations by solving the KB equation or
the ED.

We remark that quantum coherence is essential for the existence of the
oscillating relaxation regime. As an example, if we spoil the coherence by
considering time-dependent random interactions instead of static interactions,
the magnetization is expected to decay monotonically: After replacing $J_{ij}$
with Brownian variables $J_{ij}(t)$, Eq. (10) is replaced by

$\left(\partial_{t}+{\Gamma}\right)\delta
G^{K}_{\uparrow\downarrow}(t,t)=-\frac{AJ}{8}\delta
G^{K}_{\uparrow\downarrow}(t,t),$ (13)

as derived in the supplementary material [50]. This results in $m(t)\sim
e^{-(\Gamma+AJ/8)t}$ with a simple exponential decay, on contrary to the
existence of different dynamical behaviors in the static case.

We also point out that amazingly our criteria $A>0$ for the oscillation regime
matches the criteria proposed in [57] for the presence of the instability
towards the formation of wormholes with $\xi_{1}=\xi_{2}=1$. However, the
analysis in [57] focuses on the low-temperature regime, while in this work we
focus on high temperatures. This makes it difficult to establish a direct
relationship between the two theoretical analyses. It would be interesting if
there is some version of duality between the high-temperature and the low-
temperature limit. Since the wormhole phase is non-chaotic, it is also
interesting to study the out-of-time correlator or the operator size
distribution in regimes with different dynamical behaviors.

On the experimental side, the quench experiment recently performed in NMR
systems directly measured the magnetization evolution when turning off the
external field, and our results can be straightforwardly verified in this
experiment[61].

Note Added. Universal behaviors of auto-correlation function related to the
quench dynamics discussed here, including oscillatory versus non-oscillatory
behavior, have been related to Lanczos coefficients computed for determining
the Krylov complexity in Ref. [62].

Acknowledgment. We are especially grateful to the invaluable discussions with
Hui Zhai, whose advice is indispensable for the whole work. We thank Riqiang
Fu, Yuchen Li, Xinhua Peng, Xiao-Liang Qi and Ren Zhang for their helpful
discussions. PZ is partly supported by the Walter Burke Institute for
Theoretical Physics at Caltech.

## References

  * Hohenberg and Halperin [1977] P. C. Hohenberg and B. I. Halperin, Theory of dynamic critical phenomena, Rev. Mod. Phys. 49, 435 (1977).
  * Ódor [2004] G. Ódor, Universality classes in nonequilibrium lattice systems, Rev. Mod. Phys. 76, 663 (2004).
  * Zhai [2021] H. Zhai,  _Ultracold Atomic Physics_ (Cambridge University Press, 2021).
  * Pal and Huse [2010] A. Pal and D. A. Huse, Many-body localization phase transition, Phys. Rev. B 82, 174411 (2010).
  * Nandkishore and Huse [2015] R. Nandkishore and D. A. Huse, Many-Body Localization and Thermalization in Quantum Statistical Mechanics, Annu. Rev. Condens. Matter Phys. 6, 15 (2015).
  * Abanin  _et al._ [2019] D. A. Abanin, E. Altman, I. Bloch, and M. Serbyn, Colloquium: Many-body localization, thermalization, and entanglement, Rev. Mod. Phys. 91, 021001 (2019).
  * Turner  _et al._ [2018a] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Serbyn, and Z. Papić, Weak ergodicity breaking from quantum many-body scars, Nat. Phys. 14, 745 (2018a).
  * Lin and Motrunich [2019] C.-J. Lin and O. I. Motrunich, Exact Quantum Many-Body Scar States in the Rydberg-Blockaded Atom Chain, Phys. Rev. Lett. 122, 173401 (2019).
  * Iadecola  _et al._ [2019] T. Iadecola, M. Schecter, and S. Xu, Quantum many-body scars from magnon condensation, Phys. Rev. B 100, 184312 (2019).
  * Lin  _et al._ [2020] C.-J. Lin, V. Calvera, and T. H. Hsieh, Quantum many-body scar states in two-dimensional Rydberg atom arrays, Phys. Rev. B 101, 220304 (2020).
  * Choi  _et al._ [2019] S. Choi, C. J. Turner, H. Pichler, W. W. Ho, A. A. Michailidis, Z. Papić, M. Serbyn, M. D. Lukin, and D. A. Abanin, Emergent SU(2) Dynamics and Perfect Quantum Many-Body Scars, Phys. Rev. Lett. 122, 220603 (2019).
  * Serbyn  _et al._ [2021] M. Serbyn, D. A. Abanin, and Z. Papić, Quantum many-body scars and weak breaking of ergodicity, Nat. Phys. 17, 675 (2021).
  * Ho  _et al._ [2019] W. W. Ho, S. Choi, H. Pichler, and M. D. Lukin, Periodic Orbits, Entanglement, and Quantum Many-Body Scars in Constrained Models: Matrix Product State Approach, Phys. Rev. Lett. 122, 040603 (2019).
  * Turner  _et al._ [2018b] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Serbyn, and Z. Papić, Quantum scarred eigenstates in a Rydberg atom chain: Entanglement, breakdown of thermalization, and stability to perturbations, Phys. Rev. B 98, 155134 (2018b).
  * Oka and Aoki [2009] T. Oka and H. Aoki, Photovoltaic hall effect in graphene, Phys. Rev. B 79, 081406 (2009).
  * Kitagawa  _et al._ [2010] T. Kitagawa, E. Berg, M. Rudner, and E. Demler, Topological characterization of periodically driven quantum systems, Phys. Rev. B 82, 235114 (2010).
  * Rudner  _et al._ [2013] M. S. Rudner, N. H. Lindner, E. Berg, and M. Levin, Anomalous Edge States and the Bulk-Edge Correspondence for Periodically Driven Two-Dimensional Systems, Phys. Rev. X 3, 031005 (2013).
  * Nathan and Rudner [2015] F. Nathan and M. S. Rudner, Topological singularities and the general classification of Floquet–Bloch systems, New J. Phys. 17, 125014 (2015).
  * Yao  _et al._ [2017] S. Yao, Z. Yan, and Z. Wang, Topological invariants of floquet systems: General formulation, special properties, and floquet topological defects, Phys. Rev. B 96, 195303 (2017).
  * Xu  _et al._ [2022] P. Xu, W. Zheng, and H. Zhai, Topological micromotion of Floquet quantum systems, Phys. Rev. B 105, 045139 (2022).
  * [21] P. Xu and T.-S. Deng, Topologically protected boundary discrete time crystal for a solvable model, arXiv:2210.15222 .
  * Li  _et al._ [2018] Y. Li, X. Chen, and M. P. A. Fisher, Quantum zeno effect and the many-body entanglement transition, Phys. Rev. B 98, 205136 (2018).
  * Li  _et al._ [2019] Y. Li, X. Chen, and M. P. A. Fisher, Measurement-driven entanglement transition in hybrid quantum circuits, Phys. Rev. B 100, 134306 (2019).
  * Skinner  _et al._ [2019] B. Skinner, J. Ruhman, and A. Nahum, Measurement-induced phase transitions in the dynamics of entanglement, Phys. Rev. X 9, 031009 (2019).
  * Gullans and Huse [2020] M. J. Gullans and D. A. Huse, Dynamical purification phase transition induced by quantum measurements, Phys. Rev. X 10, 041020 (2020).
  * Chan  _et al._ [2019] A. Chan, R. M. Nandkishore, M. Pretko, and G. Smith, Unitary-projective entanglement dynamics, Phys. Rev. B 99, 224307 (2019).
  * Jian  _et al._ [2021] S.-K. Jian, C. Liu, X. Chen, B. Swingle, and P. Zhang, Measurement-Induced Phase Transition in the Monitored Sachdev-Ye-Kitaev Model, Phys. Rev. Lett. 127, 140601 (2021).
  * Jian  _et al._ [2020] C.-M. Jian, Y.-Z. You, R. Vasseur, and A. W. W. Ludwig, Measurement-induced criticality in random quantum circuits, Phys. Rev. B 101, 104302 (2020).
  * Zhang  _et al._ [2022] P. Zhang, C. Liu, S.-K. Jian, and X. Chen, Universal Entanglement Transitions of Free Fermions with Long-range Non-unitary Dynamics, Quantum 6, 723 (2022).
  * [30] Y.-N. Zhou, Generalized Lindblad Master Equation for Measurement-induced Phase Transition, arXiv:2204.09049 .
  * yoon Choi  _et al._ [2016] J. yoon Choi, S. Hild, J. Zeiher, P. Schauß, A. Rubio-Abadal, T. Yefsah, V. Khemani, D. A. Huse, I. Bloch, and C. Gross, Exploring the many-body localization transition in two dimensions, Science 352, 1547 (2016).
  * Smith  _et al._ [2016] J. Smith, A. Lee, P. Richerme, B. Neyenhuis, P. W. Hess, P. Hauke, M. Heyl, D. A. Huse, and C. Monroe, Many-body localization in a quantum simulator with programmable random disorder, Nat. Phys. 12, 907 (2016).
  * Bernien  _et al._ [2017] H. Bernien, S. Schwartz, A. Keesling, H. Levine, A. Omran, H. Pichler, S. Choi, A. S. Zibrov, M. Endres, M. Greiner, V. Vuletić, and M. D. Lukin, Probing many-body dynamics on a 51-atom quantum simulator, Nature 551, 579 (2017).
  * Zhang  _et al._ [2023] P. Zhang, H. Dong, Y. Gao, L. Zhao, J. Hao, J.-Y. Desaules, Q. Guo, J. Chen, J. Deng, B. Liu, W. Ren, Y. Yao, X. Zhang, S. Xu, K. Wang, F. Jin, X. Zhu, B. Zhang, H. Li, C. Song, Z. Wang, F. Liu, Z. Papić, L. Ying, H. Wang, and Y.-C. Lai, Many-body Hilbert space scarring on a superconducting processor, Nat. Phys. 19, 120 (2023).
  * Lindner  _et al._ [2011] N. H. Lindner, G. Refael, and V. Galitski, Floquet topological insulator in semiconductor quantum wells, Nature Phys 7, 490 (2011).
  * Rechtsman  _et al._ [2013] M. C. Rechtsman, J. M. Zeuner, Y. Plotnik, Y. Lumer, D. Podolsky, F. Dreisow, S. Nolte, M. Segev, and A. Szameit, Photonic Floquet topological insulators, Nature 496, 196 (2013).
  * [37] Google-Quantum-AI and Collaborators, Quantum information phases in space-time: Measurement-induced entanglement and teleportation on a noisy quantum processor, arXiv:2303.04792 .
  * Hazzard  _et al._ [2014] K. R. A. Hazzard, B. Gadway, M. Foss-Feig, B. Yan, S. A. Moses, J. P. Covey, N. Y. Yao, M. D. Lukin, J. Ye, D. S. Jin, and A. M. Rey, Many-Body Dynamics of Dipolar Molecules in an Optical Lattice, Phys. Rev. Lett. 113, 195302 (2014).
  * Yan  _et al._ [2013] B. Yan, S. A. Moses, B. Gadway, J. P. Covey, K. R. A. Hazzard, A. M. Rey, D. S. Jin, and J. Ye, Realizing a lattice spin model with polar molecules, Nature 501, 521 (2013).
  * Choi  _et al._ [2017] J. Choi, S. Choi, G. Kucsko, P. C. Maurer, B. J. Shields, H. Sumiya, S. Onoda, J. Isoya, E. Demler, F. Jelezko, N. Y. Yao, and M. D. Lukin, Depolarization dynamics in a strongly interacting solid-state spin ensemble, Phys. Rev. Lett. 118, 093601 (2017).
  * Kucsko  _et al._ [2018] G. Kucsko, S. Choi, J. Choi, P. C. Maurer, H. Zhou, R. Landig, H. Sumiya, S. Onoda, J. Isoya, F. Jelezko, E. Demler, N. Y. Yao, and M. D. Lukin, Critical thermalization of a disordered dipolar spin system in diamond, Phys. Rev. Lett. 121, 023601 (2018).
  * Smale  _et al._ [2019] S. Smale, P. He, B. A. Olsen, K. G. Jackson, H. Sharum, S. Trotzky, J. Marino, A. M. Rey, and J. H. Thywissen, Observation of a transition between dynamical phases in a quantum degenerate Fermi gas, Sci. Adv. 5, eaax1568 (2019).
  * Schauß  _et al._ [2012] P. Schauß, M. Cheneau, M. Endres, T. Fukuhara, S. Hild, A. Omran, T. Pohl, C. Gross, S. Kuhr, and I. Bloch, Observation of spatially ordered structures in a two-dimensional Rydberg gas, Nature 491, 87 (2012).
  * Signoles  _et al._ [2021] A. Signoles, T. Franz, R. Ferracini Alves, M. Gärttner, S. Whitlock, G. Zürn, and M. Weidemüller, Glassy Dynamics in a Disordered Heisenberg Quantum Spin System, Phys. Rev. X 11, 011011 (2021).
  * Gabardos  _et al._ [2020] L. Gabardos, B. Zhu, S. Lepoutre, A. M. Rey, B. Laburthe-Tolra, and L. Vernac, Relaxation of the Collective Magnetization of a Dense 3D Array of Interacting Dipolar $S=3$ Atoms, Phys. Rev. Lett. 125, 143401 (2020).
  * Krojanski and Suter [2006] H. G. Krojanski and D. Suter, Reduced decoherence in large quantum registers, Phys. Rev. Lett. 97, 150503 (2006).
  * Lovrić  _et al._ [2007] M. Lovrić, H. G. Krojanski, and D. Suter, Decoherence in large quantum registers under variable interaction with the environment, Phys. Rev. A 75, 042305 (2007).
  * Kamenev [2011] A. Kamenev,  _Field Theory of Non-Equilibrium Systems_ (Cambridge University Press, 2011).
  * Stefanucci and van Leeuwen [2013] G. Stefanucci and R. van Leeuwen,  _Nonequilibrium Many-Body Theory of Quantum Systems: A Modern Introduction_ (Cambridge University Press, 2013).
  * [50] See Supplementary Material for: (1) The formulation of the equilibrium Green’s function; (2) The symmetry of the Green’s function; (3) Details of the numerical simulation; (4) Explicit calculation of the self energy; (5) Derivation of the normal mode equation with regular interactions; (6) Derivation of the normal mode equation with Brownian interactions; (7) Susceptibility in thermal equilibrium.
  * Kitaev [2015] A. Kitaev, Talk given at kitp program: Entanglement in strongly-correlated quantum matter (2015).
  * Maldacena and Stanford [2016] J. Maldacena and D. Stanford, Remarks on the Sachdev-Ye-Kitaev model, Phys. Rev. D 94, 106002 (2016).
  * Song  _et al._ [2017] X.-Y. Song, C.-M. Jian, and L. Balents, Strongly Correlated Metal Built from Sachdev-Ye-Kitaev Models, Phys. Rev. Lett. 119, 216601 (2017).
  * Davison  _et al._ [2017] R. A. Davison, W. Fu, A. Georges, Y. Gu, K. Jensen, and S. Sachdev, Thermoelectric transport in disordered metals without quasiparticles: The Sachdev-Ye-Kitaev models and holography, Phys. Rev. B 95, 155131 (2017).
  * Guo  _et al._ [2019] H. Guo, Y. Gu, and S. Sachdev, Transport and chaos in lattice sachdev-ye-kitaev models, Phys. Rev. B 100, 045140 (2019).
  * Gu  _et al._ [2020] Y. Gu, A. Kitaev, S. Sachdev, and G. Tarnopolsky, Notes on the complex Sachdev-Ye-Kitaev model, J. High Energ. Phys. 02 (2020), 157.
  * Zhou  _et al._ [2021] T.-G. Zhou, L. Pan, Y. Chen, P. Zhang, and H. Zhai, Disconnecting a traversable wormhole: Universal quench dynamics in random spin models, Phys. Rev. Research 3, L022024 (2021).
  * Sachdev and Ye [1993] S. Sachdev and J. Ye, Gapless spin-fluid ground state in a random quantum Heisenberg magnet, Phys. Rev. Lett. 70, 3339 (1993).
  * Baldwin and Swingle [2020] C. L. Baldwin and B. Swingle, Quenched vs Annealed: Glassiness from SK to SYK, Phys. Rev. X 10, 031026 (2020).
  * Zhang  _et al._ [2021] P. Zhang, Y. Gu, and A. Kitaev, An obstacle to sub-AdS holography for SYK-like models, J. High Energ. Phys. 03 (2020), 094.
  * [61] Y. C. Li et.al, to appear.
  * [62] R. Zhang and H. Zhai, Universal Hypothesis of Autocorrelation Function from Krylov Complexity, arXiv post, the same day.

See pages 1, of supp.pdf See pages 0, of supp.pdf