Many-body localization with synthetic gauge fields in disordered Hubbard chains
MMany-body localization with synthetic gauge fields in disordered Hubbard chains
Kuldeep Suthar, Piotr Sierant,
1, 2 and Jakub Zakrzewski
1, 3 Instytut Fizyki im. Mariana Smoluchowskiego, Uniwersytet Jagiello´nski, Łojasiewicza 11, 30-348 Krak´ow, Poland ICFO- Institut de Sciences Fotoniques, The Barcelona Institute of Science and Technology, 08860 Castelldefels (Barcelona), Spain Mark Kac Complex Systems Research Center, Uniwersytet Jagiello´nski, Krak´ow, Poland. ∗ (Dated: January 28, 2020)We analyze the localization properties of the disordered Hubbard model in the presence of a synthetic mag-netic field. An analysis of level spacing ratio shows a clear transition from ergodic to many-body localizedphase. The transition shifts to larger disorder strengths with increasing magnetic flux. Study of dynamics of lo-cal correlations and entanglement entropy indicates that charge excitations remain localized whereas spin degreeof freedom gets delocalized in the presence of the synthetic flux. This residual ergodicity is enhanced by thepresence of the magnetic field with dynamical observables suggesting incomplete localization at large disorderstrengths. Furthermore, we examine the effect of quantum statistics on the local correlations and show that thelong-time spin oscillations of a hard-core boson system are destroyed as opposed to the fermionic case. I. INTRODUCTION
The phenomenon of many-body localization (MBL) has at-tracted a significant interest in condensed matter physics overpast several years, both theoretically [1–4] and experimen-tally [5–9]. The MBL is an extension of Anderson local-ization (which describes localization of single-particle eigen-states in the presence of disorder potential) to highly excitedeigenstates of interacting many-body systems. The character-istic properties of MBL phase are Poisson eigenvalue statis-tics [10–16], an absence of thermalization [17–20], a vanish-ing transport [21–23], and the logarithmic spreading of the en-tanglement entropy [24–26]. Starting from early works, manyaspects of the topic have been examined to date. The ultracoldatomic systems are ideal platform to explore the localizationphenomena due to their ability to tune dimensionality, natureof applied disorder, atomic interactions, lattice geometry andsynthetic gauge fields. The existence of MBL phase has beenconfirmed in recent quantum gas experiments using quantumsimulators such as optical lattices [5, 6, 27, 28] and trappedions [8]. Moreover, the flux dependent mobility edge of dis-ordered chain and the significance of symmetry on the local-ization and transport properties in the presence of syntheticgauge field have been explored in experiments [29, 30].The nonergodic behaviour of disorder Hubbard chain atstrong disorder and temporal evolution of its correlation func-tions have been examined theoretically in [26, 31]. The SU(2)symmetry of the model limits the full MBL as the charge de-grees of freedom are localized but spins remain delocalizedand reveal subdiffusive dynamics [32–34]. The partial MBLis due to the separation of time scale between charge and spinsectors in the presence of correlated disorder [35] and the de-cay rate of the transport strongly depends on the density ofsingly occupied sites in the initial state [36]. The implemen-tation of synthetic gauge field in ultracold atoms and recentlyon disorder systems allows one to study the effect of time-reversal symmetry breaking on the localization properties ofthe disordered fermions [29, 37]. Furthermore, it is possible ∗ [email protected] to couple internal atomic spin states which allows for an in-terpretation in terms of an additional synthetic lattice dimen-sion [38, 39]. The gauge field in synthetic dimension leadsto flux ladders where the hopping and atomic interactions be-tween different states can be controlled in experiments [40].The presence of the random gauge fields delocalizes the in-teracting system [41] and many-body states become extendedinto a single Landau level [42]. Despite numerous studieson the topic of MBL, the effects of time-reversal symmetrybreaking on localization phenomena have been sparsely stud-ied [43] and this is the subject of the present study.In particular, we examine the spectral and dynamical prop-erties of disordered Hubbard chains in the presence of thesynthetic gauge fields. At lower disorder strengths, in er-godic regime, the breaking of time-reversal invariance (TRI)by gauge field results in spectral statistics well described byGaussian Orthogonal Ensemble (GOE) of random matricesinstead of Gaussian Unitary Ensemble (GUE) as expected forbroken TRI. This is due to a residual discrete symmetry. Onlywhen the residual reflection symmetry is broken by local fieldor asymmetric tunneling rate of spin-up and down fermions,the level statistics is characterized by GUE. The time dynam-ics of charge and spin correlations for random initial statesreveal the localization of charges and a subdiffusive decay ofspin correlations. The introduction of a synthetic flux dampsthe spin oscillations and finally delocalizes them. Further-more, the entanglement entropy confirms the delocalizationof spins in the presence of the synthetic gauge fields.The paper is organized as follows. In Sec. II we introducethe Hubbard model with synthetic gauge field and disorder.In Sec. III we analyze the effect of symmetry breaking andsynthetic flux on the spectral properties of the system. Wefurther study the local charge and spin dynamics of fermionsand hard-core bosons in Sec. IV. In Sec. V we examine thebipartite entanglement entropy. The local correlations in thepresence of spin-dependent disorder are discussed in Sec.VI.Finally, we conclude in Sec. VII. a r X i v : . [ c ond - m a t . d i s - nn ] J a n FIG. 1. Schematic visualization of the system studied (1) in whichup (down)-spin corresponds to upper (lower) rung of the ladder.
II. THE MODEL
We consider interacting spin-1/2 fermions in a quasi-onedimensional lattice. The two spin components can be inter-preted as the realization of a two-leg ladder geometry, withtwo legs corresponding to two spin states. Two componentsmay be realized as in the recent experiment [5] for K atomswith spin up state being | F, m F (cid:105) = | , − (cid:105) ≡ | ↑(cid:105) and thespin-down state | , − (cid:105) ≡ | ↓(cid:105) . The two spin componentsmay be coupled by e.g. Raman coupling realizing the Hamil-tonian [39, 44] ˆ H = ˆ H + ˆ H sb , wherein ˆ H = − (cid:88) j,σ (cid:16) J ˆ c † j,σ ˆ c j +1 ,σ + Ke − iγj ˆ c † j, ↑ ˆ c j, ↓ + H . c . (cid:17) + U (cid:88) j ˆ n j, ↑ ˆ n j, ↓ + (cid:88) j,σ (cid:15) j ˆ n j,σ . (1)Here, j and σ = ↑ , ↓ are the spatial and spin indices, J is thehopping amplitude between neighbouring lattice sites on thesame leg and K is the strength of complex hopping in thesynthetic dimension, ˆ c † j,σ (ˆ c j,σ ) creates (annihilates) fermionswith spin σ at site j , and the occupation number operator ˆ n j,σ = ˆ c † j,σ ˆ c j,σ . The local “charge” density is n j = n j, ↑ + n j, ↓ while the spin magnetization is m j = n j, ↑ − n j, ↓ . The on-siteinteraction strength of two spins is assumed to be repulsivei.e. U > and (cid:15) j is the uniform distributions of a randomspin-independent on-site potential, (cid:15) j ∈ [ − W/ , W/ with W being the disorder amplitude. The parameter γ = 2 k R a determines the synthetic magnetic flux with the flux per pla-quette of the lattice being Φ = γ/ π . Here k R is the recoilwave vector of the Raman beams and a is the lattice spacing[44]. The finite value of γ leads to complex hoppings alongthe rungs of the ladder (corresponding to flipping the spins).The system length along the synthetic dimension is two, al-though it may be increased by trapping few Fermi mixturestogether. The gauge is chosen in such a way that the Peierlsphase is along the rung and not on the legs of the ladder - com-pare Fig. 1. The above model Hamiltonian without a syntheticgauge field has been realized in quantum gas experiments inoptical lattices [5] similarly physics with the synthetic dimen-sion in clean, disorder-free experiments has also been stud-ied [40].The hopping amplitude sets the unit of energy scale, J =1 and we use periodic boundary conditions throughout thiswork. We study the system at unit filling (the total number offermions N = N ↑ + N ↓ = 1 is conserved). Note that such a < r > L = 6 L = 7 L = 8 W < r > (W - W c ) L ν < r > W c ≅ ν ≅ (a)(b) FIG. 2. Mean gap ratio (cid:104) r (cid:105) around the center of the spectra (cid:15) ≈ . as a function of disorder strength W for TRI case of K = 1 , γ = 0 of Hamiltonian (1) in the absence (a) or in the presence (b) of thesymmetry breaking local term (2) for h b = 0 . . The overlappingdistinct spectra in case (a) lead to almost Poissonian level statisticsfor any W (main figure). The inset in (a) shows that diagonalizationof a single symmetry block of the Hamiltonian allow to observe thecrossover between GOE-like behavior at low W to Poisson statisticsfor larger W . The inset in panel (b) shows the finite size scaling ofthe data yielding the critical disorder strength and exponent of er-godic to MBL phase transition. The number of disorder realizationsconsidered are , , and for the system size L = 6 , , and , respectively. situation is sometimes called a half-filling in condensed mat-ter community. In the absence of the gauge field ( K = 0) , themodel Hamiltonian preserves parity, time-reversal and SU(2)spin symmetry [45], however, pseudospin SU(2) [46] andparticle-hole symmetries are broken due to the on-site disor-der potential. The parity and SU(2) spin symmetry can beremoved by adding a local weak magnetic field ( h b ) at theedge of the chain [26]. The symmetry breaking Hamiltonianis ˆ H sb = h b (ˆ n , ↑ − ˆ n , ↓ ) . (2)Additionally, the presence of gauge field in the synthetic di-mension breaks the time-reversal symmetry of the system.We use an exact diagonalization for different system sizeand flux quanta obtaining both spectral and time evolutionproperties of systems studied. While using more sophisti-cated techniques one could extend the study to slightly largersystems, we restrict ourselves in this exploratory approach tosmall L (cid:54) system amenable to exact diagonalization. Weestimate influence of finite size effects on our results by com-parison with smaller system sizes. While the question of apossible instability of MBL phase in thermodynamic limit isdebated [47–50], we believe that our results are robust on ex-perimentally relevant time scales. < r > L = 6L = 7L = 8 h b < r > ∆ J γ = 0γ = 1 γ = 0γ = 1 FIG. 3. Mean gap ratio for different system sizes in the delocalizedregime ( W = 3 ) as a function of the strength of the symmetry break-ing local magnetic field (left) or the assumed difference between tun-neling rates for spin up and down fermions (right) which also breaksthe symmetry. The top row represents TRI case, γ = 0 where pres-ence of symmetries results in a superposition of independent spec-tra yielding Poisson-like statistics. For complex fluxes, the residualreflection symmetry between spin up and down fermions (for unitfilling) combined with TRI leads to generalized time-reversal sym-metry (see discussion in the text) leading to an apparent GOE-like be-haviour. Only breaking this symmetry by a local field or differencein tunneling rates, the GUE-like statistics corresponding to brokenTRI is observed. III. SPECTRAL STATISTICS
We consider first the statistical spectral properties of themodel. A widely used signature of localization propertiesof a many-body system is the so called gap ratio statistics[10]. From eigenvalue spectrum of the many-body Hamil-tonian one computes the average gap ratio (cid:104) r (cid:105) , where r n =min( δE n , δE n +1 ) / max( δE n , δE n +1 ) with E n being the en-ergy levels and δE n = E n +1 − E n is the spacing between twoconsecutive eigenvalues. It turns out that the mean gap ratiomay be correlated with the localization properties of eigen-states. For a localized system described by the Poisson levelstatistics, the mean gap ratio is (cid:104) r Poisson (cid:105) = 2ln2 − ≈ . ; whereas the ergodic system is described by (cid:104) r (cid:105) cor-responding to GOE (cid:104) r GOE (cid:105) = 0 . for real Hamiltonianmatrix and GUE (cid:104) r GUE (cid:105) = 0 . in the absence of (the gen-eralized) time-reversal symmetry [51]. The mean gap ratiomay depends on energy (e.g. for systems with the mobilityedge) [11, 52–54] thus it should be determined in the inter-val of energies where one expects the dynamics to be similar.Thus, from now on, we average r n over the states lying in thecenter of the spectrum for which (cid:15) n ≡ ( E n − E min ) / ( E max − E min ) ≈ . where E min and E max are the smallest andlargest eigenvalues for a given disorder realization. Finallywe average the obtained r n over disorder realizations.Such a mean gap ratio as a function of disorder strength forHamiltonian (1) is shown in Fig. 2. The on-site interactionis fixed at U = 1 . Consider first γ = 0 case - then TRI in Hamiltonian (1) is preserved. Without the symmetry breakingterm, Hamiltonian (2), the individual spectra from diagonal-izations consist of independent subsets of eigenvalues due toconserved quantities. Hence, the value of the average gap ratio (cid:104) r (cid:105) is close to the Poisson value for arbitrary disorder strength W vis. [Fig. 2(a)]. To observe the transition between ergodicand MBL phases, we use the explicit forms of generators ofSU(2) symmetry S z = 12 (cid:88) j (ˆ n j, ↑ − ˆ n j, ↓ ) , S + = ( S − ) † = (cid:88) j ˆ c † j, ↑ ˆ c j, ↓ (3) S = 12 ( S + S − + S − S + ) + ( S z ) . (4)For K = 1 , the total number of up/down fermions, N ↑ /N ↓ is not conserved and [ H , S z ] (cid:54) = 0 . However, the Hamilto-nian still commutes with S x = ( S + + S − ) / and S opera-tors as can be checked by a direct calculation. Calculating theHamiltonian matrix in a basis composed of eigenstates of S and S x and performing exact diagonalization within a singleblock of the matrix, we explicitly observe the crossover be-tween ergodic and MBL regimes as demonstrated in the insetof Fig. 2(a).Alternatively, the crossover can be observed by applyingthe symmetry breaking local Hamiltonian (2) [26] indeed,as Fig. 2(b) demonstrates at weak W , (cid:104) r (cid:105) reaches the value (cid:104) r GOE (cid:105) , while at strong disorder it approaches (cid:104) r Poisson (cid:105) .This indicates the existence of two phases at γ = 0 .In order to extract the critical disorder strength W c and ex-ponent ν of transition we use the finite-size scaling technique.We collapse the curves of different system sizes L as a func-tion of ( W − W c ) L /ν as shown in the inset of the Fig. 2(b).The critical disorder strength is the value of W at which thecurves cross or merge and indicates the ergodic to MBL tran-sition. We estimate W c ∼ . and ν ∼ . . Observe that asin the spinless fermions case (equivalent to the spin Heisen-berg chain) [11] the critical exponent breaks the Harris boundwhich is for one dimensional system [55]. In the latter caseit has been suggested that the possible explanation may be re-lated to the character of the transition resembling a Kosterlitz-Thouless transition [56] with possible important logarithmiccorrections.The remaining symmetries of Hamiltonian (1) for TRI case( γ = 0 ) may be alternatively removed by making the tunnel-ing J spin dependent, i.e. J ↑ = J ↓ + ∆ J . Both cases arevisualized in Fig. 3 - top row. At W = 3 the system is char-acterized by almost Poissonian mean gap ratio in the presenceof symmetries. The finite value of h b couples different sym-metry classes and reveals the true extended character of thesystem with (cid:104) r (cid:105) corresponding to GOE (top left panel). Avery similar behavior is observed when spin dependent tun-neling is present (top right panel). In both cases the larger thesystem the smaller value of the symmetry breaking parameteris required to fully break symmetry constrains.Let us now consider the case of nonvanishing phase γ intunnelings - compare Hamiltonian (1). Non-zero γ breaks < r > L = 6 L = 7 L = 8 < r > W < r > (a) γ = 0.1 (b) γ = 0.5 (c) γ = 1 FIG. 4. Transition from GOE-like to Poissonian behavior, i.e. fromdelocalized to localized regime for Hamiltonian (1) possessing, for γ > , a generalized time reversal symmetry for different values of γ as indicated in the figure. Observe that the crossover point betweenphases does not depend on γ . simple TRI - one might naively expect in that case the GUE-like behavior in delocalized regime. Yet, as shown in the bot-tom row in Fig. 3 without h b (or spin dependent tunneling)the mean gap ratio (cid:104) r (cid:105) ≈ . points towards GOE statistics.This is explained by the fact that while standard TRI is broken,there exists a generalized TRI in our system (1), namely TRIcombined with reflection in x − y plane (i.e. change of the spindirection). That generalized TRI leads to GOE statistics (foran excellent discussion of symmetries in different universalityclasses see [57]).With the introduction of h b or spin dependent tunneling thisgeneralized TRI symmetry is broken and, once the symmetryis fully broken, the GUE statistics is fully recovered in thedelocalized regime as shown in the bottom row of Fig. 3.Consider now the crossover from the delocalized to local-ized regime for nonzero γ , i.e., in the presence of a syntheticfield. Let us discuss first the case with h b = 0 and a general-ized TRI - compare Fig. 4. Please note that the observed GOE– Poisson transition in mean gap ratio values seems not to de-pend on γ (once it is sufficiently big to mix different symme-tries). This could in principle be quantified by comparison offinite size scaling parameters obtained - we refrain from sucha procedure as the obtained values of W c must be regarded asonly rough estimates due to small system sizes . Thus we relyfor this observation on the crossing points of curves in Fig. 4for different system sizes.On the contrary, in the presence of additional local field [i.e.adding the term (2) to the Hamiltonian (1)] the transition be-tween GUE and Poisson-like behaviour becomes dependenton γ as apparent from data presented in Fig. 5. One mightthink that this shift of the transition is related to the way inwhich the symmetry is broken, namely by a local magneticfield. Such a local perturbation might differently affect thetransition to the localized phase for different value of γ . How- < r > L = 6 L = 7 L = 8 < r > < r > W < r > (a) γ = 0.01 (b) γ = 0.05 (c) γ = 0.5 (d) γ = 1 FIG. 5. Transition from GUE-like to Poissonian behavior, i.e. fromdelocalized to localized regime for Hamiltonian (1) with additionalsymmetry breaking term (2) with h b = 0 . . Contrary to h b = 0 casepresented in Fig. 4 now the crossover point between phases signifi-cantly depends on γ . ever, we have checked, by adding random fields on all of thelattice sites that this is not the case. IV. DYNAMICS AND DENSITY CORRELATIONS
To study the dynamical properties, we choose random Fockstate | ψ (0) (cid:105) as an initial state for the temporal evolutionand examine the dynamics of the local correlations and theentanglement entropy for L = 8 . The local charge andspin correlations are C ( t ) = A (cid:80) j (cid:104) ρ j ( t ) ρ j (0) (cid:105) and S ( t ) = B (cid:80) j (cid:104) m j ( t ) m j (0) (cid:105) , where ρ j = n j − ¯ n with the averagedensity ¯ n = (cid:80) j n j /L , and A and B are normalization con-stant such that C (0) = S (0) = 1 . Recall that n j ( m j ) are thesum (difference) site occupations of up and down polarizedfermions, respectively. We assume a unit filling, ¯ n = 1 . Thememory of the initial state is lost for ergodic system and thisleads to decay of the charge and spin correlations. Interestingsituation occurs for the localized phase where (in the absenceof symmetry breaking local magnetic field or other effects de-stroying the system symmetries) the charge sector appears lo-calized but spin sector does not show localization [32, 35].The absence of fully localized phase is due to the SU(2) sym-metry of the model. The choice of a different random potentialfor up and down polarized fermions, a random magnetic fieldor a weak spin asymmetry which breaks the SU(2) spin sym-metry recovers the full MBL phase [34, 58, 59]. Moreover,the subdiffusive transport of spins is due to a singular randomdistribution of effective spin exchange interactions [33] and atlong time evolution the transport is strongly suppressed [35].The particle transport rate is exponentially small in J/W , andthe rate depends strongly on the initial state of the system. Thestates occupied with only doublons or holons exhibit full MBL C ( t ) K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 t -0.500.51 S ( t ) (a)(b) FIG. 6. Time evolutions of the local charge and spin correlationsof random Fock states for U = 1 , W = 32 and h b = 0 . (a) Thecharge remains localized in the presence of the synthetic flux. (b) Inthe absence of the flux, when K = 0 , the spin delocalizes with sub-diffusive decay of correlations. Once K is finite, the spins oscillate(red curve) and the S ( t ) of K = 0 (black curve) provides an enve-lope to these oscillations. This is apparent from the inset plot wherethe oscillations in the evolution during large time are concentrated.At γ = 0 . , the oscillations are damped and S ( t ) saturates to zero.The damping of the oscillations are pronounced as the flux value in-creases, as evident for γ = 1 . The gauge field delocalizes the spins.Here, the average over disorder is performed for realizations. phase as for such states the delocalization rate is zero [36].This is the current understanding of dynamical properties ofdisordered Hubbard model. Considering the coupling of spinup and down particles allows for further studies of ergodicitybreaking in the model.We concentrate mostly on the localized side of the tran-sition. The long time evolution of density and spin correla-tions are shown in Fig. 6. We examine the behaviour of thecorrelations at several values of the synthetic flux to observeits influence on the dynamics. For γ = 0 , K = 0 (i.e. nocoupling between the up and down polarized fermions) thecharge correlation C ( t ) saturates at a finite value for an exem-plary W = 32 value. Also at zero flux ( γ = 0 ), the couplingof spin degrees of freedom through finite hopping along thesynthetic dimension ( K = 1 ) does not alter the behaviour of C ( t ) and the charge remains localized. In the presence of afinite synthetic flux γ = 0 . and , after a transient time ofthe order of /J , the charge again seems to saturate. Yet acareful analysis of C ( t ) reveals a residual very slow algebraicdecay, C ( t ) ∝ t − α with α = 0 . for K = γ = 1 as com-pared to order of magnitude less α for γ = 0 . This may bedue to the fact that at fixed W = 32 we are slightly closerto the MBL crossover for K = γ = 1 than for K = γ = 0 case. Nevertheless, α = 0 . is significantly smaller thanthe exponent governing residual decay of imbalance at MBLtransition in system of spinless fermions [60, 61] thus we be-lieve that the charge sector of the correlated fermions remainslocalized in the presence of the synthetic gauge field for a suf-ficiently strong disorder. C ( t ) K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 t -0.500.51 S ( t ) (a)(b) FIG. 7. The local correlations of the hard-core bosons in the localizedregime W = 32 . The initial states are the random Fock states and U = 1 . (a) The charge is localized in the presence of gauge field. (b)The spins are delocalized, and the spin oscillations saturate to zero atlong time evolution. In contrast to the fermionic S ( t ) [Fig. 6(b)], thespins do not remain oscillatory at K = 1 , γ = 0 and the oscillationsare damped as time evolves. Consider now the evolution of the spin correlation S ( t ) which is shown in Fig. 6(b). In the absence of the flux, thedecay of S ( t ) at long time is suppressed and eventually spinsare localized [26, 35]. For K = 0 the intermediate timedecay of S ( t ) is monotonic (modulo fluctuations resultingfrom relatively small number of disorder realizations equalto ) while the presence of non-zero K leads to oscilla-tions in S ( t ) . The frequency of oscillations is simply K asmay be verified in the inset of Fig. 6(b). Interestingly theoscillations have the envelope given by the S ( t ) curve for K = 0 . Denoting by H K the term of H proportional to K and by H H = H − H K , we verify that [ H H , H K ] = 0 as long as γ = 0 . Thus, the evolution operator factorizes: exp( − iH t ) = exp( − iH H t ) exp( − iH K t ) and average num-ber of fermions with spin σ at site j is given by n jσ ( t ) = (cid:104) ψ (0) | e iH H t e iH K t ˆ n jσ e − iH K t e − iH H t | ψ (0) (cid:105) , (5)so that the time evolution of spin correlation function is asuperposition of dynamics determined by the H H Hamilto-nian and oscillations driven by the H K term. The long-time average of the oscillations is zero, which is in agree-ment with previous study on the localization of coupled chainswith correlated disorder [62]. In the presence of a finite flux, [ H K , H H ] (cid:54) = 0 , the oscillations are damped and the dampingrate is growing with γ . Thus, the spin sector of the systemgets delocalized by the introduction of gauge field.To examine the effect of quantum statistics, we further an-alyze the time dependence of local correlations when insteadof fermions we consider hard-core bosons in the Hamiltonian(1). The correlations of hard-core bosons for different valuesof K and γ are shown in Fig. 7. The behaviour of chargecorrelation, C ( t ) , is similar to the fermionic case and its fi-nite saturation value represents the localization of particles in C ( t ) K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 t -0.500.51 S ( t ) (a)(b) FIG. 8. Same as Fig. 6 but in the additional presence of the localmagnetic field term, [Eq.(2)] breaking the generalized TRI symme-try. The influence of coupling K and synthetic flux becomes similarto that of hard-core bosons. The strength of local field h b = 0 . . the presence of the synthetic flux. The saturation value of C ( t ) decreases with (cid:54) γ (cid:54) at K = 1 . The spin sectoris localized as S ( t ) saturates to a finite value when the twohyperfine spin states are decoupled. The coupling of thesetwo states, even for γ = 0 , leads to oscillations, which incontrast to fermionic S ( t ) , are damped to zero. A finite flux γ only weakly affects the damping and spin sector becomesfully delocalized. The stark contrast of S ( t ) between hard-core bosons and fermions for K = 1 and γ = 0 is a verynice demonstration of effects induced by quantum statisticsand different commutation properties of fermions and hard-core bosons (for ground state properties this observation goesback to [63]).Finally, we return to the fermionic system and calculatetime evolution of charge and spin correlation functions asshown in Fig. 8 for non-zero value of symmetry breakingfield h b . In this case, even at γ = 0 the commutator [ H H + H sb , H K ] (cid:54) = 0 and the oscillations of S ( t ) are damped.The effects of broken generalized TRI symmetry, while alter-ing significantly spectral properties of the system, have rela-tively small effect on time dynamics of correlation functions C ( t ) and S ( t ) for γ > . V. ENTANGLEMENT ENTROPY GROWTH
To corroborate our study of spin delocalization in the pres-ence of the synthetic gauge field we investigate the entangle-ment entropy growth in the system. The two-leg ladder systemcan be partitioned in several ways with two prime options:the first version is to cut the ladder perpendicular to the di-rection of the spatial dimension in two equal (left-right) partsand the second approach is to cut the ladder perpendicular tothe direction of synthetic dimension decoupling the hyperfinestates. This allows us to study both the entanglement betweentwo left and right ladder subsystems but also between an en- S E ( t ) K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 t S E ( t ) (a)(b) FIG. 9. (a) Left-right bipartite entanglement entropy in the local-ized regime. With no coupling of spin up and spin down polarizedfermions a typical logarithmic growth of S E is observed after initialtransients. Coupling K induces oscillations of the entanglement en-tropy. Complex flux γ leads to damping of these oscillations withalso faster entropy growth. (b) Bipartite entanglement entropy dy-namics when rungs of the ladder are cut, i.e. a partition devided upand down polarized fermions. semble of up and down polarized fermions.Regardless of the splitting the von Neumann entanglemententropy is defined in a standard way as S E ( t ) = − Tr ρ A ( t ) ln ρ A ( t ) , (6)where the two subsystems are denoted as A and B with ρ A ( t ) = Tr B | ψ ( t ) (cid:105) (cid:104) ψ ( t ) | being the reduced density matrixof subsystem A . Here Tr B is the trace over degrees of free-dom of subsystem B .In the MBL regime the entropy of entanglement of initiallyseparable state should grow logarithmically in time [24, 25]after an initial transient eventually saturating for a finite sys-tem size [25]. Such a behavior is indeed observed in ourmodel for decoupled spin chains ( K = γ = 0 ) both for thestandard left-right splitting [Fig. 9(a)] and for the up-downsplitting [Fig. 9(b)].The coupling between up and down polarized fermions in-troduces strong oscillations in the entanglement entropy, withthe same period as revealed by the spin correlation function- compare Fig. 9(a). Similarly, the entropy curve for K = 0 forms an envelope (here low lying envelope) to the oscilla-tions. Again this fact is related to vanishing commutators be-tween different kinetic energy terms in the Hamiltonian (1),as described in the previous Section. The presence of com-plex flux makes these commutators non-zero and introducesthe damping of the oscillations. The resulting entropy growthis significantly faster than in the K = γ = 0 case reflectingthe delocalization of spin sector.Similar, even more spectacular oscillations of the entangle-ment entropy are observed when the partition devides up anddown polarized fermions (i.e. the rungs in the ladder are cut).The coupling between up and down fermions ( K = 1 ) makes C ( t ) K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 t S ( t ) (a)(b) FIG. 10. (a) Charge C ( t ) and (b) spin S ( t ) correlation functions foruncorrelated disorder (cid:15) jσ ∈ [ − W/ , W/ , the chemical potential issite and spin dependent. The system size is L = 8 , disorder strength W = 32 , various values of coupling K and flux γ are considered.Both C ( t ) and S ( t ) saturate signalling localization of both chargeand spin sectors. the entropy oscillate maximally, with lower envelope given by K = 0 curve as before while the upper envelope is systemsize and disorder strength dependent. Again the non-zero fluxintroduces damping of these oscillations - on the time scalediscussed the entropy growth practically saturates at upper al-lowed value [Fig. 9(b)]. VI. SPIN-DEPENDENT DISORDER
Up till now we considered the disorder identical for bothspin components. Such a situation (albeit for quasiperi-odic disorder) is realized in experiments [5]. However, itis possible (and feasible) experimentally to consider a situ-ation when disorder is different (and uncorrelated) for bothspin components with the last term in the Hamiltonian (1)changed to (cid:80) j,σ (cid:15) j,σ ˆ n j,σ where now (cid:15) j,σ are independent ran-dom variables drawn from uniform distribution in the inter-val [ − W/ , W/ . The resulting charge and spin correlationfunctions are shown in Fig. 10. We observe a clear localiza-tion of both charge and spin sectors as C ( t ) and S ( t ) saturateto constant values after initial oscillations. The fact that addi-tion of a (small amount of) disorder in spin sector, in presenceof strong charge disorder is sufficient to fully localize the sys-tem was also observed in [59]. The saturation values are lowerfor K = 1 than for K = 0 and, interestingly, are independentof the flux γ .To complete the picture we consider the case of “spin dis-order” namely, we replace the last term in Hamiltonian (1)by (cid:80) j (cid:15) j (ˆ n j, ↑ − ˆ n j, ↓ ) so the disorder affects the spin degreeof freedom. In this case, for K = 0 , the system possessespseudo-spin SU(2) symmetry that can be utilized in studies ofits localization properties [64]. Charge and spin correlationfunctions for various values of K and γ are shown in Fig. 11. C ( t ) K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 t S ( t ) (a)(b) FIG. 11. (a) Charge C ( t ) and (b) spin S ( t ) correlation functionsfor spin disorder (cid:15) jσ = − (cid:15) j ¯ σ (where ¯ σ denotes spin opposite to σ ). The system size is L = 8 , disorder strength W = 32 , variousvalues of coupling K and flux γ are considered. The saturation of S ( t ) signals localization of the spin sector whereas the charge sectorremains delocalized. IPR K = 0, γ = 0 K = 1, γ = 0 K = 1, γ = 0.2 K = 1, γ = 1 n IPR (a)(b)
FIG. 12. Averaged over disorder realizations inverse participation ra-tio (a) for Hamiltonian (1) and (b) for disorder uncorrelated betweenspin components - leading to full MBL of both charge and spin.
We observe that the spin sector is localized as S ( t ) quicklysaturate whereas the charge correlations decay suggesting de-localization of charge degrees of freedom.Note that when disorder becomes spin-dependent the timeevolution of both charge and spin correlation function is nolonger dependent on the phase of the coupling γ between spincomponents. While for individual realizations of disorder theevolution may differ substantially, the differences average outwhen taking the disorder average. Once the spin degree offreedom becomes localized the information transfer betweenthe rungs of the ladder stops regardless of γ value. It is worthmentioning that similar independence on the phase of flux isalso observed in the entanglement entropy for uncorrelatedand “spin disorder”.One may make this observation even stronger. Consider theinverse participation ratio of n -th eigenstate | n (cid:105) : I n = 1 / (cid:88) i |(cid:104) n | b i (cid:105)| , (7)where | b i (cid:105) are basis functions (we chose a natural basis ofFock states on lattice sites). We define the averaged over dis-order inverse participation ratio as IPR( n ) = I n with theoverbar denoting disorder average. For identical random dis-order for both spin components [Hamiltonian (1)] IPR definedin this way depends on K and γ - compare Fig. 12(a). Signif-icant values of IPR indicate that localization is not complete.For the disorder uncorrelated between spin components withthe same amplitude W we observe significantly smaller IPR(indicating strong localization), moreover, IPR becomes in-dependent of γ (recall that IPR(n) are obtained after disorderaveraging). Thus not only C ( t ) or S ( t ) but also other ob-servables may be expected to be γ -independent once the spincomponent is localized.In this respect we observe a clear asymmetry between a realdimension (along the chain) and the synthetic dimension (rep-resented by spin components). The localization of the latter isessential for γ -independence of disorder averaged observableswhile the charge localization plays little role. Observe that forpure “spin disorder”, as shown in Fig. 11 the correlations for K = 1 are flux independent while charge degrees of freedomare delocalized. VII. CONCLUSIONS
We have analyzed properties of the disordered chain ofspin-1/2 fermions in the presence of the synthetic magneticfield. While the spectral properties such as the average gapratio indicate the transition to many-body localized phase, thetime dynamics suggest that MBL is realized in charge (den-sity) sector only. The presence of the synthetic magnetic fielddelocalizes the spin sector as revealed by the decay of spintime correlation function. Similarly, the entropy of entangle- ment in the system grows much faster in the presence of themagnetic field flux.Interestingly the spectral properties of the system stronglydepend on the realized symmetries providing a nice exampleof the effects due to a generalized TRI (i.e. a TRI combinedwith a discrete symmetry). In effect, despite the standard timereversal symmetry being broken the spectral properties in thedelocalized regime resemble that of GOE unless additionalsymmetry breaking terms are introduced into the model.A comparison of the dynamics of correlation functions forfermions and hard-core bosons in the absence of the flux( γ = 0 ) but when the driving term K couples up and down po-larized particles shows sensitivity of the observed phenomenato quantum statistics. For fermions, due to their commutationrelations, a kinetic energy part corresponding to the transitionbetween spin up and down fermions commutes with the restof the Hamiltonian. In effect, the exact quantum dynamics isgiven by a rapidly oscillating solution whose one envelope isgiven by K = 0 (i.e. no Rabi coupling) solution. This behav-ior is absent for hard-core bosons.Last but not least, we considered different types of randomdisorder, in particular spin dependent disorder that leads tostrong localization in the spin sector. In such a case time-dependent observables become, after disorder average, inde-pendent of flux γ . We believe that the system sizes examinedin present work are sufficient to grasp robust features of con-sidered systems on time scales relevant for experiments withultracold atoms. ACKNOWLEDGMENTS
Interesting discussions with T. Chanda, D. Delande andK. ˙Zyczkowski on subjects related to this work are acknowl-edged. The support by PL-Grid Infrastructure was importantfor this work. This research has been supported by NationalScience Centre (Poland) under projects 2015/19/B/ST2/01028(P.S.), 2018/28/T/ST2/00401 (doctoral scholarship – P.S.) and2016/21/B/ST2/01086 (K.S. and J.Z.). [1] L. Fleishman and P. W. Anderson, Phys. Rev. B , 2366(1980).[2] F. Alet and N. Laflorencie, Comptes Rendus Physique , 498(2018).[3] S. A. Parameswaran and R. Vasseur, Reports on Progress inPhysics , 082501 (2018).[4] D. A. Abanin, E. Altman, I. Bloch, and M. Serbyn, Rev. Mod.Phys. , 021001 (2019).[5] M. Schreiber, S. S. Hodgman, P. Bordia, H. P. L¨uschen, M. H.Fischer, R. Vosk, E. Altman, U. Schneider, and I. Bloch, Sci-ence , 842 (2015).[6] P. Bordia, H. P. L¨uschen, S. S. Hodgman, M. Schreiber,I. Bloch, and U. Schneider, Phys. Rev. Lett. , 140401(2016).[7] P. Bordia, H. Lschen, U. Schneider, M. Knap, and I. Bloch,Nature Physics , 907 (2016). [8] J. Smith, A. Lee, P. Richerme, B. Neyenhuis, P. W. Hess,P. Hauke, M. Heyl, D. A. Huse, and C. Monroe, Nature Physics , 907 (2016).[9] K. Xu, J.-J. Chen, Y. Zeng, Y.-R. Zhang, C. Song, W. Liu,Q. Guo, P. Zhang, D. Xu, H. Deng, K. Huang, H. Wang, X. Zhu,D. Zheng, and H. Fan, Phys. Rev. Lett. , 050507 (2018).[10] V. Oganesyan and D. A. Huse, Phys. Rev. B , 155111 (2007).[11] D. J. Luitz, N. Laflorencie, and F. Alet, Phys. Rev. B ,081103 (2015).[12] M. Serbyn and J. E. Moore, Phys. Rev. B , 041424 (2016).[13] R. Vasseur, A. J. Friedman, S. A. Parameswaran, and A. C.Potter, Phys. Rev. B , 134207 (2016).[14] A. Maksymov, P. Sierant, and J. Zakrzewski, Phys. Rev. B ,224202 (2019).[15] P. Sierant and J. Zakrzewski, Phys. Rev. B , 104205 (2019).[16] P. Sierant and J. Zakrzewski, “Model of level statistics fordisordered interacting quantum many-body systems,” (2019), arXiv:1907.10336 [cond-mat.dis-nn].[17] A. Pal and D. A. Huse, Phys. Rev. B , 174411 (2010).[18] M. Serbyn, Z. Papi´c, and D. A. Abanin, Phys. Rev. Lett. ,127201 (2013).[19] D. J. Luitz, N. Laflorencie, and F. Alet, Phys. Rev. B ,060201 (2016).[20] P. Sierant, D. Delande, and J. Zakrzewski, Phys. Rev. A ,021601 (2017).[21] K. Agarwal, S. Gopalakrishnan, M. Knap, M. M¨uller, andE. Demler, Phys. Rev. Lett. , 160401 (2015).[22] Y. Bar Lev, G. Cohen, and D. R. Reichman, Phys. Rev. Lett. , 100601 (2015).[23] M. Kozarzewski, P. Prelovˇsek, and M. Mierzejewski, Phys.Rev. B , 235151 (2016).[24] J. H. Bardarson, F. Pollmann, and J. E. Moore, Phys. Rev. Lett. , 017202 (2012).[25] M. Serbyn, Z. Papi´c, and D. A. Abanin, Phys. Rev. Lett. ,260601 (2013).[26] R. Mondaini and M. Rigol, Phys. Rev. A , 041601 (2015).[27] J.-y. Choi, S. Hild, J. Zeiher, P. Schauß, A. Rubio-Abadal,T. Yefsah, V. Khemani, D. A. Huse, I. Bloch, and C. Gross,Science , 1547 (2016).[28] H. P. L¨uschen, P. Bordia, S. Scherg, F. Alet, E. Altman,U. Schneider, and I. Bloch, Phys. Rev. Lett. , 260401(2017).[29] F. A. An, E. J. Meier, and B. Gadway, Phys. Rev. X , 031045(2018).[30] C. Hainaut, I. Manai, J.-F. Clment, J. C. Garreau, P. Szriftgiser,G. Lemari, N. Cherroret, D. Delande, and R. Chicireanu, Na-ture Communications , 1382 (2018).[31] P. Prelovˇsek, Phys. Rev. B , 144204 (2016).[32] P. Prelovˇsek, O. S. Bariˇsi´c, and M. ˇZnidariˇc, Phys. Rev. B ,241104 (2016).[33] M. Kozarzewski, P. Prelovˇsek, and M. Mierzejewski, Phys.Rev. Lett. , 246602 (2018).[34] M. ´Sroda, P. Prelovˇsek, and M. Mierzejewski, Phys. Rev. B ,121110 (2019).[35] J. Zakrzewski and D. Delande, Phys. Rev. B , 014203 (2018).[36] I. V. Protopopov and D. A. Abanin, Phys. Rev. B , 115111(2019).[37] Y.-J. Lin, R. L. Compton, K. Jimnez-Garca, J. V. Porto, andI. B. Spielman, Nature , 628 (2009).[38] O. Boada, A. Celi, J. I. Latorre, and M. Lewenstein, Phys. Rev.Lett. , 133001 (2012).[39] A. Celi, P. Massignan, J. Ruseckas, N. Goldman, I. B. Spiel-man, G. Juzeli¯unas, and M. Lewenstein, Phys. Rev. Lett. ,043001 (2014).[40] L. F. Livi, G. Cappellini, M. Diem, L. Franchi, C. Clivati,M. Frittelli, F. Levi, D. Calonico, J. Catani, M. Inguscio, andL. Fallani, Phys. Rev. Lett. , 220401 (2016).[41] C. Cheng and R. Mondaini, Phys. Rev. A , 053610 (2016). [42] S. D. Geraedts and R. N. Bhatt, Phys. Rev. B , 054303 (2017).[43] J. Major, M. Płodzie´n, O. Dutta, and J. Zakrzewski, Phys. Rev.A , 033620 (2017).[44] D. Suszalski and J. Zakrzewski, Phys. Rev. A , 033602(2016).[45] F. H. L. Essler, H. Frahm, F. Ghmann, A. Klmper, and V. E. Ko-repin, The One-Dimensional Hubbard Model (Cambridge Uni-versity Press, 2005).[46] S. Zhang, Phys. Rev. Lett. , 120 (1990).[47] J. ˇSuntajs, J. Bonˇca, T. Prosen, and L. Vidmar, “Quantum chaoschallenges many-body localization,” (2019), arXiv:1905.06345[cond-mat.str-el].[48] D. A. Abanin, J. H. Bardarson, G. D. Tomasi, S. Gopalakr-ishnan, V. Khemani, S. A. Parameswaran, F. Pollmann, A. C.Potter, M. Serbyn, and R. Vasseur, “Distinguishing localiza-tion from chaos: challenges in finite-size systems,” (2019),arXiv:1911.04501 [cond-mat.str-el].[49] P. Sierant, D. Delande, and J. Zakrzewski, “Thouless timeanalysis of anderson and many-body localization transitions,”(2019), arXiv:1911.06221 [cond-mat.dis-nn].[50] R. K. Panda, A. Scardicchio, M. Schulz, S. R. Taylor, andM. ˇZnidariˇc, “Can we study the many-body localisation tran-sition?” (2019), arXiv:1911.07882 [cond-mat.dis-nn].[51] Y. Y. Atas, E. Bogomolny, O. Giraud, and G. Roux, Phys. Rev.Lett. , 084101 (2013).[52] S. Nag and A. Garg, Phys. Rev. B , 060203 (2017).[53] S.-H. Lin, B. Sbierski, F. Dorfner, C. Karrasch, and F. Heidrich-Meisner, SciPost Phys. , 002 (2018).[54] P. Sierant and J. Zakrzewski, New Journal of Physics ,043032 (2018).[55] A. Chandran, C. R. Laumann, and V. Oganesyan, (2015),arXiv:1509.04285 [cond-mat.dis-nn].[56] A. Goremykina, R. Vasseur, and M. Serbyn, Phys. Rev. Lett. , 040601 (2019).[57] F. Haake, Quantum Signatures of Chaos (Springer, 2013).[58] G. Lemut, M. Mierzejewski, and J. Bonˇca, Phys. Rev. Lett. , 246601 (2017).[59] B. Leipner-Johns and R. Wortis, Phys. Rev. B , 125132(2019).[60] E. V. H. Doggen, F. Schindler, K. S. Tikhonov, A. D. Mirlin,T. Neupert, D. G. Polyakov, and I. V. Gornyi, Phys. Rev. B ,174202 (2018).[61] T. Chanda, P. Sierant, and J. Zakrzewski, (2019),arXiv:1908.06524 [cond-mat.stat-mech].[62] Y. Zhao, S. Ahmed, and J. Sirker, Phys. Rev. B , 235152(2017).[63] F. Cr´epin, N. Laflorencie, G. Roux, and P. Simon, Phys. Rev.B , 054517 (2011).[64] X. Yu, D. Luo, and B. K. Clark, Phys. Rev. B98