Absence of dynamical localization in interacting driven systems

Using a numerically exact method we study the stability of dynamical localization to the addition of interactions in a periodically driven isolated quantum system which conserves only the total number of particles. We find that while even infinitesimally small interactions destroy dynamical localization, for weak interactions density transport is significantly suppressed and is asymptotically diffusive, with a diffusion coefficient proportional to the interaction strength. For systems tuned away from the dynamical localization point, even slightly, transport is dramatically enhanced and within the largest accessible systems sizes a diffusive regime is only pronounced for sufficiently small detunings.


Introduction
The last decade has witnessed an increasing interest in isolated out-of-equilibrium (OOE) quantum systems both in theory and experiment. In particular, periodically-driven, or Floquet many-body systems have become one of the most active areas in condensed matter theory and many-body physics, due to the exciting possibilities of finding novel OOE phases of matter, without an analog in static systems.
Driven interacting systems typically heat up, until they approach a nonequilibrium steadystate (NESS) "locally identical" to a featureless state of maximal Gibbs entropy. That is, local observables cannot differentiate between the two states [1,2,3].
This "local identity" between the NESS and the maximal Gibbs entropy state is broken when conservation laws exist, for example in noninteracting systems [4,5]. The systems still approach a NESS, but this NESS is "locally identical" to a state which maximizes the Gibbs entropy under all the constraints imposed by the conserved quantities. A large body of work has concentrated on "Floquet engineering" of such noninteracting systems, namely using periodic driving as a tool to create effective time-independent Hamiltonians of desirable but difficult to implement static Hamiltonians [6,7,8].
For interacting systems it is possible to avoid the featureless NESS by the addition of a quenched disorder. In the absence of driving, these systems are many-body localized (MBL) [9,10,11,12,13,14,15,16], and for sufficiently high driving frequencies they do not heat up indefinitely [15,16]. Such systems were also shown to host nontrivial eigenstate order, corresponding to a new OOE phase of matter, colloquially dubbed a "time-crystal" [17] and subsequently studied in both experimental and theoretical studies [18,19,20,21].
In the present work we explore a different potential avenue to suppress the indefinite heating of interacting driven systems by starting with a noninteracting translationally invariant system which is dynamically localized [22,23,24,25,26] (see also the closely related study of a disordered system [27], which appeared while our work was in preparation). The special form of the driving term in dynamically localized systems leads to a frozen stroboscopic dynamics and therefore to disorder-free spatial localization. Whether adding interactions to such disorder-free localized systems results in heating, and what is the nature of the transport are the main questions we consider in our work. We would like to stress that the mechanism behind this localization, while of quantum nature, is not due to Anderson localization, as happens in related disorder-free kicked rotor systems which unfortunately are also dubbed "dynamically localized" [28,29,30]. The interacting versions of such systems were recently studied, but are not of direct relevance to our work [31,32].
The structure of this work is as follows: We present the model in Section 2, briefly survey dynamical localization adapted to the many-body context in Section 3 and in Section 4 present our main results showing that interactions destroy dynamical localization and cause diffusive dynamics with a diffusion coefficient proportional to the interaction. We end with our conclusions, a brief discussion of the implications of our results for the locality of Floquet effective Hamiltonians and an outlook.

Model
In this work we consider a driven Hamiltonian of the form, where the static Hamiltonian is given bŷ Here L is the length of the lattice,ĉ † m creates a spinless fermion at site m,n m is the number operator, J is the hopping and U is the interaction strength, and the driving protocol is where T = 2π/ω is the period of the driving, ω is the driving frequency and A is the amplitude of the drive. This model corresponds to a constant potential gradient with a strength proportional to A, the sign of which is flipped every half period.

Dynamical localization
For ∆ = 0, the model (1) is noninteracting and exactly solvable [22,23]. For specific ratios of A/ω it exhibits dynamical localization, namely, an initially localized density excitation does not spread over time. For convenience of the reader in this Section we adapt the derivation of dynamical localization of Ref. [23] to the many-body setting. For this purpose it is instructive to perform a time-dependent unitary transformation which eliminates the drive. This can be always achieved, since setting |ψ = V (t) |φ and using the Schrödinger equation, i∂ t |ψ = H |ψ it is easy to show that |φ satisfies i∂ t |φ =H (t) |φ , with the transformed Hamiltonian, Therefore if the original Hamiltonian is of the following form, the second term can be eliminated by choosing, where The transformed Hamiltonian is then given by, While the original driving term was eliminated, it appears that nothing was gained since the transformed Hamiltonian is still time-dependent. For a spatially uniform and temporally periodic force on the system, namely,Ĥ further progress can be made. For such a drive the destruction and creation operators transform as,â and the number operators are invariant under this transformation,n m =ĉ † mĉm =â † mâm . Therefore a generic Hamiltonian of the form, will be transformed to,H with time-dependence appearing only in the hopping term and with an arbitrary function g({n i }) depending only on the number operators. For simplicity we will now restrict the discussion to a single band by setting h nm = J 2 (δ n,m−1 + δ n,m+1 ). The total current is given byĴ which can be diagonalized toĴ wheref k = n exp (ikn)â n / √ L. For g ({n i }) = 0, namely in the absence of interactions and external potentials, the Hamiltonian can be diagonalized simultaneously withĴ, and therefore also the transformed Floquet operator, can be computed. Although the current is not conserved its expectation can be readily calculated, Dynamical localization occurs when the total number of particles transported in one period vanishes, that is, when the integral of the current over one period yields zero for all k, In this case the Floquet operator (14) reduces to the identity and its spectrum collapses such that all its eigenvalues are degenerate, making the stroboscopic dynamics trivial. Concretely, for the driving we use in this work (3), A (t) = A |t| (for |t| < T /2), and we have, which vanishes for e iAT /2 = 1, or A/ω = 2n n ∈ Z. When an external potential or interactions are present the Hamiltonian cannot be diagonalized in the basis ofĴ and the above derivation does not apply since quasimomentum is not conserved. Notwithstanding, one might still expect some residual suppression of transport, since using the periodicity of exp[iA (t)] = n C n exp [iωnt] one can cast the generic Hamiltonian (11) into the form which includes a static Hamiltonian with a hopping term renormalized by C 0 . This description is however oversimplified, due to the presence of a nontrivial driving term, which is coupled to the hopping. In the following Section we will explore what happens to transport in an interacting model, where the first term is set to vanish, namely C 0 = 0.

Results
We study density transport in a driven and interacting system. For this purpose we calculate the density-density correlation function, which describes the spreading of density excitation in the system. Here Z is the dimension of the relevant Hilbert space, and the operators are written in the Heisenberg picture with respect to the driven Hamiltonian. Intuitively this describes the spreading of a density excitation in the presence of driving. Here, we fix the number of fermions to N = (L − 1) /2, and use odd sizes L, such that an excitation at the center is at the same distance from the left and right (open) boundaries. To compute the excitation profile we stroboscopically evolve two initial states over n periods, |ψ R (n) =Û (nT ) n 0 − 1 2 |ψ 0 and |ψ L (n) =Û (nT ) |ψ 0 , using a numerically exact method. The method is based on projecting the time evolution operator to an orthonormal basis of a truncated Krylov space spanned by the initial state, which is updated at each time step (for a pedagogical review of the method, see Sec. 5.1.2 of Ref. [12]). We then calculate the matrix element C x (nT ) = ψ L (n) n x − 1 2 ψ R (n) . For the driving protocol we use in this work (3) the propagatorÛ (nT ) over n periods is simplified whereĤ ± corresponds to the driven Hamiltonian (1) in the first and second halves of the period and T is the time ordering operator. The initial state |ψ 0 is sampled randomly from the Haar measure, which due to Lévy's Lemma approximates the trace in (19) with an error which is exponentially small in the system size [33,12]. This computational scheme allows us to reach systems sizes as large as L = 31 sites, which corresponds to a Hilbert space dimension of 300, 540, 195 (see Ref. [12] for more details on the numerical procedure).   (24)). The colored lines indicate system size L = 31 and the grey lines systems size L = 27, simulated with the same parameters. The interaction strength is fixed at ∆ = 0.5.
To characterize the transport we calculate the mean-square displacement (MSD) of the density excitation, and the time-dependent diffusion coefficient [34,35,12,36], We note in passing that within the linear response theory, namely when the drive is taken to be a small perturbation to the equilibrium state, this diffusion coefficient computed in the limit of t → ∞ corresponds to the value calculated from the Kubo formula. In this work we however do not consider the linear response regime; our driving amplitudes are large and take the state of the system very far from equilibrium. Therefore there is no reason to expect that the diffusion coefficient calculated from (22) and the Kubo formula coincide. In fact it is not even clear a priori if transport in such case will be diffusive. The dynamical localization of noninteracting systems, which was discussed in the previous Section, serves as an example when the diffusion coefficient (depending on the drive amplitude) is either zero or infinity and there is either no transport or the transport is ballistic. We now proceed to examine the stability of dynamical localization to the addition of interactions. To this end we set the system to the noninteracting dynamical localization point by setting the amplitude of the drive to, Throughout this work we fix the period of the drive to be T = 5, such that the frequency, ω = 2π/T ≈ 1.256 is smaller then the one-particle and many-body bandwidths, to allow effective energy distribution in the system using single and many-particle rearrangements.
While the system is in the noninteracting dynamical localization point, the left panel of Fig. 1 shows sub-ballistic spreading of the density excitation, indicating that dynamical localization is unstable under the addition of interactions. The right panel shows fixed time cuts through the excitation profile. Interestingly, even after relatively long times these profiles do not approach a Gaussian form, suggesting that the density-density correlation function C x (t) is not well described by a diffusion equation. We also added an isoline at C x (t) = 0 (the background of C x (t) is slightly negative due to the trivial correlation imposed by particle number conservation), which shows a peculiar feature at short times. We attribute this feature to short time transient dynamics which roughly corresponds to the location of the peak in the time dependent diffusion coefficient D(t) shown in Fig. 2. The asymptotic transport is only reached at later times.
To characterize the nature of transport, and to estimate the finite size effects in the system we calculate the MSD (21) for different interaction strengths and several system sizes (not all shown for readability). From the left panel of Fig. 2 it is clear that while dynamical localization is destroyed by the addition of interactions, for weak interactions transport is still significantly suppressed and signalled by a slow growth of the MSD, which allows us to reach very long times before our results are affected by finite size effects. This time becomes shorter for stronger interactions, as could be inferred from the deviation in MSD between the various system sizes. After a short time of fast transient transport the MSD appears to enter a diffusive regime, as could be seen from a linear growth of the MSD. To see it more clearly we compute the time-dependent diffusion coefficient (22) in the middle panel of Fig. 2. It shows a clear saturation to finite plateaus, which become longer for larger system sizes (while the height of the plateaus remains the same) and indicate that the departure from the plateau at later times is a finite size effect. The asymptotic diffusion coefficient is extracted from the height of the plateaus for different interaction strengths and is presented in the right panel of Fig. 2. While the overall dependence on the interaction strength is nonlinear, for sufficiently weak interaction strengths we find D ∼ ∆.
We now fix the interaction strength to ∆ = 0.5, such that we are able to reach a pronounced diffusive regime within our simulation times for A/A 0 = 1, and detune away from the dynamical localization point by slightly changing the amplitude of the drive. We characterize the strength of the detuning from the dynamical localization point by   is much faster in this regime even for a few percent detuning from the dynamical localization point, requiring very large system sizes to capture the asymptotic transport. Both the MSD and the time-dependent diffusion coefficient show diffusive behavior for small detuning δ = ±0.03, while for larger (in magnitude) detuning δ = −0.06 it appears that our results do not capture the asymptotic transport even for the largest system size L = 31.
Finally we examine the broadening of the Floquet spectrum upon addition of interactions. As explained in Sec. 3, the spectrum of the operatorÛ (T, 0) in (20) collapses to a single point at the noninteracting (∆ = 0) dynamically localized point δ = 0, resulting in a trivial stroboscopic dynamics. A perturbative reasoning suggests that a finite interaction will lift this massive degeneracy and broaden the quasienergy spectrum. We have confirmed this expectation by evaluating the quasienergies α , defined asÛ (T, 0) |α = exp (−i α T ) |α , by exactly diagonalizingÛ (T, 0) up to systems sizes L = 14. The width of spectrum was calculated from the standard deviation σ ω of the α (after the quasienergies were shifted to the center of the band to avoid spurious widening due to wrapping around the Floquet Brillouin zone). The results are shown in Fig. 4, where we see that the spectrum broadens upon addition of interactions, however the location of the minimal (but still finite) width remains at zero detuning, δ = 0 (see (24)). This is consistent with our observation that slowest transport occurs at this point (Fig. 3).

Discussion
In this work we have examined the stability of dynamical localization to the addition of interparticle interactions. We have found that while dynamical localization is destroyed by interactions, transport is significantly suppressed when the interacting system is at the noninteracting dynamically localized point (see right panel of Fig. 4 where the minimum of the width occurs at zero detuning, δ = 0 for all considered interaction strengths ∆). The slow transport allows us to find clear evidence of diffusive transport with a diffusion constant proportional to the interaction strength, indicating that scattering is dominated by interparticle collisions. We have also studied the dependence of the width of the quasienergy spectrum on the interaction strength. While the spectrum is completely degenerate for zero interactions, leading to trivial stroboscopic dynamics, for finite interaction the degeneracy is lifted and the width appears to be proportional to ∆ before it saturates to its maximum value fixed by the driving frequency. From Fig. 4 one can see that for a fixed interaction strength the width increases with system size. This leads us to speculate that in the thermodynamic limit an arbitrarily weak interaction causes the spectrum to occupy the full allowed quasienergy range. While it is tempting to conclude that the width of the spectrum corresponds to the rate with which the dynamics unfolds, this is not true in general. For example, a finite width occurs for Floquet-MBL systems while the dynamics is frozen [15,16]. We therefore can only conclude that when the spectrum broadens dynamical localization may be destroyed, but the precise nature of the underlying dynamics has to be assessed by other means (as we have done in Figs. 2 and 3).
One interesting question on which our work might have some bearing is the locality of the effective Hamiltonian,Ĥ F . This quantity is defined fromÛ (T, 0) in (20) asÛ (T, 0) = exp −iĤ F T and is the generator of the stroboscopic time evolution. It has been speculated that for generic Floquet systems which are not MBL the effective HamiltonianĤ F is nonlocal [37,1]. Since it is generally believed that spreading of correlations in systems with nonlocal static Hamiltonians will be superballistic, [38,39,40,41] one may wonder whether the diffusive transport we observe in this work points toward a local effective Hamiltonian. We would like to argue that due to the ambiguity in the definition of the effective Hamiltonian the question is not well posed. The ambiguity follows from the fact thatÛ (T, 0) and therefore all stroboscopically measured observables are invariant under the transformationĤ F =Ĥ F + ω α n α |α α|, where |α are the eigenstates ofÛ (T, 0) and n α ∈ Z. Since in general the projectors |α α| are nonlocal objects this renders the locality ofĤ F a property which cannot be inferred from physical observables. While we believe that the discussion of locality ofĤ F is irrelevant as far as physical observables are concerned, it remains an interesting question whether the ultimate fate of driven interacting systems is encoded in the locality of the corresponding Floquet operator,Û (T, 0), which is a physical quantity and whose locality can be assessed, for example using the operator entanglement entropy [42]. We leave this question for future work.