Pumping approximately integrable systems
PPumping approximately integrable systems
Florian Lange, Zala Lenarˇciˇc, ∗ and Achim Rosch Institute for Theoretical Physics, University of Cologne,Z¨ulpicher Straße 77a, D-50937 Cologne, Germany
Weak perturbations can drive an interacting many-particle system far from its initial equilibriumstate if one is able to pump into degrees of freedom approximately protected by conservation laws.This concept has for example been used to realize Bose-Einstein condensates of photons, magnons,and excitons. Integrable quantum systems, like the one-dimensional Heisenberg model, are charac-terized by an infinite set of conservation laws. Here we develop a theory of weakly driven integrablesystems and show that pumping can induce large spin or heat currents even in the presence ofintegrability breaking perturbations, since it activates local and quasi-local approximate conservedquantities. The resulting steady state is qualitatively captured by a truncated generalized Gibbsensemble with Lagrange parameters that depend on the structure but not on the overall amplitudeof perturbations nor the initial state. We suggest to use spin-chain materials driven by terahertzradiation to realize integrability-based spin and heat pumps.
A simple classical example for a weakly driven sys-tem is a well-insulated greenhouse. Due to the approxi-mate conservation of the energy within the greenhouse,even weak sunlight can lead to high temperatures in itsinterior, which can be computed from the simple rateequation for the energy transfer. Similarly, large spinaccumulation can be achieved in systems with approxi-mate spin conservation [1]. Using approximate conser-vation of the number of photons, magnons, or excitonpolaritons, one can use pumping by light to reach densi-ties which allow for the realization of Bose-Einstein con-densates [2–4]. Number-conserving collisions induce aquasi-equilibrium state in these systems, which can beefficiently described by introducing a chemical potentialwhose value is determined by balancing pumping and de-cay processes. Related theoretical approaches that de-scribe electron-phonon systems far from equilibrium areso-called two-temperature models [5]: here one uses thatthe energy of the electrons and phonons are approxi-mately separately conserved to introduce two differenttemperatures for the subsystems.Integrable many-particle systems, like the one-dimensional (1D) fermionic Hubbard model or the XXZHeisenberg model, are described by an infinite number of(local or quasi-local) conservation laws [6–11]. In closedintegrable systems those prevent the equilibration intoa simple thermal state, e.g., after a sudden change ofparameters. Instead the system can be described by ageneralized Gibbs ensemble (GGE) [12–21] ρ ∼ exp (cid:32) − (cid:88) i λ i C i (cid:33) (1)where C i are the conserved quantities and λ i the corre-sponding Lagrange parameters. It has also been shownexperimentally [22] that GGEs for a Lieb-Liniger modelcan provide highly accurate descriptions of interactingbosons in 1D. ∗ [email protected] a b FIG. 1. (a) A well-insulated greenhouse exposed to sunshinecan heat up significantly since energy within it is approxi-mately conserved. (b) As the heat current in spin chain ma-terials is approximately conserved even weak terahertz radia-tion can induce large heat current. Material candidates musthave appropriate crystal structure, schematically denoted bydashed lines indicating alternating chemical bonds.
Many materials are described with high accuracy by in-tegrable models [23], however, weak integrability break-ing terms and the coupling to thermal phonons implythat in equilibrium these systems are described by sim-ple thermal states, ρ ∼ e − βH , instead of GGEs. Theproximity to the integrable point and the presence ofapproximate conservation laws leads to enhanced spinor heat conductivities (within linear-response theory)[24–26] and also to a slow relaxation after a quantumquench (via GGE-prethermalization) towards the equi-librium state [27].We will show that – as in the greenhouse example, seeFig. 1 – such an approximately integrable system canbe driven far from its thermal equilibrium by weak per-turbations arising, e.g., from a driving periodic in timeor from coupling to a non-thermal bath. In order tobalance the constant heating due to driving the systemhas to be weakly open, e.g., by coupling to a phononbath. As we will demonstrate this mechanism can beused for example to create large spin and heat currents.Besides the quasi 1D systems considered by us, also ap-proximately many-body localized systems are character-ized by infinitely many approximate conservation lawswhich may lead to a strong response to driving [28, 29]. a r X i v : . [ c ond - m a t . s t r- e l ] J un Results
Weakly driven system.
We consider an interact-ing many-body system that is approximately describedby Hamiltonian H and characterized by a finite or in-finite number of (quasi-)local conserved quantities C i ,[ H , C i ] = 0, one of them being H . Energy and otherconservations are weakly broken by coupling to thermalor non-thermal baths and/or perturbations periodic intime. For simplicity we assume periodic boundary con-ditions and a (discreet) translational invariance. We de-scribe the system with density matrix ρ whose dynamicsis governed by the Liouvillian super-operator ˆ L ,˙ ρ = ˆ L ρ, ˆ L = ˆ L + (cid:15) ˆ L , (2)where ˆ L can be split into the dominant unitary Hamil-tonian evolution ˆ L ρ = − i [ H , ρ ] and perturbation ˆ L ofstrength (cid:15) . We are interested in the limit of small (cid:15) for t → ∞ where a unique (Floquet) steady state ρ ∞ is ob-tained. The general structure of perturbation theory inthis case has, e.g., been discussed in Refs. [30–32]. In thislimit, ρ ∞ can be approximated by ρ = lim (cid:15) → lim t →∞ ρ with ˆ L ρ = 0 according to Eq. (2). We assume and latersupport numerically that ρ is approximately describedby a GGE, see Eq. (1).Here it is essential to note that – as in the greenhouseexample discussed above – the parameters λ i are not de-termined by the initial state but by the form of the weakperturbations ˆ L . Our central goal is to compute the λ i .We first discuss the case of Lindblad dynamics, whereperturbation theory linear in (cid:15) can be used, and then fo-cus on Hamiltonian dynamics where we have to consider (cid:15) contributions. Markovian perturbation.
Within the Markovian ap-proximation one can use the Lindblad form for ˆ L [33].Note that Lindblad dynamics is considered here mainlyfor pedagogical purposes (formulas are simpler) while noLindblad approximation is used for the models studiedbelow. The coefficients λ i that fix the GGE are deter-mined from the condition that the change of the approx-imately conserved quantities has to vanish in the steadystate (cid:104) ˙ C i (cid:105) = Tr ( C i ˆ L ρ ) = Tr ( C i (cid:15) ˆ L ρ ) ! = 0 , (3)where we used that ˆ L ρ = − i [ H , ρ ] = 0. Relation(3) yields a set of coupled equations for λ i , where thenumber of equations is equal to the number of conservedquantities. We define the super-projector ˆ P onto thetangential space of GGE density matrix,ˆ P X ≡ − (cid:88) i,i (cid:48) ∂ρ ∂λ i ( χ − ) ii (cid:48) Tr ( C i (cid:48) X ) , (4)using χ ii (cid:48) = − Tr ( C i ∂ρ /∂λ i (cid:48) ). Then the conditions for ρ can be compactly written asˆ L ρ = 0 , ˆ P ( ˆ L ρ ) = 0 . (5) This equation can also be derived by considering higherorder perturbations in (cid:15) , see Methods for details. Hamiltonian perturbation.
For Hamiltonian dynam-ics (cid:15) ˆ L ρ = − i [ H , ρ ], where H may be a sum of sev-eral integrability breaking perturbations. Perturbationtheory linear in (cid:15) vanishes, Tr ( C i (cid:15) ˆ L ρ ) = 0 for all λ i .Therefore one has to expand to order (cid:15) and Eq. (5) isreplaced by ˆ L ρ = 0 , ˆ P ( ˆ L ˆ L − ˆ L ρ ) = 0 . (6)Since ˆ P ( ˆ L ρ ) = 0, ˆ L ρ is not in the kernel of ˆ L − .For periodic driving this equation has to be interpretedwithin the Floquet formalism, see Methods. Model.
As discussed in the introduction, our goal isto describe a situation which can be realized experimen-tally in spin-chain materials driven by lasers operatingin the terahertz regime. We assume that spin chains areapproximately described by a spin-1/2 XXZ Heisenbergmodel, possibly in the presence of an external magneticfield B , H = (cid:88) j J S + j S − j +1 + S − j S + j +1 ) + ∆ S zj S zj +1 − BS zj . (7)The system is driven out of equilibrium by a weak(integrability-breaking) time-dependent perturbation H d = (cid:15) d J (cid:88) j (cid:16) ( − j +1 S j · S j +1 sin( ωt ) + ( − j S zj cos( ωt ) (cid:17) , (8)with driving frequency ω . This specific term has beenchosen because it can induce heat and spin currents (ascan be shown by a symmetry analysis), and becauseit can be realized experimentally. Such staggered ex-change couplings and staggered magnetic fields arise nat-urally in certain compounds with (at least) two magneticatoms per unit cell when coupled to uniform electric andmagnetic fields, respectively [34–37]. See Fig. 1 for aschematic drawing of such a compound and Methods forconcrete experimental suggestions. Therefore H d can berealized by shining a laser (typically at terahertz frequen-cies) onto the sample. In this case (cid:15) is proportional tothe laser power. Note that for T = 0 and B = 0 inthe adiabatic limit, ω →
0, Eqs. (7,8) realize an adia-batic Thouless pump, where per pumping cycle one spinis transported by one unit cell [38]. We will be interestedin the opposite regime of large ω and large (effective)temperatures.Formally the periodic perturbation H d would drive thesystem to infinite temperature [39–42] (up to remainingconservation laws [43], possibly through a prethermal-like regime [44]). In a solid state experiment this is pro-hibited by the coupling to phonons and, ultimately, tothe thermal environment of the experimental setup. Wemimic this effect by coupling the spin system to a bathof Einstein phonons, H ph0 = ω ph (cid:80) j a † j a j + . . . , wheredots stand for the couplings to further reservoirs whichguarantee that the phonon system is kept at fixed tem-perature T ph , ρ ph ∼ e − H ph0 /T ph . See Methods for detailson finite size calculation using a broadened distributionof phonon energies. The (weak) coupling to the spin sys-tem is described by H ph = (cid:15) ph J (cid:88) j (cid:16) S j · S j +1 ( a j + a † j ) (9)+ γ m ( S xj S zj +1 + S zj S xj +1 )( a j + a † j ) (cid:17) . To obtain a unique steady state it is essential to break allsymmetries, including the S z conservation. Relativisticeffects which relax S z are mimicked by γ m in our ap-proach. We expect γ m (cid:28) γ m = 1 within our numer-ics as this is found to minimize finite size effects, withouta qualitative influence on the results. Besides phononsalso other integrability breaking perturbations exist inreal materials, including defects, which typically domi-nate at the lowest temperatures. For high temperaturesof the order of J (relevant for the considered setup) it isrealistic to assume that phonon coupling dominates.In the presence of a periodic perturbation, Eq. (8), inthe long-time limit the density matrix is changing peri-odically, ρ ( t → ∞ ) = (cid:80) n e − iωnt ρ ( n ) with ρ ( n ) † = ρ ( − n ) , n ∈ Z . Within the Floquet formalism one therefore pro-motes the steady-state density matrix to a vector andLiouville operator to a matrix, see Methods. For weakdriving, (cid:15) d →
0, only the n = 0 sector remains and theGGE ansatz, Eq. (1), simply reads ρ ( n )0 = ( ρ ⊗ ρ ph ) δ n, where we included also the phonon density matrix, seeabove. Steady state.
We will use two different approaches todetermine an approximate solution for the steady statedensity matrix. First, we will parametrize ρ , Eq. (1),with a small number of (quasi-)local conserved quanti-ties, C i , i = 1 , . . . , N C . In an alternative approach, fea-sible for small systems, we take all conserved quantitiesinto account: local and non-local, commuting and non-commuting. While the second approach is formally exactin the limit (cid:15) d , (cid:15) ph →
0, the first one is, perhaps, moreintuitive and can be computed for larger system sizes.For the XXZ Heisenberg model an infinite set of mutu-ally commuting local conserved quantities C i is known,see Methods. C is the total spin C = (cid:80) i S zi and C = H XXZ . Importantly, C is the heat current [45], C = J H ( B = 0). In addition there also exist (in-finite) sets of quasi-local commuting conserved quanti-ties [8–10]. As shown in [8, 46] the spin-reversal parity-odd family has an overlap with the spin current J S at∆ < J . Therefore both heat and spin current could showa large response to a weak perturbation. For our anal-ysis we choose three or five ( N C = 4 , N C = 6) most - - - ( d / ph ) = - - - ( d / ph ) = ba FIG. 2. Effective force F in the space of Lagrange parame-ters ( β, λ ) using e − βH − λ C as an ansatz for the generalizedGibbs ensemble. Parameters: J = ∆ = − B = ω = ω ph = T ph . Lagrange parameters ( β, λ ) are plotted in units 1 /J and1 /J , respectively. (a) In the absence of an external driving, (cid:15) d = 0, the stable fixed point (red dot) is given by the thermalensemble, β = 1 /T ph , λ = 0. (b) When the system is drivenby H d ( (cid:15) d = (cid:15) ph ), it heats up and λ becomes finite. local conserved quantities C i , i = 1 , ..., N C −
1. From thequasi-local sets we include as a single (effective) operatorthe conserved part of spin current J cS , computed numer-ically [25, 47]. For details see Methods. In the presenceof an external magnetic field, Eq. (7), the heat currentalso has, in addition to C , a spin current component, J H = C − BJ S .For the visualization of our results it is useful to definegeneralized forces F i in the space of Lagrange parametersby rewriting ˆ P ˙ ρ = (cid:80) i ∂ρ ∂λ i F i such that ˙ λ i ≈ F i , F i = (cid:88) i (cid:48) ( χ − ) ii (cid:48) Tr ( C i (cid:48) (cid:15) ˆ L ˆ L − (cid:15) ˆ L ρ ) (10)computed using exact diagonalization, see Methods. Thevector F is a function of the Lagrange parameters λ i which points into the direction of the steady state stablefixed point obtained from F i = 0. In the absence ofdriving (Fig. 2a) one obtains the expected thermal statewith T = T ph while all other Lagrange parameters λ i vanish. For finite driving the GGE is activated and the λ i become finite (Fig. 2b). To obtain the steady state,we solve χ F = 0 using Newton’s method.For the second approach, performed on small N -sitesystems, we first numerically construct a basis in the setof all (local and non-local) conserved operators, Q = {| n (cid:105)(cid:104) m | with E m = E n } , where H | n (cid:105) = E n | n (cid:105) . Dueto degeneracies we find (for finite B and ∆ (cid:54) = J ) about2 · N elements Q i ∈ Q . In the limit (cid:15) d , (cid:15) ph → ρ ∞ has to fulfill ˆ L ρ ∞ = 0and therefore can be exactly written as a linear combina-tion of the Q i , ρ ∞ = (cid:80) α i Q i . Using Eq. (6), we thereforefind that the steady state density matrix for (cid:15) d , (cid:15) ph → L Q mn = − Tr ( Q † m (cid:15) ˆ L ˆ L − (cid:15) ˆ L Q n ) , (11)where ˆ L , ˆ L are Floquet matrices, see Methods. Notethat only the relative (cid:15) d /(cid:15) ph and not the absolutestrength of perturbations determine ρ , as can be seenby dividing the equations χ F = 0 or L Q ρ = 0 by (cid:15) .In Fig. 3 we show the expectation value of the energy N N N N
160 5 10 15 20 25 300.0.0020.004 C N H N a a Ε d Ε ph ( d / ph ) N C N C Ε d Ε ph J H N b ( d / ph ) FIG. 3. Expectation values of (a) energy and (b) heat cur-rent densities for a weakly driven spin chain, (cid:15) d , (cid:15) ph → (cid:15) d and phononcoupling (cid:15) ph . Red solid lines: exact result taking into ac-count all 7969 conservation laws of a system of N = 12sites. (a) For the energy accurate results are already ob-tained with a GGE ensemble based on N C = 4 (dot-dashedlines) or N C = 6 (dashed lines) conserved quantities. (b)Also the heat current J H = C − BJ S is qualitatively welldescribed by the GGE ensemble but quantitative deviationsare larger. Inset: Finite size analysis for (local) C based onGGE ensemble with N C = 6 conserved quantities. Parame-ters: J = 1 , ∆ = 0 . , B = − . , ω = 1 . ω ph , ω ph = T ph = 1. ( d / ph ) J S / N N C = 4 N C = 6exact FIG. 4. For vanishing magnetic field a spin current (but noheat current) is generated within our model for finite ratios of (cid:15) d /(cid:15) ph . The expectation value of spin current density is againmaximal for (cid:15) d /(cid:15) ph ≈
1. Parameters: J = 1 , ∆ = 0 . , ω =1 . ω ph , ω ph = T ph = 1 , N = 12. and of the heat current densities as functions of (cid:15) d /(cid:15) ph taking into account N C = 4, N C = 6, and all conservedquantities. The energy density expectation value is al-ready obtained with good accuracy for N C = 4 and evenbetter for N C = 6. The heat current vanishes both inthermal equilibrium, (cid:15) d →
0, and for (cid:15) ph →
0, where thesystem is described by an infinite temperature state withfinite magnetization, ρ ∼ e − λ S z and (cid:104) H (cid:105) = − B (cid:104) S z (cid:105) .It takes its largest value for (cid:15) ph ∼ (cid:15) d . For the currents adescription in terms of N C = 4 or 6 is qualitatively butnot quantitatively accurate. Our study strongly suggeststhat further quasi-local conserved quantities contribute,as discussed in quench protocols [15–17], see also Ref.[25]. For the chosen parameters our results depend onlyweakly on the system size N , see inset of Fig. 3. Systemsize analysis is performed for N C = 6 since the solutionbased on all conservations cannot be obtained for largersystems.Our setup can also be used to create spin currents.Whilst, by symmetry (bond-centered rotation in real andspin space by π around y axis), a finite external field B is needed to obtain a finite heat current, this is not thecase for the spin current. Fig. 4 displays the spin currentdensity as a function of (cid:15) d /(cid:15) ph for B = 0. Qualitativelyone obtains a behavior rather similar to the results for N C B J H N a B J N J S C b FIG. 5. (a) Heat current J H , (b) spin current J S , and C den-sities as a function of external magnetic field B obtained froma GGE ensemble with N C = 6 conserved quantities (dashed)or from an exact calculation (solid) including all conserva-tions. Parameters: ( (cid:15) d /(cid:15) ph ) = 2 . J = 1 , ∆ = 0 . , ω =1 . ω ph , ω ph = T ph = 1 , N = 12. the heat current shown in Fig. 3 with a maximum in thespin current for (cid:15) ph ∼ (cid:15) d .The external magnetic field B is a parameter whichcan easily be tuned experimentally. Fig. 5 shows heatand spin current densities as a function of external mag-netic field B for ( (cid:15) d /(cid:15) ph ) = 2 .
5. Note that the sign ofthe magnetic field determines the sign of the heat cur-rent (cid:104) J H (cid:105) = (cid:104) C (cid:105) − B (cid:104) J S (cid:105) . All main features of theB-dependence are semi-quantitatively reproduced by thetruncated GGE with N C = 6. For very large mag-netic fields the convergence to the steady state fixedpoint becomes slow as transitions rates connecting sec-tors with different magnetization are strongly suppressed,see Methods for further details. Discussion
We have demonstrated that driving approximately inte-grable systems activates and pumps into approximatelyconserved quantities. Perhaps the most simple experi-mental setup to measure the pumping effect predictedin this work, is to use a terahertz laser that excites aspin chain material like Cu-benzoate where by symmetrystaggered terms of the form (8) are expected [34, 35]. Asa consequence of the induced heat currents it is antici-pated that the system cools down on one side while itheats up on the other. The direction of the effect can becontrolled either by changing the direction of the laserbeam or the sign of the external magnetic field B .For the chosen parameters, the spin and heat cur-rents expressed in dimensionless units appear to be rathersmall of the order of 10 − . While these values can def-initely be increased by tuning parameters, for examplethe external magnetic field, it is important to note thatthe currents are actually quite large compared to the typ-ical heat or spin currents obtained in bulk materials. Tocreate a heat current of similar size in a good heat con-ductor like Cu (assuming J ∼ k B ·
100 K, 5 ˚A for the dis-tance of the spin chains, and κ Cu ≈
400 Wm − K − ) onewould need a temperature gradient of several 10 Km − .Similarly, to create a (transversal) spin current of compa-rable size in a heavy element like Pt using the spin-Halleffect (assuming ρ Pt ≈ µ Ω cm and α Pts ≈
10 % for thespin Hall angle [48]) one needs electric fields of the or-der of 10 Vm − or sizable current densities of the orderof 10 Am − . These numbers are even more remarkablewhen one takes into account that the electron densitiesin Cu or Pt are at least an order of magnitude higherthan the spin density for spin-chains with a distance of5 ˚A.While our study has focused on the steady state, itis instructive to discuss the relevant time scales for itsbuildup. For this argument we consider a quench whereat time t = 0 an initial state is perturbed both by theintegrable part of the Hamiltonian and by small non-integrable perturbations. At short times of the order ofseveral 1 /J the initial state will prethermalize [27, 49–52]into a GGE where the values of the conserved quanti-ties, (cid:104) C i (cid:105) , are set by the initial conditions (with small corrections from the perturbations [52, 53]). Furthertime evolution can be approximately described by a GGEwith time-dependent Lagrange parameters. Their time-dependence is determined by perturbations which assertforces F i ∼ (cid:15) , such that dλ i /dt ≈ F i . Governed bythe perturbations the system will loose the memory ofits initial condition on a time scale of order 1 /(cid:15) and re-lax to the steady state (obtained from F i = 0) which is,in general, completely unrelated to the prethermalizedstate. Note that the same approach predicts ordinarythermalization in the absence of external driving.Our results suggest that the concept of generalizedGibbs ensembles has a much broader range of applica-tion than previously anticipated, now extended to opensystems where symmetries are not exact and integrabilityis weakly broken. A truncated GGE proved to be use-ful for qualitative description, however, it showed quan-titative discrepancies most probably due to disregardedquasi-local conserved quantities, as observed already inquench protocols [15, 16]. We are planning a future studytailored to address this issue systematically. It would beinteresting to develop integrability-based methods simi-lar to the quench-action approach [15, 16, 54, 55] to treatsuch situations.Most important for applications is that the integrabil-ity is not required to be realized exactly but only ap-proximately. Efficient pumping requires only that thepumping rates are of the same order of magnitude asthe loss rates arising from integrability breaking terms.Especially the integrability based creation of large spincurrents could find its application in future spintronicsdevices. Methods
Perturbing around ρ . The central equations (5) or(6), used to determine the density matrix ρ in the limit (cid:15) →
0, have to be consistent and can also even be derivedby considering perturbations around ρ , ρ ∞ = ρ + δρ .First, the leading δρ correction to (cid:104) ˙ C i (cid:105) , Eq. (3), aris-ing from Tr ( C i ˆ L δρ ) which is nominally of the same or-der as Tr ( C i (cid:15) ˆ L ρ ) vanishes trivially as Tr ( C i [ H , δρ ]) =Tr ( δρ [ C i , H ]) = 0.For arbitrary ρ , δρ is exactly given by δρ = − ˆ L − (cid:15) ˆ L ρ , where ˆ L − is a short-hand notation forlim η → ( ˆ L − η ˆ1) − with the infinitesimal regularizer η .The correct expansion point ρ is found if lim (cid:15) → δρ =0. Below we show that for the projection operator ˆ P ,Eq. (4), ˆ L − ˆ P ∼ O ( (cid:15) − ) , (12)which would yield ˆ L − ˆ P (cid:15) ˆ L ρ ∼ O (1). This contradictsour perturbative approach unless ˆ P ˆ L ρ = 0, as set byour condition Eq. (5).Eq. (12) is a consequence of the fact that ˆ P projectsonto the tangential space to GGE density matrix. Inthis space ˆ L vanishes by definition, ˆ L ( ∂ρ /∂λ i ) = 0,and ˆ L = ˆ L + (cid:15) ˆ L is therefore of order (cid:15) . Technically,this can be seen by using the general relation( ˆ X + ˆ Y ) − − ˆ X − = − ( ˆ X + ˆ Y ) − ˆ Y ˆ X − (13)for ˆ X = ˆ P (cid:15) ˆ L ˆ P , ˆ Y = ˆ L + ˆ Q(cid:15) ˆ L ˆ Q + ˆ P (cid:15) ˆ L ˆ Q + ˆ Q(cid:15) ˆ L ˆ P , (14)with ˆ Q = ˆ1 − ˆ P and ˆ X + ˆ Y = ˆ L . Thenˆ L − ˆ P = ( ˆ X + ˆ Y ) − ˆ P = ˆ X − ˆ P − ( ˆ X + ˆ Y ) − ˆ Q ˆ Y ˆ P ˆ X − ˆ P ∼ O ( (cid:15) − ) + O (1) (15)The second term is O (1) as ˆ L ˆ P = 0 and thereforeˆ Y ˆ P ∼ O ( (cid:15) ). The divergence of ˆ L − ˆ P for (cid:15) → Staggered hopping and magnetic field modula-tion.
Sizable staggered g-tensors leading to staggeredB-fields have been observed in a number of different com-pounds [34–37]. Similarly an external electric field willdistort the crystalline structure in these materials, lead-ing to staggered exchange couplings linear in homoge-neous electric fields. An example of such a material isCu-benzoate [34] with the above modulations allowed bysymmetry for electric (magnetic) fields applied in the 010(001) crystallographic direction. In this system the stag-gered g-tensor has been measured to be approximately0.08 [35], the size of the staggered exchange coupling isunknown. For simplicity, we assume in Eq. (8) that thetwo staggered terms are of the same size.
Conservation laws of the XXZ Heisenberg model.
An infinite set of local conserved quantities C i of theHeisenberg model H XXZ = H ( B = 0) can be ob-tained using the boost operator O b = − i (cid:80) j jh j,j +1 (where H XXZ = (cid:80) j h j,j +1 ) from the recursion relation[ O b , C i ] = C i +1 for i > C = (cid:80) j S zj , C = H XXZ [7]. In general, C i are operators involving maximally i neighboring sites. Importantly, C in the absence of ex-ternal magnetic field equals the heat current J H ( B = 0) = C = J (cid:88) j ( S (cid:48) j × S (cid:48)(cid:48) j +1 ) · S (cid:48) j +2 (16)with rescaled spin operators S (cid:48) aj = √ λ a S aj , S (cid:48)(cid:48) aj = (cid:112) λ z /λ a S aj for λ z = ∆ /J , λ x = λ y = 1 . In the pres-ence of external magnetic field, Eq. (7), heat current has in addition to C also a spin current component, J H = J (cid:88) j ( S (cid:48) j × S (cid:48)(cid:48) j +1 ) · S (cid:48) j +2 − BJ S . (17)As understood recently there also exist families ofquasi-local conserved quantities [8–10], which are mostlydisregarded in our study with the exception of a spin-reversal parity-odd operator, J cS . The latter is con-structed as the conserved part of the spin current op-erator J S , J S = i J (cid:88) j ( S + j S − j +1 − S − j S + j +1 ) (18) J c S = (cid:88) ˜ n | ˜ n (cid:105)(cid:104) ˜ n | J S | ˜ n (cid:105)(cid:104) ˜ n | where | ˜ n (cid:105) are simultaneous eigenstates of the C i . Sinceit is known that the spin current has an overlap with thequasi-local family [46] for ∆ < J , the conserved J cS contains quasi-local components (and, possibly, non-local components not contributing in the thermodynamiclimit). Floquet formulation.
For a periodically driven systemdescribed by ˙ ρ = ˆ L ( t ) ρ with ˆ L ( t + T ) = ˆ L ( t ) the den-sity matrix changes periodically in the long-time limit.Therefore it is useful to split it into Floquet components, ρ = (cid:88) n e − inωt ρ ( n ) , n ∈ Z (19)with ρ ( − n ) = ρ ( n ) † and ω = 2 π/T . The Flo-quet components are combined into the vector ρ =( . . . ρ ( − , ρ (0) , ρ (1) , . . . ). The Liouvillian is promoted toa (static) matrix ˆ L nm = inωδ nm + ˆ L n − m with ˆ L n − m = T (cid:82) T ˆ L ( t ) e iω ( n − m ) t dt . Using this notation, all resultsobtained for static Liouvillian super-operators directlytranslate to the time-periodic case. Within our setup, H , all approximate conservation laws C i and the GGEdensity matrix ρ are static and therefore the projectionoperator ˆ P , Eq. (4), projects onto the n = 0 Floquet sec-tor only. The steady state condition, Eq. (6), thus meansthat the approximately conserved quantities do not growafter averaging over an oscillation period. To second or-der in (cid:15) d only transitions from the n = 0 to the n = ± L n = 0 for | n | > F (d) i = 2 πN (cid:15) (cid:88) i (cid:48) ( χ − ) ii (cid:48) (cid:88) m,k ρ m ( C i (cid:48) ,m − C i (cid:48) ,k ) ×× (cid:110) |(cid:104) k | H (+)d | m (cid:105)| δ ( E k − E m − ω )+ |(cid:104) k | H ( − )d | m (cid:105)| δ ( E k − E m + ω ) (cid:111) (20)where we used H eigenstates | m (cid:105) with H | m (cid:105) = E m | m (cid:105) ,matrix elements ρ m = (cid:104) m | ρ | m (cid:105) , C i,m = (cid:104) m | C i | m (cid:105) , andthe notation H d = (cid:15) d (cid:16) e iωt H ( − )d + e − iωt H (+)d (cid:17) . Notethat Eq. (20) contains – as expected – transition rateswell-known from Fermi’s golden rule. Eq. (20) is eval-uated for finite systems of size N by replacing the δ function by a Lorentzian (1 /π ) η/ ( ω + η ) ( η = 0 . J for N = 12).Eq. (20) is only valid for situations where all con-servation laws commute with each other, with C i = (cid:80) m | m (cid:105) C i,m (cid:104) m | , see below for a brief discussion of thenon-commuting case. Phonon coupling.
As written in the main text, weassume that the phonon system always remains at equi-librium, ρ ph ∼ e − H ph0 /T ph . Using Eq. (10), after tracingover phonons, we obtain for the generalized force F (ph) i = 2 π(cid:15) (cid:88) i (cid:48) ( χ − ) ii (cid:48) (cid:88) m,k ρ m ( C i (cid:48) ,m − C i (cid:48) ,k ) × J (cid:0) |(cid:104) k | S j · S j +1 | m (cid:105)| + γ |(cid:104) k | S xj S zj +1 + S zj S xj +1 | m (cid:105)| (cid:1) × (cid:0) ( n B ( E m − E k ) + 1) A (ph) ( E m − E k )+ n B ( E k − E m ) A (ph) ( E k − E m ) (cid:1) (21)where n B ( E ) = 1 / ( e E/T ph −
1) is the equilibrium Bose dis-tribution evaluated at the temperature T ph and A (ph) ( ω )is the phonon spectral function. For our finite size cal-culation we broaden the spectral function of the Einsteinphonons using A (ph) ( ω ) = Θ( ω ) ωω ph η √ π e − ( ω − ω ph ) /η .This choice of broadening ensures detailed balance rela-tions (necessary to obtain a thermal state in the absenceof driving) and the positivity of phonon frequencies (nec-essary for stability). For all plots we use η = 0 . J . How-ever, we have checked that similar results are obtained,e.g., for η = 0 . J for magnetic fields up to | B | = 2 J .For larger fields η = 0 . J does not provide a sufficientamount of relaxation between sectors with different mag-netization and convergence becomes slow and unstable.For η = 0 . J larger fields, | B | (cid:46) J , can be reached. Implementation of non-commuting conservationlaws.
As discussed in the main text, a complete basisof all non-local commuting or non-commuting conservedquantities is given by Q = {| n (cid:105)(cid:104) m | with E m = E n } whichsolve the equation ˆ L Q i = 0 for Q i ∈ Q . Using the ex-act eigenstates of H it is straightforward to evaluateEq. (11) where we use for our finite size calculations thebroadening procedures described above. As a technicaldetail we note that, when one follows this procedure, onehas to evaluate in the phonon sectors integrals of thetype (cid:82) A (ph) ( ω (cid:48) ) ω − ω (cid:48) n B ( ω (cid:48) ) dω (cid:48) numerically. For efficient eval-uations we use interpolating functions for these integrals. GGE estimation for other conserved quantities.
To provide further support for our claim that truncatedGGEs give a semi-quantitative description of our weakly N C B H N a B C N b FIG. 6. (a) The energy density and (b) the expectation valueof another conserved quantity C (4-spin operator) as a func-tion of magnetic field B , obtained from calculation using allconserved quantities (solid) and a GGE with N C = 6 (quasi-)local conserved quantities (dashed). open system we show in Fig. 6 additional comparisonof the (cid:104) H (cid:105) and (cid:104) C (cid:105) as a function of magnetic field B at ( (cid:15) d /(cid:15) ph ) = 2 .
5, comparing as in the main text theexact calculation including all conserved quantities andthe truncated GGE with N C = 6 (quasi-)local conservedquantities. The GGE ansatz captures the right magni-tude and the correct behaviour in the dependence on B also for more complicated 4-spin operators like C . Weuse same parameters as for the Fig. 5 in the main text:( (cid:15) d /(cid:15) ph ) = 2 . , J = 1 , ∆ = 0 . , ω = 1 . ω ph , ω ph = T ph = 1 , N = 12. Acknowledgements
We acknowledge useful discussions with S. Diehl, F. H.L. Essler, M. Fagotti, E. Ilievski, M. Mierzejewski, J.De Nardis, T. Prosen, and M. C. Rudner, H. F. Leggfor reading the manuscript, and financial support of theGerman Science Foundation under CRC 1238 (projectC04) and CRC TR 183 (project A01).
Author contributions
A. R. and Z. L. designed the study, Z. L. and F. L. per-formed analytical calculations and F. L. implemented thenumerical codes, all authors analyzed the results and con-tributed to the manuscript. [1] Kikkawa, J. M. & Awschalom, D. D. Resonant spin am-plification in n -type GaAs. Phys. Rev. Lett. , 4313–4316 (1998).[2] Klaers, J., Schmitt, J., Vewinger, F. & Weitz, M. Bose-Einstein condensation of photons in an optical microcav-ity. Nature , 545–548 (2010).[3] Demokritov, S. et al.
Bose–Einstein condensation ofquasi-equilibrium magnons at room temperature underpumping.
Nature , 430–433 (2006).[4] Kasprzak, J. et al.
Bose–Einstein condensation of excitonpolaritons.
Nature , 409–414 (2006).[5] Allen, P. B. Theory of thermal relaxation of electrons inmetals.
Phys. Rev. Lett. , 1460–1463 (1987).[6] Faddeev, L. Algebraic aspects of the Bethe ansatz. In-ternational Journal of Modern Physics A , 1845–1878(1995).[7] Grabowski, M. P. & Mathieu, P. Structure of the con-servation laws in integrable spin chains with short rangeinteractions. Ann. Phys. , 299–371 (1995).[8] Prosen, T. & Ilievski, E. Families of quasilocal conserva-tion laws and quantum spin transport.
Phys. Rev. Lett. , 057203 (2013).[9] Mierzejewski, M., Prelovˇsek, P. & Prosen, T. Identify-ing local and quasilocal conserved quantities in integrablesystems.
Phys. Rev. Lett. , 140601 (2015).[10] Ilievski, E., Medenjak, M. & Prosen, T. Quasilocalconserved operators in the isotropic Heisenberg spin-1 / Phys. Rev. Lett. , 120601 (2015).[11] Ilievski, E., Medenjak, M., Prosen, T. & Zadnik, L.Quasilocal charges in integrable lattice systems.
Journalof Statistical Mechanics: Theory and Experiment ,064008.[12] Rigol, M., Dunjko, V., Yurovsky, V. & Olshanii, M. Re-laxation in a completely integrable many-body quantumsystem: An
Ab Initio study of the dynamics of the highlyexcited states of 1d lattice hard-core bosons.
Phys. Rev.Lett. , 050405 (2007).[13] Pozsgay, B. The generalized Gibbs ensemble for Heisen-berg spin chains. Journal of Statistical Mechanics: The-ory and Experiment , P07003.[14] Fagotti, M. & Essler, F. H. Stationary behaviour of ob-servables after a quantum quench in the spin-1/2 Heisen-berg
XXZ chain.
Journal of Statistical Mechanics: The-ory and Experiment , P07012.[15] Wouters, B. et al.
Quenching the anisotropic Heisenbergchain: Exact solution and generalized Gibbs ensemblepredictions.
Phys. Rev. Lett. , 117202 (2014).[16] Pozsgay, B. et al.
Correlations after quantum quenchesin the
XXZ spin chain: Failure of the generalized Gibbsensemble.
Phys. Rev. Lett. , 117203 (2014).[17] Ilievski, E. et al.
Complete generalized Gibbs ensemblesin an interacting theory.
Phys. Rev. Lett. , 157201(2015).[18] Vidmar, L. & Rigol, M. Generalized Gibbs ensemble inintegrable lattice models.
Journal of Statistical Mechan-ics: Theory and Experiment , 064007.[19] Essler, F. H. L. & Fagotti, M. Quench dynamics andrelaxation in isolated integrable quantum spin chains.
Journal of Statistical Mechanics: Theory and Experiment , 064002.[20] Ilievski, E., Quinn, E. & Caux, J.-S. From interacting particles to equilibrium statistical ensembles.
Phys. Rev.B , 115128 (2017).[21] De Luca, A., Collura, M. & De Nardis, J. Non-equilibrium spin transport in the XXZ chain: steadyspin currents and emergence of magnetic domains. https://arxiv.org/abs/1612.07265 (2016).[22] Langen, T. et al.
Experimental observation of a general-ized Gibbs ensemble.
Science , 207–211 (2015).[23] Mourigal, M. et al.
Fractional spinon excitations in thequantum Heisenberg antiferromagnetic chain.
NaturePhysics , 435–441 (2013).[24] Jung, P., Helmes, R. W. & Rosch, A. Transport in almostintegrable models: Perturbed Heisenberg chains. Phys.Rev. Lett. , 067202 (2006).[25] Jung, P. & Rosch, A. Spin conductivity in almost inte-grable spin chains. Phys. Rev. B , 245108 (2007).[26] Jung, P. & Rosch, A. Lower bounds for the conductivitiesof correlated quantum systems. Phys. Rev. B , 245104(2007).[27] Bertini, B., Essler, F. H. L., Groha, S. & Robinson,N. J. Prethermalization and thermalization in modelswith weak integrability breaking. Phys. Rev. Lett. ,180601 (2015).[28] De Luca, A. & Rosso, A. Dynamic nuclear polarizationand the paradox of quantum thermalization.
Phys. Rev.Lett. , 080401 (2015).[29] De Luca, A., Rodr´ıguez-Arias, I., M¨uller, M. & Rosso, A.Thermalization and many-body localization in systemsunder dynamic nuclear polarization.
Phys. Rev. B ,014203 (2016).[30] Cirac, J. I., Blatt, R., Zoller, P. & Phillips, W. D. Lasercooling of trapped ions in a standing wave. Phys. Rev. A , 2668–2681 (1992).[31] Benatti, F., Nagy, A. & Narnhofer, H. Asymptotic en-tanglement and Lindblad dynamics: a perturbative ap-proach. Journal of Physics A: Mathematical and Theo-retical , 155303 (2011).[32] Li, A. C. Y., Petruccione, F. & Koch, J. Perturbativeapproach to Markovian open quantum systems. ScientificReports , 4887 (2014).[33] Petruccione, F. & Breuer, H.-P. The theory of open quan-tum systems (Oxford Univ. Press, 2002).[34] Affleck, I. & Oshikawa, M. Field-induced gap in Cu ben-zoate and other s = antiferromagnetic chains. Phys.Rev. B , 1038–1056 (1999).[35] Nojiri, H., Ajiro, Y., Asano, T. & Boucher, J. Magneticexcitation of s = 1 / New Journal of Physics , 218 (2006).[36] Kimura, S. et al. Collapse of magnetic order of the quasione-dimensional ising-like antiferromagnet BaCo V O in transverse fields. Journal of the Physical Society ofJapan , 033706 (2013).[37] Niesen, S. K. et al. Substitution effects on the tem-perature versus magnetic field phase diagrams of thequasi-one-dimensional effective Ising spin- chain systembaco v o . Phys. Rev. B , 104419 (2014).[38] Shindou, R. Quantum spin pump in s = 1 / Journal of the Physical Society of Japan , 1214–1223 (2005). [39] Genske, M. & Rosch, A. Floquet-Boltzmann equationfor periodically driven Fermi systems. Phys. Rev. A ,062108 (2015).[40] D’Alessio, L. & Rigol, M. Long-time behavior of isolatedperiodically driven interacting lattice systems. Phys. Rev.X , 041048 (2014).[41] Lazarides, A., Das, A. & Moessner, R. Equilibrium statesof generic quantum systems subject to periodic driving. Phys. Rev. E , 012110 (2014).[42] Ponte, P., Chandran, A., Papi´c, Z. & Abanin, D. A. Pe-riodically driven ergodic and many-body localized quan-tum systems. Annals of Physics , 196–204 (2015).[43] Lazarides, A., Das, A. & Moessner, R. Periodic Thermo-dynamics of Isolated Quantum Systems.
Phys. Rev. Lett. , 150401 (2014).[44] Canovi, E., Kollar, M. & Eckstein, M. Strobo-scopic prethermalization in weakly interacting periodi-cally driven systems.
Phys. Rev. E , 012130 (2016).[45] Zotos, X., Naef, F. & Prelovˇsek, P. Transport and con-servation laws. Phys. Rev. B , 11029–11032 (1997).[46] Prosen, T. Open XXZ spin chain: Nonequilibriumsteady state and a strict bound on ballistic transport.
Phys. Rev. Lett. , 217206 (2011).[47] Mierzejewski, M., Prelovˇsek, P. & Prosen, T. Breakdownof the generalized Gibbs ensemble for current-generatingquenches.
Phys. Rev. Lett. , 020602 (2014).[48] Sinova, J., Valenzuela, S. O., Wunderlich, J., Back, C. H.& Jungwirth, T. Spin hall effects.
Rev. Mod. Phys. , 1213–1260 (2015).[49] Berges, J., Bors´anyi, S. & Wetterich, C. Prethermaliza-tion. Phys. Rev. Lett. , 142002 (2004).[50] Moeckel, M. & Kehrein, S. Interaction quench in thehubbard model. Phys. Rev. Lett. , 175702 (2008).[51] Kollar, M., Wolf, F. A. & Eckstein, M. Generalized Gibbsensemble prediction of prethermalization plateaus andtheir relation to nonthermal steady states in integrablesystems.
Physical Review B , 054304 (2011).[52] Essler, F. H. L., Kehrein, S., Manmana, S. R. & Robin-son, N. J. Quench dynamics in a model with tuneableintegrability breaking. Phys. Rev. B , 165104 (2014).[53] Mierzejewski, M., Prosen, T. & Prelovˇsek, P. Approx-imate conservation laws in perturbed integrable latticemodels. Phys. Rev. B , 195121 (2015).[54] Caux, J.-S. & Essler, F. H. L. Time evolution of localobservables after quenching to an integrable model. Phys.Rev. Lett. , 257203 (2013).[55] Caux, J.-S. The quench action.