Realistic many-body quantum systems vs full random matrices: static and dynamical properties

We study the static and dynamical properties of isolated many-body quantum systems and compare them with the results for full random matrices. In doing so, we link concepts from quantum information theory with those from quantum chaos. In particular, we relate the von Neumann entanglement entropy with the Shannon information entropy and discuss their relevance for the analysis of the degree of complexity of the eigenstates, the behavior of the system at different time scales and the conditions for thermalization. A main advantage of full random matrices is that they enable the derivation of analytical expressions that agree extremely well with the numerics and provide bounds for realistic many-body quantum systems.


I. INTRODUCTION
Advances in quantum information science and many-body quantum physics have been closely intertwined. Quantum computers, for example, are many-body quantum systems. To build the first, one needs a better understanding of the latter. At the same time, one of the main reasons for developing quantum computers is to simulate many-body quantum systems [1]. New numerical methods, such as the density matrix renormalization group employed in the studies of many-body quantum systems [2,3], rely on notions of entanglement. In particular, the accuracy of these methods requires a limited growth of entanglement with time and system size [4]. An important task in quantum information processing is the faithful transfer of quantum states. To accomplish this, one needs a precise characterization of the dynamics of many-body quantum systems at different time scales. In turn, studies of nonequilibrium quantum dynamics often touch upon the old quest of deriving statistical mechanics and thermodynamics from first principles.
Thermalization in isolated many-body quantum systems is caused by strong interactions between their particles (or quasiparticles) and is intimately related to the notion of quantum chaos [5][6][7]. Quantum chaos refers to signatures that one finds at the quantum level, such as level repulsion, that tell us whether or not the classical counterpart of the system is chaotic [8]. A main feature of a classically-chaotic system is the extreme sensitivity of its dynamics to the initial conditions. At the quantum level, one can no longer talk about phase-space trajectories, but one can still expect to find quantum signatures of classical chaos, such as those associated with spectrum properties [9]. This notion of quantum chaos has, however, been extended to systems that may not even have a classical limit. In addition to the properties of the spectrum, quantum chaos is also directly connected with the emergence of chaotic eigenstates, that is highly delocalized states that are similar to (pseudo)-random vectors [5,6]. The level of delocalization of the eigenstates is measured with quantities, such as the participation ratio and the Shannon information entropy [5,6,9].
In this work, we explore the relationship between measures of entanglement and measures of delocalization [10][11][12][13][14][15]. We show that the von Neumann entanglement entropy and the Shannon information entropy provide similar information about the system. The first requires the partial trace of the density matrix, being computationally more expensive than the second, which deals with the states of the entire system.
We also analyze the evolution of many-body quantum systems from very short times until the moment they reach equilibrium. Our analysis is based on the behavior of the survival probability (probability of finding the initial state later in time), the Shannon entropy and the entanglement entropy. We show that the short-and intermediate-time evolutions are generic when the system is taken far from equilibrium, while the long-time dynamics depends on the level of chaoticity of the initial state. At intermediate times, the decay of the survival probability may be faster than exponential [16][17][18][19][20][21][22][23][24]. At long times, the survival probability necessarily shows a power law decay [25]. From the value of the power law exponent, one may anticipate whether or not the initial state will thermalize.
We start our studies in Section II using full random matrices. These are matrices filled with random numbers. Early connections between the properties of the spectrum of quantum systems and classical chaos were made in the context of full random matrices and became known as the Bohigas-Giannoni-Schmit conjecture [26]. Full random matrices are not the most suitable models for many-body quantum systems, because they imply that all of the particles interact at the same time. Yet, they allow for the derivations of analytical results that can serve as references and bounds for the analysis of realistic models. In Section III, we then investigate the static properties and the dynamics of realistic many-body quantum systems described by one-dimensional spin-1/2 models. We consider both integrable and chaotic limits. These models are similar to those employed in experiments that use cold atoms to study nonequilibrium dynamics [27,28] and thermalization [29,30].

II. FULL RANDOM MATRICES AND THERMALIZATION
Full random matrices are matrices filled with random numbers whose only constraint is to satisfy the symmetries of the system they are trying to describe. They were used extensively by Wigner to model the spectra of heavy nuclei, which are very complex systems [31]. In this approach, interactions are treated statistically, and their details are overlooked. This simple idea led to results that agreed very well with data from real nuclei and was soon employed in the analysis of other complex systems, such as atoms, molecules and quantum dots [8]. We discuss some of the static properties of these matrices in Sections II A and II B, and some of their dynamic properties in Section II C. In Section II D, we discuss why thermalization is trivially achieved for full random matrices.

A. Eigenvalues: Density of States and Level Repulsion
There are different kinds of ensembles of full random matrices defined according to the symmetries that the matrices satisfy [8,32,33]. When modeling systems with time reversal symmetry, random matrices of the Gaussian orthogonal ensemble (GOE) are used. These are D × D real and symmetric matrices with entries from a Gaussian distribution with mean zero, In practice, GOE full random matrices can be obtained by generating a matrix with D 2 random numbers and then adding it to its transpose. The density of states of full random matrices follows the standard semicircle distribution [34], where 2E is the length of the spectrum, that is −E ≤ E ≤ E.
A key property of full random matrices and a main feature of quantum chaos is the strong repulsion between neighboring levels, as captured, for example, by the nearest-neighbor level spacing distribution P (s), where s is the spacing between neighboring rescaled energies. The unfolding procedure [9] guarantees that the mean level spacing of the rescaled eigenvalues is one. For the unfolded spectrum of GOE matrices, one finds: This contrasts with the level spacing distribution of a sequence of uncorrelated eigenvalues, where the levels are not prohibited from crossing and the distribution is Poisson, P (s) = exp(−s). Level repulsion causes the rigidity of the spectrum. As a result, the variance Σ 2 (ℓ) of the number of unfolded eigenvalues in an interval of length ℓ grows logarithmically with ℓ. For the GOE, one has: where γ e = 0.5772 . . . is Euler's constant. The level number variance of full random matrices is between the linear dependence Σ 2 (ℓ) = ℓ found for uncorrelated eigenvalues and Σ 2 (ℓ) = 0 reached by the completely rigid spectrum of the harmonic oscillator. P (s) and Σ 2 (ℓ) are complementary. The former characterizes the short-range fluctuations of the spectrum, and the latter characterizes the long-range fluctuations. Both are shown in Figures 1a and 1b, respectively.

B. Eigenstates: Delocalization and Entanglement Measures
There are various ways to quantify how much a state spreads out in a certain basis. One of them is the participation ratio P R. Given an eigenstate |ψ α = k C α k |φ k written in a basis |φ k , P R (α) = 1/ k |C α k | 4 . A comparable quantity is the Shannon information entropy, defined as: The values of P R (α) and S

(α)
Sh depend on the chosen basis. In the case of full random matrices, where one deals with ensembles of matrices filled with random numbers, the notion of basis is not well defined, but we can still associate values with the delocalization measures.
All eigenstates of full random matrices are (pseudo)-random vectors. In the particular case of GOEs, the coefficients are real random numbers from a Gaussian distribution satisfying the normalization condition. All eigenstates are therefore equivalent and lead to approximately the same values of the participation ratio, P R GOE ∼ D/3, and of the Shannon entropy: The latter is shown in Figure 1c. This result can be derived as follows. The sum is approximated by an integral, where P (C) = D/2πe −DC 2 /2 is the Gaussian probability distribution of the components. This approach guarantees that the average of the coefficients is C = 0, the average of their squares is C 2 = 1/D, and k |C α k | 2 → D ∞ −∞ C 2 P (C)dC = 1. The entropy is therefore: Notice that the Shannon entropy of random vectors is smaller than the maximum value S max Sh = ln D obtained when all weights |C α k | 2 = 1/D. Similarly, there are various methods for quantifying the amount of entanglement in a state. We focus on the von Neumann entanglement entropy, S vN [35]. Its computation requires the bipartition of the system in subsystems A and B and the partial trace of one of the two. S vN (ρ A ) is the von Neumann entropy of the reduced density matrix ρ A = Tr B ρ, where ρ is the density matrix associated with the total system. It gives: Maximum entanglement occurs when the reduced density matrix is maximally mixed, that is when it is the normalized identity matrix. In this case, S max vN ∼ ln D A , where D A is the dimension of ρ A . Studies have shown that this limit is nearly achieved by random pure states [36,37]. In particular, Page [37] Figure 1d) is slightly above our numerical results for the entanglement entropy. The values for D and D A that we select are motivated by the comparison that we later make in Section III with spin systems that have matrices of these same dimensions. We also find numerically that the results for the entanglement entropy are slightly larger than ln(0.48D A ) (solid line Figure 1d). However, since these differences are minor, when D A is large, we chose to write, in analogy with S GOE Sh , that:

C. Time Evolution: Entropy Growth and Survival Probability
To study the evolution under a GOE full random matrix, we assume that a fictitious basis of product states |φ k was used to construct the matrix and that the initial state |Ψ(0) = |φ ini is the basis vector in the middle of the matrix, k = ini = D/2. The entanglement entropy at t = 0 is therefore S vN (0) = 0, but it grows as the state evolves.
Written in the energy eigenbasis, the initial state is: Since the eigenstates |ψ α are random vectors, so is |Ψ(0) . The Shannon entropy of |Ψ(0) in the energy eigenbasis (also known as diagonal entropy; see Equation (24)) gives values very similar to what we find in Figure 1 We note, however, that to study dynamics, we consider the spreading in time of |φ ini over the other basis vectors. In Figure 2a, we show the evolution of the Shannon entropy given by: where the probability W k (t) for the initial state to be found in state |φ k at time t is: and Planck's constant is one. We can obtain W k (t) by taking the Fourier transform of the distribution P k,ini (E), if the latter is known. Alternatively, a simpler approach is to separate the contributions to W k (t) into those coming from k = ini and those from k = ini, as we do next.

Survival Probability and Power Law Decays
is the survival probability of the initial state. It corresponds to the probability of finding the initial state |Ψ(0) = |φ k=ini later in time, is the energy distribution of the initial state weighted by the components |C α ini | 2 . It is often referred to as the local density of states (LDOS) or strength function. For full random matrices, this distribution is a semicircle [17,18], as the density of states in Equation (2). This is seen in Figure 2a. The Fourier transform of a semicircle gives: where J 1 is the Bessel function of the first kind, is the variance of the energy distribution of the initial state and: is the energy of |Ψ(0) . For the initial state considered here, σ GOE ini ∼ E/2 and E GOE ini ∼ 0. The agreement between the numerical results for W ini (t) and Equation (14) is excellent, as seen in Figure 2b.
The asymptotic behavior of W GOE ini (t) for t ≫ σ −1 ini is a power law and ∝ t −3 [22,25]. This can be derived from a theorem proven in [38]. For the semicircle, . The theorem says that if 0 < ξ < 1 and η(E) is N times continuously differentiable for −E ≤ E < E, then for t → ∞, we have [38,39]: where η The dominant term for n = 0 leads to: This inverse power law decay is a manifestation of the Khalfin effect [40], which refers to the slow decay of the survival probability at long times due to the unavoidable bounds (at least the ever present lower bound) in the spectrum of any quantum system. Khalfin showed that the usual exponential decay of unstable states could not persist for long times. Later studies showed that the decay should become a power law [25,39,41,42] with the exponent depending on the properties of P ini,ini (E) [23][24][25].

Entropy Growth
We now come to Equation (11) and rewrite it as [43][44][45]: where the use of 1 − W ini (t) in the second line is motivated by the normalization condition D k=1 W k (t) = 1 and the ratio [1 − W ini (t)]/N pc by the fact that at very long times, W ini (t) should be very small and S Sh → ln N pc . N pc refers to the number of states that contribute to the evolution of the initial state. In the case of full random matrices, we would expect N GOE pc ∼ 0.48D. However, by doing an average of the numerical values of exp[S Sh (t)] at long times, after the saturation of the evolution, we actually find that N pc is slightly larger than 0.48D. Why this is so still needs to be understood. The semi-analytical expression above, with the result for N pc ∼ exp[S Sh (t)] , agrees extremely well with the numerical results for the Shannon entropy (see Figure 2c). The agreement is so good that the curves are indistinguishable. The same occurs for the entanglement entropy in Figure 2d using The entropies in Figure 2c,d initially grow nearly quadratically, as can be verified by expanding Equation (19). This gives: which matches the numerical results well for very short times. Subsequently, the data are very well fitted with a linear curve, as shown with circles in Figure 2c,d. This makes a connection with the linear behavior, S Sh,vN ∝ t, studied theoretically in realistic systems [4,[43][44][45][46][47] and also observed experimentally [30].

D. Relaxation and Thermalization
We say that an isolated finite quantum system equilibrated if, after a transient time, the initial state simply fluctuates around a steady state, remaining very close to it for most times (revivals may occur, but they take exceedingly long times to happen for large system sizes). Consider, for example, a few-body observable O evolving according to the equation, where O βα = ψ β |O|ψ α and O αα is the eigenstate expectation value of O. If the system equilibrates, there will be only fluctuations around the infinite time average O. In the case of systems without degeneracies, as in GOE matrices, this average is given by the first term in Equation (21), with the second term leading to the fluctuations. The size of the fluctuations is crucial to determine whether the system indeed reaches equilibrium [48][49][50][51][52][53][54][55][56]. For this, the fluctuations need to be small and decrease with D. Small fluctuations certainly occur when the energy spacings are not zero and the products C β * ini C α ini and the off-diagonal elements O αβ are small. This is evidently the case for full random matrices, but also for realistic many-body quantum systems with interactions, where the fluctuations decrease exponentially with system size [55]. We have found that the exponent of this decay depends on the level of delocalization of the initial state in the energy eigenbasis, with maximum exponents being found for full random matrices [55].
One can then talk about equilibration (or relaxation) without the need to invoke an environment. The loss of the initial coherence is, for all practical purposes, irreversible, because the recurrence time for systems with chaotic eigenstates is extremely long and increases with system size [57][58][59].

Infinite-Time Averages
The evolution of the Shannon entropy (11) can be written in terms of the energy eigenbasis as: is approximately: For GOE full random matrices, S Sh is close to S GOE Sh (6). It is also similar to the diagonal entropy S d , defined as [60][61][62], Equivalently, the infinite-time average of the entanglement entropy is close to S GOE vN (9). The infinite-time average of the survival probability is the inverse of the participation ratio, For GOE full random matrices, W GOE ini which is reached in Figure 2b. The size of the fluctuations decreases as 1/D. This is obtained from: On average, this term cancels out, except for Eα − E β = E δ − Eγ. Since there are no degeneracies of the energies or of the energy spacings in full random matrices, this equality requires that Eα = E δ and E β = Eγ , which gives: The standard deviation of the temporal fluctuations of Wini coincides with W ini.

Thermalization in Full Random Matrices
After relaxation, the observable is said to have thermalized if its infinite-time average coincides with a thermodynamic average. Equality can only happen in the thermodynamic limit, but thermalization is already suggested if the two averages are close and further approach each other as the system size increases. The comparison between the two averages is made explicit with the equation, where OME is the thermodynamic (microcanonical) average and N E ini ,δE is the number of energy eigenbases in the window δE taken around Eini. Equation (29) is valid when Oαα for eigenstates close in energy agree with the microcanonical average, an idea that is at the heart of statistical mechanics and has recently become known as the eigenstate thermalization hypothesis (ETH) [7,21,[63][64][65][66][67][68][69]. Notice that the ETH is not a condition for thermalization, but a statement of what one means by it. Equation (29) is trivially satisfied for the GOE full random matrices. Since the eigenstates are random vectors, Oαα is approximately the same for any eigenstate |ψα , so it can be taken out of the sum and each Oαα ∼ OME. The condition for thermalization is therefore quantum chaos, not only in the sense of level repulsion, but also in the sense of chaotic states [5,6,44,45,61,68,[70][71][72][73].
We stress that thermalization does not require an equal distribution of all probabilities, that is we do not need to have |C α k | 2 = 1/D, S (α) Sh = ln D and S (α) vN = ln DA for all eigenstates, but the eigenstates that take part in the evolution of the initial state should be statistically close to these limits. We emphasize also that, contrary to full random matrices, not all eigenstates of real systems are chaotic, even when the Hamiltonian shows level repulsion. This is further discussed in Section III. Therefore, the prerequisite for thermalization in real systems is to have equivalent values of Oαα in the energy window sampled by the initial state. This is satisfied if the eigenstates in that window are very similar to random vectors, closely fulfilling Equations (6) and (9). These two elements, energy window sampled by |Ψ(0) and chaotic states, make evident the key role of the interplay between initial state and Hamiltonian in the analysis of thermalization [68,69]. The third important element is the observable, as discussed also in [74]. The observables considered in the studies of thermalization in realistic systems are few-body observables.
As Dyson had already anticipated, the development of random matrix theory marked the beginning of a "new kind of statistical mechanics" [75], but at those early stages, the link between equilibrium and nonequilibrium dynamics, as put forward by Equation (29), was not yet well established. Current studies analyze Equation (29) for specific few-body observables, take into account the properties of the initial state and consider realistic systems that may not even be disordered, but in which strong interactions lead to stochastic behavior.
In the cases of the particular quantities studied here, quantum chaos guarantees that in realistic systems, S Sh ∼ ln(aD) and SvN ∼ ln(aDA), with a being a constant. The fact that these entropies saturate to values that follow a volume law imply that they approach thermodynamic entropies as the system size increases [61,62]. For chaotic initial states, we also have W ini ∼ D −1 , which assures thermalization.

III. REALISTIC INTEGRABLE AND CHAOTIC MODELS
Full random matrices do not model realistic systems. The Hamiltonian matrices describing realistic systems are usually sparse and banded, due to the presence of few-body (often only two-body) interactions, and in various cases, they do not even have random elements. Attempts to bring random matrix theory closer to realistic systems started already with Wigner, with the introduction of the so-called Wigner banded random matrices [34]. In these matrices, only the elements within a bandwidth around the diagonal are nonzero. Other similar approaches include two-body random ensembles and power law banded random matrices [76][77][78][79].
To bring the discussion of Section II down to earth, we consider a one-dimensional spin-1/2 model, which is similar to the systems studied experimentally with cold atoms, ion traps and nuclear magnetic resonance platforms. It is described by the Hamiltonian: where: Above, = 1, S x,y,z n are the spin operators on site n; L is the total even number of sites in the chain; J is the coupling strength; and ∆ is the anisotropy parameter. The Hamiltonian contains nearest-neighbor (NN) couplings and next-nearest-neighbor (NNN) couplings if λ = 0. We always include a small defect ǫ1 = 0.1 (Zeeman splitting 0.1 larger than that of the other sites) on the first site to break parity and spin reversal symmetries. A defect of amplitude d may also be included in the middle of the chain, on site n = ⌊L/2⌋. We consider open boundary conditions. The total spin in the z-direction, S z , is conserved. We study the largest subspace, S z = 0, so D = L!/(L/2)! 2 . In all of the figures below, except for Figure 7, we consider L = 16, so D = 12870 and DA = 256. This explains the choices of D and DA made for the figures for full random matrices in Section II.
All of the parameters are taken as positive, and J = 1 sets the energy scale. When λ, d = 0 and ∆ = 0 (we fix ∆ = 0.48), the Hamiltonian is integrable and known as the XXZ model. [XXZ refers to models where the coupling strengths in the x and y directions are the same and different from that in the z direction. The isotropic version is known as XXX.] The inclusion of NNN couplings or of a defect away from the borders of the chain breaks the integrability of the XXZ Hamiltonian [21,80]. Level repulsion occurs for the parameters considered here: λ = 1 (d = 0) and d = 0.9 (λ = 0). We refer to these two cases as the NNN model and the defect model, respectively. In common with the XXZ Hamiltonian, the first is clean, and the second has only NN couplings.
For any value of the parameters in H (30), the density of states is Gaussian, as is typical of many-body quantum systems with two-body interactions. The Gaussian shape is a consequence of combinatorics and the central limit theorem. This is the first crucial difference with respect to full random matrices that is worth emphasizing. The distributions are shown in Figure 3a-c for the three considered spin models. The Gaussian energy distribution is symmetric when ǫ1, d, ∆, λ = 0, while defects, anisotropy and λ make it asymmetric. Level spacing distribution and level number variance are shown in the middle and bottom panels of Figure 3, clearly distinguishing the integrable model from the chaotic ones. The XXZ model presents additional symmetries when the anisotropy reaches the roots of unity, ∆ = cos(π/µ), where µ ≥ 2 is an integer. The largest amount of degeneracies occur for µ = 2, but significant amounts exist also for µ = 3 [55]. This explains why we chose ∆ = 0.48 instead of 1/2, even though the presence of the small defect at the border, ǫ1, also prevents too many degeneracies.
The initial impression is that there is not much difference between Figure 3e,f and Figure 1a. However, visibly, the spectrum of the GOE full random matrix is more rigid than those for the chaotic spin models, as evident by comparing Figure 1b with Figure 3h,i. This reinforces the importance of considering different signatures of quantum chaos when comparing models.

A. Delocalization and Entanglement Measures: Basis Dependence
The Gaussian density of states implies that there are more states in the middle of the spectrum than at the edges. The states away from the borders are then expected to be more delocalized. This energy dependence on the level of delocalization of |ψα is another crucial difference between realistic systems and full random matrices.
To compute S (α) Sh , we need to choose a basis representation. This choice is made according to the problem we are studying. If one is interested in localization in real space, the site basis (known in quantum information as computational basis) is the natural choice. It corresponds to states where on each site, the spin points either up or down in the z-direction. However, if the interest is in the level of chaoticity of the states, the mean-field basis is the most appropriate one. The mean-field basis is defined by the integrable (regular) part of the Hamiltonian. In the case of the NNN and defect models, a good choice for the mean-field basis is to use the eigenstates of the XXZ model. We also verify that the differences obtained when we consider the eigenstates of the XX model are not significant. For the XXZ model, we use the eigenstates of the XX model [44,45] as a way to analyze how the level of complexity of |ψα increases with the Ising interaction.
In Figure 4a-c, we show the Shannon entropy for the XXZ, defect and NNN models in the site and mean-field basis. Independent of the basis, larger fluctuations are observed for the integrable model [70][71][72], where the eigenstates are further from random vectors than in chaotic models. The values of S (α) Sh for the XXZ Hamiltonian are also smaller overall than those for the chaotic models.
FIG. 4: Top panels: Shannon entropy in the site basis (black symbols) and in the mean-field basis (red symbols) for the XXZ (a), defect (b) and NNN (c) models. Bottom panels: normalized entanglement entropy (black symbols) and normalized Shannon entropy in the mean-field basis (red symbols) for the same models as in (a,b,c). Parameters for all panels are the same as in Figure 3.
In the mean-field basis, one clearly sees that the eigenstates of the chaotic models close to the middle of the spectrum are very complex, with S (α) Sh approaching the full random matrix value S GOE Sh (6). As we move towards the edges of the spectrum, the states become more localized. This is why thermalization is not expected to take place close to the borders of the spectrum [68,[70][71][72]. We note that the somewhat large values of S (α) Sh obtained for the NNN model at low energies occurs only when λ is large, but they are small when λ ∼ 0.5 [21,45]. This may suggest that in that region of the spectrum and for λ → 1, the basis considered is not the best mean-field basis.
In the site basis, the energy dependence is not so obvious, and S (α) Sh even surpasses S GOE Sh . Values above S GOE Sh indicate that the chosen basis is not a good mean-field basis. For the models considered, the site basis is therefore not the correct one to evaluate the proximity of the eigenstates to chaotic vectors.
In Figure 4d-f, we compare the normalized entanglement entropy S Sh /S GOE Sh for the same three models. The behavior of the two entropies is quite similar, so either one could be used in the analysis of the complexity of the states. Each entropy has its own advantages and disadvantages. In contrast to S (α) vN , the basis needs to be carefully chosen when computing S (α) Sh , but the latter is computationally less expensive, since no partial trace needs to be performed.

B. Dynamics at Intermediate Times: Generic Behaviors
The dynamics of the spin model is initiated after its preparation in an eigenstate of a certain Hamiltonian Hini. The state is then left to evolve according to H (30), where H = Hini + νV and ν is the perturbation strength.
As discussed in Section II C 1, the behavior of the survival probability is determined by the LDOS, Pini,ini(E). As ν increases from zero, the LDOS broadens [18,19,21]. When the couplings (elements of νV ) become larger than the mean level spacing, the envelope of the LDOS acquires a Lorentzian shape, and its Fourier transform leads to the exponential decay of Wini(t). By further increasing ν, the LDOS eventually becomes Gaussian, causing the Gaussian decay of Wini(t) [17][18][19][20][21][22][23][24][25]. This is not the mere quadratic decay observed for very short times for any perturbation strength, which emerges from the expansion: but a true Gaussian decay holding for larger times (see also [16] and the references therein). The Gaussian LDOS reflects the Gaussian density of states and corresponds to the maximum spread in energy of the initial state. This limiting scenario is reached by strongly perturbing the many-body quantum system, and it does not matter whether or not there is level repulsion. Several examples of Gaussian LDOS and Gaussian decay of Wini(t) have been shown for the clean XXZ and NNN models in [17][18][19][20][21][22] and also for disordered models [23,24]. In Figure 5, we compare the integrable XXZ model with the chaotic defect model. As the initial state, |Ψ(0) = |φini , we consider the Néel state, | ↓↑↓↑↓↑↓↑ . . . , where the z-polarization of the spins alternates from one site to the other. This state is often prepared in experiments with cold atoms. The evolution under both models is equivalent to having had a very strong perturbation, where ∆ is changed from ∞ to 0.48.
The envelopes of the LDOS in Figure 5a,c agree very well with a Gaussian of width: centered at: with d being zero for the XXZ model and 0.9 for the defect model. For both models, the chosen initial state is further from the middle of the spectrum and, therefore, less filled than for the NNN model, where: Comparing Figure 5a and 5c, we also see that the LDOS is better filled for the defect model than for the XXZ model, since the first is chaotic and has Eini slightly closer to the middle of the spectrum. We also note that the width of the initial state is sub-extensive, σini ∝ √ L, so the initial state is narrow in energy. This is a necessary condition for thermalization in real systems [64,65].
In Figure 5b,d, we show the evolution of the survival probability. The numerical results are very close to the analytical Gaussian decay, Wini(t) = exp(−σ 2 ini t 2 ), until the curves cross the saturation line 1/P Rini (25). In this limit of strong perturbation and initial states with energy Eini away from the edges of the spectrum, the behavior of the survival probability is general and equivalent for both integrable and chaotic many-body quantum models. For the intermediate time scales considered in Figure 5, level repulsion does not play any important role.
In Figure 6, we compare the evolution of the Shannon entropy (top panels) and the entanglement entropy (bottom panels) for the Néel state under the XXZ (a,d), defect (b,e) and NNN (c,f) models. At very short times, t ≪ σ −1 ini , the entropies increase nearly quadratically [18]. Later, as in full random matrices, the entropy growth becomes linear in time and remains as such until close to saturation. This linear behavior is expected for initial states far from the edges of the spectrum and seems independent of the presence or absence of level repulsion. What one observes, however, is a dependence on the energy of the initial state [20] and on the connectivity of the Hamiltonian [18]. For the NNN model, σini is still the same as in Equation (34), but since Eini from Equation (36) is closer to the middle of the spectrum, where there are more eigenstates and they are more delocalized than those from the XXZ and defect model (cf. Figure 4), the slope of the linear entropy increase is larger in Figure 6c,f than in Figure 6a,b,d,e.
The saturation values of the entropies for the GOE full random matrices are indicated in Figure 6 with horizontal dashed lines. For the XXZ model, the entropies saturate to values smaller than S GOE Sh,vN . For the defect model, S Sh surpasses S GOE

Sh
. For the NNN model, both S Sh and SvN are above the results for full random matrices. This is in contrast with previous studies for the evolutions of the Shannon entropy under the NNN model that start with mean-field basis vectors and analyze the spreading over other mean-field basis [44,45]. In these cases, the saturation values are not larger than S GOE Sh,vN , but very close to it. Another difference between the results for site basis initial states and mean-field initial states is the agreement with Equation (19), which occurs for the latter [44,45], but not for the former [18]. The approximation in Equation (19) seems appropriate for initial states that are directly coupled with many other basis vectors, that is their connectivity is ∝ D. This is not the case of site basis vectors. The Néel state, for instance, is directly coupled with only L − 1 states.

C. Dynamics at Long Times
At long times, the evolution of Wini(t) slows down. The decay of the survival probability becomes necessarily a power law, It is the tails of the LDOS that determine the value of the power law exponent γ. Good filling implies that, despite the  discreteness of the spectrum, the LDOS can be treated as a nearly continuous function. The behavior of Wini(t) at long times can then be obtained from the Fourier transform of the Gaussian Pini,ini(E), but taking into account also the unavoidable energy bounds E low and Eup of the spectrum [25], The above power law decay with exponent γ = 2 can be seen in Figure 7a and also in other chaotic initial states studied in [25]. The exponent is smaller than that for full random matrices, where γ = 3 (see Equation (18)), but larger than what we obtain for systems undergoing localization, where the LDOS is very sparse and γ < 1 [23][24][25]. In many-body quantum systems, exponents γ ≥ 2 indicate that the LDOS is very well filled and that the components C α ini are close to random numbers from a Gaussian distribution. In this scenario, we should expect thermalization to occur [25].
In Figure 7b, we show the evolution of the Shannon entropy. The symbols represent a fitted linear curve. The log-log plot makes it very clear that the linear growth holds only for times that are not too short. The dot-dashed line corresponds to Equation (20). It is parallel to the numerical data nearly up to the point where the linear behavior develops. For the curve of Equation (20) to actually coincide with the numerical data, we would need to substitute Npc by a small value, ∼ 20. This happens when the initial state is a site basis vector, not when it is a mean-field basis vector.
In analogy with the onset of the power law decay for Wini(t), we might wonder whether at long times, Figure 7b actually detects a behavior different from the linear growth. Right before saturation, between t ∼ 1 and t ∼ 3, a small oscillation is visible. This is also noticeable in Figure 6, especially in Panel (c). This time interval coincides with that for the survival collapse observed for the survival probability. The collapse [25,81,82] corresponds to an interference between the Gaussian decay and the emergence of the power law decay, which pushes Wini(t) very much below the saturation line (see Figure 7a and also Figure 5b,d). Thus, even though the small oscillation requires further studies, it may be an indication of an interference between two different behaviors, as in the case of the survival collapse.

IV. DISCUSSION
We compared the static and dynamical properties of GOE full random matrices and finite isolated many-body quantum systems described by one-dimensional spin-1/2 models with two-body interactions. This comparison is useful, because for full random matrices, analytical expressions can be derived and then used as references and bounds for the analysis of realistic models. Our main findings are itemized below: 1. The results for the von Neumann entanglement entropy SvN , which is a concept employed in quantum information science, and for the Shannon information entropy S Sh , which is generally used as a measurement of the degree of delocalization of quantum states, were very similar. Thus, either one can be used to measure the level of complexity of the eigenstates. The advantage of the Shannon entropy is that it is computationally less expensive than the entanglement entropy. The disadvantage is that it is strongly dependent on the basis chosen.
2. For full random matrices, all eigenstates are pseudo-random vectors and therefore lead to the same values of S (α) Sh , but the results for realistic systems depend on the region of the spectrum and on the basis selected.
3. An analytical expression was given for full random matrices for the time evolution of both entropies. It agrees extremely well with the numerical results. For the spin systems; this expression gives an upper bound for S Sh (t) and SvN (t).
4. At short times, S Sh (t) and SvN (t) show a nearly quadratic behavior. It is only at longer times that the linear increase, S Sh,vN (t) ∝ t, develops. These two behaviors seem to be independent of the presence or absence of level repulsion.
We also reviewed some of our previous studies and discussions. They include: 1. In realistic chaotic models, the spectrum is not as rigid as that of full random matrices. When comparing different chaotic models, it is appropriate to compare different signatures of chaos, such as those that detect short-range and also long-range correlations.
3. At long times, the decay of the survival probability becomes a power law, Wini(t) ∝ t −γ , with γ = 3 for full random matrices and γ max = 2 for the spin systems. The emergence of a power law decay at long times should have interesting consequences for problems associated with quantum information science and foundations of quantum mechanics. One should expect, for example, that external actions on the system, such as measurements, performed at long times may change the power law decay and recover Gaussian or exponential decays. This idea was explored in [83] for a one-body system interacting with an environment. It would be worth extending it also to many-body quantum systems. 4. Equilibration and thermalization are trivially reached under full random matrices. In realistic models, the absence of degeneracies and the presence of chaotic states in the energy window sampled by the initial state are both key elements for achieving thermal equilibrium.

V. MATERIALS AND METHODS
The numerical methods used include exact diagonalization and Expokit [84,85]. Exact diagonalization was used for Hamiltonian matrices with D < 20,000, and Expokit was employed for the dynamics of systems with larger D.