Persistent currents in Bose-Bose mixtures after an interspecies interaction quench
aa r X i v : . [ qu a n t - ph ] F e b Persistent currents in Bose-Bose mixtures after an interspeciesinteraction quench
D. Spehner , L. Morales-Molina , and S. A. Reyes Departamento de Ingenier´ıa Matem´atica, Universidad de Concepci´on, Concepci´on, Chile Univ. Grenoble Alpes, CNRS, Institut Fourier and LPMMC, F-38000 Grenoble, France Facultad de F´ısica, Pontificia Universidad Cat´olica de Chile, Casilla 306, Santiago, Chile
February 17, 2021
Abstract
We study the persistent currents and interspecies entanglement generation in a Bose-Bose mixtureformed by two atomic gases (hereafter labelled by the letters A and B ) trapped in a one-dimensionalring lattice potential with an artificial gauge field after a sudden quench from zero to strong interactionsbetween the two gases. Assuming that the strength of these interactions is much larger than the singlespecies energies and that the gas A is initially in the Mott-insulator regime, we show that the currentof the gas B is reduced with respect to its value prior the interaction quench. Averaging fast oscillationsout, the relative decrease of this current is independent of the initial visibility and Peierls phase of thegas B and behaves quadratically with the visibility of the gas A . The second R´enyi entropy of thereduced state measuring the amount of entanglement between the two gases is found to scale linearlywith the number of sites (volume law) and to be proportional to the relative decrease of the current. The manifestations of quantum coherence and entanglement in many-body systems is one of the mostchallenging problems in condensed matter physics and quantum technology. The recent experimentalrealizations with trapped ultracold atoms of analogs of electronic mesoscopic systems such as supercon-ducting quantum interference devices [1, 2, 3] has opened new perspectives in the study of matter-waveinterferences. Atomtronics focuses on the design of such atomic quantum devices, characterized by tunableparameters and low decoherence, and their applications to fundamental research and technology [4, 5].One of the striking manifestation of quantum interferences is the Aharonov-Bohm effect. In a ringpierced by a magnetic flux, this gives rise to persistent currents varying periodically with the flux, whichwere observed long ago in superconductors and normal metals [6, 7, 8, 9]. In atomic Bose-Einstein con-densates (BECs) trapped in a ring-shaped potential, an artificial gauge field can be induced by laserfields or by the Coriolis force in the presence of a rotating potential barrier [10]. This system provides anovel platform for studying persistent currents [11, 12, 13, 14]. It encompasses a rich variety of quantumphenomena, including the formation of macroscopic superpositions of clockwise and anticlockwise flowingstates at specific values of the magnetic flux when rotation invariance is broken by a localized potentialbarrier [15, 16, 17, 18, 19]. The influence of such a barrier on the current amplitude has been studied inRefs. [13, 14] for all strengths of the repulsive atomic interactions and in Ref. [20] for attractive interactions.The quantum interference effects become more complex when two atomic species are involved. It isknown that mixtures of two condensates can lead to the generation of exotic phases and to the formationof quantum droplets, among others phenomena [21, 22, 23, 24, 25, 26, 27, 28, 29]. It has been conjectured1hat the phase coherences of the two atomic gases in the mixing process may play a crucial role in theseemerging phenomena [24, 25]. A number of theoretical and experimental works have focused on quantumphases in optical lattices, in particular because of the analogy with condensed matter systems [30, 31, 32,33, 34, 35, 36, 37, 38]. For repulsive interspecies interactions, dynamical phase separation effects have beeninvestigated thoroughly (see e.g. [39] and references therein for a many-body approach). More recently,Bose-Bose mixtures consisting of a single atom interacting with a gas of atoms of a different species haveattracted a lot of attention due to the formation of polarons in such systems [40, 41, 42, 43, 44].A particular emphasis in the study of multi-component quantum gases concerns out-of-equilibriumdynamics following a sudden quench. Such quench dynamics can be investigated experimentally thanksto the high level of control on the trapping potential and the ability to tune (and even change the signof) the intra- and inter-species interactions using Feshbach resonances [45, 46, 47]. On the theoreticalside, special attention has been devoted to interaction quenches (inter-species interactions are suddenlyswitched on) [48, 49, 39, 43]. In these works, the generation of entanglement in the quench dynamics andits role in the observed phenomena have been investigated.Coming back to persistent currents in BECs, it is natural to ask about the effect on the current of thefirst atomic species of the presence of a second atomic gas trapped in the same ring. In order to change thiscurrent, the two quantum gases must interact and be entangled with each other. One may then wonderwhether the variation of the current provides in some way a measure of the amount of entanglementbetween the two gases. Another natural question concerns the dependence of the persistent current onthe phase coherence of the two gases before they start to interact.In this article, we investigate these questions by considering a binary mixture containing two bosonicatomic species A and B (which may either correspond to different types of atoms or to identical atoms indifferent internal states) trapped in a 1D-ring lattice potential in the presence of an artificial gauge field.Assuming that the two species are initially decoupled and that the gas A is in the Mott-insulator (MI)regime, we calculate analytically and numerically the time evolution of the current of B -atoms and thegeneration of entanglement between the A and B -condensates after a sudden quench from zero to stronginter-species interactions. We quantify the amount of entanglement using the Schmidt number [50] shiftedby one, K AB = tr A [(tr B | ψ AB ih ψ AB | ) ] − −
1, where | ψ AB i is the wavefunction of the binary mixture afterthe interaction quench. For small K AB this quantity is nearly equal to the 2-R´enyi entropy of entanglement S (2) AB = ln(1 + K AB ) [51], which is a measure of entanglement having similar properties as the entanglementof formation (von Neumann entropy of the reduced state). It has been shown recently that S (2) AB can bedetermined from statistical correlations between random measurements [52, 53, 54].We show that the B -current after the quench is reduced due to interactions with the gas A . Averagingout the fast oscillations with frequency equal to the inter-species interaction energy, we prove that therelative variation hJ B i of the B -current before and after the quench and the shifted Schmidt number hK AB i are proportional to each other and follow some universal laws. In particular, hJ B i is independentof the visibility of the gas B before the quench, being the same when the latter is in the MI or in thesuperfluid (SF) regimes, and does also not depend on the Peierls phase of the B -atoms. On the otherhand, hJ B i behaves quadratically with the initial visibility of the gas A . Furthermore, hK AB i and S (2) AB are linear in the number of lattice sites L (volume law of entanglement) and hK AB i /L = β hJ B i , with aproportionality factor β depending only on the filling factor ν B , the Peierls phase, and the initial visibilityof the gas B when this gas is in the MI regime, whereas when it is in the SF regime β depends on thefilling factor ν B only. We show that the generation of entanglement in the quench dynamics and itsuniversal relation with the reduction of B -current comes from particle-hole excitations in the gas A whichslow down the flow of B -atoms and are coupled to site-dependent wavefunctions of the gas B . We arguethat these results could be used for determining the amount of entanglement between the two gases fromquantities (atomic current and visibility of interference fringes) that can be measured experimentally (seee.g. [55, 56, 57, 12] for the observation of the supercurrent of a single-component gas in ring-shape trapping2otentials).Let us comment on the methods used to obtain the aforementioned results. Our analytical calculationsrely on (i) a perturbative expansion of the initial ground state and (ii) small and intermediate timeapproximations for the time propagator, valid when the inter-species interaction strength is much largerthan the tunneling and intra-species interaction energies. We point out that the method employed todetermine the Schmidt number in Sec. 7 is original and could be useful in other contexts. It consists inexpressing K AB ( t ) when the gas A is in the MI regime in terms of an effective propagator acting on theother gas (see Appendix B); by diagonalizing perturbatively the Hamiltonian in this propagator for largeinter-species interactions, we are able to determine K AB ( t ) at times of the order of the inverse tunnelingand intra-species interaction energies of the gas B . Our numerical simulations, in turn, rely on exactdiagonalization. Although they are carried out for small lattice sizes, these simulations corroborate theanalytical results which apply to much larger atom numbers and lattice sizes. This implies that finitesize effects are not relevant for the physical effects we are interested in, thus enabling their applications tomesoscopic systems of arbitrary size.The paper is organized as follows. We introduce our model of a Bose-Bose mixture in a 1D-ring latticein Sect. 2. Section 3 is devoted to a brief discussion on the visibility of a single species before the interactionquench and its behavior as function of the Peierls phase for finite lattice sizes. The calculation of the B -current after the interaction quench is performed in Sec. 4. We determine in Sec. 5 the entanglementgeneration between the two gases and show there that the time-averaged shifted Schmidt number hK AB i is proportional to L and to hJ B i , assuming an averaging time much larger than the inverse inter-speciesinteraction energy and much smaller than the inverse single species energies. In Sec. 6 we show thatquantum superpositions in the post-quench total wavefunction involving particle-hole excitations in thegas A are at the origin of the entanglement. In Sec. 7, we evaluate the Schmidt number at larger timesand show that its average hK AB i t remains constants over a large time period, suggesting a convergenceof hK AB i t at large times t to universal values which are determined analytically. Concluding remarks aregiven in Sec. 8. Some technical details on the analytical calculations are presented in the two appendices. Let us consider two atomic gases A and B trapped in the same 1D-ring lattice potential with L sites. Wedenote by ˆ a † j , ˆ a j , and ˆ n Aj = ˆ a † j ˆ a j the creation, annihilation, and number operators at site j for atoms ofthe gas A , and by ˆ b † j , ˆ b j , and ˆ n Bj the corresponding operators for the gas B . The two gases are describedby Bose-Hubbard Hamiltonians, e.g. for the gas A ˆ H A = − J A X j (e i φ A ˆ a † j +1 ˆ a j + h . c . ) + U A X j ˆ n Aj (ˆ n Aj −
1) = ˆ K A + ˆ H int A , (1)where J A and U A are the tunneling and interaction energy strengths and φ A is the Peierls phase associatedto an artificial gauge field. The sums in (1) run over all lattice sites j = 0, . . . , L − H B of the gas B is given similarly in terms of J B , U B , and φ B . We consider repulsive intra-species interactions, U A , U B >
0, and assume fixed totalatom numbers N A , N B for each species with integer filling factors ν A = N A /L , ν B = N B /L . The twogases are initially uncorrelated and in their ground states (GSs) | ψ A i and | ψ B i . At time t = 0, attractiveinteractions between the two species are suddenly switched on (see the left panel in Fig. 1). The mixturesubsequently evolves according to the Hamiltonian ( V > H AB = ˆ H A ⊗ ˆ B + ˆ A ⊗ ˆ H B + ˆ H int AB , ˆ H int AB = − V X j ˆ n Aj ˆ n Bj . (2)3e assume that inter-species interactions are much larger than the single species energies and that thegas A is initially in the MI regime, i.e., V ≫ U A , U B , J B , λ A = J A U A ≪ . (3)The gas B can be either in the MI or the SF regime; in fact both cases λ B ≪ λ B & B -atoms (called in whatfollows the B -current) and of the amount of A - B entanglement after the interaction quench. The B -current is given at time t ≥ I B ( t ) = tr B [ˆ ρ B ( t ) ˆ I B ] with ˆ ρ B ( t ) = tr A | ψ AB ( t ) ih ψ AB ( t ) | the reduceddensity matrix of the gas B , | ψ AB ( t ) i the time-evolved state of the mixture with the Hamiltonian (2), andˆ I B the current operator ˆ I B = 12 LJ B ∂ ˆ H B ∂φ B = 12i L X j (cid:0) ˆ b † j +1 ˆ b j e i φ B − h . c . (cid:1) (4)(here ~ = 1). We will focus in the following on the relative variation of B -current and its average in thetime interval [0 , t ], defined by J B ( t ) = I B (0) − I B ( t ) I B (0) , hJ B i t = 1 t Z t d t ′ J B ( t ′ ) . (5)When the tunelling energy of the gas A vanishes (i.e. when J A = 0), the B -current is time-independentand thus J B ( t ) = 0. Actually, then the initial state | ψ A i| ψ B i of the mixture has ν A atoms A per site andthe coupling ˆ H int AB acts on it as − V ν A N B , where the total number N B = P j ˆ n Bj of B -atoms is a c -number;hence | ψ AB ( t ) i is equal to | ψ A i| ψ B i at all times t up to irrelevant phases. In contrast, when J A > A has particle-hole excitations [58], thus the coupling entangles the two gases andthe B -current is modified, J B ( t ) = 0.The time evolution of the relative B -current variation, obtained numerically from exact diagonalizationsof ˆ H A , ˆ H B , and ˆ H AB , is shown in the right panel of Fig. 1. We use a second-order Suzuki-Trotterdecomposition method with time steps as small as ∆ t = 0 . J B ( t ) presents fastoscillations of period 2 π/V superimposed on a complex pattern with oscillations with larger periods. Thetime-averaged version of J B ( t ) displays two plateaus: the first one at times V − ≪ t ≪ U − A , U − B , J − B ,the second one at times of the order of the single species inverse energies. We will show analytically inSec. 4 that the value of the first plateau does not depend on the initial state of the gas B . Note that thisuniversal value holds for the relative variation of the B -current defined in (5); in contrast, the B -currentstrongly depends on the initial state of the gas B , being much smaller when B is in the MI regime thanwhen it is in the SF regime.Before calculating the B -current, we study the visibility of a single atomic species in a ring lattice witha gauge field before the interaction quench. A good indicator of the degree of phase coherence of a single species (say B ) is the visibility. The latteris estimated experimentally by measuring the interference pattern after a free expansion of the gas. Itis defined as V B = ( S max − S min ) / ( S max + S min ), S max and S min being the maximum and minimum ofthe momentum distribution S ( q ) = P i,j e i q ( i − j ) h ψ B | ˆ b † i ˆ b j | ψ B i . For J B = 0, the GS of the gas B is thephase-incoherent MI state | ψ B i = | ψ MI i having ν B atoms per site [58]. Then V B = 0 (no interferencefringes) and the B -current h ψ B | ˆ I B | ψ B i vanishes. Increasing the energy ratio λ B = J B /U B , the visibilityand current increase, with V B reaching its maximum V B = 1 in the SF limit λ B ≫ time J B / λ A Figure 1:
Left panel:
Sketch of moving atoms on a 1D-ring lattice with two interacting species A and B represented by green and red balls, in the presence of Peierls phases φ A and φ B . Right panel:
Timeevolution after the interaction quench of the relative B -current variation J B ( t ) (red curve) and its time-average (black curve), Eq. (5), divided by λ A for a lattice with L = 4 sites, N A = N B = 4 atoms of eachspecies, V = 200 U B , J B = U A = U B , J A = 0 . U B , φ A = φ B = π/
10 (from numerical calculations). Timeis in units of U − B . Dashed horizontal segment: universal value of Eq. (19). Inset: amplification of theblue box shown in the figure with formula (18) displayed in dashed line.It has been shown both theoretically and experimentally [58, 59] that the visibility of a single BECtrapped in an infinite gauge-free lattice potential behaves linearly with λ B when λ B .
1. However, wewill need in the sequel the visibility for finite lattice sizes in the presence of a gauge field and cannot relyon the results of these references. In the MI regime λ B ≪
1, one obtains by expanding perturbatively theGS | ψ B i and the momentum distribution up to second order in λ B that (see Appendix A) V B = 4( ν B + 1) λ B v L ( φ B ) (cid:2) − (4 ν B + 1) λ B w L ( φ B ) (cid:3) + O ( λ B ) (6)with v L ( φ B ) = cos( φ B − ℓφ ) if L is even (cid:0) cos( φ B − ℓφ ) + cos (cid:0) k φ − φ B (cid:1)(cid:1) if L is odd (7) w L ( φ B ) = L is even − cos( φ B − ℓφ ) + cos (cid:0) k φ − φ B (cid:1) if L is odd, (8)where φ = 2 π/L is the lattice flux, ℓ = E [ φ B /φ + 1 /
2] the angular momentum of the SF state, and k = E [ φ B /φ ] (here E is the integer part). Note that Gauge invariance implies that V B is periodic in thePeierls phase φ B with period φ . In the limit L → ∞ one finds that v L ( φ B ) , w L ( φ B ) → φ B ,so that V B becomes phase independent and one recovers the result of Refs. [58, 59]. It is easy to showfrom (6)-(8) that for finite lattice sizes L , the visibility V B reaches its minimum when φ B is equal to ahalf-integer value of φ (see Appendix A).In the opposite limit of weak interactions U B ≪ J B (SF regime λ B ≫ φ B is not close to ahalf-integer value of φ the GS of the Bose-Hubbard Hamiltonian can be approximated by the SF state(GS of the tunelling Hamiltonian). Since the visibility in the latter state is equal to 1 for all φ B ’s, itfollows that V B = 1 + O ( λ − B ). The phases φ B = ( ℓ + 1 / φ , ℓ = 0 , ± , . . . , are the points where theparabola giving the energies of the SF states with angular momenta ℓ and ℓ + 1 intersect. For such phases,quantum fluctuations of angular momentum produce fluctuations in the speed of rotation of the gas whichare expected to blur the interference pattern, henceforth reducing the visibility as compared to its value5 φ/π V L=3L=4L=5 -0,5 0 0,5 φ/π -0,200,2 C u rr e n t L=3L=4L=5 1/5 1/4 1/3-1/5-1/4 0-1/3
Figure 2:
Left panel:
Visibility of a Bose gas with a single species B trapped in a 1D-ring lattice potentialas function of the Peierls phase φ B for λ B = J B /U B = 0 . N B = L = 3 ,
4, and 5 (from numericalcalculations).
Right panel: B -current I B (0) (in arbitrary units) of the gas B prior the interaction quenchas function of φ B for the same parameters. The value of λ B is close to the transition between the MI andSF regimes.when φ B is not close to a half-integer value of φ . This is confirmed by numerical simulations for L = 3 , V B is minimum atPeierls phases equal to half-integer values of φ for any energy ratio λ B . This is supported by numericalcalculations. For instance, one observes in the left panel of Fig. 2 that V B presents a well pronouncedminimum at φ B = φ / π/L when L = 3 ,
4, and 5 for λ B = 0 .
2. Note that the phase points ( ℓ + 1 / φ are singled out by the behavior of the GS persistent current I B (0) = h ψ B | ˆ I B | ψ B i of the Bose gas as functionof φ B . The latter is periodic with period φ and drops rapidly and changes sign at these phase points,independently of strength of the repulsive interactions (i.e., of λ B ) [15, 20]. For reference, the behavior of I B (0) as function of φ B is shown in the right panel in Fig. 2.Note that in Fig. 2 the depth of the minima of V B become less pronounced as L is increased, inagreement with the expectation that the visibility is phase independent for an infinite ring. Let us stressthat the numerical simulations carried out in this work consider small lattice sizes. This enables us toobserve the dependence of physical quantities like the visibility on the Peierls phase, whose effect vanishesat large sizes. B -current We can calculate the B -current at times t ≪ U − A , U − B , J − B as follows. Using the invariance of thereduced density matrix ˆ ρ B ( t ) with respect to translations by one lattice site and introducing the eigenstates | n B i = | n B , . . . , n BL − i of ˆ n B = (ˆ n B , . . . , ˆ n BL − ) (Fock states), one has I B ( t ) = Im X n B q n B ( n B + 1) e i φ B h n B | ˆ ρ B ( t ) | n B + 1 − i , (9)where 1 i denotes the vector with components δ ij , j = 0 , . . . , L −
1, and the sum runs over all n B ∈ N L , P j n Bj = N B . Due to our small time hypothesis, hopping and intra-species interactions can be neglectedin the Hamiltonian (2) and the dynamics is solely governed by inter-species interactions. Thus | ψ AB ( t ) i =e i tV ˆ n A · ˆ n B | ψ A i| ψ B i . Plugging this expression into (9), the current can be cast as I B ( t ) = I B (0) h e − i tV ∆ˆ n A i ψ A , ∆ˆ n A = ˆ n A − ˆ n A , (10)6 V A < J B > / V A φ A = φ B = π/20, J B /U B = 0.1 φ A = π/20, φ B = π/3, J B /U B = 10 φ A = - φ B = π/3, J B /U B = 10 φ A = φ B = π/3 J B /U B = 0.1 (a) λ V V A < J B > / λ A L=3L=4L=5 (b)
Figure 3: (a)
Time averaged current variation hJ B i t divided by the square visibility V A vs. V A fordifferent values of λ B = J B /U B and Peierls phases, from numerical calculations with L = N A = N B = 4.The average is up to time t = 0 . /U B . Triangles downward: V = 100 U B , U A = U B , J B = 0 . U B , φ A = φ B = π/
20. Squares: V = 1000 U B , U A = J B = 10 U B , φ A = π/ φ B = π/
3. Triangles upward: V = 1000 U B , U A = J B = 10 U B , φ A = − φ B = π/
3. Diamonds: V = 100 U B , U A = U B , J B = 0 . U B , φ A = φ B = π/
3. In all cases J A is varied in order to change V A . The horizontal segments correspond tothe values in the RHS of Eq.(19). (b) hJ B i t divided by λ A for V = 100 U B , U A = U B , J B = 0 . U B , J A variable, φ A = φ B = π/
20, and different site numbers L = N A = N B = 3 , Inset: single species visibility vs λ = J/U for the samevalues of L and N , with dashed lines showing the linear term in Eq. (6).where I B (0) is the current in the GS of ˆ H B (obtained by substituting ˆ ρ B ( t ) by | ψ B ih ψ B | in (9)), h·i ψ A = h ψ A | · | ψ A i is the quantum expectation in | ψ A i , and we have used that h e − i tV ∆ˆ n A i ψ A is real, as will beshown below. The relative current variation reads J B ( t ) = 1 − (cid:10) e − i tV ∆ˆ n A (cid:11) ψ A , t ≪ U − A , U − B , J − B . (11)It follows from the Cauchy-Schwarz inequality that J B ( t ) ≥
0. Hence the effect of the inter-speciescoupling at time t > J B ( t ) does not depend on the initial GS of the gas B , being the same when B is initially in theSF or MI regimes, and is also independent of the Peierls phase φ B of the gas B . This independence of therelative variation of B -current after the interaction quench from φ B and λ B contrasts with the behavior ofthe B -current I B (0) before the quench. The latter oscillates with φ B (see the right panel in Fig. 2) withan amplitude which strongly depends on the visibility V B , i.e., on λ B (see the discussion at the beginningof Sec. 3). The independence of J B ( t ) from λ B and φ B is confirmed by the numerical results displayedin Fig. 3(a). Indeed, one sees in this figure that the time-averaged values hJ B i t of J B ( t ) for λ B = 10, φ B = π/ λ B = 0 . φ B = π/
20 are very close from each other, for all values of the visibility V A of the gas A . As we shall explain below, the observed dependence of hJ B i t / V A on the Peierls phase φ A ofthe gas A for small V A originates from the variation of V A with φ A predicted by Eq. (6); this is thus afinite-size effect.We next calculate the A -expectation in (11), assuming that the gas A is in the MI regime ( λ A ≪ λ A . (Let us note that although this is true forall values of λ A when φ A = 0 by symmetry of ˆ H A under j L − j mod ( L ), this is not the case for φ A = 0 and large λ A .) The GS of ˆ H A can be evaluated by treating the tunneling Hamiltonian ˆ K A in(1) perturbatively. The unperturbed Hamiltonian ˆ H int A has eigenvectors | n A i , eigenenergies E (0) A ( n A ), and7S | ψ MI i = | ν A i with ν A = ( ν A , . . . , ν A ). Noting that the energy to create a particle-hole excitation is E (0) A ( ν A ± i +1 ∓ i ) − E (0) A ( ν A ) = U A , one finds | ψ A i = | ψ MI i + λ A | ψ (1) A i + λ A | ψ (2) A i + O ( λ A ) , (12)where the first-order correction is | ψ (1) A i = − J A ˆ K A | ψ MI i = X j (cid:0) e i φ A ˆ a † j +1 ˆ a j + h . c . (cid:1) | ψ MI i (13)and the second-order correction satisfies (see appendix A) h ψ MI | ψ (2) A i = − α A L with α A = ν A ( ν A + 1) . (14)A simple calculation yields (cid:10) e − i tV ∆ˆ n A (cid:11) ψ (1) A ≡ h ψ (1) A | e − i tV ∆ˆ n A | ψ (1) A i = 2 α A X j cos( tV ∆ j,j +101 ) , (15)where ± ∆ j,j +101 stands for the eigenvalue of ∆ˆ n A in the particle-hole excitation state | ϕ A,j, ± i = | ν A ± j +1 ∓ j i . (16)Since ∆ j,j +101 = δ j +1 , − δ j, − δ j +1 , + δ j, equals 2 if j = 0, ∆ j,j +101 = − j = 0 or L −
1, and ∆ j,j +101 = 0otherwise, one finds (cid:10) e − i tV ∆ˆ n A (cid:11) ψ (1) A = 2 α A [ L − tV ) + cos(2 tV ))] . (17)Using (11), (14), (17), e − i tV ∆ˆ n A | ψ MI i = | ψ MI i , and h ψ (1) A | ψ MI i = 0, this gives J B ( t ) = 2 α A λ A [3 − tV ) − cos(2 tV )] + O ( λ A ) . (18)A good agreement between Eq. (18) and the numerical result at times t ≪ /U A , /U B , /J B is seen inFig. 1 (inset, right panel). The cosines disappear upon averaging up to a time t much larger than 1 /V ,with t ≪ /U A , /U B , /J B . Plugging the expression of the visibility V A analogous to (6), one gets hJ B i t = 6 α A λ A = 3 ν A ν A + 1) V A v L ( φ A ) (19)with errors of order λ A . We see that hJ B i t is independent of the initial coherence, number of atoms, andPeierls phase of the gas B and of the size L of the ring . Moreover, it depends quadratically on the initialvisibility of the gas A . Recalling that v L ( φ A ) → L → ∞ , the dependence on the Peierls phase φ A ofthe expression in the last member of Eq. (19) disappears in the infinite size limit.Comparing Eq. (19) with the numerical results shown in Fig. 3(a), which display hJ B i t / V A as functionof V A for two different values of λ B corresponding to the gas B being initially in the MI and SF regimesand for different phases, a good agreement is observed when V A ≪
1, the values calculated numericallybeing about 5% above the prediction of Eq. (19). Clear deviations from Eq. (19) show up for V A & . λ A ≪ B -current when the two gases are coupled at time t > B -atoms of particle-hole excitations in the gas A (Notethat the MI state of the gas A has no effect on the B -current since the attractive potential produced oneach B -atom by its coupling with the gas A in the MI state is site-independent and equal to − V ν A ; in8ther words, if the two gases are in the state | ψ MI i| ψ B i then their coupling energy is independent of thedistribution of the B -atoms on the lattice and is equal to − V ν A N B , as noted in Sec. 2.) Consider e.g. agas A in the particle-hole state | ϕ A,j, + i . Then the transfer of a B -atom from site j + 1 to site j has acoupling energy cost of 2 V (in fact, the coupling energy increases by V ( ν A + 1) when removing the B -atomfrom site j + 1 and decreases by V ( ν A −
1) when adding it on site j ). Similarly, if the gas A is in theparticle-hole state | ϕ A,j, − i , the transfer of a B -atom from site j to site j + 1 has an energy cost of 2 V .We can infer that particle-hole excitations in the gas A slow down the flow of B -atoms in the ring lattice.We shall come back to this effect in Sect. 6 after having studied the entanglement generation between thetwo gases. B -current The reduction of current is due to quantum correlations between the A - and B -atoms induced by theinterspecies coupling. We estimate the amount of entanglement between the two gases using the shiftedSchmidt number K AB ( t ) = (cid:0) tr A [ˆ ρ A ( t )] (cid:1) − − (cid:0) tr B [ˆ ρ B ( t )] (cid:1) − − , (20)where ˆ ρ A ( t ) = tr B | ψ AB ( t ) ih ψ AB ( t ) | is the reduced density matrix of the gas A at time t . The shift byone in our definition of the Schmidt number, as compared with the usual definition, makes sure that K AB ( t ) = 0 when the two gases are unentangled. Note that K AB ( t ) is symmetric under the exchange of A and B .Assuming as before that the gas A is in the MI regime and considering times t much smaller than U − A ,a simple calculation using the perturbative expansion (12) and neglecting the hopping of atoms A in thedynamics yields (see in Appendix B) K AB ( t ) = 4 Lα A λ A (cid:16) − (cid:12)(cid:12)(cid:10) e − i t ˆ W (01) B (cid:11) ψ B (cid:12)(cid:12) (cid:17) + O ( λ A ) (21)with the Hamiltonian ˆ W (01) B = ˆ H B − V ∆ˆ n B , ∆ˆ n B = ˆ n B − ˆ n B . (22)In this section, like for the calculation of the B -current, we restrict our analysis to times t much smallerthan all the inverse single species energies. We refer to this condition as defining the short time regime .The behavior of the Schmidt number at larger times will be discussed in Sec. 7. By applying (21), we firstdetermine K AB ( t ) in the short time regime in two limits: λ B ≪ B in the MI regime) in Subsect. 5.1and λ B ≫ B in the SF regime) in Subsect. 5.2. We then discuss the general case in Subsect. 5.3. A and B in the MI regime Let us first assume that both gases are in the MI regime ( λ A , λ B ≪ t ≪ U − A , U − B , theHamiltonian ˆ H B can be replaced by ˆ H int B in Eq. (21) and thus be dropped out (recall that [ ˆ H int B , ∆ˆ n B ] =0 and ˆ H B | ψ B i = E B GS | ψ B i ). The quantum expectation in this equation becomes the same as in thecalculation of the B -current, replacing A by B . This yields K AB ( t ) = 16 Lα A α B λ A λ B (cid:2) − tV ) − cos(2 tV ) + O ( λ A ) + O ( λ B ) (cid:3) (23)with α B = ν B ( ν B + 1). Averaging up to time t with V − ≪ t ≪ U − A , U − B one gets hK AB i t = 48 Lα A α B λ A λ B = 3 L ν A ν B ( ν A + 1)( ν B + 1) V A V B v L ( φ A ) v L ( φ B ) , (24)where we have used in the last equality the linear approximation of the visibility, see (6).9n Eqs. (23) and (24) we can observe that the shifted Schmidt number scales with the system size L .Note that the 2-R´enyi entropy of entanglement, which is a meaningful entanglement measure analogousto the entanglement of formation [51], is related to the Schmidt number defined in (20) by S (2) AB ( t ) = ln( K AB ( t ) + 1) ≃ K AB ( t ) , (25)where the second equality holds for K AB ( t ) ≪
1. The latter condition is satisfied for not too large L ’s in Eq. (24) since we assume λ A , λ B ≪
1. Therefore, the quench dynamics follows a volume law ofentanglement, that is, the entanglement entropy scales linearly with the system size L .Furthermore, by comparing (19) and (23) one finds that K AB ( t ) (and thus S (2) AB ( t )) is proportional tothe relative B -current reduction, S (2) AB ( t ) /L = K AB ( t ) /L = β ν B ,λ B J B ( t ) (26)with a proportionality factor β ν B ,λ B = 8 α B λ B depending on the filling factor ν B and the energy ratio λ B of the B -gas only. Using (6) again, one has β ν B ,λ B = 8 α B λ B = ν B ν B + 2 V B v L ( φ B ) , (27)showing that K AB ( t ) /L is proportional to V B J B ( t ).It is worth noting that the time-averaged Schmidt number per site has a universal value 48 α A α B λ A λ B ,see (24), which is independent of the Peierls phases and, more remarkably, of the strength V of theinter-species interactions. Moreover, the amplitude of the time oscillations of K AB ( t ) has also a universalcharacter, see (23); in particular it does not depend on V which only sets in the period of oscillations.Since both the B -current and the visibilities are measurable quantities, Eqs. (24) and (26) provide away to evaluate the amount of A - B entanglement. This can be done either by measuring the interferencepatterns of the two gases prior to the interaction quench (Eq. (24)) or by performing measurements priorand after the quench on the gas B only (Eq. (26)). Note that the usual technique to estimate entan-glement experimentally in quantum many-body systems relies on quantum state tomography (see [62, 63]and references therein) and can be applied in practice to small systems only (see also alternative tech-niques proposed in [53, 54, 39]). In contrast, estimating entanglement from observations of the interferencepatterns and measurement of the atomic current can be done without any restriction on the system size.The numerical results displayed in Fig. 4 show good agreement with Eqs. (24) and (26) for values of λ A corresponding to V A . .
2. For higher visibilities V A , it is seen in the left panel that hK AB i t dependson the Peierls phases. Note that for large L , (26 and (27) give hK AB i t /L = ( ν B / (2 ν B + 2)) V B hJ B i t , hencethe slopes of the straight lines in the right panel become φ A -independent and equal to 1 / A and B in the MI and SF regime Let us now assume that the gases A and B are in the MI and SF regimes, respectively ( λ A ≪ ≪ λ B ).We will show that, after averaging over time, the shifted Schmidt number and 2-R´enyi entropy ofentanglement show again a volume law and a proportionality to the relative (time-averaged) B -currentreduction, h S (2) AB i t /L = hK AB i t /L = β ν B , ∞ hJ B i t (28)with β ν B , ∞ depending on the filling factor ν B and being almost independent of L . For instance, for ν B = 1one has β ν B , ∞ = 0 . ± .
01 with a difference between the values for L = 3 ,
4, 5, or L → ∞ smaller than10 V A φ A = φ B = π/20φ A = φ B = π/3 time < K > λ B λ A K (t) λ B λ A (a) < J B > < K > / L V B φ A = φ B = π/20 L=3 φ A = φ B = π/20 L=4 φ A = φ B = π/20 L=5 (b)
Figure 4: (a):
Time averaged Schmidt number hK AB i t divided by λ A λ B vs. initial visibility of the gas A for a gas B in the MI regime. The parameters are L = N A = N B = 4, V = 100 U B , J B = 0 . U B , U A = 0 . U B , J A variable, φ A = φ B = π/
20 (red dots) and π/ t = 0 . /U B .Horizontal dashed segment: value of Eq. (24). Inset:
Schmidt number vs time (in units of U − B ) for thesame parameters and J A = 0 . U B . Red and blue line: K AB ( t ) and hK AB i t from numerical calculations.Dashed black lines: formula (64). (b): hK AB i t /L divided by the squared visibility of the gas B vs. relativecurrent reduction hJ B i t , for the same energy parameters, φ A = φ B = π/
20, and L = N A = N B = 3 , .
01 (see Table 1). To prove (28), we assume that the Peierls phase φ B is not close to a half-integer valueof φ and approximate the GS | ψ B i of the gas B by the SF state | ψ SF i = 1 p L N B N B ! (cid:18) L − X j =0 e i φ ℓj b † j (cid:19) N B | i = r N B ! L N B X n e i φ ℓ P j jn j p n ! · · · n L − ! | n i (29)(here | i denotes the vacuum state), making a small error O ( λ − B ). As mentioned above, in the shorttime regime t ≪ J − B , U − A one may replace ˆ W (01) B by − V ∆ˆ n B in the quantum expectation in (21), whichbecomes (cid:10) e − i t ˆ W (01) B (cid:11) ψ B = X n e i tV ( n − n ) (cid:12)(cid:12) h n | ψ SF i (cid:12)(cid:12) + O ( λ − B )= N B ! L N B N B X σ =0 (cid:18) X n + n = σ e i tV n e − i tV n n ! n ! (cid:19)X n ′ ( σ ) n ′ ! · · · n ′ L − ! + O ( λ − B )= N B X σ =0 L σ (cid:16) tV ) (cid:17) σ N B ! σ !( N B − σ )! (cid:18) L − L (cid:19) N B − σ + O ( λ − B ) . (30)Here, the sum over n ′ in the second line runs over all ( L − n ′ = ( n ′ , · · · , n ′ L − ) ∈ N L − such that P L − j =2 n ′ j = N B − σ , and in the last line we have used the multinomial formula (cid:18) X j ∈J x j (cid:19) N = N ! X P j n j = N,n j ≥ Y j ∈J x n j j n j ! if J ⊂ N is a finite subset of indices and x j ∈ R .11 ν B β ν B , ∞ (small time regime) β ′ ν B , ∞ (intermediate time regime)3 1 ≃ . ≃ . ≃ . ∞ ≃ . ≃ . ≃ . ∞ ≃ . ≃ . β ν B , ∞ of Eq. (28) and β ′ ν B , ∞ of Eq. (74) when the gas B is in the SFregime for some values of the site number L and filling factor ν B .Plugging (30) into (21), this yields K AB ( t ) L = 4 α A λ A (cid:26) − (cid:18) − L (cid:0) − cos( tV ) (cid:1)(cid:19) N B + O ( λ A ) + O ( λ − B ) (cid:27) . (31)Averaging over a time interval [0 , t ] with t ≫ V − and using Eq. (19), one obtains Eq. (28) above withthe proportionality factor β ν B , ∞ = 23 (cid:26) − π Z π d s (cid:18) − L (cid:0) − cos s (cid:1)(cid:19) N B (cid:27) ∼ (cid:26) − π Z π d s e − ν B (1 − cos s ) (cid:27) , (32)where the last expression corresponds to the thermodynamical limit N B , L ≫ ν B = N B /L .The values of β ν B , ∞ for N B = 3 ,
4, and 5 and ν B = 1 or 2 are given in Table 1. It is noteworthy that theycoincide with the value at the thermodynamical limit up to the second decimal.The time evolution of the Schmidt number determined numerically when the gas B is in the SF regimeis shown in Figs. 5(a) and (b). One sees that K AB ( t ) presents, in addition to oscillations of period 2 π/V ,a complex behavior which is likely to result from the superposition of the different frequencies associatedto the energy parameters U A , U B , and J B entering the Bose-Hubbard Hamiltonians. Nonetheless, thetime-averaged Schmidt number hK AB i t has a much simpler evolution and basically displays two plateaus.The first plateau at small time t . J − B agrees well with the predicted value obtained by replacing hJ B i t by6 α A λ A in Eq. (28). We shall derive in Sec. 7 an analytical formula giving the value of the second plateauat times t & J − B in the limit L → ∞ (see Eq. (75) below). The fact that this formula reproduces well thevalue of the second plateau for L = 4 in Fig. 5(a) is an indication that the time-averaged Schmidt numberis not sensitive to finite size effects. In Fig. 5(c) one observes that hK AB i /L is proportional to hJ B i forsmall values of hJ B i , i.e., small values of λ A , with a proportionality factor in very good agreement withthe predicted value β ν B , ∞ in Eq. (28). Interestingly, although hK AB i /L is seen to behave non-linearlywith hJ B i at larger hJ B i ’s, the curves obtained for L = 3, L = 4, and L = 5 are superposed on eachothers, suggesting once again that finite size effects are not relevant.Note that our calculations of the B -current and Schmidt number in the previous and in this sectionare restricted to times much smaller than the inverse single species energies (small time regime). In thethermodynamical limit, since these energies scale linearly with the atom numbers N A and N B , to insure thevalidity of our results up to finite times t > H A and ˆ H B as U A → U A /N A , J A → J A /N A , U B → U B /N B and J B → J B /N B .We will extend in Sec.7 the above results on the Schmidt number in both the MI regime λ B ≪ λ B ≫ U − B and J − B , respectively, keeping the assumption t much smaller than the inverse interaction energy of the A -atoms. It will be shown that the time-averagedSchmidt number is still given at these larger times by Eqs. (24) in the MI regime and by Eq. (28) with aslighly different prefactor β ν B , ∞ in the SF regime. 12 time K ( t ) / ( L λ A ) time K ( t ) / ( L λ A ) simulationtheory (a)(b) < J B > < K > / L L=4, φ A = φ B = π/20 L=3, φ A = φ B = π/20 L=4, φ A = φ B = π/3 L=5, φ A = φ B = π/20 (c) Figure 5: (a):
Schmidt number K AB ( t ) divided by Lλ A vs. time (black plain line) and its time average hK AB i t / ( Lλ A ) (red plain line) for a gas B in the SF regime. Time is in units of U − B . The parametersare L = N A = N B = 4, V = 500 U B , J B = 5 U B , U A = 0 . U B , J A = 0 . U B , and φ A = φ B = π/ (b): same on the time interval [0 , . /U B ] with formula (73)shown in green dashed lines. (c): Time-averaged Schmidt number per site hK AB i t /L vs. relative currentreduction hJ B i t , for N A = N B = L and φ A = φ B = φ with L = 3, φ = π/
20 (blue dots), L = 4, φ = π/ L = 4, φ = π/ L = 5, φ = π/
20 (green diamonds). The energyparameters are U A = 0 . U B , J B = U B , J A variable, averaging time t = 0 . /U B . Dashed line: values fromEqs. (28) with β ν B , ∞ = 0 .
54 (see Table 1). A in the MI regime, λ B arbitrary In the more general case λ B arbitrary and λ A ≪
1, one deduces from (21) and (19) that hK AB i t /L = β ν B ,λ B ,L hJ B i t , V − ≪ t ≪ U − A , U − B , J − B , (33)with a proportionality factor given by β ν B ,λ B ,L = 23 (cid:18) − π Z π d s (cid:12)(cid:12)(cid:10) e i s ∆ˆ n B (cid:11) ψ B (cid:12)(cid:12) (cid:19) = 23 (cid:18) − X n,m,n − n = m − m |h n | ψ B i| |h m | ψ B i| (cid:19) . (34)Note that | ψ B i and thus β ν B ,λ B ,L depend on L , so that we cannot claim that Eq. (33) shows a volume lawof entanglement for finite lattices and λ B ≈
1, although one may expect that β ν B ,λ B ,L converges rapidlyto its large L limiting value as in the case λ B ≫ B appear in the RHS of Eq. (33), thisequation provides a way to estimate the amount of entanglement between the two gases when λ A ≪ performing local measurements on the species B . In order to understand better the entanglement generation process in the quench dynamics, it is enlight-ening to look at the superpositions at the origin of this entanglement in the time-evolved wavefunction ofthe two gases. With this aim, we determine this wavefunction | ψ AB ( t ) i to second order in λ A .13isregarding as before the tunneling Hamiltonian ˆ K A and the commutator [ ˆ H int AB , ˆ H B ] in the dynamicalevolution up to times t ≪ U − A , U − B , J − B and recalling that ˆ H B | ψ B i = E B GS | ψ B i , we find from (2) and (12) | ψ AB ( t ) i = e − i t ˆ H int AB e − i t ˆ H int A e − i tE B GS (cid:0) | ψ MI i + λ A | ψ (1) A i + λ A | ψ (2) A i + O ( λ A ) (cid:1) | ψ B i . (35)Let us expand the initial GS of B in the Fock basis as | ψ B i = P n B c n B | n B i and separate the second-ordercorrection to the GS of A as the sum of h ψ MI | ψ (2) A i| ψ MI i and Π ⊥ MI | ψ (2) A i , where Π ⊥ MI = 1 − | ψ MI ih ψ MI | isthe projector onto the orthogonal subspace to the MI state. Using (13) and (14) and recalling that theenergy to create a particle-hole excitation is equal to U A , we obtain | ψ AB ( t ) i = e − i t ( V ν A N B + E (0) A ( ν A )+ E B GS ) (cid:26) (1 − λ A α A L ) | ψ MI i| ψ B i + λ A √ α A e − i tU A X j X n B c n B (cid:16) e i φ A e − i tV ( n Bj +1 − n Bj ) | ϕ A,j, + i + e − i φ A e i tV ( n Bj +1 − n Bj ) | ϕ A,j, − i (cid:17) | n B i (cid:27) + λ A e − i tE B GS e − i t ˆ H int AB e − i t ˆ H int A Π ⊥ MI | ψ (2) A i| ψ B i + O ( λ A ) , (36)where | ϕ A,j, ± i are the particle-hole states, see 16. It is easy to convince oneself that the contribution ofthe term in the last line of (36) to the reduced density matrix ˆ ρ B ( t ) = tr A | ψ AB ( t ) ih ψ AB ( t ) | of the gas B is of order λ A or higher. Since we are interested in the Schmidt number and the B -current, which bothdepend on ˆ ρ B ( t ) only (see (20) and (9)), this term can be dropped out. Disregarding also the irrelevantphase factor, the wavefunction of the two gases at time t is given by(1 − λ A α A L ) | ψ MI i| ψ B i + λ A √ α A e − i tU A X j, ± e ± i φ A | ϕ A,j, ± i| ψ B,j, ± ( t ) i + O ( λ A ) (37)with | ψ B,j, ± ( t ) i = X n B c n B e ∓ i tV ( n Bj +1 − n Bj ) | n B i . (38)The corrections of order λ A to the separable state | ψ MI i| ψ B i in (37) are at the origin of the entanglementbetween the two gases. Thus entanglement comes from the coupling of particle-hole excitations | ϕ A,j, ± i in the gas A with the states | ψ B,j, ± ( t ) i of the gas B at each lattice site j . Let us note that the latterstates are periodic in time with period 2 π/V ; in particular, they come back to their initial j -independentvalue | ψ B i at times t m = 2 mπ/V , m = 1 , , . . . , leading to a disappearance of entanglement (this is validwhen t m ≪ U − A , U − B , J − B ; see the next section for larger times). This explains the observed oscillationsof K AB ( t ) with period 2 π/V in the inset of Fig. 4(a) and in Fig. 5(b).We now look at the effect of the aforementioned superpositions on the B -current I B ( t ) = tr B [ˆ ρ B ( t ) ˆ I B ].Plugging into this formula the reduced density matrix ˆ ρ B ( t ) obtained from the wavefunction 37, one gets I B ( t ) = (1 − λ A α A L ) I B (0) + λ A α A X j, ± (cid:10) ˆ I B (cid:11) ψ B,j, ± ( t ) + O ( λ A ) . (39)The reduction of the B -current comes from the fact that the currents h ˆ I B i ψ B,j, ± ( t ) in the states | ψ B,j, ± ( t ) i are smaller than I B (0). In fact, a simple calculation shows that these currents averaged in the time interval[0 , t ] with t ≫ /V are equal to (cid:10) h ˆ I B (cid:11) ψ B,j, ± (cid:11) t = 12i L X i/ ∈{ j − ,j,j +1 } h ψ B | (cid:0) e i φ B ˆ b † i +1 ˆ b i − h . c . (cid:1) | ψ B i . (40)This expression is nothing but the current for the initial state | ψ B i in the ring after having removed thethree sites i = j − j , and j + 1. Replacing h ˆ I B (cid:11) ψ B,j, ± ( t ) by ( L − I B (0) /L in (39), one gets back formula1419) for the relative current variation. We conclude that the universal reduction of B -current is due tothe coupling of particle-hole excitations in the gas A with the states (38) of the B -gas. The latter have atime-averaged current reduced by a factor ( L − /L , in qualitative agreement with the intuitive energeticargument given in Sec. 4.Let us stress that as a consequence of the energy scale separation (3), all states appearing in thewavefunction (37) have energies much higher than the GS energy E GS AB ≈ − V N A N B for the post-quenchHamiltonian ˆ H AB . More precisely, (cid:10) ˆ H AB (cid:11) ψ MI ⊗ ψ B ≡ − V ν A N B + E (0) A ( ν A ) + E B GS + O ( J A ) (41) (cid:10) ˆ H AB (cid:11) ϕ A,j, ± ⊗ ψ B,j, ± ( t ) = − V ν A N B − V (cid:10) ˆ n Bj +1 − ˆ n Bj (cid:11) ψ B + E (0) A ( ν A ) + U A + (cid:10) ˆ H B (cid:11) ψ B,j, ± ( t ) + O ( J A ) . In particular, the GS of the post-quench Hamiltonian ˆ H AB does not play any role in the quench dynamicsat the times we are considering. This can be explained from the fact that the initial state | ψ A i| ψ B i hasan energy lying in the highly excited part of the spectrum of ˆ H AB . This is illustrated in the numericalresults of Fig. 6 displaying the probabilities p i = |h ψ AB ( t ) | Φ i i| of finding the wavefunction of the twogases in the eigenstates | Φ i i of ˆ H AB (note that the p i are time independent). The peak with value p ∗ i closeto 1 corresponds to the eigenstate of ˆ H AB with closest eigenvalue to E (0) ∗ AB = − V ν A N B + E (0) A ( ν A ) + E B GS .Actually, the state | ψ MI i| ψ B i in (37) is an eigenstate of ˆ H (0) AB = ˆ H int A + ˆ H B + ˆ H int AB with eigenvalue E (0) ∗ AB .Thus by perturbation theory ˆ H AB = ˆ H (0) AB + ˆ K A has an eigenvector | Φ ∗ AB i differing from this state by O ( J A /V ) with eigenenergy E ∗ AB = E (0) ∗ AB + O ( J A ). The corresponding probability p ∗ i is the square of thecoefficient multiplying | ψ MI i| ψ B i in (37), p ∗ i = (cid:12)(cid:12) h ψ MI |h ψ B | ψ AB ( t ) i + O ( J A /V ) (cid:12)(cid:12) = 1 − λ A α A L + O ( λ A ) + O ( J A /V ) . (42)The small non-zero probabilities shown in the two insets are the contribution of the linear terms in (37),which are at the origin of the entanglement. The peak in the upper inset comes from the time-independentcontribution in these linear terms obtained by forgetting all Fock states with n Bj = n Bj +1 in the definition(38) of | ψ B,j, ± ( t ) i ; such states have inter-species interaction energy − V ν A N B and contribute to the time-averaged current. In this section, we determine the Schmidt number at times t that can be of order U − B or J − B but mustbe small with respect to the inverse interaction energy of the gas A , i.e.,0 ≤ t ≪ ( N A U A ) − . (43)Having in mind that the single species energy of the gas A can be tuned to be smaller than U B and J B (this can be done experimentally using Feshbach resonances or by varying the atom number N A ), we definean intermediate time regime by U − B , J − B . t ≪ ( N A U A ) − . Note that one may exchange the role of A and B in the previous inequalities, using the symmetry property of K AB ( t ). The main results establishedbelow are:(i) Formula (64) generalizes Eq. (23) to larger times, giving the time evolution of K AB ( t ) when the gas B is in the MI regime.(ii) Under the same condition, the time-averaged Schmidt number hK AB i t is still given by Eq. (24) inthe intermediate time regime. In the numerical simulations shown in Fig. 7(a), one indeed observesthat hK AB i t stays constant over a wide range of time, even at times larger than U − A .15 Eigenvalues / U B | < ψ ( t ) | φ i > | -1000 -500 0 -800 E AB* / U B Figure 6: Probabilities p i = |h ψ AB ( t ) | Φ i i| of finding the wavefunction of the two gases at time t in aneigenstate | Φ i i of the post-quench Hamiltonian ˆ H AB , for N A = N B = L = 4, V = 200 U B , J B = U A = U B , J A = 0 . U B , φ A = φ B = π/
10, and t = 1 /U B (from numerical exact diagonalization). We checked thatthe numerically evaluated values of p i do not change when changing t , as it should be. Each red pointcorresponds to a given eigenstate of ˆ H AB with an eigenenergy E iAB represented in the horizontal axis (inunits of U B ) and the associated probability p i represented in the vertical axis. The eigenenergy E ∗ AB atwhich p i displays a maximum is defined in the text. The insets show amplifications in particular regionsof the spectrum. For energies E iAB . − U B the values of p i are smaller than 2 × − .(iii) Formula (73) generalizes Eq. (31) to larger times, giving the time evolution of K AB ( t ) when the gas B is in the SF regime.(iv) Under the same condition, Eq. (75) shows that the time-averaged shifted Schmidt number hK AB i t takes in the intermediate time regime a constant value slightly above the short-time value; we findnumerically that hK AB i t is constant even at times larger than U − A and U − B , see Fig. 7(b).The analytical calculations leading to these results are technically more involved than those of the previoussections. Actually, to be able to estimate K AB ( t ) in the intermediate time regime one needs to determinethe spectrum and eigenvectors of the Hamiltonian ˆ W (01) B appearing in Eq. (21). This is done in the twolimiting cases λ B ≪ B in the MI regime) and λ B ≫ B in the SF regime) in the two nextsubsections. A and B in the MI regime As before, we consider a gas A initially in the MI regime and inter-species interaction strengths V muchlarger than the single species energies, see (3). We assume in this subsection that, in addition, the gas B is also in the MI regime, i.e. λ B ≪
1. We determine the shifted Schmidt number in these limits byrelying on formula (21), which gives a good approximation of K AB ( t ) at times t satisfying (43), as shownin Appendix B. Using the spectral decomposition of e − i t ˆ W (01) B , this gives K AB ( t ) = 4 Lα A λ A (cid:20) − X n,n ′ cos (cid:0) t ( w n − w n ′ ) (cid:1) |h φ n | ψ B i| |h φ n ′ | ψ B i| (cid:21) , (44)where | φ n i and w n are the eigenstates and eigenenergies of the Hamiltonian ˆ W (01) B . In order to estimatethese eigenstates and energies, we diagonalize ˆ W (01) B through a two-fold application of time-independentperturbation theory. 16
20 40 60 80 100 time φ A = φ B = π/20φ A = φ B = π/3 L λ A < K > (a) φ A = φ B = π/20φ A = φ B = π/3φ A = φ B = π/4 time < K > t / ( L λ A ) (b) Figure 7: (a):
Time-averaged Schmidt number per site hK AB i t /L divided by λ A as function of time t (inunits of U − B ) for a gas B in the MI regime up to t = 100 /U B (from numerical calculations). Parameters: L = N A = N B = 4, V = 100 U B , J B = 0 . U B , U A = 0 . U B , and J A = 0 . U B , φ A = φ B = π/
20 (blackline) and π/ (b): Same for a gas B inthe SF regime. Parameters: L = N A = N B = 4, V = 500 U B , J B = 5 U B , U A = 0 . U B , J A = 0 . U B , φ A = φ B = π/
20 (black line), π/ π/ U − B shown in the inset. The horizontal blue and green dashed lines display the values from Eq. (75) inthe short and intermediate time regimes, respectively. Note that the phase π/ φ / H B of the gas B as a perturbation of − V ∆ˆ n B . Thelatter unperturbed Hamiltonian has degenerated eigenvalues w (0)∆ = − V ∆ and eigenprojectorsˆΠ ∆ = X n, (∆ n ) =∆ | n ih n | (45)with ∆ ∈ {− N B , . . . , N B } , where we have set (∆ n ) ≡ n − n . Hereafter, we omit the label B referringto the gas B to simplify notations. To lowest order in U B /V , the eigenstates of ˆ W (01) are given by theeigenstates of the projected Bose-Hubbard HamiltonianˆΠ ∆ ˆ H ˆΠ ∆ = ˆ H int ˆΠ ∆ + ˆΠ ∆ ˆ K ˆΠ ∆ (46)and the corresponding eigenvalues are the first-order corrections to the unperturbed energies w (0)∆ . Recallthat the kinetic Hamiltonian reads ˆ K = − J B X e i φ ij ˆ b † j ˆ b i , (47)where < i, j > refers to nearest neighbor sites i, j on the ring lattice and φ ij = ± φ B if j = i ± λ B ≪
1, one may diagonalize the Hamiltonian (46) using once again perturbation theory, now withthe intraspecies interaction term ˆ H int ˆΠ ∆ as the unperturbed part and the projected kinetic HamiltonianˆΠ ∆ ˆ K ˆΠ ∆ as the perturbation. Denote as before by E (0) ( n ) and by E (0)0 the eigenenergies and GS energyof ˆ H int . For any ∆ ∈ {− N B , . . . , N B } and any energy ǫ in the spectrum of ˆ H int ˆΠ ∆ , let us introducethe set S ∆ ,ǫ = { n ∈ N L ; (∆ n ) = ∆ , E (0) ( n ) = ǫ, P j n j = N B } and let ( c m,n ) m ∈S ∆ ,ǫ be the normalizedeigenvectors of the matrix ( h p | ˆ K | m i ) p,m ∈S ∆ ,ǫ with eigenvalues E (1) ( n ), n ∈ S ∆ ,ǫ . The eigenenergies and17igenstates of ˆ W (01) are given by w n = − V (∆ n ) + E (0) ( n ) + E (1) ( n ) + O ( λ B J B ) + O (cid:16) U B V (cid:17) (48) | φ n i = | φ (0) n i + λ B | φ (1) n i + λ B | φ (2) n i + O ( λ B ) + O (cid:16) U B V (cid:17) (49)with (∆ n ) = ∆ and E (0) ( n ) = ǫ , where | φ (0) n i = X m ∈S ∆ ,ǫ c m,n | m i , λ B | φ (1) n i = X p, (∆ p ) =∆ ,E (0) ( p ) = ǫ | p ih p | ˆ K | φ (0) n i E (0) ( n ) − E (0) ( p ) . (50)In particular, since the GS of ˆ H int is non-degenerated and given by the MI state | ψ MI i = | ν i , the set S ∆ ,ǫ reduces to a single element ν when ∆ = 0 and ǫ = E (0)0 . Then (compare with (13)) | φ (0) ν i = | ψ MI i , | φ (1) ν i = √ α B X , ∆ ij =0 e i φ Bij | ν + 1 j − i i , (51)where we have set∆ ij = − ∆ ji = δ ,j − δ ,i − δ ,j + δ ,i = i = 0, j = 11 if ( i = 0, j / ∈ { , } ) or ( i / ∈ { , } , j = 1)0 if i, j / ∈ { , } . (52)Furthermore, the second-order correction is given when ∆ = 0 and ǫ = E (0)0 by λ B | φ (2) ν i = X m,p = ν, (∆ m ) =(∆ p ) =0 h m | ˆ K | p ih p | ˆ K | ψ MI i ( E (0)0 − E (0) ( m ))( E (0)0 − E (0) ( p )) | m i− X m = ν, (∆ m ) =0 |h m | ˆ K | ψ MI i| ( E (0)0 − E (0) ( m )) | ψ MI i . (53)The next step is to calculate the scalar products h φ n | ψ i appearing in (44) up to second order in λ B .To this end, we combine (49) with the perturbative expansion of the GS | ψ i of the gas B analog to (12).Since ˆ K only couples the MI state to Fock states with a single particle-hole excitation and the energy tocreate such an excitation is equal to U B , one has h m | ˆ K | ψ MI i 6 = 0 if and only if E (0) ( m ) = E (0)0 + U B .Consequently, by (50), λ B h φ (1) n | ψ MI i = U B h φ (0) n | ˆ K | ψ MI i if (∆ n ) = 0 and E (0) ( n ) = E (0)0 + U B − λ B h φ (0) n | ψ (1) i if (∆ n ) = 00 if (∆ n ) = 0, (54)where we have used (13) (with A replaced by B ) in the second equality. Expanding | ψ i and | φ (0) n i inpowers of λ B and using h φ (0) n | ψ MI i = δ n,ν , h φ (1) ν | ψ MI i = h ψ (1) | ψ MI i = 0, and (54), we obtain (cid:12)(cid:12) h φ n | ψ i (cid:12)(cid:12) = λ B Re (cid:8) h φ (2) ν | ψ MI i + h φ (1) ν | ψ (1) i + h ψ MI | ψ (2) i (cid:9) + O ( λ B )+ O ( λ B U B V ) if n = νλ B |h φ (0) n | ψ (1) i| + O ( λ B ) + O ( λ B U B V ) if n = ν , (∆ n ) = 0 O ( λ B ) + O ( λ B U B V ) + O ( λ B U B V ) if n = ν , (∆ n ) = 0.(55)18e proceed to evaluate the scalar products appearing in the RHS of (55). One easily shows using (53)that, in analogy with (14), h ψ MI | φ (2) ν i = − X , ∆ ji =0 α B = − α B ( L − , (56)where we have taken advantage of (52) in the last equality. Furthermore, h φ (1) ν | ψ (1) i is given by the sameexpression multiplied by −
2. Replacing these values into (55), this yields (cid:12)(cid:12) h φ ν | ψ i (cid:12)(cid:12) = 1 − λ B α B + O ( λ B ) + O ( λ B U B V ) . (57)Moreover, one has in view of (13) h φ (0) n | ψ (1) i = √ α B X , ∆ ij =(∆ n ) c ν +1 j − i ,n e i φ ij if E (0) ( n ) = E (0)0 + U B E (0) ( n ) = E (0)0 + U B have a single particle-hole excitation, that is, there are of theform | n i = | ν + 1 l − k i with k = l not necessarily nearest neighbors, one deduces that K AB ( t ) = 8 Lα A α B λ A λ B (cid:20) − X k = l,k or l ∈{ , } A kl cos (cid:16) t (cid:0) − ∆ kl V + U B + E (1) B ( ν B + 1 l − k ) (cid:1)(cid:17)(cid:21) (59)with A kl ≡ α B (cid:12)(cid:12) h φ (0) ν +1 l − k | ψ (1) i (cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) X , ∆ ij =∆ kl c ν B +1 j − i ,ν B +1 l − k e − i φ ij (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) . (60)The error terms inside the square bracket in (59) are of order λ B , U B / ( λ B V ), λ B J B t , and U B t/V .Therefore, by inspection on (52) and (59), K AB ( t ) presents oscillations with frequencies 2 V ± U B + O ( J B )and V ± U B + O ( J B ) around the mean value hK AB i t = 48 Lα A α B λ A λ B . In particular, Eq. (24) remainsvalid when the Schmidt number is averaged over time intervals [0 , t ] satisfying V − ≪ t ≪ U − A , ( λ B J B ) − , V U − B (61)(recall that Eq. (24) was shown in Sect. 5.1 to hold in the small time regime t ≪ U − B , U − A only). Notethat for large atom numbers N A , N B , all energy parameters in these inequalities must be multiplied bythe corresponding atom number, see the discussion at the end of Sec. 5.2. The averaged Schmidt number hK AB i t evaluated numerically is displayed as function of t in Fig. 7(a). One observes that hK AB i t isconstant in the time range [100 V − , J − B ]. Note that a good agreement with the predicted value fromEq. (24) is obtained even for times t not satisfying the condition t ≪ U − A .Formula (59) is valid in the intermediate time regime t ≈ J − B (with J B ≫ U A , U B /V ). We now derivea simpler and more explicit expression for times 0 ≤ t . U − B ≪ J − B . Then the first-order correction E (1) B to the eigenenergies (48), which is of order J B , can be neglected and one gets K AB ( t ) = hK AB i t − Lα A λ A λ B X ∆= ± , ± cos (cid:0) t ( − V ∆ + U B ) (cid:1) X m, (∆ m ) =∆ ,E (0) B ( m )= E (0)0 ,B + U B (cid:12)(cid:12) h φ (0) m | ψ (1) B i (cid:12)(cid:12) . (62)19ut for any eigenenergy ǫ of ˆ H int B one has X m ∈S ∆ ,ǫ (cid:12)(cid:12) h φ (0) m | ψ (1) B i (cid:12)(cid:12) = X m ∈S ∆ ,ǫ (cid:12)(cid:12) h m | ψ (1) B i (cid:12)(cid:12) . (63)Furthermore, in view of Eq. (13), for any m = ν B it holds h m | ψ (1) B i = √ α B e i φ Bij if m = ν B + 1 j − i with i, j nearest neighbors and h m | ψ (1) B i = 0 otherwise. This yields, using (52) again,1 L K AB ( t ) = 8 α A α B λ A λ B h − tV ) cos( tU B ) − tV ) cos( tU B ) i if t ≪ U − A , J − B , V U − B . (64)Note that this formula agrees with (23) in the short time regime t ≪ U − B , U − A . The time evolution of K AB ( t ) calculated numerically by exact diagonalization is compared with our analytical result (64) in theinset of Fig. 4(a), showing a good agreement. The presence of the frequency U B in Eq. (64), which is theparticle-hole excitation energy for the species B, indicates that the entanglement does not only depend onthe inter-species interactions but also on the post-quench intra-species interactions.It is worth observing that the Schmidt number K AB ( t ) being symmetric in A and B , the wholecalculation above remains valid by exchanging the role of A and B . Therefore, the time-averaged Schmidtnumber is also given by (24) for averages up to time t satisfying V − ≪ t ≪ U − B , ( λ A J A ) − , V U − A , insteadof (61). Similarly, (64) is correct with U B replaced by U A when t ≪ U − B , J − A , V U − A . Furthermore, asshown in Appendix B, Eqs. (24) and (64) can be easily extended to Bose-Bose mixtures in arbitrary finitelattices Λ in dimension D = 1, 2, or 3 with periodic boundary conditions, see Eqs. (B.12) and (B.13). A and B in the MI and SF regimes Let us now determine the Schmidt number in the intermediate time regime when λ B ≫ B in the SFregime). As in the previous subsection, we rely on the separation of energy scales (3) and determine thespectrum and eigenvectors of ˆ W (01) using perturbation theory, the perturbation being the Bose-HubbardHamiltonian ˆ H of the gas B . One must diagonalize the projected HamiltonianˆΠ ∆ ˆ H Π ∆ = Π ∆ ˆ K Π ∆ + O ( U B ) = Π ∆ ˆ K ′ + O ( U B ) , (65)where ˆΠ ∆ is the projector (45). In the last equality in Eq. (65) we have introduced the kinetic Hamiltonianˆ K ′ of an open chain Λ ′ = { , . . . , L − } with ( L −
2) sites,ˆ K ′ = − J B L − X j =2 ( e i φ B ˆ b † j +1 ˆ b j + h . c . ) . (66)The equality Π ∆ ˆ K Π ∆ = Π ∆ ˆ K ′ follows from the fact that ∆ ij = 0 if and only if i, j / ∈ { , } , see (52).The eigenvectors of Π ∆ ˆ K ′ are given by | n , n i| φ ′ σ,k i with n − n = ∆ and n + n = σ , | φ ′ σ,k i being the eigenvectors of ˆ K ′ for an atomic gas with N σ ≡ N B − σ atoms trapped in the chain Λ ′ , i.e.,ˆ K ′ | φ ′ σ,k i = w ′ σ,k | φ ′ σ,k i . Thus the perturbative eigenvalues and eigenvectors of ˆ W (01) are w n ,n ,k = − V ( n − n ) + w ′ σ,k + O ( U B ) + O (cid:16) J B V (cid:17) (67) | φ n ,n ,k i = | n , n i| φ ′ σ,k i + O ( λ − B ) + O (cid:16) J B V (cid:17) . (68)Using the spectral decomposition e − i t ˆ K ′ = P k e − i tw ′ σ,k | φ ′ σ,k ih φ ′ σ,k | , we find h ψ | e − i t ˆ W (01) | ψ i = X n ,n e i tV ( n − n ) h ψ | n , n i e − i t ˆ K ′ h n , n | ψ i (69)20p to errors of order λ − B , J B /V , U B t , and J B t/V . Note that Eq. (69) applies to times satisfying t ≪ U − B ,in addition to (43), because we have neglected the interaction term Π ∆ ˆ H int in (65); it is, however, justifiedto use (69) for times t ≈ J − B . To simplify we assume that − φ / ≤ φ B < φ /
2, so that the gas B haszero angular momentum ( ℓ = 0). Replacing | ψ i by the SF state | ψ SF i , a calculation similar to the oneperformed in Sect. 5.2 leads to h ψ | e − i t ˆ W (01) | ψ i = N B X σ =0 (cid:18) tV ) L (cid:19) σ N B ! σ !( N B − σ )! (cid:18) L − L (cid:19) N B − σ R σ ( t ) (70)with R σ ( t ) ≡ h ψ ′ SF ,σ | e − i t ˆ K ′ | ψ ′ SF ,σ i . (71)Here, | ψ ′ SF ,σ i is the SF state of a gas with N σ atoms and angular momentum ℓ = 0 trapped in the ringlattice Λ ′ ring = { , · · · , L − } with ( L −
2) sites. Note that the identity (70) does not hold when ℓ = 0 dueto the different values of φ in the lattices Λ and Λ ′ ring .We now show that in the thermodynamical limit the factor R σ ( t ) is equal to e − i tE SF ,σ , where E SF ,σ = − N σ J B cos φ B is the GS energy of a non-interacting Bose gas with N σ atoms in the ring lattice Λ ′ ring .Actually, let us denote by ˆ K ′ the first quantization version of (66) and by | ϕ ring0 i the GS of a single atomin the ring Λ ′ ring . The kinetic Hamiltonian on Λ ′ ring differs from ˆ K ′ by the presence of hopping termsbetween sites 2 and L −
1, namely by ˆ Q = − J B ( e i φ B | ih L − | +h . c . ). We can use second-order perturbationtheory to expand | ϕ ring0 i in terms of the eigenstates | ϕ ′ k i of ˆ K ′ and to evaluate h ϕ ring0 | e − i t ˆ K ′ | ϕ ring0 i .Observing that |h ϕ ′ k | ˆ Q | ϕ ′ i| is typically much smaller than L − for large L ’s (the GS and first excitedstates of the chain Λ ′ are delocalized on the whole chain), this gives h ϕ ring0 | e − i t ˆ K ′ | ϕ ring0 i ≃ e − i tE ′ + o ( L − ),where E ′ is the GS energy of a single atom in the chain Λ ′ . Since the SF state is the N σ -fold tensor productof the single atom state | ϕ ring0 i and noting that o ( L − ) = o ( N − σ ), one infers that for N σ ≫ R σ ( t ) = h ϕ ring0 | e − i t ˆ K ′ | ϕ ring0 i N σ ∼ e − i tN σ E ′ ∼ e tN σ J B cos φ B , (72)where we have approximated E ′ by the GS energy E ring0 = − J B cos φ B of a single atom in the ring Λ ′ ring .Replacing (70) and (72) into (21), one obtains in the thermodynamical limit N B , L ≫ K AB ( t ) L ≃ α A λ A (cid:16) − exp (cid:8) − ν B (cid:0) − cos( tV ) cos(2 tJ B cos φ B ) (cid:1)(cid:9)(cid:17) . (73)This formula is valid for Peierls phases φ B satisfying − φ / < φ B < φ / t much smaller than U − B , U − A , and V J − B . Notice that one recovers the result of Eq. (31) when t ≪ J − B . On the other hand,averaging K AB ( t ) up to times t with J − B ≪ t ≪ U − B , U − A , V J − B , one is led to hK AB i t L ≃ α A λ A β ′ ν B , ∞ , β ′ ν B , ∞ = 23 (cid:18) − lim k →∞ Z πk d s πk e − ν B (1 − cos( s ) cos( s/k )) (cid:19) . (74)This result easily follows by (i) making the change of variables s = V t in the integral between 0 and t of the exponential in (73); (ii) approximating the second cosine by cos( s/k ) with k = E [ V / (2 J B cos φ B )](this is justified since t ≪ V J − B ); (iii) using the periodicity of the resulting integrand; (iv) taking k = O ( V /J B ) → ∞ .We may summarize the results of this subsection and those of Sec. 5.2 by the following formula hK AB i t L = 6 α A λ A ( β ν B , ∞ if V − ≪ t ≪ J − B β ′ ν B , ∞ if J − B ≪ t ≪ U − B , U − A , V J − B , (75)21here β ν B , ∞ and β ′ ν B , ∞ are given by (32) and (74); their numerical values are given in Table 1 for ν B = 1and 2.We compare in Fig. 5(b) formula (73) with the Schmidt number K AB ( t ) evaluated numerically. Onesees that while Eq. (73) describes well the time evolution of K AB ( t ) in the time range [0 , J − B ], at latertimes discrepancies appear; this is likely to be due to contributions of larger frequencies 1 /U B and 1 /U A that have been neglected in (73). In spite of this, one observes in Fig. 7(b) that the averaged Schmidtnumber hK AB i t evaluated numerically for L = N A = N B = 4 and λ B = 5 is constant over a wide rangeof time [ J − B , U − B ], in analogy with what happens when B is in the MI regime. This shows that theaforementioned contributions of the frequencies 1 /U A and 1 /U B average out to zero. For Peierls phases φ A = φ B = π/
20 and π/
3, the constant value agrees reasonably well with the value predicted by Eq. (74)even for times t not satisfying the conditions t ≪ U − B and t ≪ U − A , even though a slight dependenceon the phase is observed for time t larger than J − B (in contrast, as seen in the inset, the two curves forthe phases π/
20 and π/ t ≪ J − B ). On the other hand, hK AB i t is smallerfor φ A = φ B = φ / π/
4, both in the short time regime t . J − B and in the intermediate and largetime regimes J − B . t ≤ U − B . As pointed out above, Eqs. (28), (74), and (75) are not valid at thisphase value, since the GS of the gas B is a superposition of SF states with angular momenta ℓ = 0 and ℓ = 1. This explains the small disagreement between our numerical and analytical results. We observeon Fig. 7(b) that the values of hK AB i t /L for φ A = φ B = φ / φ A = φ B = π/
20 or π/ . λ A .The fact that Eqs. (73) and (74)), which have been derived assuming a large ring, describe well thebehavior of K AB ( t ) at times t . J − B and its averaged value at times t ≫ J − B for a ring with L = 4 sitesis an indication that finite size effects do not play an important role. We have studied the dynamics of persistent currents in a binary mixture in the presence of an interactionquench between the species. Assuming that one of the gas components is in the MI regime, we founduniversal relations between the persistent current of the other species, the visibility before the quench,and the amount of inter-species entanglement. In particular, for the short time regime, the analyticalexpressions of the time average of the current show that it is proportional to the 2-R´enyi entropy ofentanglement and scales quadratically with the visibility before the quench.Our results are in principle amenable to experimental realizations. In fact, 1D-ring lattice trappingpotentials have been already achieved [60] and artificial gauge fields can be created in various ways [10].Interaction quench dynamics can be generated using the ability to tune atomic interactions with Feshbachresonances [45, 46, 47]. Finally, time-of-flight techniques for obtaining information on the atomic currentare available. In particular, the winding number and current-phase relationship can be extracted from thespiral interference pattern obtained after releasing the trap in the plane of the ring and letting the BECinterfere with another BEC confined initially near its center [55, 56, 57, 12].Some open problems require further studies. Firstly, even though our numerical simulations indicatea saturation of the Schmidt number at large times, at least after having performed a local time averagingto suppress oscillations with frequency equal to the inter-species interaction strength (see Fig. 7), wecan not give a definite answer to the question of a possible convergence of the quench dynamics to astationary state. Such a stationary state, if it exists, will certainly not be a thermal state in view of thestate probability distribution in the eigenbasis of the post-quench Hamiltonian shown in Fig. 6. In orderto get more insight on the convergence properties of the out-of-equilibrium dynamics, one would need toperform simulations for larger system sizes, which is difficult with our exact diagonalization approach; ouranalytical approach is also of limited use for this problem because it applies to small or intermediate times.Secondly, it would be of interest to study the current and entanglement evolutions when one of the gascomponents is in the SF regime, instead of the MI regime. One then expects that the gas component could22ave a stronger effect on the current of the other component. This will be investigated in a forthcomingwork. Thirdly, it would be worth adapting our methods to the case of a single impurity atom interactingwith a gas of atoms of a different species, taking N A = 1 and N B ≫ Acknowledgments
D.S. acknowledges support from the Fondecyt project N A Visibility of a BEC trapped in a 1D-ring lattice potential
In this appendix we determine the visibility of the interference fringes of a single BEC trapped in a 1D-ringlattice potential in the presence of an artificial gauge field. As mentioned in the main text, the interferencepattern is obtained after releasing the potential trap and letting the atomic gas expand freely, thus it isdirectly related to the momentum distribution S φ ( q ) = L − X i,j =0 e i q ( i − j ) h ψ φ | ˆ b † i ˆ b j | ψ φ i , (A.1)where | ψ φ i is the ground state (GS) of the Bose-Hubbard Hamiltonian ˆ H φ of a single species (Eq. (1) inthe main text), φ is the Peierls phase, ˆ b † i and ˆ b j are the creation and annihilation operators of a bosonicatom at sites i and j , and L is the number of lattice sites. The visibility is defined as V φ = S max − S min S max + S min , (A.2) S max and S min being the maximum and minimum of S φ ( q ). By gauge invariance, V φ is periodic in φ with period φ = 2 π/L . Actually, the Hamiltonians ˆ H φ + φ and ˆ H φ and their GSs are related by a gaugetransformation, e.g. | ψ φ + φ i = e i φ ˆ X | ψ φ i with ˆ X = P j j ˆ b † j ˆ b j , implying that S φ + φ ( q ) = S φ ( q − φ ).We determine analytically the momentum distribution and visibility when the Bose gas is in the Mott-insulator (MI) regime, i.e., for small energy ratio λ = J/U . We rely on standard perturbation theory, theperturbation being the kinetic term ˆ K φ = − J L − X j =0 (e i φ ˆ b † j +1 ˆ b j + h . c . ) . (A.3)As in the main text, the Fock states diagonalizing ˆ H int = ˆ H φ − ˆ K φ are denoted by | n i = | n , . . . , n L − i , E (0) ( n ) are the corresponding eigenenergies, the filling factor ν = N/L is assumed to be an integer, and | ψ MI i = | ν i = | ν, . . . , ν i is the MI state (GS of ˆ H int with energy E (0)0 = E ( ν )). For small λ , the GS of theBEC reads | ψ φ i = | ψ MI i + λ | ψ (1) φ i + λ | ψ (2) φ i + O ( λ ) (A.4)with the first and second-order corrections λ | ψ (1) φ i = − U ˆ K φ | ψ MI i = λ L − X j =0 (cid:0) e i φ ˆ b † j +1 ˆ b j + h . c . (cid:1) | ψ MI i (A.5) λ | ψ (2) φ i = X m,p = ν h m | ˆ K φ | p ih p | ˆ K φ | ψ MI i ( E (0)0 − E (0) ( m ))( E (0)0 − E (0) ( p )) | m i − X m = ν |h m | ˆ K φ | ψ MI i| ( E (0)0 − E (0) ( m )) | ψ MI i . (A.6)23ere, we have used h ψ MI | ˆ K φ | ψ MI i = 0 and the fact that the energy to create a particle-hole excitation ∝ ˆ b † i ˆ b j | ψ MI i is equal to E (0) ( ν + 1 i − j ) − E (0)0 = U , i = j , (A.7)where 1 j = (0 , . . . , , , . . . ,
0) is the vector having a single nonzero component equal to unity at site j .Noting that h m | ˆ K φ | ψ MI i = 0 if E (0) ( m ) = E (0)0 + U , one deduces from (A.6) that λ h ψ MI | ψ (2) φ i = − U h ψ MI | ˆ K φ | ψ MI i = − λ (cid:13)(cid:13) ψ (1) φ (cid:13)(cid:13) = − ν ( ν + 1) λ L . (A.8)Substituting (A.4) into (A.1), the momentum distribution reads S φ ( q ) = L − X i,j =0 e i q ( i − j ) h ψ MI | ˆ b † i ˆ b j | ψ MI i + 2 λ Re (cid:26) L − X i,j =0 e i q ( i − j ) h ψ MI | ˆ b † i ˆ b j | ψ (1) φ i (cid:27) +2 λ Re (cid:26) L − X i,j =0 e i q ( i − j ) h ψ MI | ˆ b † i ˆ b j | ψ (2) φ i (cid:27) + λ L − X i,j =0 e i q ( i − j ) h ψ (1) φ | ˆ b † i ˆ b j | ψ (1) φ i + O ( λ ) . (A.9)A simple calculation yields h ψ MI | ˆ b † i ˆ b j | ψ MI i = νδ ij h ψ MI | ˆ b † i ˆ b j | ψ (1) φ i = ν ( ν + 1)(1 − δ i,j ) L − X l =0 (cid:16) e i φ δ i,l δ j,l +1 + e − i φ δ i,l +1 δ j,l (cid:17) h ψ MI | ˆ b † i ˆ b j | ψ (2) φ i = ν ( ν + 1) (cid:20) (1 − δ i,j ) L − X l =0 (cid:16) ( ν + 1)e φ δ i,l δ j,l +2 + ν e φ δ i,l − δ j,l +1 (A.10)+( ν + 1)e − φ δ i,l +1 δ j,l − + ν e − φ δ i,l +2 δ j,l (cid:17) − νLδ i,j (cid:21) h ψ (1) φ | ˆ b † i ˆ b j | ψ (1) φ i = h ψ MI | ˆ b † i ˆ b j | ψ (2) φ i + 3 ν ( ν + 1) Lδ i,j , where, according to the periodic boundary conditions, δ j,L must be identified with δ j, , etc... Note thatthe expressions (A.10) for i = j immediately follow from b † i b i | ψ MI i = ν | ψ MI i , h ψ MI | ψ (1) φ i = 0, and (A.8).Plugging (A.10) into (A.9), one finds S φ ( q ) = Lν (cid:26) ν + 1) λ h(cid:16) − L (cid:17) cos( q − φ ) + 1 L cos (cid:0) q ( L −
1) + φ (cid:1)i +6( ν + 1)(2 ν + 1) λ h(cid:16) − L (cid:17) cos (cid:0) q − φ ) (cid:1) + 2 L cos (cid:0) q ( L −
2) + 2 φ (cid:1)i + O ( λ ) (cid:27) . (A.11)One easily checks that this expression fulfills S φ + φ ( q ) = S φ ( q − φ ) (gauge invariance) and S − φ ( q ) = S φ ( − q ). As a consequence, the visibility (A.2) satisfies V φ ∓ φ = V φ and it is enough to let φ variesbetween 0 and φ / π/L .A non-perturbative expression of S ( q ) has been determined in Refs. [58, 59] for infinite lattices in theabsence of gauge field. Expanding the result of the first reference for small λ ’s, one recovers Eq.(A.11) inthe special case φ = 0 and for an infinite ring ( L → ∞ ), up to an irrelevant constant term. However, asstressed in the main text, we are interested in this article in finite rings with a non-zero gauge field.We now determine the extrema of S φ ( q ). We assume L ≥ λ ≥
0. One has ∂S φ ∂q = − λLν ( ν + 1) (cid:26)(cid:16) − L (cid:17)h sin( q − φ ) + sin (cid:0) q ( L −
1) + φ (cid:1)i +3 λ (2 ν + 1) (cid:16) − L (cid:17)h sin (cid:0) q − φ ) (cid:1) + sin (cid:0) q ( L −
2) + 2 φ (cid:1)i + O ( λ ) (cid:27) . (A.12)24e first obtain the extrema of S φ to lowest order in λ . Neglecting the terms of order λ , ∂S φ /∂q vanisheswhen q − φ = − ( q ( L −
1) + φ ) + 2 mπ or q − φ = q ( L −
1) + φ − (2 m + 1) π , i.e., q = q m = 2 mπL or q = q ′ m = (2 m + 1) π − φL − m an arbitrary integer. At these extremal points, S φ ( q ) takes the values S φ ( q m ) = Lν n ν + 1) λ cos (cid:0) q m − φ (cid:1) + O ( λ ) o S φ ( q ′ m ) = Lν n ν + 1) (cid:16) − L (cid:17) λ cos (cid:0) q ′ m − φ (cid:1) + O ( λ ) o . (A.14)To determine which of the extremal points q m , q ′ m corresponds to the global maximum (minimum), we areleft with the task of maximizing (minimizing) C m ( φ ) = cos( q m − φ ) and C ′ m ( φ ) = cos( q ′ m − φ ) over all m ∈ { , . . . , L − } and m ∈ { , . . . , L − } , respectively.Under the assumption 0 ≤ φ ≤ π/L , one has max m { C m ( φ ) } = C ( φ ) = cos φ . Furthermore, it followsfrom the bound 1 − /L < cos( π/L ) that for any m , − cos (cid:16) πL (cid:17) < (cid:16) − L (cid:17) C ′ m ( φ ) < cos (cid:16) πL (cid:17) ≤ cos φ . (A.15)Thus, to lowest order in λ , the maximum of S φ ( q ) is reached for q = q = 0 modulo 2 π . Since the term oforder λ in the derivative in Eq.(A.12) also vanishes for q = q , the exact momenta at which S φ ( q ) reachesits maximum differ from q , q ± π, . . . by at most O ( λ ). The value of the maximum is S max = Lν n ν + 1) λ cos φ + 6( ν + 1)(2 ν + 1) λ cos(2 φ ) + O ( λ ) o . (A.16)Similarly, min m { C m ( φ ) } = C L/ ( φ ) = − cos φ if the number of sites L is even and min m { C m ( φ ) } = C ( L +1) / ( φ ) = − cos( π/L − φ ) if L is odd. By (A.15), the minimum of S φ ( q ) is, to lowest order in λ ,reached for q = q L/ = π modulo 2 π in the first case and for q = q ( L +1) / = π + π/L modulo 2 π in thelast case. Again, the term of order λ in the derivative in (A.12) vanishes at these momenta, thus theseformulas give the locations of the momenta minimizing S φ ( q ) up to corrections of order O ( λ ). The valueof the minimum is S min = Lν n − ν + 1) λ cos φ + 6( ν + 1)(2 ν + 1) λ cos(2 φ ) + O ( λ ) o if L is even Lν n − ν + 1) λ cos (cid:0) πL − φ (cid:1) + 6( ν + 1)(2 ν + 1) λ cos (cid:0) πL − φ (cid:1) + O ( λ ) o if L is odd.(A.17)We note that when φ = π/L (respectively φ = 0 and L is odd), S φ ( q ) reaches its maximum (minimum) attwo phase points in the interval [0 , π ], namely at q = 0 and q = 2 π/L (respectively at q ( L +1) / = π + π/L and q ( L − / = π − π/L ), up to errors O ( λ ).Even though we are considering a gas in the MI regime, it is convenient to introduce the angularmomentum of the superfluid (SF) state, ℓ = E [ φ/φ + 1 / E denotes the integer part, and similarlylet k = E [ φ/φ ]. The introduction of the integers ℓ and k helps to write the result in a concise form forarbitrary phases φ ∈ R , taking into account the periodicity and symmetry properties V ± φ + φ = V φ . Weconclude from (A.16) and (A.17) that if L is even, L ≥
4, the visibility is given by V φ = 4( ν + 1) λ cos (cid:0) φ − ℓφ (cid:1) + O ( λ ) , L even. (A.18)If L is odd, L ≥
3, it is given by V φ = 2( ν + 1) λ h cos (cid:16) φ − φ + kφ (cid:17) + cos (cid:0) φ − ℓφ (cid:1)i × (cid:26) − (4 ν + 1) λ h cos (cid:16) φ − φ + kφ (cid:17) − cos (cid:0) φ − ℓφ (cid:1)i(cid:27) + O ( λ ) , L odd. (A.19)25 A / π I n t e n s it y π/10π/2π/3 -2 -1 0 1 2 q/ π I n t e n s it y φ=π/2φ=0φ=π/10φ=π/3 φ V φ V Figure 8:
Upper panels: momentum distribution S φ ( q ) (in arbitrary units) as function of q for a singleBEC in the MI regime (left panel, λ = J/U = 0 .
05) and in the SF regime (right panel, λ = 1), fromnumerical calculations; N = 4 atoms are trapped in a 1D-ring lattice potential with L = 4 sites andPeierls phases φ = 0 , π/ , π/
3, and π/ Lower panels: visibility of a Bose gas in the MI regime asfunction of the Peierls phase for N = L = 4 (left panel) and N = L = 5 (right panel), with λ = J/U = 0 . L coincideswith the visibility of an infinite chain, V φ =0 = 4( ν + 1) λ + O ( λ ). The latter has been determined inRef. [58]. In contrast, for finite odd numbers of sites V φ =0 = 2( ν + 1)(1 + cos( π/L )) λ + O ( λ ) differs fromthe aforementioned gauge-free result of Ref. [58], which holds in the large L limit only.It is easy to show from Eqs. (A.18) and (A.19) that V φ is minimum at half-integer values of φ . If L is even, this follows from the fact that cos( φ − ℓφ ) is minimum at such phase values. If L is odd, thefirst-order contribution to the visibility in (A.19) is minimum for φ = mφ and φ = ( m + 1 / φ , with m ∈ Z . By comparing the second-order corrections, one sees that V φ reaches its minimum for the latterphase point.The visibility is displayed as function of φ in Fig. 8 (lower panels) for ring lattices with L = 4 and 5 sitesand λ = 0 .
01. The values obtained from Eqs. (A.18)-(A.19) agree well with those obtained by locatingnumerically the extrema of the momentum distribution. The latter is obtained numerically via an exactdiagonalization of the Bose-Hubbard Hamiltonian. We have also calculated numerically the momentumdistribution and visibility for a gas in the SF regime ( λ = 1) and in the transition regime ( λ ≃ . S φ ( q ) are the same in both regimes. Forinstance, it is seen in Fig. 8 that S φ ( q ) reaches its maximum at q = ℓφ modulo 2 π in both regimes.26 Schmidt number as an expectation of a local time-evolution operator
In this appendix we derive formula (21) for the Schmidt number K AB ( t ) quantifying the entanglementbetween the two gases after interspecies interactions have been switched on.Our analysis applies here to trapping potentials forming an arbitrary finite lattice Λ with L sites indimension D = 1, 2, or 3, with periodic boundary conditions. We denote by z the number of nearestneighbors of a given site i ∈ Λ. For instance, z = 2 for a 1 D ring lattice and z = 4 for a 2 D square lattice.As we will show below, the Schmidt number turns out to depend on the geometry of the lattice through L and z only. The kinetic part ˆ K B of the Bose-Hubbard Hamiltonian of the gas B is given by Eq. (47)in the main text, where the sum runs over all pairs of nearest neighbor sites i, j ∈ Λ and the phases φ ij ∈ { φ B , − φ B } are translation-invariant, invariant under rotation over the origin, and satisfy φ ji = − φ ij for any neighboring sites i and j , in order to insure that ˆ K B be self-adjoint. A similar expression holdsfor the kinetic Hamiltonian ˆ K A of the gas A .The wavefunction of the Bose mixture at time t is | ψ AB ( t ) i = e − i t ˆ H AB | ψ A i ⊗ | ψ B i , where ˆ H AB is theHamiltonian after interactions have been turned on [see Eq. (2) in the main text] and | ψ A i and | ψ B i arethe GSs of the gases A and B . We assume that the gas A is in the MI regime ( λ A ≪ | ψ A i hasthe form given in (A.4). Expanding in powers of λ A , the reduced density matrix of the gas A readsˆ ρ A ( t ) = tr B [ | ψ AB ( t ) ih ψ AB ( t ) | ] = ˆ ρ (0) A ( t ) + λ A ˆ ρ (1) A ( t ) + λ A ˆ ρ (2) A ( t ) + O ( λ A ) (B.1)with ˆ ρ (0) A ( t ) = tr B (cid:2) e − i t ˆ H AB | ψ MI ih ψ MI | ⊗ | ψ B ih ψ B | e i t ˆ H AB (cid:3) ˆ ρ (1) A ( t ) = tr B (cid:2) e − i t ˆ H AB (cid:0) | ψ (1) A ih ψ MI | + | ψ MI ih ψ (1) A | (cid:1) ⊗ | ψ B ih ψ B | e i t ˆ H AB (cid:3) (B.2)ˆ ρ (2) A ( t ) = tr B (cid:2) e − i t ˆ H AB (cid:0) | ψ (2) A ih ψ MI | + | ψ MI ih ψ (2) A | + | ψ (1) A ih ψ (1) A | (cid:1) ⊗ | ψ B ih ψ B | e i t ˆ H AB (cid:3) . For times t satisfying t ≪ U − A , the kinetic Hamiltonian ˆ K A of the gas A can be dropped out, so thatˆ H AB ≃ ˆ H int A + ˆ H int B + ˆ K B | {z } = ˆ H B + ˆ H int AB , ˆ H int AB = − V ˆ n A · ˆ n B = − V X j ∈ Λ ˆ n Aj ˆ n Bj (B.3)(for simplicity we do not write the identity operators explicitly, e.g., ˆ H int A stands for ˆ H int A ⊗ ˆ B ). Makingthis approximation in (B.2) and using ˆ H int A | ψ MI i = E (0)0 ,A | ψ MI i and ˆ H int AB | ψ MI i = − V ν A N B | ψ MI i , one easilyfinds that ˆ ρ (0) A ( t ) is independent of time and given byˆ ρ (0) A ( t ) = | ψ MI ih ψ MI | . (B.4)Consequently, plugging (B.1) and (B.4) into Eq. (20) in the main text,( K AB ( t ) + 1) − = 1 + 2 λ A h ψ MI | ˆ ρ (1) A ( t ) | ψ MI i + 2 λ A h ψ MI | ˆ ρ (2) A ( t ) | ψ MI i + λ A tr (cid:2) ˆ ρ (1) A ( t ) (cid:3)(cid:1) + O ( λ A ) . (B.5)The linear term in λ A vanishes. In fact, h ψ MI | ˆ ρ (1) A ( t ) | ψ MI i = 2Re h ψ MI | ψ (1) A i by (B.2) and (B.3), and h ψ MI | ψ (1) A i = 0, see (A.5). Generalizing (A.8) to the case of a general lattice Λ, one has λ A h ψ MI | ˆ ρ (2) A ( t ) | ψ MI i = 2 λ A Re h ψ MI | ψ (2) A i = − λ A α A Lz (B.6)with α A = ν A ( ν A + 1). 27herefore, the only time-dependent contribution to K AB ( t ) is due to the trace of ˆ ρ (1) A ( t ) in (B.5), whichcan be estimated as follows. Let us set ˆ A (1) = | ψ (1) A ih ψ MI | + | ψ MI ih ψ (1) A | . Since | ψ (1) A i = − J − A ˆ K A | ψ MI i ,see A.5, it holds (cid:12)(cid:12) h n A | ˆ A (1) | m A i (cid:12)(cid:12) = α A X δ n A ,ν A +1 j − i δ m A ,ν A + ( n A ↔ m A ) (B.7)with | ν A + 1 j − i i the Fock state differing from | ψ MI i by the presence of an extra particle at site j and ahole at site i . By evaluating the trace in the Fock basis {| n A i} diagonalizing both ˆ n A and ˆ H int A , one getsfrom (B.2) and (B.3)tr A (cid:2) ˆ ρ (1) A ( t ) (cid:3) = X n A ,m A (cid:12)(cid:12) h n A | ˆ ρ (1) A ( t ) | m A i (cid:12)(cid:12) = X n A ,m A (cid:12)(cid:12) h n A | ˆ A (1) | m A ih ψ B | e i t ( ˆ H B − V m A · ˆ n B ) e − i t ( ˆ H B − V n A · ˆ n B ) | ψ B i (cid:12)(cid:12) = 2 α A X (cid:12)(cid:12) h ψ B | e i t ( ˆ H B − V ν A ˆ N B ) e − i t ( ˆ H B − V ν A ˆ N B − V (ˆ n Bj − ˆ n Bi )) | ψ B i (cid:12)(cid:12) . (B.8)Thanks to the commutation of the Bose-Hubbard Hamiltonian ˆ H B with the total number operator ˆ N B ofthe B -atoms and since | ψ B i is an eigenstate of ˆ H B , it holds (cid:12)(cid:12) h ψ B | e i t ( ˆ H B − V ν A ˆ N B ) e − i t ( ˆ H B − V ν A ˆ N B − V (ˆ n Bj − ˆ n Bi )) | ψ B i (cid:12)(cid:12) = (cid:12)(cid:12) h ψ B | e − i t ( ˆ H B − V (ˆ n Bj − ˆ n Bi )) | ψ B i (cid:12)(cid:12) . (B.9)Now, it follows from the translation-invariance of ˆ H B and its invariance under rotations around theorigin that h ψ B | e − i t ( ˆ H B − V (ˆ n Bj − ˆ n Bi )) | ψ B i = h ψ B | e − i t ( ˆ H B − V (ˆ n B − ˆ n B )) | ψ B i if i, j are nearest neighbors, (B.10)where 1 denotes an arbitrary nearest neighbor site to the origin 0. Collecting the above results, one finds K AB ( t ) = 2 Lzα A λ A (cid:16) − (cid:12)(cid:12) h ψ B | e − i t ( ˆ H B − V (ˆ n B − ˆ n B )) | ψ B i (cid:12)(cid:12) (cid:17) + O ( λ A ) . (B.11)For a 1D-ring lattice, this formula reduces to Eq. (21 ) in the main text.Recall that Eq. (B.11) has been obtained by neglecting the kinetic Hamiltonian ˆ K A of the gas A inthe dynamics following the interaction quench. Since K AB ( t ) is of order λ A , this is justified a posteriori provided that the product J A t is much smaller than λ A . In fact, by including in our calculations thecorrections due to ˆ K A estimated from time-dependent perturbation theory, one finds that such correctionsare of order ( J A t ) and λ A ( J A t ) (all terms proportional to J A vanish due to h ψ MI | ˆ K A | ψ MI i = 0). Hencethe expression in the RHS of (B.11) approximates well K AB ( t ) for times satisfying t ≪ U − A , as one maysuspect from the fact that neither J A nor U A appear in this expression.Let us note that it is straightforward to extend the arguments of Sect. 7.1 to the case of an arbitrarylattice Λ. One finds that when the gas B is in the MI regime, the time-averaged Schmidt number reads hK AB i t = 8 Lα A α B z (2 z − λ A λ B (B.12)for times satisfying (61). Moreover, Eq. (64) reads for an arbitrary lattice1 L K AB ( t ) = 4 α A α B zλ A λ B h z − − z −
1) cos( tV ) cos( tU B ) − tV ) cos( tU B ) i (B.13)if t ≪ U − A , J − B , V U − B . eferences [1] K. C. Wright, R. B. Blakestad, C. J. Lobb, W. D. Phillips, and G. K. Campbell, Phys. Rev. Lett. , 025302 (2013).[2] C. Ryu, P. W. Blackburn, A. A. Blinova, and M. G. Boshier, Phys. Rev. Lett. , 205301 (2013).[3] S. Eckel, J. G. Lee, F. Jendrzejewski, M. Murray, C. W. Clark, C. J. Lobb, W. D. Phillips, M.Edwards, and G. K. Campbell, Nature (London) , 200 (2014).[4] B. T. Seaman, M. Kr¨amer, D. Z. Anderson, and M. J. Holland, Phys. Rev. A , 023615 (2007)[5] L. Amico and A. M. Boshier, in: R. Dumke, Roadmap on Quantum Optical Systems , Vol. (Journal of Optics, 2016), p. 093001[6] D. S. Deaver and W. M. Fairbank, Phys. Rev. Lett. , 43 (1961)[7] N. Byers and C. N. Yang, Phys. Rev. Lett. , 46 (1961)[8] L. Onsager, Phys. Rev. Lett. , 50 (1961)[9] L. P. L´evy, G. Dolan, J. Dunsmuir, and H. Bouchiat, Phys. Rev. Lett. , 2074 (1990)[10] J. Dalibard, F. Gerbier, G. Juzeliunas, and P. ¨Ohberg, Rev. Mod. Phys. , 1523 (2011).[11] Y. Kagan, N. V. Prokofev, and B. V. Svistunov, Phys. Rev. A , 045601 (2000)[12] T. Haug, J. Tan, M. Theng, R. Dumke, L-C. Kwek, and L. Amico, Phys. Rev. A , 013633 (2018)[13] M. Cominotti, D. Rossini, M. Rizzi, F. Hekking, and A. Minguzzi, Phys. Rev. Lett. , 025301(2014)[14] D. Aghamalyan, M. Cominotti, M. Rizzi, D. Rossini, F. Hekking, A. Minguzzi, L. C. Kwek, R.Dumke, and L. Amico, New. J. Phys. , 045023 (2015)[15] A. J. Leggett, Prog. Theoret. Phys. Suppl. , 80 (1980)[16] D. W. Hallwood, T. Ernst, and J. Brand, Phys. Rev. A , 063623 (2010)[17] D. Solenov and D. Mozyrsky, Phys. Rev. A , 061601(R) (2010)[18] A. Nunnenkamp, A. M. Rey, and K. Burnett, Phys. Rev. A , 053604 (2011)[19] C. Schenke, A. Minguzzi, and F. W. J. Hekking, Phys. Rev. A , 053636 (2011)[20] P. Naldesi, J. P. Gomez, V. Dunjko, H. Perrin, M. Olshanii, L. Amico, and A. Minguzzi,arXiv:1901.09398[21] G. Modugno, M. Modugno, F. Riboli, G. Roati, and M. Inguscio, Phys. Rev. Lett. , 190404(2002).[22] K. G¨unter, T. St¨oferle, Henning Moritz, M. K¨ohl, and T. Esslinger, Phys. Rev. Lett. , 180402(2006).[23] S. B. Papp, J. M. Pino, and C. E. Wieman, Phys. Rev. Lett. , 040402 (2008).[24] J. Catani, L. De Sarlo, G. Barontini, F. Minardi, and M. Inguscio, Phys. Rev. A , 011603(R)(2008). 2925] G. Thalhammer, G. Barontini, L. De Sarlo, J. Catani, F. Minardi, and M. Inguscio, Phys. Rev.Lett. , 210402 (2008).[26] S. Sugawa, R. Yamazaki, S. Taie, Y. Takahashi, Phys. Rev. A , 011610(R) (2011).[27] D. J. McCarron, H. W. Cho, D. L. Jenkin, M. P. K¨oppinger, and S. L. Cornish, Phys. Rev. A ,011603(R) (2011).[28] I. Ferrier-Barbut, M. Delehaye, S. Laurent, A. T. Grier, M. Pierce, B. S. Rem, F. Chevy, and C.Salomon, Science , 1035 (2014).[29] C. R. Cabrera, L. Tanzi, J. Sanz, B. Naylor, P. Thomas, P. Cheiney, and L. Tarruel, Science ,301 (2017).[30] A. Kuklov and B. Svistunov, Phys. Rev. Lett. , 100401 (2003).[31] A. Kuklov, N. Prokof’ev, and B. Svistunov, Phys. Rev. Lett. , 030403 (2004).[32] A. Kuklov, N. Prokof’ev, and B. Svistunov, Phys. Rev. Lett. , 050402 (2004).[33] C. Menotti and S. Stringari, Phys. Rev. A , 045604 (2010).[34] Y. Li, M. R. Bakhtiari, L. He, and W. Hofstetter, Phys. Rev. B , 144411 (2011).[35] Y. Li, L. He, and W. Hosfstetter, New J. Phys. , 093028 (2013).[36] L. Morales-Molina, S. A. Reyes, and E. Arevalo, EPL , 36004 (2016).[37] V. Penna and A. Richaud, J. of Phys: Conf. Series , 012011 (2019).[38] I. Morera, G. E. Astrakharchik, A. Polls, and B. Juli´a-D´ıaz, Phys. Rev. Research , 022008(R)(2020).[39] S. I Mistakidis, G. C. Katsimiga, P. G. Kevrekidis, and P. Schmelcher, New J. Phys. , 043052(2018)[40] F. Meinert, M. Knap, E. Kirilov, K. Jag-Lauber, M. B. Zvonarev, E. Demler, and H.-C. N¨agerl,Science , 945 (2017)[41] F. Grusdt, G. E. Astrakharchik, and E. Demler, New J. Phys. , 103035 (2017)[42] A. G. Volosniev and M.-W. Hammer, Phys. Rev. A , 031601(R) (2017)[43] S. I. Mistakidis, G. C. Katsimiga, G. M. Koutenkakis, Th. Busch, and P. Schmelcher, Phys. Rev.Lett. , 183001 (2019)[44] F. Theel, K. Keiler, S. I. Mistakidis, and P. Schmelcher, New J. Phys. , 023027 (2020)[45] F. Meinert, M. J. Mark, E. Kirilov, K. Lauber, P. Weinmann, A. J. Daley, and H.-C. N¨agerl, Phys.Rev. Lett. , 053003 (2013).[46] F. Meinert, M. J. Mark, K. Lauber, A. J. Daley, and H.-C. N¨agerl, Phys. Rev. Lett. , 205301(2016).[47] S. Greschner, G. Sun, D. Poletti, and L. Santos, Phys. Rev. Lett. , 215303 (2014).[48] S. Campbell, M.A. Garc´ıa-March, T. Fogarty, and T. Busch, Phys. Rev. A , 013617 (2014).3049] L. Morales-Molina, S. A. Reyes, and M. Orszag, Phys. Rev. A , 033629 (2012).[50] D. Spehner, J. Math. Phys. , 075211 (2014).[51] R. Horodecki, P. Horodecki, M. Horodecki, and K. Horodecki, Rev. Mod. Phys. , 865 (2009).[52] S.J. van Enk and C.W.J. Beenakker, Phys. Rev. Lett. , 110503 (2012).[53] A. Elben, B. Vermersch, C.F. Roos, and P. Zoller, Phys. Rev. A , 052323 (2019).[54] T. Brydges et al. , Science , 260 (2019)[55] L. Corman, L. Chomaz, T. Bienaim´e, R. Desbuquois, C. Weitenberg, S. Nascimb`ene, J. Dalibard,and J. Beugnon, Phys. Rev. Lett. , 135302 (2014)[56] S. Eckel, F. Jendrzejewski, A. Kumar, C. J. Lobb, and G. K. Campbell, Phys. Rev.X , 031052(2014)[57] R. Mathew, A. Kumar, S. Eckel, F. Jendrzejewski, G. K. Campbell, M. Edwards, and E. Tiesinga,Phys. Rev. A , 033602 (2015)[58] F. Gerbier, A. Widera, S. F¨olling, O. Mandel, T. Gericke, and I. Bloch, Phys. Rev. Lett. , 050404(2005); ibid , Phys. Rev. A , 053606 (2005).[59] K. Sengupta, N. Dupuis, Phys. Rev. A , 033629 (2005)[60] L. Amico, D. Aghamalyan, F. Auksztol, H. Crepaz, R. Dumke, and L. C. Kwek, Sci. Rep. , 4298(2014)[61] L. Morales-Molina, S. A. Reyes, and E. Arevalo, EPL , 36001 (2020).[62] R. Schmied and P. Treutlein, New J. Phys. , 065019 (2011)[63] B. P. Lanyon et al. , Nature Physics13