Amos Chan

Saumya Shivam

David A Huse

Andrea De Luca

Many-Body Quantum Chaos and Space-time Translational Invariance

We study the consequences of having translational invariance in space and in time in many-body quantum chaotic systems. We consider an ensemble of random quantum circuits, composed of singlesite random unitaries and nearest neighbour couplings, as a minimal model of translational invariant many-body quantum chaotic systems. We evaluate the spectral form factor (SFF) as a sum over many-body Feynman diagrams, which simplifies in the limit of large local Hilbert space dimension q. At sufficiently large t, diagrams corresponding to rigid translations dominate, reproducing the chaotic behavior of random matrix theory (RMT). At finite t, we show that translational invariance introduces additional mechanisms via two novel Feynman diagrams, known as the crossed and deranged diagrams, which delay the emergence of RMT. Our analytics suggests the existence of exact scaling forms which describe the approach to RMT behavior in the scaling limit where both t and L are large while the ratio between L and L Th (t), the many-body Thouless length, is fixed. We numerically demonstrate, with simulations of two distinct circuit models, that in such a scaling limit, most microscopic details become unimportant, and the resulting scaling functions are largely universal, remarkably being only dependent on a few global properties of the system like the spatial dimensionality, and the space-time symmetries.

Understanding the chaotic properties of quantum systems is a notoriously hard problem. A fruitful direction has been opened by the combination of two ingredients: First, fingerprints of an underlying chaotic dynamics are visible in the Hamiltonian spectrum of quantum systems [START_REF] Bohigas | Characterization of chaotic quantum spectra and universality of level fluctuation laws[END_REF]; second, spectral properties are best discussed in statistical terms [START_REF] Altshuler | Repulsion of energy levels and conductivity of small metal samples[END_REF]. This approach eliminates dependence on the microscopic details of the studied systems and brings out the universal characteristics of an ensemble of statistically similar Hamiltonians, which are captured by the random matrix theory (RMT) contrained only by symmetries [START_REF] Brody | Random-matrix physics: spectrum and strength fluctuations[END_REF][START_REF] Mehta | Random Matrices[END_REF]. RMT provides a prototype of thermalising dynamics for which the eigenstate thermalisation hypothesis [START_REF] Deutsch | Quantum statistical mechanics in a closed system[END_REF][START_REF] Srednicki | Chaos and quantum thermalization[END_REF][START_REF] Deutsch | Eigenstate thermalization hypothesis[END_REF] is confirmed [START_REF] Rigol | Thermalization and its mechanism for generic isolated quantum systems[END_REF].

However, RMT fails to reproduce the local structure of interactions of many-body quantum systems which results in a complex geometry and correlation in the Fock space [START_REF] Basko | Metal-insulator transition in a weakly interacting many-electron system with localized single-particle states[END_REF][START_REF] De Luca | Ergodicity breaking in a model showing many-body localization[END_REF][START_REF] De Luca | Anderson localization on the bethe lattice: Nonergodicity of extended states[END_REF][START_REF] Biroli | Delocalized glassy dynamics and many-body localization[END_REF][START_REF] Roy | Fock-space correlations and the origins of many-body localization[END_REF][START_REF] Roy | Percolation in fock space as a proxy for many-body localization[END_REF]. For this reason, random unitary circuits (RUC) have been proposed as toy models which utilize RMT while incorporating a notion of locality and dimensionality. In the simplest formulation, time evolution of RUC is performed by acting with randomly generated unitary gates on pairs of nearest neighbours in a spin lattice (Fig. 1a) [START_REF] Nahum | Quantum entanglement growth under random unitary dynamics[END_REF][START_REF] Keyserlingk | Operator hydrodynamics, otocs, and entanglement growth in systems without conservation laws[END_REF]. These models have proven fruitful in developing a unifying picture of the out-of-equilibrium dynamics of generic many-body systems with predictions for the entanglement growth [START_REF] Nahum | Quantum entanglement growth under random unitary dynamics[END_REF][START_REF] Zhou | Emergent statistical mechanics of entanglement in random unitary circuits[END_REF][START_REF] Li | Quantum zeno effect and the many-body entanglement transition[END_REF][START_REF] Skinner | Measurementinduced phase transitions in the dynamics of entanglement[END_REF][START_REF] Li | Measurement-driven entanglement transition in hybrid quantum circuits[END_REF][START_REF] Chan | Unitary-projective entanglement dynamics[END_REF][START_REF] Gullans | Dynamical purification phase transition induced by quantum measurements[END_REF][START_REF] Bao | Theory of the phase transition in random unitary circuits with measurements[END_REF][START_REF] Jian | Measurement-induced criticality in random quantum circuits[END_REF][START_REF] Zabalo | Critical properties of the measurement-induced transition in random quantum circuits[END_REF], and the outof-time-ordered correlators [START_REF] Keyserlingk | Operator hydrodynamics, otocs, and entanglement growth in systems without conservation laws[END_REF][START_REF] Nahum | Operator spreading in random unitary circuits[END_REF][START_REF] Rakovszky | Diffusive hydrodynamics of out-of-time-ordered correlators with charge conservation[END_REF][START_REF] Khemani | Operator spreading and the emergence of dissipation in unitary dynamics with conservation laws[END_REF]. More recently, Floquet random unitary circuits (FRUC) have been introduced by applying repeatedly the same set of random gates (Fig. 1b) [START_REF] Chan | Solution of a minimal model for many-body quantum chaos[END_REF][START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF][START_REF] Chan | Eigenstate correlations, thermalization, and the butterfly effect[END_REF][START_REF] Friedman | Spectral statistics and many-body quantum chaos with conserved charge[END_REF][START_REF] Chan | Spectral lyapunov exponents in chaotic and localized many-body quantum systems[END_REF][START_REF] Moudgalya | Spectral statistics in constrained many-body quantum chaotic systems[END_REF]. FRUC have given access to the study of non-trivial spectral properties in extended many-body systems. In particular, for the spectral form factor (SFF) [START_REF] Chan | Solution of a minimal model for many-body quantum chaos[END_REF][START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF][START_REF] Friedman | Spectral statistics and many-body quantum chaos with conserved charge[END_REF][START_REF] Chan | Spectral lyapunov exponents in chaotic and localized many-body quantum systems[END_REF][START_REF] Moudgalya | Spectral statistics in constrained many-body quantum chaotic systems[END_REF][START_REF] Haake | Quantum Signatures of Chaos[END_REF][START_REF] Bertini | Exact spectral form factor in a minimal model of many-body quantum chaos[END_REF][START_REF] Bertini | Random matrix spectral form factor of dual-unitary quantum circuits[END_REF][START_REF] Flack | Statistics of the spectral form factor in the self-dual kicked ising model[END_REF][START_REF] Kos | Many-body quantum chaos: Analytic connection to random matrix theory[END_REF][START_REF] Cotler | Black holes and random matrices[END_REF][START_REF] Cotler | Chaos, complexity, and random matrices[END_REF][START_REF] Saad | A semiclassical ramp in syk and in gravity[END_REF][START_REF] Gharibyan | Onset of random matrix behavior in scrambling systems[END_REF][START_REF] Li | Spectral statistics of non-hermitian matrices and dissipative quantum chaos[END_REF][START_REF] Winer | Hydrodynamic theory of the connected spectral form factor[END_REF][START_REF] Winer | Spontaneous symmetry breaking, spectral statistics, and the ramp[END_REF][START_REF] Vasilyev | Monitoring quantum simulators via quantum nondemolition couplings to atomic clock qubits[END_REF][START_REF] Joshi | Probing many-body quantum chaos with quantum simulators[END_REF][START_REF] Kos | Chaos and ergodicity in extended quantum systems with noisy driving[END_REF],

K(t, L) ≡ Tr[W (t)] Tr[W † (t)] (1)
where W (t) is the time evolution operator for time t, L is the system size and . . . indicates the ensemble average, it has been argued that the RMT behavior is recovered only for t > t Th (L), with t Th (L) the SFF Thouless time.

t Th (L) is an intrinsic time scale which generally grows unbounded with the system size L (with the exception of the dual-unitary circuits [START_REF] Bertini | Exact spectral form factor in a minimal model of many-body quantum chaos[END_REF][START_REF] Bertini | Random matrix spectral form factor of dual-unitary quantum circuits[END_REF][START_REF] Flack | Statistics of the spectral form factor in the self-dual kicked ising model[END_REF][START_REF] Akila | Particletime duality in the kicked ising spin chain[END_REF]). Its origin traces back to the existence of domain walls separating growing chaotic subregions [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF][START_REF] Garratt | Many-body quantum chaos and the local pairing of feynman histories[END_REF][START_REF] Garratt | Many-body delocalisation as symmetry-breaking[END_REF]. In this article, we consider the effect of translational invariance in space and time on the SFF. We introduce a spatially translational invariant (TI) version of the random phase model (RPM) [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF] on a d-dimensional lattice of length L, which can also be time-periodic (Floquet) or not. We will refer to the four setups resulting from the combination of TI and time-periodicity as cases (a), (b), (c) and (d) as illustrated in Fig. 1. We show that the SFF is exactly computable in the limit of large local Hilbert space dimension q via a diagrammatic expansion made up of contractions between Tr[W (t)] and Tr[W † (t)] (respectively top and bottom layer in Fig. 3b) of [START_REF] Bohigas | Characterization of chaotic quantum spectra and universality of level fluctuation laws[END_REF]. Before providing the explicit derivation, we outline the main results. At large t t Th (L) (but still with t t Heis (L), the Heisenberg time which is exponentially large in the system size), only ladder diagrams, corresponding to rigid translations of the top layer Tr[W (t)] w.r.t. the bottom layer Tr[W † (t)], contribute (see Fig. 3 top). This reproduces exactly the RMT predictions, i.e. K(t, L) ∼ K RMT (t, L), with

K RMT (t, L) ≡          1 , w/o symm. -(a) t , Floquet -(b) L d , TI -(c) tL d , TI + Floq. -(d)
.

(

) 2
The first two lines are standard and result from replacing the time evolution with a random matrix drawn from the circular unitary ensemble (CUE), either re-drawn at every time step (a), or repeated in time (b) [START_REF] Mehta | Random Matrices[END_REF]. The remaining lines of Eq. ( 2) can be understood observing that TI on a square lattice leads to L d momentum sectors, modeled as independent unitary blocks (still drawn from the CUE). This justifies the factors of L d for cases (c) and (d).

Additionally, we characterise the corrections to RMT at large but finite t t Th (L). Our diagrammatic calculations at infinite q indicate the existence of scaling forms in the scaling limit where both time t and the system size L are large but the ratio L/L Th (t) is kept fixed. Here, L Th (t) denotes the Thouless length, defined as the inverse function of the aforementioned Thouless time, i.e. L Th (t Th (L)) = L. Specifically, we obtain for the relevant cases (b), (c) and (d), the scaling forms lim

L,t→∞ L/L Th (t)=x K F (t, L) -t = κ F (x) , lim L,t→∞ L/L Th (t)=x L -1 K TI = κ TI (x) , (3) lim
L,t→∞ L/L Th (t)=x L -1 K TIF -t = κ TIF (x) .
for d = 1, and for general d below in Sec. III. Remarkably, not only the validity of the scaling forms (3) obtained at q → ∞ is confirmed by our numerics at finite q, but we also have evidence that the value of the scaling functions κ(x) in each case is universal, being independent of the microscopic details of the model, which only affect the non-universal L Th (t).

The paper is organized as follows. In Sec. I, we introduce the random circuit models in the preference of space-time translational invariance. In Sec. II, we discuss the numerical simulations at finite q confirming the existence of the scaling limit and its universality. In Sec. III, we derive the scaling forms (3) using the large-q diagrammatical computations. Lastly, in Sec. IV, we discuss the discrepencies between the finite-q numerics and large-q analytics, and the possible physical mechanism behind them.

I. MODELS

The random phase model (RPM) is defined by a quantum circuit which is a matrix product W (t) = t t =1 w(t ) where w(t) = w 2 (t) w 1 (t) is a q N × q N operator, with N the number of sites. The model can be defined on arbitrary lattices, but here we focus on integer lattices L in d-dimension with periodic boundary conditions and length L, so that L ≡ Z d L and N = L d . w 1 (t) = r∈L u(r, t) generates transformations at each site, with q ×q unitary matrices u(r, t) chosen from the CUE; w 2 (t) couples neighbouring sites and is diagonal in the basis of site orbitals with matrix elements [w 2 (t)] a1,...a N ;a1,...a N = exp ı r,r ϕ ar,a r (r, r , t) , where r, r are the nearest neighbours in L and a r ∈ {1, . . . , q}. We take each coefficient ϕ ar,a r (r, r , t) to be a Gaussian random variable with mean zero and variance , which effectively controls the coupling between neighbouring spins.

For the temporally and spatially random RPM (Fig. 1a), all unitaries u(r, t) and phases ϕ(r, r , t) are drawn independently. Correlation exists only between unitaries in Tr[W (t)] and their unique conjugates in Tr[W † (t)], which gives K(t, L) = 1 for all q. In [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF], the Floquet RPM (Fig. 1b) was considered where all gates are drawn independently in space but are constant in t.

Here, we also consider the TI RPM (Fig. 1c), where the gates are such that u(r, t) = u(t), φ(r, r , t) = φ (µ) (t), whenever rr = e µ , with e µ the unit vector in the µ spatial direction. The Floquet TI RPM (Fig. 1d), where the gates are also constant in time, is arguably the most realistic set-up.

II. NUMERICS

To test the universality of the scaling forms in (3), we numerically evaluate the SFF for RPM and additionally the brick wall model (BWM) defined in the SM [53] at finite q = 2, 3, t ≤ 8, and large L. Because of the computational cost, higher d are currently out of reach and we restrict to d = 1. The comparison between two different models is shown in Fig. 2, together with the infinite-q analytical predictions shown in green (see [START_REF] Deutsch | Eigenstate thermalization hypothesis[END_REF][START_REF] De Luca | Ergodicity breaking in a model showing many-body localization[END_REF][START_REF] Roy | Fock-space correlations and the origins of many-body localization[END_REF] below). In all cases, we see an excellent collapse among the different models and times, consistent with the scaling form and pointing at the existence of a universal scaling function. We stress that the only free parameter in this procedure is the Thouless length L Th (t), which rescales the horizontal axis for each t. In our procedure, we fix it by imposing that the numerical data at different t's all cross at a reference value x 0 = L/L Th (t) and equal the infinite-q expression (see Appendix F 1 for details). L/L Th (t) ; in all panels the RPM has q = 3, while the BWM has q = 2 in the left and center panels, and q = 3 in the right panel. The infinite-q scaling forms are plotted in green. The main panel of the right plot shows t = 2, 3, 4 from light to dark red. In the inset, L Th is plotted against t for TIF BWM (red) and TIF RPM (blue). The TIF RPM is not shown in the main panel because the apparent L Th (t) is too small (as seen in the inset) to be reliably estimated [53].

10 1 10 4 K F (t, L) -t d = 1 F (x) t = 2, F-BWM t = 3, F-BWM t = 4, F-BWM t = 5, F-BWM t = 6, F-BWM t = 4, F-RPM t = 5, F-RPM t = 6, F-RPM t = 7, F-RPM t = 8, F-RPM
For Floquet circuits (Fig. 2 left), using (3), we obtain an exceptionally good collapse for both models with the analytic infinite-q calculation (see [START_REF] De Luca | Ergodicity breaking in a model showing many-body localization[END_REF] below). Note that while the SFF for the Floquet RPM had already been computed in [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF], the universality of the corresponding scaling function in the scaling limit had not been observed before. However, for both TI (non-Floquet) and TI Floquet circuits (Fig. 2 middle and right), the scaling functions which emerge from numerical data are not well described by those computed at infinite q. The physical mechanism behind the discrepancies are discussed further in Sec. IV.

III. ANALYTICS

Now we sketch the exact large-q analytics of SFF for the TI, Floquet, and TI FLoquet cases systematically. The analytics allows us to derive the emergence of the RMT behaviour (2) in large-t, and the existence of scaling forms (3) describing the approach towards such emergence. Remarkably, as demonstrated in Sec. II the scaling forms are largely universal and depend only on q, spatial dimensionality, and the space-time symmetries. The dependence of universality classes in the spatial dimensionality was first observed in [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF] for Floquet systems and is even more striking for TI ones: In d = 1, corrections are controlled by crossed diagrams where sub-intervals in the top layer are rigidly contracted with those in the bottom layer (Fig. 3 middle); instead, in d ≥ 2, corrections are generated by deranged defect diagrams, where confinement forces excitations (on top of ladder diagrams) to be dilute (Fig. 3

bottom).

A. Translational invariant case

Here we compute the SFF of the TI-RPM in the limit q → ∞. To compute the SFF, we perform the ensemble averages in two steps: (i) ensemble average over the CUEs u(r, t ); (ii) ensemble averages over the random phases. Within a fixed time slice t , there are N copies of u(r, t ) / u † (r , t ) on the top/bottom layer. Following [START_REF] Chan | Solution of a minimal model for many-body quantum chaos[END_REF] (see also [53]), at the leading order in large q [54], the ensemble average is expressed as a sum over permutations, σ ∈ S(N ), pairing u(r, t ) with u † (σ(r), t ). Additionally, at leading order in large q, one is forced to take the same permutation on all t time slices, i.e. SFF is a single sum over σ ∈ S(N ) [53]. We now turn to the average over the random phases and we will see that it is natural to interpret it as a cost function associated to each σ. Expanding the orbital sum from phases at time slice t , we have a1,...,a N exp ı r,r ϕ ar,a r (t ) -ϕ a σ(r) ,a σ(r) (t ) where the sum is over nearest neighbor pair of sites. We see that, in large q, cancellations of the phases are only possible whenever σ maps nearest neighbour sites onto nearest neighbours (preserving the orientation). Using that e ıϕ = e -/2 and that all time slices contribute equally, we arrive at the expression

K TI = σ∈S N e -t(dN -N pb (σ)) (4)
where N pb (σ) = #{r, µ | σ(r + e µ ) -σ(r) = e µ } is the number of preserved bonds in any direction by the permutation σ (# denotes the cardinality of a set). The sum in Eq. ( 4) can be reorganized by grouping all σ ∈ S N with the same N pb (σ). We observe that N pb (σ) = N pb (τ σ) for all translations τ ∈ S N . Since the subgroup of translations is isomorphic to the lattice itself L, we arrive at

K TI = N dN n=0 A d (N , n)e -n t , (5)
where

A d (N , n) = #{σ ∈ S N /L | N pb (σ) = dN -n} and n = dN -N pb (σ)
is the number of broken bonds. Computing exactly the A d (N , n) poses a non-trivial combi-natorial problem, which nevertheless simplifies in d = 1 or for large N , as we show below. However, it is easy to see that A d (N , n = 0) = 1, corresponding to the identity equivalence class. Therefore in the limit t → ∞, we recover the RMT result for this case in Eq. ( 2). Generalising this construction, one can see that K converges to the dimension of the group of spatial symmetries. As example, for a two-dimensional TI circuit on a square lattice with rotational symmetry by angles of π/2, K TI = N = 4 in the limit of large t.

For d = 1, σ ∈ S L /Z L can be represented as cyclic permutations of L elements: Represent σ ∈ S L /Z L as σ ≡ (σ(1), σ(2), . . . ) and define an associated cyclic permutation with cycle σ . Then, σ with a fixed number of broken bonds n can be obtained as follows: We first partition L ≡ Z L = I 1 ∪ I 2 . . . I n into n adjacent nonempty intervals. Then we take any cyclic permutation σ of n elements such that σ(i) = σ(i + 1) + 1 (mod n) and rigidly map I i → I σ(i) . Clearly, the resulting mapping breaks exactly n bonds and all possible mappings can be uniquely constructed in this way. As an example, see Fig. 3e where L = 8, n = 3 and σ = (132). This leads to

A 1 (L, n) = L n a n (6)
where the binomial factor counts the partitionings of L and a n the possible σ. Although an explicit expression for the a n is not available, they correspond to a wellstudied sequence whose exponential generating function is known [START_REF] Sloane | The OEIS Foundation Inc[END_REF][START_REF] Stanley | Enumerative combinatorics volume 1 second edition[END_REF]. For both L, t large with fixed x = L/L Th (t) < 1, L Th (t) = e t , we can take L n ∼ L n /n! using the dominated convergence theorem [START_REF] Bartle | The Elements of Integration and Lebesgue Measure[END_REF] and sum over n to obtain the scaling form lim

L,t→∞ L/L Th (t)=x L -1 K d=1 TI = e -x (1 -ln(1 -x)) ≡ κ d=1 TI (x) .
(7) L Th (t) denotes the Thouless length. Note also that throughout the article, the limit of large q is always taken before the limits of large t and L.

For d > 1, evaluating the multiplicities A d (N , n) is a much harder task as they result from the interplay between permutations and the geometry of L. Nevertheless, the problem simplifies in the limit of large L at fixed n, as it corresponds to a dilute regime where a fixed number of bonds is broken in a very large system. Consider first a transposition which exchanges two sites. This will generally break 4d bonds (2d neighbours for each site), and therefore A d>1 (N , 4d) ∼ N 2 /2. More generally, we first pick the positions of k well-separated spatial defects and then we consider the possible ways of permuting them without leaving fixed points, so that precisely n = 2dk bonds are broken (Fig. 3g). These deranged defect diagrams lead to the asymptotic expansion

A d (N , n = 2dk) N →∞ = N k k! d k , d > 1 ( 8
)
where d k are the number of derangements [START_REF] Stanley | Enumerative combinatorics volume 1 second edition[END_REF][START_REF] Graham | Concrete mathematics: a foundation for computer science[END_REF] (i.e. permutation with no fixed points) of k elements [START_REF] Sloane | The OEIS Foundation Inc[END_REF]. The error we make in ( 8) is related to situations where the defects are close to one another thus breaking less bonds, but these are sub-leading in N . Once again, we consider the limit t, L → ∞ at fixed x = N /N Th (t), with N Th (t) = e 2d t the Thouless volume, and obtain lim

L,t→∞ N /N Th (t)=x N -1 K d>1 TI = e -x 1 -x ≡ κ d>1 TI (x) (9)
We stress that the apparent difference between the scaling functions ( 7) and ( 9) has a fundamental origin: In d = 1, extended intervals can be rigidly exchanged paying a cost only at their boundary; instead, in d > 1, exchanging two extended domains has a cost which grows with their boundary, i.e. extended defects are suppressed by confinement. Therefore, the leading contribution at large L and t is given by well-separated single-site excitations.

B. Floquet case

Before analysing the effect of combining time periodicity and translation invariance, we review the calculation of K F (t, L) for the Floquet RPM (Fig. 1a) [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF]. Here the single-site unitary gates are constant in time but are random in space u(r, t) = u(r). Thus, in the diagrammatic expansion of the SFF, we can choose any permutation σ ∈ S t to pair the t-copies of u(r) in the top layer with those of u † (r) in the bottom one. In the limit of large q, only time translations contribute, i.e. σ t (k) = k + t (mod t) with t = 0, . . . , t -1 [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF]. Therefore, we get a many-body diagram by choosing a configuration t(r) (color) for each site r ∈ L. After averaging over the random phases, it was shown that K F (t, L) = Z Potts , with Z Potts the partition function of a t-state Potts model with a Boltzmann weight across all bonds W t,t = e -t(1-δ t,t ) . At large t, the partition function is dominated by the t ferromagnetic groundstates where all sites have the same color, leading to the RMT prediction K F ∼ t. As t approaches t Th from above, excitations from the t ferromagnetic groundstates -the lowest-lying excitation being the domain wall states [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF][START_REF] Garratt | Many-body quantum chaos and the local pairing of feynman histories[END_REF] -become important. To access such corrections, in 1D, one makes use of the transfer matrix to write Z Potts = Tr[W L ]. Computing the spectrum of W and evaluating Z Potts in the scaling limit, this leads to the scaling form [53] lim

L,t→∞ L/L Th (t)=x K d=1 F -t = e x -x -1 ≡ κ d=1 F (x) (10)
with x = L/L Th (t) and L Th (t) = e t /t. In higher dimension, Z Potts cannot be easily computed for finite L and t. Nevertheless, in the scaling limit where L, t are both large, corrections to the RMT SFF are associated with diluted excitations where the color is changed with respect to the ground state's one. The position of the excitation can be chosen in ∼ N n /n! ways and each of them can be assigned any of the t -1 remaining color, with a cost e -2n dt . As a consequence, setting x = N /N Th (t) and N Th (t) = e 2d t /t, lim

L,t→∞ N /N Th (t)=x t -1 K d>1 F = e x ≡ κ d>1 F (x) . (11)
Intriguingly, note that, in contrast with Eqs. ( 7) and ( 9) which are divergent for any x ≥ 1, Eqs. ( 10) and ( 11) remain always smooth for finite x. This can be understood observing that infinite-q TI systems are mapped onto stat-mech models with non-local interactions, so that the scaling function is associated with the exchange of distant domains (d = 1) or defects (d > 1); on the contrary for Floquet systems, the resulting Potts model has purely local interactions.

C. Translational invariant Floquet case

We can now turn to the TI Floquet case. The same matrix CUE matrix u and u † appear tN times respectively in the top and bottom layer. However, at large q only the subgroup S N × Z N t ⊂ S tN , corresponding to arbitrary spatial permutations σ and time translations t(r) at each site. In d = 1, as explained in the TI case, the permutation σ corresponds to crossed diagrams, where spatial intervals in the top layer Tr[W (t)] are mapped onto intervals in the bottom one Tr[W † (t)] (e.g. Fig. 3e). Then, the cost associated to the average over the phases depends on the choices of t(r): within the same interval, the cost is given by the Boltzmann weights W as in the Floquet case; instead, between different intervals, the cost is always e -t irrespectively of the choice of t's at the interface. To account for the resulting combinatorics, it is useful to introduce the partition function Z(ω) = Tr((W + ωR) L ) = n ω n Z n (t, L), where R t,t = e -t is a rank 1 matrix with constant coefficients. In words, n counts the number of intervals and the factors Z n (t, L) contain the sum over all possible the colors with a n intervals. This leads to the explicit formula valid for arbirtrary t and L

K d=1 TIF (t, L) = ∞ n=0 a n Z n (t, L) . (12)
In the scaling limit, we obtain the behavior [53] lim

L,t→∞ L/L Th (t)=x L -1 K d=1 TIF -t = ln e -x 1 -x ≡ κ d=1 TIF (x) . ( 13
)
In d > 1, the large-t dominant contribution corresponds to ladder diagrams in space and a ferromagnetic groundstate in the color, leading as expected to K d>1 TIF (t, L) ∼ tL d . Corrections at large L, t are once again obtained by diluted excitations which can have two different origins: derangements as in [START_REF] Basko | Metal-insulator transition in a weakly interacting many-electron system with localized single-particle states[END_REF], which are now deranged defect diagrams in space-time; or color changes as in [START_REF] De Luca | Anderson localization on the bethe lattice: Nonergodicity of extended states[END_REF]. The two effects combine multiplicatively [53] lim

L,t→∞ N /N Th (t)=x (N t) -1 K d>1 TIF = 1 1 -x ≡ κ d>1 TIF (x) , (14)
with N Th (t) = e 2d t /t.

IV. DISCUSSION

The numerics shown in Fig. 2 (middle and right) shows a discrepancy between the infinite-q analytics and the finite-q numerics in the presence of space translation invariance. There can be different justifications behind this discrepancy: One possibility is that finite t corrections decay very slowly for TI systems. This is qualitatively confirmed in the infinite-q case, inspecting how the limit in ( 7) is approached increasing t [53]. A slow convergence is also to be expected due to the presence of singularities at finite x in [START_REF] Deutsch | Eigenstate thermalization hypothesis[END_REF][START_REF] Roy | Fock-space correlations and the origins of many-body localization[END_REF]. More probably, we have indications that Eq. ( 10) is more robust than [START_REF] Deutsch | Eigenstate thermalization hypothesis[END_REF]: By looking at the Floquet/TI RPM with p-site unit cell at infinite q, we find that the scaling function [START_REF] De Luca | Ergodicity breaking in a model showing many-body localization[END_REF] is independent of p while [START_REF] Deutsch | Eigenstate thermalization hypothesis[END_REF] is not [53]. Still, it might appear puzzling that finite-q numerics shows a collapse to a scaling function indepedent of q, which is nonetheless not in agreement with the infinite-q analytics. To elucidate this aspect, we propose a simple qualitative scenario. First, we observe that in the RPM, the two parameters and q control respectively the coupling in the space and time directions and the RMT behavior emerges when long-range order is established in both directions. In d = 1, corrections to RMT are then controlled by dilute excitations which break ordered domains, either in space or time, i.e. in the leading order, we have K(t, L) ∼ K RMT (t, L) + g 0 (L/L Th,0 (t)) + g 1 (L/L Th,1 (t)); where the subscript 0 and 1 refer respectively to the time and space directions, with the corresponding correlation lengths L Th,0/1 (t). The functions g 0 (x) and g 1 (x) tend to zero as x → 0, and are expected to be model independent and only dependent on the symmetries. In the large q limit at fixed L and t, the coupling in the time direction becomes infinitely strong with L Th,0 (t) q→∞ -→ ∞ so that only the function g 1 (x) survives in the decomposition above. For the Floquet RPM, we expect the relevant length scale to be L Th (t) ∼ L Th,1 (t)

L Th,0 (t), so that the scaling function is dominated by spatial domain walls, i.e. κ d=1 F (x) ∼ g 1,F (x). However, for TI systems, L Th (t) ∼ L Th,0 (t) L Th,1 (t), and therefore κ d=1 TI (x) ∼ g 0,TI (x). However, as pointed out above, this scaling function is not accessible if q → ∞ before L and t, thus explaining the observed discrepancy between the numerics and the analytics. In practice, the universal scaling function g 0,TI (x) observed in Fig. 2 results from temporal domain walls where contractions in different time slices take different permutation values. We stress that, because of unitarity, the SFF cannot diverge exponentially in t, therefore temporal domain walls must contribute both positively and negatively, in a distinct contrast to the spatial domain walls discussed in [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF][START_REF] Garratt | Many-body quantum chaos and the local pairing of feynman histories[END_REF]. This will be discussed in an upcoming work [START_REF] Chan | [END_REF].

A few additional comments are in order. Firstly, it would be beneficial to justify the universality which emerges from our work by means of a well-defined renormalization procedure. The main difficulty in this direction are the lack of locality and positivity of the resulting stat-mech model. Secondly, it is natural to expect that the scaling regime we identified is also visible in other quantities, as time-dependent correlation functions like out-of-time-ordered correlators. Thirdly, quasimomentum is conserved in TI lattice systems and affects its spectral properties but does not lead to the transport of an extensive conserved quantity because of Umklapp scattering; it is therefore interesting to explore its interplay with U (1) conserved charges. This will be discussed in an upcoming work [START_REF] Chan | [END_REF].

V. ACKNOWLEDGEMENTS.

Supplementary Material Many-Body Quantum Chaos and Space-time Translational Invariance

In this supplementary material we provide additional details about: A) SFF for translational invariant RPM 1. Derivation of Eq. 4 To derive Eq. 4, we compute the ensemble average of the TI-RPM in two steps: the ensemble average of (i) the unitaries drawn from the CUE in w 1 (t ), and (ii) the random phases in w 2 (t ). To carry out (i), we recall that the ensemble averages of unitaries u i,j drawn the CUE can be evaluated as

u i1j1 u i2j2 . . . u i j u * i 1 ,j 1 u * i 2 ,j 2 . . . u * i ,j = σ,τ ∈S Wg(σ -1 τ, q) k=1 δ i σ(k) ,i k δ j τ (k) ,j k
where S is the symmetric group of elements, and Wg(σ -1 τ, q) is the Weingarten function [START_REF] Samuel | u(n) Integrals, 1/n, and the Dewit-'t Hooft Anomalies[END_REF][START_REF] Mello | Averages on the unitary group and applications to the problem of disordered conductors[END_REF]. Wg(σ -1 τ, q) can be expressed as a polynomial in 1/q, and is of the leading order in q when σ = τ . We apply the above formula to average over unitaries in a fixed time slice t . In the large-q limit, the sum over permutation is dominated by the term with σ = τ ∈ S N , where N is the number of sites in the lattice. Together with all t time slices, the ensemble average is a sum over a vector of permutations, σ = (σ 1 , σ 2 , . . . , σ t ), where σ i ∈ S N . Observe that the average over the random phases can only maintain or decrease the order in 1/q of a given contribution. This implies that in the large-q limit, the leading contribution must have

σ 1 = σ 2 = • • • = σ t = σ ∈ S N ,
so that the number of sums over the site orbitals is maximized. Next, we perform the ensemble averages of the random phases in (ii). Expanding the orbital sum from phases at time slice t , we have and Tr[W † (t)] respectively. The delta functions in (SA.1) will lead to cancellations of phases in the exponent, when σ(r + e µ ) -σ(r) = e µ . In the large-q limit, each remaining phase ϕ arar+e µ will give rise to a factor of e ıϕ = e -/2 . Furthermore, since all time slices take the same contraction σ, the factors from (SA.1) are raised to the t-th power. This gives us the expression Eq. (4).

One dimensional case

Consider first of all a possible partitioning of L into n intervals I 1 , . . . , I n . We define the set

M n = {σ ∈ S n /Z n |∀i , σ(i + 1) = σ(i) + 1 (mod n)} . (SA.2)
As explained in the main text, we denote as a n its cardinality:

a n = #M n . The sequence a n satisfies n a n z n n! = e -z (1 -ln(1 -z)) , (SA.3)
which converges for all |z| < 1. Therefore, taking L large in Eq. ( 6), we have lim

L,t→∞ L/L Th (t)=x L -1 K d=1 TI = ∞ n=0 a n x n n! = e -x (1 -ln(1 -x)) ≡ κ d=1 TI (x) , (SA.4)
where x = L/L Th (t) and L Th (t) = e t , as given in Eq. ( 7) of the main text.

x ≥ 1 regime

According to the analysis in the previous section, we saw that the scaling function κ d=1 TI (x) diverges when x ≥ 1. We can investigate how does it diverge at finite but large L. In order to do so, we consider the exact expression for finite L at infinite q and d = 1, in Eqs. ( 4) and [START_REF] Srednicki | Chaos and quantum thermalization[END_REF].

κ d=1 TI (x; L) ≡ L n=0 L n a n x n L -n , (SA.5)
where we fixed x = L/L Th . Since for x ≥ 1, the expression in (SA.5) diverges at large L, its behavior is dominated by the large n expansion of the a n coefficients which can be deduced from (SA.3) and reads

a n n! e(n + 1) ≡ ãn (SA.6)
We can thus split the sum (SA.5) as

κ d=1 TI (x; L) L0 n=0 L n (a n -ãn )x n L -n + L n=0 L n ãn x n L -n , (SA.7)
where L 0 is large but finite. The first term in (SA.7) converges to the finite contribution

L0 n=0 L n (a n -ãn )x n L -n L→∞ -→ L0 n=0 (a n -ãn )x n n! (SA.8)
To evaluate the second term in (SA.7), we use the integral representation of the factorial as

n! = ∞ 0 dt e -t t n (SA.9) which leads to L n=0 L n ãn x n L -n = L ex(L + 1) ∞ 0 dt t e -Lt (1 + tx) L+1 -1 (SA.10)
For large L and x > 1, this last integral can be estimated via saddle point at t * = (x -1)/x, which leads to

κ d=1 TI (x; L) ∝ exp[L(ln x -1 + 1/x)] (SA.11)
For x = 1, the integral is dominated by small t * ∝ 1/L, which leads eventually to

κ d=1 TI (x; L) ln L 2e + O(1)
(SA.12)

Higher dimensional case

For higher dimension, we cannot obtain an explicit expression for K d>1 TI at finite N and t. However, we are still interested in the limit of large L and t. In this case, once again we are interested in the dilute limit, where a permutation only exchanges a fixed number of sites without leaving fixed points. For the number of derangements d n , we have the exponential generating function

∞ n=0 d n z n n! = e -x 1 -x . (SA.13)
Then, using ( 5) and ( 8) in the main text, we arrive at lim

L,t→∞ N /N Th (t)=x N -1 K d>1 TI (t, L) = ∞ n=0 x n d n n! = e -x 1 -x ≡ κ d>1 TI (x) , (SA.14)
where x = N /N Th (t) and N Th (t) = e 2d t , as given in Eq. ( 9) of the main text.

Appendix B: SFF for the Floquet RPM

As reviewed in the main text, the calculation of the SFF in the Floquet RPM [START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF] at q → ∞ can always be mapped into the partition function of a corresponding Potts model

K F (t, L) = Z Potts , (SB.1)
whose Boltzmann weights are given by the matrix W t,t = e -t(1-δ t,t ) and t, t = 0, . . . , t -1.

One dimensional case

In 1d, we can simply express the partition function by using the Transfer matrix

K d=1 F (t, L) = Tr[W L ] . (SB.2)
The spectrum of the matrix W contains two eigenvalues λ + > λ -: λ -= 1 -exp(-t) with degeneracy t -1 and λ + = 1 + e -t (t -1). Consequently,

K d=1 F (t, L) = (t -1)λ -(t) L + λ + (t) L . (SB.3)
In this case, setting L Th (t) = e t /t and x = L/L Th (t), we have lim

L,t→∞ L/L Th (t)=x K d=1 F (t, L) -t = e x -x -1 ≡ κ d=1 F (x) , (SB.4)
as given in Eq. ( 10) of the main text.

Higher dimensional case

For d > 1, the calculation of the Potts partition function poses a non-trivial problem. For d = 2, integrability can be used. However, here we are interested in the scaling limit where both the number of sites N and the time t are large. Let's set z = e -t and consider the small z (large ) expansion. We focus on d = 2 for simplicity, but the procedure can be easily extended to any d > 1. This corresponds to a low-temperature expansion in the ferromagnetic phase. At the zero-th order, Z Potts is simply given by the t possible groundstates. The leading correction is obtained changing the color of one site. This can be done on any site, choosing any among the t -1 remaining colors and will break 4 bonds; the corresponding contribution to the partition function is therefore N (t -1)z 4 . Higher orders in z are obtaining changing colors at more sites. For instance, considering two site changes, we have two cases, according to whether the two sites are nearest neighbours or not. The first case gives a contribution 4(t -1)N z 6 + 4(t -1)(t -2)N z 7 . The second case gives instead N (N -5)(t -1) 2 z 8 . In the limit, L, t → ∞ with fixed x = N tz 4 = N /N Th (t), we clearly see that only the last contribution survives. This corresponds to a dilute limit, which once accounting for the permutations among defects, takes in any d > 1 the form lim

L,t→∞ N /N Th (t)=x t -1 K d>1 F (t, L) = ∞ n=0 x n n! = e x ≡ κ d>1 F (x) , (SB.5)
with x = N /N Th N Th (t) = z -2d /t = e 2d t /t, as given in Eq. ( 11) of the main text.

Appendix C: SFF for translational invariant Floquet RPM

In the TI Floquet RPM, the diagrams appearing in the expansion of the SFF can be put in correspondence of i) spatial permutations σ; ii) time translations at each site. We can thus face the problem in two steps, first fixing the permutation σ and then summing over the possible choices of time translations at each site with fixed σ. We analyse separately the 1d from the d > 1 case.

One dimensional case

As explained for the temporal random case, fixing the permutation σ is equivalent to partitioning the L sites into intervals I 1 , . . . , I n and then mapping the intervals on the top layer onto intervals on the bottom layer using cyclic permutations σ ∈ M n in (SA.2). For a fixed choice of the interval and of σ, the sum over the possible choices of the time translations can be still written in terms of a Potts-like partition function with modified Boltzmann weights: bonds inside the same interval are given by the matrix W; bonds at the interface between two different intervals always give the trivial Boltzmann weight e -t . Since choosing a partition I 1 , . . . , I n can be done choosing the bonds where the edges of the intervals are, we can rewrite the SFF introducing the generalised partition function

Z(ω) = Tr[(W + ωR) L ] = ∞ n=0 ω n Z n , (SC.1)
where we introduced the rank 1 matrix R t,t = e -t . The coefficients Z n in the power series expansion of Z(ω) contain all the configurations where n intervals are present and thus the trivial Boltzmann weight R t,t is used. We thus have

K d=1 TIF (t, L) = L L n=0 a n Z n , (SC.2)
where, as before, a n is defined as the cardinality of the set in (SA.2). The coefficients Z n can be computed explicitly from the spectrum of W + ωR. Similarly to (SB.3), we have

Z(ω) = (t -1)λ -(t, ω) L + λ + (t, ω) L (SC.3)
where λ -(t, ω) = λ -(t) = 1 -e -t and λ + (t, ω) = 1 + e -t ((t -1) + ωt). It follows that

Z 0 = K d=1 F (t, L) , Z n = λ + (t, ω = 0) L L n t e t + t -1 n , n > 1 .
(SC.4) Taking once again the scaling limit, we obtain lim

L,t→∞ L/L Th (t)=x L -1 K d=1 TIF (t, L) -t = κ d=1 F (x) + e x (κ d=1 TI (x) -1) = -x -ln(1 -x) ≡ κ d=1 TIF (x) (SC.5)
with x = L/L Th (t) and L Th (t) = e t /t, as given in Eq. ( 13) of the main text.

Higher dimensional case

An explicit formula for the SFF at finite t and L is hard to derive in this case, because it would require the computation of the Potts model partition function in d-dimensions where some bonds have been removed. We focus on the scaling limit L, t → ∞ but with x = N /N Th (t) kept constant. As explained above this correspond to a dilute limit. We thus have the expansion

K d>1 TIF (t, L) ∼ tN ∞ n=0 (tN ) n n! d n e -2d nt ∞ m=0 N m m! (t -1) n e -2d mt . (SC.6)
The origins of the expression are explained as follows: In the first square bracket, the sum over n accounts for the deranged defects on top of the rigid spatial translation which can be chosen in N n ∼ N n /n! ways. The factor of t accounts for t choices of the t variables for these defects, and the factor d n counts their derangements. In the second square bracket, the sum over m accounts for the defects in the sum over the t variables. We can choose m defects in N m ∼ N m /m! ways and each of them can be independently changed with the t -1 other colours different from the background. Taking again the scaling limit lim

L,t→∞ N /N th (t)=x (tN ) -1 K d>1 TIF (t, L) = κ d>1 F (x) × κ d>1 TI (x) = 1 1 -x ≡ κ d>1 TIF (x) , (SC.7)
with N Th (t) = e 2d t /t and x = N /N Th (t), as given in Eq. ( 14) of the main text.

Appendix D: Generalised unit cell in large-q

In this appendix, we demonstrate the robustness of the scaling forms. To this end, we compute the scaling forms in the large-q limit of the modified circuit models, namely TI RPM (Floquet RPM) that are invariant under p-site (pdiscrete-time) translation in space (time). The 2-site and 2-discrete-time translational invariant model is illustrated in Fig. S1a andb respectively. A summary of the results is as follows. For TI RPM with p-site translational invariance, we show that for d = 1, the scaling form κ TI ) for TI RPM with p-site translational invariance. We will first present the derivation for p = 2 then generalise it to any p ∈ N + . Following the derivation of Eq. ( 5), we have for p = 2, a sum over two permutations, σ 1 and σ 2

K (2) TI (t, L eff ) = σ1∈S odd σ2∈Seven e -t(L-N pb (σ[σ1,σ2])) , (SD.1)
where L eff = L/2 is the number of unit cells, S odd / S even are permutations restricted to the odd/even sublattices and N pb (σ) is the number of preserved bonds as defined in the main text. The permutation σ is defined composing

σ 1 and σ 2 σ[σ 1 , σ 2 ](i) = σ 1 (i) if i is odd σ 2 (i) if i is even (SD.2)
It is difficult to evaluate (SD.1) exactly, but a lower bound can be identified. Observe that one can always restrict the broken bonds to be placed in between site i and i + 1 for even i. This means K

TI > L eff

L eff n=0 A 1 (L eff , n)
, where A 1 is given in Eq. ( 6) in the main text. In addition to the above configurations, one can find another set of diagrams where the broken bonds are placed between site i and i + 1 with odd i. Putting the two cases together and accounting for the overcounting of the cases where no bonds are broken we have

K (2) TI (t, L eff )/L eff ≥ L eff n1=0 L eff n 1 a n1 e -n1 t + L eff n2=0 L eff n 2 a n2 e -n2 t -1 (SD.3)
This expression is easily generalised to arbitrary p to

K (p) TI (t, L eff )/L eff ≥ p L eff n=0 L eff n a n e -n t -(p -1) (SD.4)
with L eff = L/p in general.

In the scaling limit where t, L → ∞ with fixed x = L eff /L Th = L eff e -t , we define

κ (p) TI (x) ≡ lim L,t→∞ L eff /L Th (t)=x L -1 eff K (p) TI (t, L) . (SD.5)
Note that with this definition κ

TI (x = 0) = 1, which corresponds to the RMT behavior at times. From (SD.4), we obtain the lower bound for the scaling function

κ (p) TI (x) ≥ pκ (1) TI (x) -(p -1) , (SD.6)
Note that this equation implies the existence of a divergence at x = 1 for all p's.

Higher-dimensional systems

For d > 1, with the foresight that we will use the scaling limit where only dilute deranged defects matter, we define a generic higher-dimensional TI RPM as follows: Consider a TI RPM made up of super-sites of p qudits, each of which evolve under independently drawn CUE, and couple with z other qudits under the (independently drawn) random phase gates defined in the main text. Additionally, we demand that the translational invariance of the supersites is preserved.

Consider the scaling limit where t, L → ∞ at fixed x = N eff /N Th with N eff = N /p and N Th = e z t . In this limit, SFF is a sum over deranged defect diagrams labelled by (k 1 , k 2 , . . . , k p ) where k i is the number of deranged dilute defects of the i-th type of qudits in the supersite. Following the derivation of Eq. ( 5), ( 8) and (9), in the limit of large N , we write

N -1 eff K (p) TI ∼ p i=1 N eff ki=0 d ki N ki eff k i ! N ki Th = N -1 eff K (1) TI p , (SD.7)
where again d k is the number of derangements of k elements. Using the definitions

κ (p) TI (x) ≡ lim N ,t→∞ N eff /N Th (t)=x N -1 eff K (p) TI (t, L) , (SD.8)
we arrive a more compact statement in the scaling limit

κ (p) TI = κ (1) TI (x) p .
(SD.9)

Again, with this convention we have κ

TI (x = 0) = 1, which corresponds to the RMT behavior at large times.

SFF for Floquet RPM with p-discrete-time translational invariance

The calculation of SFF for RPM with p-discrete-time translation invariance closely follows the one for RPM with 1-discrete-time translation invariance (see main text and [START_REF] Chan | Solution of a minimal model for many-body quantum chaos[END_REF][START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF]). For all dimensions d, upon averaging over the CUE in the large-q limit following [START_REF] Chan | Solution of a minimal model for many-body quantum chaos[END_REF][START_REF] Chan | Spectral statistics in spatially extended chaotic quantum many-body systems[END_REF], SFF is mapped to a Potts model with DOF t i at each site i. t i can take t eff = t/p number of possible states, corresponding to the t eff possible ladder diagrams. The averages over the random phases give an effective Boltzmann weight W t,t = e -t eff (1-δ t,t ) eff with eff = /p for nearest neighbour pairs of state t and t . Therefore, the SFF K (p) F (t eff , L) can be written as the partition function of a t eff -state Potts model with Boltzmann weight parametrized by eff .

One-dimensional systems

In d = 1, we define in the scaling limit,

κ (p) F (x) ≡ lim L,t→∞ N /N Th (t)=x K (p) F (t eff , L) -t eff .
(SD.10)

where L Th (t eff ) = e t eff eff /t eff = pe t /t. κ (p) F (x) can now be related to κ

F (x) as

κ (p) F (x) = e x -x -1 = κ (1)
F (x) .

(SD.11)

In other words, after properly defining the Thouless length, the scaling function is the same independently of p. Note that at large time, we have κ F (x = 0) = 0, and the RMT result is reproduced.

Higher-dimensional systems

For d ≥ 1, we define N Th = e z eff t eff /t eff with t eff = t/p and eff = p . Furthermore, we define

κ (p) F (x) ≡ lim L,t→∞ N /N Th (t)=x t -1 eff K (p) F (t eff , L) .
(SD.12)

In the scaling limit, κ

F coincides with κ

F as

κ (p) F (x) = κ (1)
F (x) , (SD. [START_REF] Roy | Fock-space correlations and the origins of many-body localization[END_REF] i.e. the scaling function is independent of p. Note again that at large time, we have κ F (x = 0) = 0, and the RMT result is reproduced.

where we have added subscripts to the traces to emphasize the Hilbert spaces on which matrices W and V act. The dual tensors for TI BWM, with a computational basis b = {b 1 , b 2 . . . , b 2t }, can be similarly written as

[v 1 ] b,b = t µ=1 ũb 2µ ,b 2µ+1 (µ) , (SF.3a) [v 2 ] b,b = t µ=1 ũb 2µ-1 ,b 2µ (µ) , (SF.3b)
where ũij,kl = u jl,ik , u ij,kl being a CUE random matrix acting on two sites, and ũ is the dual of the unitary gate, which is non-unitary in general. The dual Floquet operator is given by V = v 2 v 1 and V (L) = V L/2 .

SFF for translational invariant circuits

We elaborate on the middle panel of Figure 2, and provide further evidence of consistency with the predicted scaling form in Eq. 7. For both RPM and BWM, K TI (t, L) is computed by exact diagonalizing the dual Floquet operator V . In figure S2, we re-plot the data for RPM (left), BWM with q = 2 (middle), and also add the data for BWM with q = 3 in the right panel (whose importance is spelled out in the discussion on TIF circuits), for which only three system sizes were accessible. Averaging was done over 12000 -15000 realizations of V for RPM, 8000 -15000 realizations for BWM with q = 2, and over 4000 -10000 realizations for BWM with q = 3.

t = 2 t = 3 t = 4
FIG. S2. K(t, L)/L vs. L/L Th for TI-RPM at q = 3 (left), TI-BWM at q = 2 (middle), and TI-BWM at q = 3 (right).

Now we describe the details of how we obtain L Th (t) and plot the collapse in Figure 2. The horizontal axes in 2 are all scaled such that the finite-q numerics for K TI (t, L)/L is equal to the infinite-q scaling function κ d=1 TI (x) (Eq. 7) at a specific point x 0 < 1. We choose x 0 = 0.95 and determine L(t) for a given t such that K TI (t, L)/ L = κ d=1 TI (x 0 ). Since x 0 = L(t)/L Th (t) by definition, we get the corresponding Thouless time L Th (t) = L(t)/x 0 . Lastly, we rescale the horizontal axis to be x = L/L Th (t), so the numerics and the scaling function can be compared directly. We visually show this procedure in Figure S3, where the left panel shows the horizontal lines for a few x 0 that we draw to extract L(t), and the right panel shows L Th (t) vs t for those x 0 . We find little difference between the behaviour at different x 0 , and exponential fits to L Th (t) are consistent with the expected scaling form in Eq. 7.

L Th (t) for BWM (both with q = 2 and q = 3) were similarly calculated and is shown in the middle and the right panel of Figure S7.

As additional checks for consistency, we consider time direction simulations, using sparse matrix multiplication and trace evaluation, in Figure S4, where the left panel shows K TI (t, L)/L vs t for different L, averaged over 6000 -9000 realizations of W (t). Using this data, t Th (L) can be obtained analogous to L Th (t) by looking at the intersection of K TI (t, L)/L with a constant δ = κ d=1 TI (x 0 ) for some x 0 , as plotted in the middle panel. Note that here we are only interested in the scaling form of t Th (L) hence we ignore the normalization arising due to different δ. Good logarithmic fits indicate consistency with the form of L Th (t) = Le -eff t , where the fitted eff = 2 is due to finite q effects. The L Th (t)

x 0 = 0.85 x 0 = 0.9 x 0 = 0.95 2.51e 0.17t 2.68e 0.17t 2.99e 0.16t FIG. S3. Left : KTI(t, L)/L vs L for different L using space direction simulations of TI RPM, with the intersections L(t) with constant horizontal lines κ d=1 TI (x0) for a few x0 used to compute L Th (t) = L/x0 . Right : L Th (t) vs t for different choices of x0. The corresponding exponential fits (solid lines) show reasonable agreement with the expected form from Eq. 7 (L Th = e eff t ), with some effective eff = 2 due to finite-q effects. scaling form is then further corroborated in the right panel, showing a collapse with 0.2t -log(L), with eff being consistent across all plots. A logarithmic fit is consistent with the expectation that L Th (t) = Le -eff t , with an effective eff = 2 due to finite-q corrections, and also consistent with the fit in Fig. S3. Right : KTI(t, L)/L vs 0.2t -log(L).

SFF for Floquet circuits

Here, we elaborate on the left panel of Figure 2. We first show in Figure S5 the collapse of K F (t, L) -t vs L/L Th (t) separately for RPM with q = 3, = 2 (left), BWM with q = 2 (middle), and additionally BWM with q = 3(right); and we find that for all these circuits, the infinite-q scaling function in Eq. 10 is in excellent agreement with finite-q numerics for these circuits. Note that these are all space direction simulations performed through sparse matrix multiplication and trace computation and not through ED. The data for RPM has been averaged over 4000 -10000 realizations, and for BWM with q = 2 over 3000 -10000 realizations, while for BWM with q = 3, we averaged over 4000 -10000 realizations. L Th (t) were obtained in a similar fashion as for TI RPM (Figure S3), fixing x 0 = 3, and are shown in the left panel of Figure S7.

3 K F (t, L) -t d = 1 F (x) t = 2 t = 3 t = 4
FIG. S5. K(t, L) -t vs. L/L Th for F-RPM at q = 3 (left), F-BWM at q = 2 (middle), and F-BWM at q = 3 (right).

SFF for translational invariant Floquet circuits

In this part of the Appendix, we focus on transnational invariant Floquet circuits, and elaborate on the right panel of Figure 2. First, we separately show in Figure S6, K TIF /L -t vs L/L Th (t) for RPM with q = 3, = 2 (left), BWM with q = 2 (middle), and BWM with q = 3. The data is averaged over 8000 -15000 realizations of the dual Floquet operator V for RPM, 7000 -1000 realizations for BWM with q = 2, and 4000 -10000 realizations for BWM with q = 3. As can be seen in the figure, there aren't enough points for x = L/L Th (t) < 1 for RPM with q = 3 and BWM with q = 2, in order to see a collapse for x < 1, hence cannot be compared reliably for x > 1 either. This can be seen as a consequence of small L Th (t) for those two models (as shown in the right panel of Figure S7). TIF-RPM with q = 3, = 2, in particular, has the additional problem where the extracted L Th (t) in fact are probing the region beyond the RMT ramp, when t > t Hei ≈ q L /L. Note that this is not an issue for F-RPM, since t Hei = q L for that model. These problems are not present in the data for BWM with q = 3 (right panel of Figure S7), because L Th (t) are sufficiently large for the accessible values of t, allowing data points for x < 1. Note that L Th (t) for t = 5 is obtained by stochastic sampling of the trace through multiplication with complex random vectors, and is not as accurate as the points for t ≤ 4.

To further investigate on the dynamics of TIF-RPM with q = 3, = 2, we look at K TIF /L vs t using time direction simulations, and show that t Th is of the order of t Hei for accessible L. The left panel of Figure S8 plots K TIF /L vs t, averaged over 8000 -10000 realizations of the Floquet operator W , along with the RMT behavior K RMT , given by K RMT (t, L) = tL t < t Hei q L t ≥ t Hei .

(SF.4)

It can be observed in the middle panel that the K TIF /L approaches K RMT /L at times larger than t Hei . The middle panel shows a rolling average of the data from the plot in the left panel. To quantify t Th , one could draw a horizontal line near zero and look at the intersection of the difference between K TIF /L and K RMT /L. A large t Th is consistent with the small L Th (t) we see in Figure S7, and is also reflected in the ratio of consecutive level spacings r =min(δ n , δ n+1 )/max(δ n , δ n+1 ), where δ n = |φ n -φ n-1 | is the difference between consecutive eigenphases {φ n }, which is plotted in the right panel, averaging over 100 realizations of the Floquet operator for each model. The plot shows that the spacing ratio for TIF-RPM hasn't converged to the GUE value for the accessible system sizes, unlike the other two models we have studied. L/L Th for TIF-RPM at q = 3 (left), TIF-BWM at q = 2 (middle), and TIF-BWM at q = 3 (right). F-RPM (q = 3) F-BWM (q = 2) F-BWM (q = 3) TI-RPM (q = 3) TI-BWM (q = 2) TI-BWM (q = 3) TIF-RPM (q = 3) TIF-BWM (q = 2) TIF-BWM (q = 3) As mentioned in the main text, for TI and TIF models, the differences between the scaling collapse of the finite-q numerics and the infinite-q solution (Fig. 2), can be due to (i) a genuinely different scaling function for finite-q, or (ii) a slow convergence in t to the same scaling function in infinite-q. To investigate the latter possibility, in Figure S9, we compare the finite -L, -t, infinite-q results for the SFF with the infinie-q scaling functions reached at large L and t, for Floquet, TI and TI Floquet models. we find that finite t corrections decay very slowly for TI systems compared to the Floquet systems. L/L Th (t) S9. For all three plots, the green curves are the infinite-t, infinite-L, infinite-q scaling functions given in Eq. 10, 7, and 13 (from left to right respectively). Left: KF(t, L) -t vs L/L Th (t) for different t where the red curves denote the finite-t, finite-L, infinite-q expression for KF(t, L)-t (Eq. SB.3), with L Th for F-RPM (from Figure S7). Middle: KTI(t, L)/L vs L/L Th different t, alongside the finite-t, finite-L, infinite-q expression for KTI(t, L)/L (Eq. SA.5) in reds, with L Th for TI-BWM with q = 2. Right: KTIF(t, L)/L -t vs L/L Th different t, alongside the finite-t, finite-L, infinite-q expression for KTIF(t, L)/L -t (Eq. SC.2) in reds, with L Th for TI-BWM with q = 3.

FIG. 1 .

 1 FIG. 1. Illustrations of the different types of RUC for the RPM: (a) Temporally and spatially random RPM; (b) Floquet (and spatially random) RPM; (c) TI (and temporally random) RPM, and (d) TI Floquet RPM. For each case, gates of the same colors are identical.

2. One dimensional case 3 .

 3 Higher dimensional case B) SFF for Floquet RPM. 1. One dimensional case 2. Higher dimensional case C) SFF for translational invariant Floquet RPM 1. One dimensional case 2. Higher dimensional case D) Generalised unit cell in large-q 1. SFF for TI RPM with p-site translational invariance 2. SFF for Floquet RPM with p-discrete-time translational invariance E) Brick wall model (BWM) F) Numerical methods and results 1. SFF for translational invariant circuits 2. SFF for Floquet circuits 3. SFF for translational invariant Floquet circuits 4. Comparison between scaling forms and finite-t, finite-L, and infinite-q solutions Appendix A: SFF for translational invariant RPM 1. Derivation of Eq. 4

2 αβ 2 /

 22 ϕa r a r+eµ -ϕ a r a r+eµ δ ar,a σ(r) , (SA.[START_REF] Bohigas | Characterization of chaotic quantum spectra and universality of level fluctuation laws[END_REF] where Dϕ αβ = exp -ϕ √ 2π dϕ αβ , and where {a} and {a } labels the orbital degrees of freedom in Tr[W (t)]

F

  FIG.S1. Illustrations of (a) Floquet (and spatially random) RPM which is invariant under 2-discrete-time translation; and (b) TI (and temporally random) RPM which is invariant under 2-site translation. For each case, gates of the same colors are identical.

5 - 4 . 8 FIG

 548 FIG. S4. Left : KTI(t, L)/L vs t for different L obtained using time direction simulations of TI RPM. Middle : t Th (L) vs L, computed by numerically solving KTI(t Th , L)/L = δ. A logarithmic fit is consistent with the expectation that L Th (t) = Le -eff t , with an effective eff = 2 due to finite-q corrections, and also consistent with the fit in Fig.S3. Right : KTI(t, L)/L vs 0.2t -log(L).

10

 10

FIG. S7 .

 S7 FIG.S7. L Th (t) vs. t for Floquet circuits (left), TI circuits (middle), and TI Floquet circuits (right)

  FIG.S8. Left : KTIF(t, L)/L vs. t for TIF-RPM with q = 3, = 2, for L = 4, 5, 6, 7, 8 (from light to dark blue), and the dashed lines denote the corresponding RMT behavior defined in Eq. SF.4. Middle : KTIF(t, L)/L -KRMT(t, L)/L vs. t where the vertical lines correspond to tHei for the L corresponding to the same color. The intersection points with a horizontal line near zero give an estimate of t Th which appear larger than tHei. Right: Nearest neighbour level spacing ratio r in the zero momentum sector for the three TIF models studied in this paper.

  FIG.2. Scaled SFF vs. x = L/L Th at finite q for Floquet RPM and BWM (left), TI RPM and BWM (middle), and TI Floquet BWM (right)

	K TI (t, L)/L eff	10 4 10 7	d = 1 TI (x) t = 3, TI-BWM t = 4, TI-BWM t = 5, TI-BWM t = 6, TI-BWM t = 6, TI-RPM t = 7, TI-RPM t = 8, TI-RPM t = 5, TI-RPM		10 6 10 12 -t TIF (t, L)/L eff	L Th (t)	30 0 10 20	t 2 4 6 8
					K
		10 1			10 0
			10 -1	10 0	10 1			10 -1	10 0
			L/L Th (t)				L/L Th (t)

  FIG. S6. K(t, L)/L eff -t vs.

			d = 1 TIF (x)			d = 1 TIF (x)			d = 1 TIF (x)
	10 12	t = 4	10 12	t = 2	10 12	t = 2
	K TIF (t, L)/L -t	10 6	t = 5 t = 6 t = 7 t = 8	2K TIF (t, L)/L -t	10 6	t = 3 t = 4 t = 5 t = 6	2K TIF (t, L)/L -t	10 6	t = 3 t = 4
		10 0			10 0			10 0
		10 -1	10 0	10 1	10 -1	10 0	10 1	10 -1	10 0	10 1
			L/L Th (t)			L/L Th (t)			L/L Th (t)

AC and ADL warmly thank John Chalker for his guidance in related projects. DAH thanks Grace Sommers and Michael Gullans for a related collaboration. DAH is supported in part by NSF QLCI grant OMA-2120757. AC is supported by fellowships from the Croucher foundation and the PCTS at Princeton University.

Appendix E: Brick wall model (BWM)

The one-dimensional BWM is defined by a quantum circuit which is a matrix product

where w(t) = w 2 (t) w 1 (t) is a q L × q L operator.

is a tensor product of q 2 × q 2 unitary matrices u 2n-1,2n (t) chosen from the circular unitary ensemble (CUE) and acting on site 2n -1 and 2n.

is again a tensor product of u 2n,2n+1 (t) drawn from CUE except that the unitary gate acts on site 2n and 2n + 1.

The model is defined with periodic boundary condition with u L,L+1 (t) acting on site L and 1.

For the temporally and spatially random BWM, each unitary gate u n,n+1 (t) is drawn independently. For Floquet BWM, we take unitary gates acting on different pairs of sites to be independently drawn, while gates acting on the same pair of sites at different discrete time t to be identically drawn, i.e. w(t) = w(t ) for t = t and W (t) = w t . For TI BWM, gates in w 1 (t) (and separately in w 2 (t)) acting on different pairs of sites at the same time are identically drawn, but gates acting on different discrete times are independently drawn. For Floquet TI BWM, u n,n+1 (t) = u n ,n +1 (t ) for even integer n and for all t. The equation separately holds true for odd integer n.

Appendix F: Numerical methods and results

We simulate the RPM with q = 3, = 2 and d = 1, and the BWM with q = 2, 3, both in d = 1 and with periodic boundary conditions. (We find that RPM with q = 2 and large is not fully chaotic in the sense that the linear ramp does not appear in late t.) We consider two types of numerical simulations: "Time-direction" simulations and "space-(or dual-)direction" simulations, which involve performing calculations by acting matrices in the Hilbert space H = L C q and the dual Hilbert space H = t C q respectively. Although Figure 2 uses only dual direction numerics, here we discuss both the directions.

For the time-direction simulations with translational invariant (but not Floquet) circuits, the computation of K(t, L) vs t for different L involves sparse matrix multiplication of locally-supported random unitaries, and taking the trace. For translational invariant Floquet circuits, we compute K(t, L) for a fixed L with exact diagonalizations (ED).

For the space-direction simulations, we define a dual Floquet operator V , which we explicitly construct for RPM. V for BWM can be similar constructed following the procedure below. Informally, V is taken to be the tensor product of the first column of tensors in Fig. 1b and1c. More precisely, we introduce the computational basis b = {b 1 , . . . , b t } with each b µ = 1, . . . , q. The dual tensor for RPM can now be written as