Out-of-equilibrium phase transitions induced by Floquet resonances in a periodically quench-driven XY spin chain

We consider the dynamics of an XY spin chain subjected to an external transverse field which is periodically quenched between two values. By deriving an exact expression of the Floquet Hamiltonian for this out-of-equilibrium protocol with arbitrary driving frequencies, we show how, after an unfolding of the Floquet spectrum, the parameter space of the system is characterized by alternations between local and non-local regions, corresponding respectively to the absence and presence of Floquet resonances. The boundary lines between regions are obtained analytically from avoided crossings in the Floquet quasi-energies and are observable as phase transitions in the synchronized state. The transient behaviour of dynamical averages of local observables similarly undergoes a transition, showing either a rapid convergence towards the synchronized state in the local regime, or a rather slow one exhibiting persistent oscillations in the non-local regime, where explicit decay coefficients are presented.


Introduction
The study of out-of-equilibrium systems has recently moved to the forefront of research in many-body physics, motivated largely by the desire to create interesting states of matter with properties beyond those achievable in equilibrium. In periodically-driven systems, which exist in nature but can also be created in the laboratory, driving can be used as a probe or to generate non-trivial phases [1][2][3]. Floquet theory [4] is then used to describe these situations. Whereas all relevant properties of the periodic dynamics are encoded in the Floquet Hamiltonian, finding exact solutions remains a severe challenge, and these are usually limited to few-body systems. This poses severe limitations to the application of Floquet theory when a large number of components are involved, as in condensed matter physics. However, in freefermion models and related spin chain systems, the Floquet Hamiltonian sometimes takes a simpler form due to their underlying algebraic structure [5]. These models can then shed light on periodically-driven systems, while describing experimentally-relevant situations [6][7][8], and being relevant in various fields of physics [1,9,10]. Previous works have explored some of these scenarios [11,12], and numerical studies have also been performed for non-integrable systems [13][14][15][16], where many-body resonances are of major interest [17][18][19][20][21][22][23].
Yet, exact expressions for Floquet Hamiltonians at arbitrary driving frequencies and explicit connections between the dynamics and the quasi-energy spectrum, from which much could be learned on equilibration and synchronization in driven systems, are lacking in most cases. In this paper, we consider an antiferromagnetic spin-1/2 XY chain subject to periodic quenches of a transverse magnetic field. Although the equilibrium properties of this model are well-known [24][25][26], and despite studies of its out-of-equilibrium behaviour in certain scenarios (e.g., [27][28][29][30][31][32][33]), an exact description of its Floquet dynamics valid in the entire parameter space is currently missing. The XY and related models have been fundamental for understanding quantum phase transitions [32,34,35] and quantum quenches [36,37], and can certainly play a similar role for the physics of periodic driving. Its Floquet Hamiltonian is obtained exactly, and we investigate its local/non-local nature for different driving parameters. A link is estab-lished between the locality of the Floquet Hamiltonian and the presence of avoided crossings between its eigenstates, which can be interpreted as Floquet resonances. These resonances are generally not observable in the thermodynamic limit, but the free-fermionic nature of the model makes them visible in the single-particle spectrum, leading to a slow convergence of observables towards a synchronized state (i.e., the steady state of the stroboscopic time evolution [11,12]) with persistent oscillations in the non-local regime. These results are expected to be generic for any system mappable to a free fermion model. This study also allows to naturally connect Floquet resonances [19] with the appearance of long-range correlations as presented in Refs. [38,39]. The phenomena analyzed are also expected to be observable in string order parameters, characterizing similar transitions in the Floquet ground state of driven Ising models [40,41]. This paper is organised as follows. In Sec. 2, we start by providing a brief overview of Floquet theory. This is subsequently applied to a periodically driven XY spin chain, whose corresponding Floquet Hamiltonian is explicitly derived (see Secs. 3 and 4) and its effective interactions thoroughly discussed (see Sec. 5). In Sec. 6, a detailed account is presented on the correspondence between the locality, or lack thereof, of the Floquet Hamiltonian and the dynamics of observables. The last section is reserved for conclusions and possible future avenues of research.

Floquet Theory
Analogously to the Bloch theorem for crystalline solids, given a time periodic Hamiltonian H(t) =Ĥ(t + T ) with T = 2π/ω 0 , one can construct a complete set of solutions of the Schrödinger equation, the so-called Floquet states [1,4,11], which are periodic in time up to a phase factor, |Ψ α (t) = e −i [T ] α t |Φ α (t) .
Here the states |Φ α (t) are called Floquet modes and are periodic |Φ α (t + T ) = |Φ α (t) , with α representing a set of quantum numbers. The real quantities [T ] α are known as the Floquet quasienergies. Due to the periodicity, these quasienergies can be translated to the first Brillouin zone [−π/T, π/T ], just like the quasimomenta of Bloch waves. To derive expressions for the Floquet modes and quasienergies, we first notice that, by definition, a Floquet state obeys the Schrödinger equation with evolution operatorÛ (t 2 , t 1 ) = T t exp −i t 2 t 1Ĥ (t)dt . One can easily show that, for timeperiodic Hamiltonians, this evolution operator obeys that U (t 0 + nT, t 0 ) = U (t 0 + T, t 0 ) n . In terms of the Floquet modes the preceding evolution equation becomeŝ Eq. (3) tells us, on the one hand, how Floquet modes evolve and, on the other, that for stroboscopic times t = t 0 +nT , the Floquet modes are eigenfunctions of the evolution operator Moreover, since a unitary operator can always be written as the exponential of a Hermitian operator, this allows us to introduce the so-called Floquet HamiltonianĤ The Floquet Hamiltonian is time-independent, but depends both on the driving period, as well as on the arbitrary choice of the initial time. Since Floquet Hamiltonians for different initial times are related by a gauge transformation, t 0 is referred to as the "Floquet gauge". Thus, the properties of the system at stroboscopic times t 0 + nT can be determined by studying the static Floquet Hamiltonian. Finally, it can be shown that the evolution of a Floquet mode Φ α (t) (where the superscript T now explicitly shows its dependence on the period) from t 0 to t = t 0 +nT +δt is given by 3 Driving the XY spin chain In the following, Floquet theory will be applied to the XY spin chain, where we pay special attention to the structure of the Floquet Hamiltonian. The instantaneous XY Hamiltonian is given by with superexchange parameter J > 0, anisotropy γ ∈ [0, 1], and h(t) a time-periodic external magnetic field with non-negative average value. The spin-1/2 operators S α j , α ∈ {x, y, z}, are taken to obey periodic boundary conditions. Performing a Jordan-Wigner transformation and a Fourier transform to spinless fermion operators {c † k , c k }, and taking care of the boundary conditions [24] (see Appendix A for details of the derivation), this instantaneous Hamiltonian can be mapped to which (assuming even filling in subsectors k ∪ −k) is a direct sum of operators acting on twodimensional subspaces spanned by {|0 , c † k c † −k |0 }. Observing the system stroboscopically at times t 0 + nT for n ∈ N, the evolution operator can be written as The so-called Floquet Hamiltonian H can be obtained by diagonalizing the evolution operator for a single driving cycle. In each two-dimensional momentum subspace, this operator can be repesented as with Pauli matrices σ α k . Solving the eigenvalue problem provides the Floquet modes and the quasienergies. After some algebra one obtains the Floquet modes with normalization factor and corresponding eigenvalues λ σ = e −iσ T (k)T with σ = ±1. Here T (k) are the quasi-energies given by This allows us to write the following pseudospin representation of the Floquet Hamiltonian (up to an overall constant) with It is, however, more illuminating to write the resulting Floquet Hamiltonian in terms of the Jordan-Wigner fermions in real space, leading to the following result: where we can see that the driving may induce non-local interactions A and B . So far the results are general, as they only rely on the algebraic structure of the system. The problem is now that obtaining explicitly exact, analytical expressions for Euler parameters u k,α is not a trivial matter. One protocol that allows for such solutions consists on periodically quenching between two values of the time dependent external field h(t), In this case, one can use the algebra of the Pauli matrices to determine the effective parameters u k,α . In the following, we take t 0 = 0 for simplicity. For this Floquet gauge, we obtain: where E i (k) are the eigenvalues of the XY Hamiltonian with field h i , ∆ k = θ 2,k − θ 1,k , and θ i,k is the corresponding Bogoliubov angle for that same external field. This is defined as Detailed expressions of the parameters u k,α for arbitrary Floquet gauge choices can be found in Appendix B.
Before analyzing the influence of these possibly non-local interactions, let us discuss their origin by considering the Floquet-Magnus expansion [13,42]. In the infinite driving frequency limit T → 0, the Floquet Hamiltonian reduces to the time-averaged Hamiltonian H av . Eigenstates of H av whose single-particle energy difference (within a single momentum sector) equals an integer multiple of the driving frequency 2π/T are said to be quasi-degenerate or resonant and, in the same way that a small perturbation on a stationary Hamiltonian can strongly couple degenerate eigenstates of the unperturbed Hamiltonian, deviations of e −iH [t 0 ,T ] F T from e −iHavT can strongly couple quasi-degenerate eigenstates [17][18][19][20][21][22][23]. These so-called Floquet resonances give rise to avoided crossings that result into the Floquet eigenstates being linear combinations of two such (near-) resonant states. As this mechanism cannot be represented by local operators, the Floquet Hamiltonian contains long-range interactions if resonances are present (note that other mechanisms can also lead to non-local terms in certain circumstances [43,44]). In such situations, the Magnus expansion is not expected to converge [13,[45][46][47][48], hence exact solutions are important for the description of these resonances.
In Fig. 1 we show the quasi-energy spectrum as a function of T , where the intervals within which avoided crossings occur can be clearly distinguished. Since quasi-energies are only defined up to shifts 2πm/T, m ∈ N, in Fig. 1 the quasi-energies are restricted to the first Brillouin Zone [−π/T, π/T ], which is known as the 'folding' of the spectrum [13].

Unfolding of the quasi-energy spectrum
By construction the Floquet Hamiltonian possesses the same eigenstates as the evolution operator, and its eigenvalues ± T (k) = ∓arccos(u k,0 )/T follow from the logarithm of the corresponding eigenvalues of U k (t 0 + T, t 0 ). As such, these quasi-energies are defined up to shifts 2πm/T , with m ∈ N, whose branch choice is seen as a discrete gauge freedom [49]. While the gauge choice for the quasi-energies is immaterial for the time evolution of the system, leading to a Floquet Hamiltonian that is similarly only defined up to a gauge freedom, its effects on the latter are, however, much more involved. As will be argued in the following, this folding might introduce unphysical interactions in the local Hamiltonian, and the spectrum needs to be explicitly unfolded in order to highlight the effects of Floquet resonances. Although there exists no general way of doing so [19], in the current model we have the advantage that, due to the free-fermionic structure, the Floquet Hamiltonian factorises into two-dimensional subspaces and has a symmetric spectrum. Therefore, the single-particle quasi-energy spectrum can be explicitly unfolded by either adding or subtracting the corresponding multiple of 2π/T whenever a resonance condition is met (see Fig. 2 for a graphical representation of the unfolding).

Avoided Crossings
Here, we present two equivalent ways of explicitly unfolding the Floquet spectrum. The first one takes advantage of the symmetry of the spectrum of the time-averaged and Floquet Hamiltonians, which implies that an avoided crossing takes place whenever the value of a quasienergy approaches the edge of the Brillouin zone. This is precisely the mechanism that keeps the spectrum found in (12) inside the interval [−π/T, π/T ]. For a given value of k and as the period increases, the quasi-energy T (k) will have a first resonance at T = π/E av (k). For periods above this value, we can unfold the quasienergy out of the first Brillouin zone by shifting it by 2π/T . The spectrum will now be confined inside the second Brillouin zone, and in order to further unfold it we need a shift of −2 × 2π/T for periods larger than 2π/E av (k), where a second resonance takes place. By repeating this process iteratively, the full unfolding can be written as with Θ(x) the Heaviside function. An alternative way of unfolding the spectrum is by making the connection with the γ → 0 limit of driven XX models. In this limit, the driving becomes trivial since the XX Hamiltonians commute at all values of h and J, leading to a Floquet operator whose eigenstates are exactly those of the time-averaged Hamiltonian H av . The spectrum can now be unfolded by demanding that the Floquet Hamiltonian reduces to H av at all values of the driving period T if we take this limit. Namely, for γ = 0 the evolution operator in the pseudo-spin representation is with E av (k) = |J cos(k) − h av | and h av = αh 1 + (1 − α)h 2 , while the Floquet Hamiltonian reduces to It is clear that this Hamiltonian shares the eigenstates of the time-averaged Hamiltonian, where the eigenvalues have acquired an appropriate shift compared to E av because of the term arccos [cos(E av (k)T )] /T . The shifts introduced by the arccos function can be explicitly undone by defining the unfolded quasienergies as In the limit γ → 0, at each value of the driving period the Floquet Hamiltonian now reduces to H av , which has a guaranteed unfolded spectrum. It is easy to show that Eqs. (17) and (20)  = H av at all driving frequencies, as returned by the proposed unfolding. Another example is the static case, which can be achieved by either h 1 = h 2 , α = 0 or α = 1, for which the Floquet Hamiltonian reduces to the static one.

Induced Interactions in the Floquet Hamiltonian
Now that the spectrum has been unfolded, we can return to the couplings A and B in the Floquet Hamiltonian given by Eq. (15). In particular, we are interested in whether or not the driving introduces significant couplings at ranges larger than the nearest-neighbour interactions present in the instantaneous Hamiltonian. With our choice of driving protocol, these couplings are given by To visualise the non-locality of the the Floquet Hamiltonian in the parameter space we introduce the measure This quantity measures the average strength of interactions up to range L, excluding those of range 1 as well as the effective on-site potential. Therefore, it allows us to probe only interaction ranges which were not present in the static Hamiltonians and, as such, the nonlocal long-range behaviour of the Floquet Hamiltonian. As mentioned before, as T goes to zero, one retrieves the time-averaged Hamiltonian, and ρ = 0. When T increases, one would naively expect an enhancement of longer range interactions, reflected in an increase of ρ.  Fig  3, where it can be seen that inside the local regime (highlighted by the pale blue background) only the local couplings are non-neglibible, whereas in the non-local regime all the couplings play a role in the Floquet Hamiltonian. These local to non-local transitions are induced by the appearance of avoided crossings, as indicated in Fig. 1, since the latter are the origin of the long-range interactions. As such, the boundary lines separating these two kinds of regions, drawn as red lines in Fig. 3, can be analytically derived by looking for these resonance conditions.
To do so, we start by noticing that, due to the symmetry of the spectrum of the timeaveraged Hamiltonian, resonances begin to occur whenever the quasienergies get close to the first Brillouin zone (as can also be observed from Fig. 1). If we were to start at small values of T we expect the quasienergies to correspond to the eigenvalues of the time-averaged Hamiltonian, given by with h av = αh 1 + (1 − α)h 2 the average magnetic field, which we take to be non-negative without loss of generality. This dispersion has its maximum for k = π, so we expect that this would be the first momentum mode to reach the boundary of the Brillouin zone. The last mode to reach it corresponds with the minimum k of the dispersion, given by This will define bands within which there are avoided crossings, delimited by the conditions with m ∈ N, giving the following equations and We can now use the set of equations (27) to find the boundary lines that appear in the density plots of Fig. 3. Consider, for instance, the left plot of Fig. 3 in the h av ≥ J(1 − γ 2 ) region, and let T where the set of h c is the family of values of average magnetic field where two boundaries intersect. In particular, in the left panel of Fig. 3 the end points of the non-local region in the (T, h 1 )-plane correspond to the collection of points (T c , h 1c ) = (π/2, 2(2m + 1) − 1/2) for m = 1, 2, 3, . . .. In the same way, one can derive the corresponding expressions for these intersections in the regimes h av < J(1 − γ 2 ) and J(1 − γ 2 ) ≤ h av < J.
To form a more explicit criterion for the locality of the interactions, we can look at the asymptotic decay of the couplings with increasing distance. We consider the couplings to behave as for large and N . Here, f ( ) is the dominant decaying term in the thermodynamic limit and g( ) the leading finite size correction to this decay. Let us first find it in the nonlocal regions, where the arguments of the sums in Eq. (21) present discontinuities due to the unfolding of the quasi-energy spectrum at points {k i } n i=1 for which the Floquet states are resonant. In the limit of large N all summations can be exchanged by integration, and defining a k = T (k)u k,z 1 − u 2 k,0 and integrating by parts then returns Thus, the couplings decay as f ( ) ∝ −1 in the non-local regions, as can be seen in Fig. 5. The deviations from the predicted decay for large are due to finite size corrections, which can similarly be obtained, as shown in Appendix C.
As for the local regions, there are no discontinuities in the integrand, and further analysis is necessary: We can see that all the prefactors in this asymptotic series are proportional to odd derivatives of a k . However, it can be verified that this function is symmetric both around k = 0 and k = π, and thus its odd derivatives at these points vanish. Therefore, this series shows that the couplings must decay faster than any polynomial of −1 , as can be seen in Fig. 6, where we can see an exponential decay of the couplings. This faster than algebraic decay of the function f ( ) is to be seen as a precise definition of "locality" in the present context.

Influence on the dynamics
A natural question to ask is how these features impact the evolution of observables towards the (synchronized) steady-state regime. It was previously observed that the presence of a limited number of Floquet resonances in finite-size systems has a major influence on the dynamics, leading to long-lived temporal fluctuations [19]. Similarly, in Ref. [38] open driven systems are considered where analogous regions to the ones in Fig. 1 are defined by the quick or slow spatial decay of correlation functions.
To explore this, suppose we prepare the system at an initial time in the ground state |GS of a static Hamiltonian (with e.g. h(t) = h 1 ) and analyse its stroboscopic properties after introducing driving. For a given observable O, the stroboscopic time evolution is given by where |Ψ(nT ) = U (nT, 0) |GS , with n ∈ N. This can be separated into two contributions: a decaying part tending to zero as n → ∞ and a (synchronized) steady-state [11] determined by the diagonal ensemble in the basis of Floquet eigenstates. Let us work these expression explicitly for different observables. Consider the following operators, which can be written in pseudo-spin representation as Each of these operators can be expanded as O = k>0 O k . Similarly, in all performed calculations the system is initially prepared in the ground state |GS of an XY Hamiltonian (for instance with h = h 1 ), which can also be expanded in momentum and Floquet modes as Following the previous sections, the time evolution can be explicitly calculated. Taking U k (t = nT + δt, 0) = U k (δt, 0)[U k (T, 0)] n with δt ∈ [0, T ) leads to The dynamical expectation values of the above operators are then given by These can be naturally separated in a steady-state (synchronized) term depending only on δt and a remaining time-dependent decaying term, Let us analyse each of these distinct behaviors.

Decay to the synchronized state
The dynamical behaviour of observables can be shown to differ greatly depending on the different phases of the system (see Fig. 7). In the non-local regime, convergence towards the steady-state (synchronized) regime is robust, happening as t −1/2 . In the local regime on the other hand, convergence generally behaves as t −3/2 , although specific driving protocols also allow for a slower decay as t −3/4 or t −1/2 . This behaviour is independent of the specific observable, and here it will be shown to follow from the dispersion relation of the quasienergies as function of momentum k. Provided O can be decomposed in momentum subsectors, in the thermodynamic limit the decaying part of the observable can be written as k,+ (0) GS k . The latter contains both an observable-dependent part with the off-diagonal matrix elements of O k and an observable-independent part containing the overlaps of the Floquet modes with the initial state.
The long-time dynamics of such expressions then follows from the stationary-phase approximation, where the dominant contributions to the integral will follow from stationary points where the phase θ T (k) becomes extremal, i.e. θ T (k) = 0. Each one of these points contributes to the dominant term in the long-time dynamics. First of all, it can be easily checked that the dispersion θ T (k) always has extremal points at the boundaries of the Brillouin zone, where θ T (0) = θ T (π) = 0. However, the crucial observation is that, at k = 0 and k = π, the prefactor of these contributions also vanishes, f (0) = f (π) = 0. This behaviour is independent of the observable, only following from the overlaps in f (k).
Indeed, let us look at the k = 0, π sectors of the XY Hamiltonian. For these modes, the off-diagonal elements in Eq. (7) vanish, and the Hamiltonian is already diagonal in the basis {|0 k , c † k |0 k } k=0,π , reflecting the fact that these modes cannot be doubly occupied (c † 0 c † 0 = c † π c † −π = 0). In order to respect the separation in sectors of parity of N + M of the Hilbert space, we only consider initial states with either a |0 k or c † k |0 k wavefunction component for these momenta. The same holds for the Floquet modes, which are by construction restrained to one of these parity subspaces, and therefore one of the overlaps in f (k) will always vanish at these points. From the stationary phase method, the contributions to the integral of these two stationary points can then be approximated as Note that this does not necessitate f (k 0 ) = 0, since this first-order contribution vanishes by the symmetry of the integral. If k = 0 and k = π are the only stationary points, these are the only two contributions, both scaling as n −3/2 . This is precisely what happens in the local regime, where the dispersion remains qualitatively similar to that of the time-averaged Hamiltonian, leading to a decay of t −3/2 . It can also easily be checked that additional stationary points arise in the non-local regime, corresponding to the avoided crossings in the Floquet spectrum. At these points, f (k) is not expected to vanish, leading to contributions which can be approximated as Once such an additional stationary point is present, as is the case in the non-local regime, this n −1/2 term will dominate the time evolution, leading to a slow decay of t −1/2 . Within the local regime, there is still the possibility of obtaining a t −1/2 decay, albeit with smaller deviations from the steady state, if h av < J(1 − γ 2 ), when the dispersion has one extremum in the (0, π) interval. Remarkably, it is also possible to obtain a t −3/4 decay by fine-tuning the system parameters. In all derivations so far, it was assumed that θ T (k 0 ) = 0 at k 0 = 0, π. However, in some specific cases the second and the third order derivatives vanish at either k = 0 or k = π (an explicit expression can be straightforwardly obtained, although it is rather unwieldy). Then, the contributions from these two stationary points are given by If no other stationary points are present, the n −3/4 terms will dominate the decay. These three decay regimes are illustrated in Fig. 8. This shows that in the non-local regimes, the system generally exhibits the slowest decaying exponent. Furthermore, it exhibits persistent oscillations, the frequency of which corresponds to the quasi-energy difference between the resonant states in the avoided crossings. These can be interpreted as Rabi oscillations where the resonant states repeatedly emit and reabsorb energy to/from the driving. This change of behaviour as we transition between local to nonlocal regions can also be connected to crossovers between area-law to volume-law scaling in the entanglement entropy of the steady state, which was similarly shown to be reflected in dynamical phase transitions [39].

Synchronized state
In previous studies, it has been shown that the synchronized state can be described by a Gibbs ensemble. In most instances the temperature associated to that description turns out to be infinite [14][15][16]. In our case, however, due to the restriction of Floquet resonances to the single-particle sector, the system does not heat up to an infinite-temperature state, The value of γ in the second case (green dots) results from imposing that the second derivative of the quasi-energy is zero at k = 0. In all cases we plot the absolute value of the observable so that the algebraic decay (shown by solid lines with corresponding matching colors) can be appreciated in a log-log scale. In the insets we show the total expectation value for each observable in the regime corresponding to a decay ∼ n −1/2 .  Fig. 1. The transitions from local to non-local regions are accompanied by the presence of nonanalyticities. The vertical green lines correspond to non-analyticities inside non-local regions. but is instead described by a non-trivial Periodic Generalized Gibbs Ensemble (PGGE) [12]. This implies that the resulting steady-state values of physically relevant observables will be highly dependent on the driving protocol, and possibly prone to display signals of phase transitions. This behaviour is indeed observed in Fig. 9, where the steady-state observables present non-analyticities precisely at the edges of local and non-local regimes. Furthermore, non-analyticities can be found well within the non-local regions, as indicated by the vertical pale green lines in the same figure. These correspond to the change in the number of resonant modes and it can be shown that their location in the parameter space also obeys Eqs. (27). Note that, even though the dynamics are independent of the unfolding of the spectrum, the match between the transitions present in the Floquet Hamiltonian and those in the expectation values of the observables further justifies our choice of unfolding. This approach also circumvents the use of an extended Hilbert space [50] for which exact solutions are out of reach. Finally, from Fig. 9 we can see that by choosing certain driving parameters one can target specific combinations of values for the magnetization and correlation functions as a means of engineering systems with particular properties.

Conclusions
In this work, we considered an anisotropic XY spin chain driven by a periodically-quenched external transverse field. The algebraic structure of the theory, together with the quench-based choice of protocol, have allowed us to obtain an exact expression for the Floquet Hamiltonian and its eigenstates. While the choice of Floquet Hamiltonian for a fixed t 0 is unique only up to discrete gauge choices, we have shown that a simple unfolding procedure allows to clearly highlight the different non-equilibrium phases. We have established the presence of regions with local and non-local interactions depending on the driving parameters, and present analytic expressions for the phase diagram as determined by the structure of the avoided crossings of the quasi-energies representing resonances. In the local regions this Floquet Hamiltonian is well approximated by the time-averaged one, whereas the non-local region is characterized by the presence of resonances in the Floquet single-particle spectrum, necessitating long-range interactions. For these systems, there exists a maximal value of the driving period above which the interactions remain non-local. The different phases can be observed in both the dynamics and the steady-state observables. Non-local regions are distinguished by a slow decay combined with persistent oscillations originating from Floquet resonances, whereas in local regions a quicker decay to the synchronized state can be expected, though this depends on the parameters of the driving. In future work, we will return to the questions of how more general interactions influence these effects, and how they manifest themselves for other types of observables such as longer-range and/or dynamical correlations.

A The Antiferromagnetic XY Model
We consider the spin-1/2 antiferromagnetic chain [24] described by the Hamiltonian where J > 0 is a superexchange coupling, γ ∈ [0, 1] is the anisotropy parameter, h is an external magnetic field, S α j for α ∈ {x, y, z} are spin-1/2 operators and we consider periodic boundary conditions S α j+N = S α j . Recall that the raising and lowering spin operators are defined as S ± j = S x j ± iS y j . After performing a Jordan-Wigner transformation where c and c † are operators for spinless fermions with the usual anticommunitation relations {c j , c † j } = δ jj , we obtain the following Hamiltonian where P ≡ e iπ N =1 c † c ≡ e iπM is the operator measuring the parity of the total number of fermions, which is related to the total magnetization alongẑ. The last term in Eq. 39 corresponds to the transformation of the boundary term under periodic boundary conditions. Since [H XY , P ] = 0 the Hilbert space splits into two subspaces with even or odd N + M . Moreover, as the periodic boundary conditions for the spin operators are, for c-operators, transformed into the condition (−1) N e iπ N =1 c † c c N +1 = c 1 , this implies that in the N + M even sector we have c N +1 = c 1 , while for the odd sector we obtain, in turn, c N +1 = −c 1 . All in all, the Hamiltonian takes the form where H + obeys periodic boundary conditions while H − obeys anti-periodic ones. Finally, we write the c-operators in Fourier modes where k takes values in a Brillouin Zone that depends on N and M , as shown in Table 1.
This Hamiltonian can then be diagonalised by a Bogoliubov transformation with a Bogoliubov angle defined by sin(θ k ) = −λ k / ξ 2 k + λ 2 k and cos(θ k ) = −ξ k / ξ 2 k + λ 2 k , where λ k = γ sin(k) and ξ k = cos(k) − h/J. In terms of Bogoliubov fermions, the Hamiltonian becomes with excitation spectrum The model undergoes a quantum phase transition when h = J, with an antiferromagnetic ground state for h < J, and a paramagnetic one for h > J. The dispersion relation (45) has a minimum at k = 0 when h ≥ J(1 − γ 2 ), and at k = arccos h/J(1 − γ 2 ) when h < J(1 − γ 2 ). In the latter case, k = 0 becomes a local maximum. In both situations, there is a maximum for k = π. Another approach to the problem comes from noticing that the Hamiltonian (42) allows a pseudospin 1/2 representation since it connects the Fock vacuum |0 to the state c † k c † −k |0 for each k. We thus introduce the following representation Thus operators acting on this space, and expressible in terms of sum of operators acting on k-labeled subspaces (e.g. O = k>0 O k ) can be written, in this representation, as sums of 2 × 2 matrices whose entries are (O k ) ij = i|O k |j for i, j ∈ {1, 2}. For the Hamiltonian (42) this allows to write it as where H k = E(k)n k ·σ k with σ k = (σ x k , σ y k , σ z k ) being Pauli matrices, and n k = (0, − sin(θ k ), cos(θ k )). This representation is related to Eq. (44) through a particle-hole transformation of the negative momentum Bogoliubov modes. Thus, the diagonalized Hamiltonian presents two energy branches E σ (k) = σE(k) for σ ∈ {−1, 1} with eigenvectors The pseudospin representation is particularly useful as it allows to express time evolution as rotations on the Bloch sphere around the axis n k , that is,

B Other Floquet gauge choices and continuous time evolution
In the main text we consider focus on the case where t 0 = 0 an consider mostly the stroboscopic dynamics. Suppose, however, that we start the evolution of the system at an arbitrary Floquet gauge t 0 and we want to evolve it up to a time t 0 + δt, δt ∈ [0, T ]. Then the evolution operator in a given k sector is given by where H (i) k is the Hamiltonian with field h i , with eigenvalues E i (k) and Bogoliubov angles θ (i) k . Using an Euler angle representation, we will write this evolution operator as where the coefficients u k,l (δt) are given by k . The coefficients are given in a piecewise manner for δt in the three intervals [0, αT − t 0 ], [αT − t 0 , T − t 0 ], and [T − t 0 , T ], respectively. This result will allow us to do any arbitrary time evolution and to obtain the Floquet modes by diagonalising it for δt = T . We denote the Euler parameters of the stroboscopic time evolution operator U k (t 0 + T, t 0 ) as u k,p ≡ u k,p (δt = T ), p ∈ {0, x, y, z}, whose expressions turn out to be

C Finite-size corrections
To derive expressions describing the finite size effects, an easy way to obtain the correction function g( ) is to replace the boundaries of the integral by the lowest and highest momenta SciPost Physics Submission allowed by the quantization: A ≈ 1 π Re π−π/N π/N a k e ik dk = 1 π (a π−π/N sin((π − π/N ) ) − a π/N sin((π/N ) ) + O( −2 ) with sinc(x) = sin(x)/x. This finite size correction is the same in both the local and non-local regime. The same analysis can be done for B yielding the same asymptotic behaviour in both the local and non-local regimes. The finite size corrections however are slightly different, and are given by