Quantum boomerang effect for interacting particles
Jakub Janarek, Dominique Delande, Nicolas Cherroret, Jakub Zakrzewski
QQuantum boomerang effect for interacting particles
Jakub Janarek,
1, 2, ∗ Dominique Delande, † Nicolas Cherroret, ‡ and Jakub Zakrzewski
1, 3, § Institute of Theoretical Physics, Jagiellonian University in Krakow, (cid:32)Lojasiewicza 11, 30-348, Krak´ow, Poland Laboratoire Kastler Brossel, Sorbonne Universit´e, CNRS,ENS-Universit´e PSL, Coll`ege de France; 4 Place Jussieu, 75004 Paris, France Mark Kac Complex Systems Centre, Jagiellonian University in Krakow, (cid:32)Lojasiewicza 11, 30-348, Krak´ow, Poland (Dated: June 9, 2020)When a quantum particle is launched with a finite velocity in a disordered potential, it may sur-prisingly come back to its initial position at long times and remain there forever. This phenomenon,dubbed “quantum boomerang effect”, was introduced in [Phys. Rev. A , 023629 (2019)]. In-teractions between particles, treated within the mean-field approximation, are shown to partiallydestroy the boomerang effect: the center of mass of the wave packet makes a U-turn, but does notcompletely come back to its initial position. We show that this phenomenon can be quantitativelyinterpreted using a single parameter, the average interaction energy. I. INTRODUCTION
Anderson Localization (AL) [1], i.e. inhibition of trans-port in disordered media, has been the source of vari-ous, often counter-intuitive phenomena discovered overthe last 60 years [2]. Already in one dimension (1D)the fact that even the tiniest random disorder generi-cally leads to a full localization of eigenstates is totallyagainst a classical way of thinking. This effect is a clearmanifestation of the inherently interferometric nature ofAL, typically explained as the effect of quantum waveinterference. The phenomenon was observed in many ex-perimental setups including light [3, 4], sound waves [5]as well as matter waves [6–10]. The closely related phe-nomenon of Aubry-Andr´e localization in quasi-periodicpotentials [11] was also observed in a cold atomic setting[12]. Last years have also led to a number of studies ofthe many-body counterpart of AL [13, 14], the so calledmany-body localization (MBL) (for recent reviews see[15, 16]). It has already been observed in cold atomicexperiments in quasi-periodic potentials [17, 18]. Whilestudies of MBL have been very extensive, even its veryexistence has been questioned recently [19], provoking avivid debate [20–22].The physics of a single particle in a random poten-tial, particularly in 1D, has much stronger foundations,while still bringing novel features such as studies of ran-dom (or quasi-random) potentials revealing the presenceof mobility edge [23, 24], i.e. situations where local-ized and extended states appear at different energies.Even for a pure random standard case, one may findnew counter-intuitive phenomena as exemplified by the quantum boomerang effect [25]. As a classical boomerangreturns to the initial position from which it was launched,the center of mass of a wave packet with a nonzero ini-tial velocity returns to its origin due to AL. The effect ∗ [email protected] † [email protected] ‡ [email protected] § [email protected] is quite general and occurs in Anderson localized multi-dimensional systems [25] whose Hamiltonians preservetime-reversal invariance.The aim of this work is to investigate how interac-tions between particles affect the boomerang effect. Ourstudy focuses on the limit of weak interactions, a regimewhere it was previously shown that AL of wave pack-ets is replaced by a subdiffusive evolution at very longtime [26–30] - see however [31, 32]. Computing the tem-poral evolution of a many-body wave packet in a disor-dered potential is in general a formidable task, even in1D. When interactions are sufficiently weak however, onemay use a mean-field approximation and describe the dy-namics with a one-dimensional Gross-Pitaevskii equation(GPE). Within this formalism, we provide an in depthnumerical analysis of the quantum boomerang effect. Forsimplicity, we study the weak-disorder case, for which ananalytic description of the boomerang effect in the ab-sence of interaction is available. At short times, of theorder of a few disorder scattering times, we observe thatthe dynamics of the wave packet is essentially not affectedby interactions, with an initial ballistic flight followed bya U-turn of the center of mass, slowly returning towardsits initial position. The main effect in this regime is asmall modification of the scattering time and scatteringmean free path due to interactions. At longer time, weobserve that the center of mass, instead of returning to itsinitial position, stops at a finite distance from it, whichincreases with the interaction strength. We show thatthis phenomenon can be interpreted in terms of a breaktime, the time scale beyond which interference effects aretypically destroyed by interactions. We finally show thatin the presence of interactions, the boomerang effect canbe quantitatively described in terms of a single scalingparameter, a variant of the interaction energy, computedat the break time.The paper is organized as follows. After formulatingthe problem and introducing the main parameters in Sec.II, we numerically analyze the influence of interactionson the boomerang effect in Sec. III. In Sec. IV, we thenperform several numerical studies which establish thatthe dynamics can be described by a single parameter. a r X i v : . [ c ond - m a t . d i s - nn ] J un The appendix shows that the short-time dynamics canalso be understood using the same parameter. We finallyconclude and briefly discuss open questions.
II. THE MODEL
To study the boomerang effect in the presence of in-teractions, we use the one-dimensional, time dependentGross-Pitaevskii equation (GPE) [33]: i (cid:126) ∂ψ ( x, t ) ∂t = (cid:20) p m + V ( x ) + g | ψ ( x, t ) | (cid:21) ψ ( x, t ) , (1)where V ( x ) is a disordered potential, while g representsthe interaction strength. Wave functions are normalizedto unity, (cid:82) | ψ ( x, t ) | d x = 1. Throughout this work weassume that the disordered potential is a Gaussian un-correlated random variable, i.e. V ( x ) = 0 , V ( x ) V ( x (cid:48) ) = γδ ( x − x (cid:48) ) , (2)where the overbar denotes the average over disorder real-izations. The parameter γ measures the disorder strengthand determines the characteristic time and length scalesfor scattering: the mean scattering time τ (the typi-cal time for a particle to be scattered by the disorder)and the corresponding mean scattering length (cid:96) , usu-ally called mean free path. The uncorrelated disordermodel is sufficient to capture the main features of theboomerang effect, which depends only on a small set ofwell defined parameters, τ , (cid:96) and time t itself. In theBorn approximation, τ and (cid:96) are given by [34]: τ = (cid:126) k mγ , (cid:96) = (cid:126) k m τ = (cid:126) k m γ . (3)Non-interacting 1D systems remain always strongly lo-calized, i.e. their eigenstates decay exponentially overthe localization length ξ loc = 2 (cid:96) [35].Following [25], as the initial state for time evolution wetake a Gaussian wave packet with mean velocity (cid:126) k /m : ψ ( x, t = 0) = (cid:18) πσ (cid:19) / exp( − x / σ + ik x ) , (4)where σ and k are chosen such that the initial wave func-tion is sharply peaked around (cid:126) k in momentum space,i.e. k σ (cid:29)
1. Moreover, we focus on the weak-disorderregime where the mean free path is much longer thanthe de Broglie wavelength, that is k (cid:96) (cid:29) (cid:104) x ( t ) (cid:105) = (cid:90) x | ψ ( x, t ) | d x , (5) as: (cid:104) x ( t ) (cid:105) (cid:96) = f (cid:18) tτ (cid:19) , (6)where f is a universal function whose Taylor expansionat short time and asymptotic behavior at long time areexactly known [25]. In particular (cid:104) x ( t ) (cid:105) ≈ t(cid:96) /τ at shorttime t (cid:28) τ , meaning that the initial motion is ballistic,and (cid:104) x ( t ) (cid:105) ≈ (cid:96) ln( t/ τ ) τ t (7)at long time t (cid:29) τ . The two assumptions of weak dis-order, k (cid:96) (cid:29) , and narrow wave packet in momentumspace, k σ (cid:29) , make it possible to have a well con-trolled non-interacting limit [37]. If they are not valid,the wave packet will contain a broad energy distribution,and consequently a distribution of scattering time (whichdepends on energy), making the analysis more tedious.However, the phenomena described below are expectedto be very similar.Another important ingredient of the model are inter-actions. Our interest is only in small values of g , whichcorresponds to a weak interaction regime [38–40].In the following, we study numerically the propaga-tion of the initial wave function ψ ( x, t = 0) for differ-ent disorder realizations. It yields the averaged den-sity profile | ψ ( x, t ) | , from which we compute the CMP,Eq. (5). The numerical technique is as follows. Theone-dimensional configuration space is discretized on aregular grid, over which the wave function is computed.The temporal propagation is performed using a split-step algorithm of step ∆ t , alternating propagations ofthe linear part of the GPE, exp (cid:2) − i ( p / m + V )∆ t/ (cid:126) (cid:3) ,and of the nonlinear part, which is simply a phase factorexp( − ig | ψ ( x, t ) | ∆ t/ (cid:126) ) . The linear part of the evolutionoperator is expanded in a series of Chebyshev polynomi-als, as described in [41–44].Throughout our work, we express lengths in units of1 /k , times in m/ (cid:126) k , and energies in (cid:126) k /m . The in-teraction strength is expressed in units of (cid:126) k /m , andthe disorder strength γ in units of (cid:126) k /m . Numerical results have been obtained on a discretizedgrid of size 4000 /k (sufficiently large for the wave func-tion to be vanishingly small at the edges; open boundaryconditions have been used) divided into 20000 points, sothat discretization effects are negligible. We have used adisorder strength γ = 0 . (cid:126) k /m , so that, in the Bornapproximation (3), one has k (cid:96) = 5 , i.e. the disorder isweak. Note that, while the true mean free path and scat-tering time may slightly differ from their expressions (3)at the lowest order Born approximation, the dynamicsof the non-interacting quantum boomerang effect is stillgiven by Eq. (6), provided corrected values of τ and (cid:96) are used. We have performed calculations for variouswidths σ of the initial wave packet ψ ( x, t = 0) , as well asvarious values of the interaction strength g . t h x i g = 0 . g = 0 . g = 1 . g = 1 . g = 2 . g = 2 . FIG. 1. CMP time evolution for an initial wave packet with σ = 10 /k , for different values of the interaction strength g .Here, (cid:104) x (cid:105) is in units of 1 /k and time t in units of m/ (cid:126) k .The legend indicates curves from bottom to top. All curveshave been averaged over more than 5 × disorder realiza-tions. The short-time behavior remains almost unchanged,whereas the long-time evolution clearly depends on the inter-action strength. The error bars represent statistical averageerrors. Center-of-mass trajectories among different disorderrealizations are normally distributed, such that we use thestandard error of the mean as estimator of the errors. g . . . . . h x i ∞ FIG. 2. Long-time average (cid:104) x (cid:105) ∞ vs. interaction strength g for a wave packet with initial width σ = 40 /k . (cid:104) x (cid:105) ∞ is inunits of 1 /k , and g is in units of (cid:126) k /m . Data have beenaveraged over more than 5 × disorder realizations. III. THE BOOMERANG EFFECT WITHINTERACTIONSA. Role of the interaction strength
In Fig. 1 we present results obtained for a wavepacket with initial width σ = 10 /k for various inter- action strengths g . Similarly to the non-interacting case,after the initial ballistic motion of the center of mass,we observe a subsequent reflection towards the origin,that is a boomerang effect. However, the long-time be-havior is affected by interactions. For non-zero valuesof g, the center of mass does not return to the originbut saturates at some finite value. In other words, theboomerang effect is only partial [45]. This can be under-stood as follows. For g = 0 , the disorder is static andfull Anderson localization sets in; the complicated in-terference between multiply scattered paths leads to fulllocalization at infinite time, and also to full return of thecenter of mass to its initial position. For non-zero butsmall g , the nonlinear term in Eq. (1), g | ψ ( r, t ) | , playsthe role of a small additional effective potential which is time-dependent , thereby adding a fluctuating phase alongeach scattering path. This breaks the interference be-tween multiple scattering paths and thus destroys bothAnderson localization [26–30] and the full boomerang ef-fect at long time. For all studied widths of the wavepacket, namely k σ = 5 , ,
20 and 40 we observe asimilar saturation effect. Of course, the phase scram-bling progressively develops over time. The characteris-tic break time over which it kills coherent transport andboomerang effect is discussed in detail in the sequel ofthis paper. Note that the interpretation of the effect ofinteractions in terms of a decoherence mechanism may bequestioned at very long time, where thermalization comesinto play and affects the momentum distribution as wellas the dynamics of the system [46, 47]. In our case, bothdisorder and interactions are small, so that thermaliza-tion takes place at times significantly longer than the onesconsidered here. Neglecting thermalization neverthelessrestricts our analysis to the regime where the long-timeCMP (cid:104) x (cid:105) ∞ (cid:28) (cid:96) (see Eq. (8) below), which constrainsthe maximum value of g .To study in detail the long-time evolution, we run nu-merical simulations up to time t max ≈ τ . Fromthese simulations we calculate the long-time average ofthe CMP, (cid:104) x (cid:105) ∞ , defined in the following way: (cid:104) x (cid:105) ∞ = 1 t − t (cid:90) t t (cid:104) x ( t ) (cid:105) d t , (8)where we choose t ≈ τ , t ≈ τ . The re-sults are essentially independent of these bounds, pro-vided they are much longer than τ . This definition pro-vides us with a very good estimate of the infinite-timeCMP. Figure 2 shows the dependence of (cid:104) x (cid:105) ∞ on the in-teraction strength g . For small values of g, the CMPdependence is quadratic in g , (cid:104) x (cid:105) ∞ ∝ g , and becomesapproximately linear for larger g . In the following, wewill mostly concentrate on the quadratic regime of smallinteractions, leaving the more difficult case of larger g forfuture studies. h x i a) no interaction σ = 10 σ = 20 σ = 40 t h x i b) g = 2 . FIG. 3. a) Comparison of non-interacting CMP (cid:104) x ( t ) (cid:105) forwave packets of widths σ = 10 /k , /k , /k . All threecurves overlap, indicating that (cid:104) x ( t ) (cid:105) is σ -independent in thenon-interacting limit. b) Same as a) but for non-zero inter-action strength g = 2 . (cid:126) k /m (the legend indicates curvesfrom top to bottom). Here, the saturation point of (cid:104) x ( t ) (cid:105) ishigher for initially narrower wave packets. The center of mass (cid:104) x ( t ) (cid:105) is in units of 1 /k and time t in units of m/ (cid:126) k . Errorbars are only shown in panel a) for σ = 10 /k to indicatetheir order of magnitude. B. Role of the wave-packet width
Another important parameter is the width σ of theinitial wave packet. In Fig. 3, we show (cid:104) x ( t ) (cid:105) for dif-ferent σ values. While, for g = 0 , it follows the an-alytic prediction (7) independently of σ, for interactingparticles the long-time behavior strongly depends on σ. A simple qualitative explanation is that the destructionof the boomerang effect is controlled by the nonlinearterm g | ψ ( x ) | in Eq. (1). This term is larger for a spa-tially narrow wave packet, so that interference betweenscattered waves are suppressed at shorter time, giving ahigher (cid:104) x (cid:105) ∞ value.Although the boomerang effect is affected by a changeof either g or σ in the interacting case, one may guessthat the CMP is not a function of these two independentparameters. Indeed, closer investigation reveals that sim-ilar “trajectories” of (cid:104) x ( t ) (cid:105) can be achieved for differentcombinations of g and σ . In particular, the same valuesof (cid:104) x (cid:105) ∞ can be obtained from different initial states, pro-vided g is properly adjusted. This property is illustratedin Fig. 4, where we have computed (cid:104) x ( t ) (cid:105) for various val-ues of σ and have adjusted g so to that the curves fall ontop of each other. t h x i σ = 5 , g = 0 . σ = 10 , g = 1 . σ = 20 , g = 1 . σ = 40 , g = 2 . FIG. 4. CMP (cid:104) x ( t ) (cid:105) vs. time for different initial states chosensuch that all of them saturate around the value (cid:104) x (cid:105) = 0 . /k .All the curves overlap. Results have been averaged over 16000disorder realizations. The CMP is shown in units of 1 /k andtime is in units of m/ (cid:126) k . Error bars are shown only for σ = 5 /k . IV. UNIVERSAL SCALING OF THEINTERACTING BOOMERANG EFFECT
From the results of the previous section, it is natural toask whether or not the interacting boomerang effect, andmore specifically its long-time average (cid:104) x (cid:105) ∞ – Eq. (8) –can be described in terms of a single parameter, in thespirit of scaling approaches well-known in the context ofAnderson localization of non-interacting particles [48]. A. Break time
Before attempting to rescale the CMP, let us introducean important parameter that will turn useful in the fol-lowing. We recall that in the non-interacting limit, theCMP at long time is given by Eq. (7). If we neglect thelogarithmic part, (cid:104) x ( t ) (cid:105) decays as t − . It suggests thatone may identify a characteristic time scale connectedwith weak interactions which is inversely proportional to g . In our analysis we call this time scale the break time t b and define it by the relation: (cid:104) x ( t b ) (cid:105) g =0 = (cid:104) x (cid:105) ∞ ( g ) , (9)where for the left-hand-side we use the analytical predic-tion of the non-interacting theory, Eq. (7). A given timescale in quantum mechanics corresponds to some energyscale, here the break energy : E b = 2 π (cid:126) t b , (10)which will turn out to be a key parameter in our rescalingof the CMP. .
00 0 .
01 0 .
02 0 .
03 0 . E int ( t = 0) . . . . . h x i ∞ σ = 10 σ = 40 σ = ∞ FIG. 5. Dependence of the long-time average (cid:104) x (cid:105) ∞ , Eq.(8), on the initial interaction energy. For wave packets E int ( t = 0) = g/ (2 √ πσ ), while E int ( t = 0) = gρ / σ = ∞ . The interaction energy is expressedin units of (cid:126) k /m , and the CMP in units of 1 /k . Error barsrepresent standard deviation of the averaged points. .
00 0 .
01 0 .
02 0 .
03 0 . E int ( t = 0) . . . . . . E b σ = 10 σ = 40 σ = ∞ FIG. 6. Break energy E b vs. initial interaction energy E int ( t = 0). The break energy is defined by E b = 2 π (cid:126) /t b ,where the break time t b , calculated according to Eq. (9), isthe characteristic time beyond which the boomerang effectdisappears. Energies are in units of (cid:126) k /m . Error bars onbreak energy points represent error on break time calculationfrom long time average of CMP. B. Rescaling of the CMP: first attempt
A closer inspection of Fig. 4 reveals that, when tryingto superimpose the CMP curves, broader wave packetsrequire higher interaction strengths. A first natural can-didate to characterize the CMP is thus the interaction energy, E int ( t ) = g (cid:90) | ψ ( x, t ) | d x . (11)We recall that the total energy, conserved by the GPE, isthe sum of the non-interacting part and the interactionenergy: E tot = (cid:104) p (cid:105) / m + (cid:104) V (cid:105) + E int . A related im-portant quantity is the interaction energy at initial time, E int ( t = 0) = g/ (2 √ πσ ). Notice that while the CMPcurves in Fig. 4 are obtained from quite different val-ues of g and σ , they are all associated with comparableinteraction energies at time t = 0. This is a very clearhint that E int ( t = 0) is a crucial parameter for describingthe impact of interactions on the boomerang effect. Wecan thus try to rescale the results by plotting them vs.the initial interaction energy. This is done for (cid:104) x (cid:105) ∞ inFig. 5, and for the break energy – see Eq. (10) – in Fig.6. In these plots we also show, for comparison, the for-mal limit σ → ∞ of infinitely large wave packets, where ψ ( x ) → √ ρ exp( ik x ) reduces to a plane wave with E int ( t = 0) = gρ /
2. In this limit, the sole parameter gρ controls the boomerang effect. Note that despite theflatness of the density profile in the plane-wave limit, inpractice it is still possible to study the boomerang effectby measuring an effective CMP defined as: (cid:104) x ( t ) (cid:105) σ = ∞ ≡ m (cid:90) t (cid:104) p ( t (cid:48) ) (cid:105) d t (cid:48) . (12)We have verified numerically, in particular, that in thenon-interacting case the definition (12) agrees with re-sults for wave packets of finite width, thus with the the-oretical prediction (7).The curves in Figs. 5 and 6, which correspond to dif-ferent values of σ , are qualitatively similar, despite thefact that the g values are widely different. This sug-gests that the interaction energy is indeed an importantparameter. Moreover, note that the break energy E b iscomparable (within a factor 4) to the initial interactionenergy. In particular, because E int ( t = 0) is proportionalto the interaction strength g for all initial states, at smallvalues of E int ( t = 0) we see the expected linear depen-dence of E b with E int ( t = 0). Nevertheless, it is clearthat plots based on E int ( t = 0) do no fall on the sameuniversal curve: (cid:104) x (cid:105) ∞ and E b deviate from each other as σ is changed, approaching the upper limit curve σ = ∞ as σ increases (this limit is represented by green squaresymbols in Figs. 5 and 6). C. Nonlinear energy at break time
The reason why E int ( t = 0) does not allow for a uni-versal rescaling of the CMP stems from the fact that theinteraction energy E int ( t ) varies significantly from t = 0onwards. This evolution is shown in Fig. 7 for two finitevalues of σ and for the plane-wave limit σ = ∞ .The figure reveals that the time evolution of the in-teraction energy generally consists of two stages. In a t . . . . . . . . E i n t ( t ) σ = 10 σ = 40 σ = ∞ FIG. 7. Interaction energy vs. time for different initial statesand interaction strengths. The plot shows data for wave pack-ets with [ σ = 10 /k , g = 0 . (cid:126) k /m ], [ σ = 40 /k , g =1 . (cid:126) k /m ], and [ σ = ∞ , g = 45 . (cid:126) k /m, ρ = 0 . k ]. E int is in units of (cid:126) k /m and time t in units of m/ (cid:126) k . Thelegend indicates curves from bottom to top. Initial states arechosen to clearly emphasize the randomization of the wavefunction amplitude at long time. It is, however, genericallypresent for any interaction strength. In the plane-wave limit,randomization doubles E int over a short time scale compara-ble to the scattering time τ , which then remains stationaryat long time. For finite σ , randomization is is visible at shorttime, but followed by a slow decay at long time, because ofwave-packet spreading. From t = t b onwards, however, thedecay is very slow (the location of t b is indicated by arrows). first stage, which takes place on a time scale of a fewscattering times, the interaction energy roughly doubles.This can be understood by noticing that E int ( t ) dependson the fourth moment of the field, see Eq. (11), whichobeys: | ψ ( x, t ) | = | ψ ( x, t ) | + Var (cid:16) | ψ ( x, t ) | (cid:17) . (13)In the plane-wave limit where the profile is flat, the ini-tial density variance is zero. During the temporal evolu-tion however, the density develops fluctuations depend-ing on the realization of the disorder, which makes theinteraction energy increase. The factor 2 enhancement isobtained by assuming that, after a few scattering times, ψ ( x, t ) is a complex Gaussian random variable. The vari-ance in Eq. (13) is then | ψ ( x, t ) | , so that | ψ ( x, t ) | is doubled. This randomization of the complex wavefunction amplitude is similar to the appearance of op-tical speckles in scattering media [49]. It implies that E int ( t (cid:29) τ ) = 2 E int ( t = 0) in the limit σ = ∞ . Forwave packets of finite size σ, the effect is also present,albeit slightly smaller [50].In the plane-wave limit σ = ∞ , the interaction en-ergy remains constant once the randomization processhas ended. For finite σ on the contrary, a second stage occurs, where E int decreases in time, see Fig. 7. This de-crease is due to the spreading of the wave packet, whichbecomes more and more dilute. The spreading is initiallyfast, and then quickly slows down.The time evolution of E int makes a detailed rescal-ing analysis of the boomerang effect very complex. Thecurves in Fig. 7, however, suggest the simple idea ofusing as scaling parameter the interaction energy at thebreak time t b , instead of the initial interaction energy.Although at finite σ such a rescaling can only be approx-imate, since wave packets keep evolving slowly in timebeyond t b (indicated as arrows in Fig. 7), we show belowthat it provides satisfactory results.Before applying this strategy, a final adjustment mustbe performed. The quantum boomerang being a dynam-ical effect governed by the GPE (1), its evolution is notstrictly governed by the interaction energy, g | ψ ( x, t ) | / nonlinear energy , defined as E NL ( t ) = 2 E int ( t ) = g (cid:90) | ψ ( x, t ) | d x . (14)In the appendix, we show that the E NL is related to thenonlinear self energy and discuss the dynamical behaviorof the system at short time. We show analytically andnumerically that it is indeed E NL , rather than E int , whichgoverns the evolution. D. Rescaling of the boomerang effect
We can now re-analyze the boomerang effect using thenonlinear energy at the break time, E NL ( t = t b ), as acontrol parameter of interactions. We show in Fig. 8the long-time average (cid:104) x (cid:105) ∞ of the CMP as a functionof the nonlinear energy calculated at the break time. Incontrast with Fig. 5, now all points collapse on a sin-gle universal curve. As expected, in the regime of smallnonlinear energy, (cid:104) x (cid:105) ∞ shows a quadratic dependence.In Fig. 9 we also compare the break energy E b =2 π (cid:126) /t b with the nonlinear energy at the break time t b for wave packets of size σ = 10 /k , /k and σ = ∞ ,for various values of g . There is a compelling evidencethat these quantities are very similar. A small differenceis observed in the plane-wave limit, which we attributeto the slow residual decay of (cid:104) x (cid:105) due to early stage ofthermalization, which leads to an underestimation of thebreak energy. Such a good agreement shows that a simplemodel involving a single parameter, the nonlinear energyat the break time, captures quantitatively the essentialfeatures of a complex dynamical process like the quantumboomerang effect for interacting particles. As discussedin the appendix, the same parameter turns out to alsocontrol the short-time dynamics and the correction ofthe scattering time due to interactions. .
00 0 .
01 0 .
02 0 .
03 0 .
04 0 . E NL ( t = t b ) . . . . . h x i ∞ σ = 10 σ = 40 σ = ∞ FIG. 8. Long-time averages (cid:104) x (cid:105) ∞ of the CMP for differentinitial states vs. the nonlinear energy computed at the breaktime, E NL ( t = t b ). (cid:104) x (cid:105) ∞ is expressed in units of 1 /k , whileenergy is in units of (cid:126) k /m . All results have been averagedover more than 5 × disorder realizations. In this repre-sentation, data fall on the same single curve. . . . . . . . . a) E NL ( t = t b ) E b . . . e n e r g y b) g . . . c) FIG. 9. Comparison of the break energy, E b , and the nonlin-ear energy at break time, E NL ( t = t b ), for increasing values ofthe interaction strength. Results are shown for wave packetsof size a) σ = 10 /k , b) σ = 40 /k and c) σ = ∞ . En-ergy is in units of (cid:126) k /m , interaction strength g in units of (cid:126) k /m . All results have been averaged over 5 × disorderrealizations. V. CONCLUSION
We have analyzed the quantum boomerang phe-nomenon in the presence of interactions on the basis ofthe one-dimensional Gross-Pitaevskii equation. We havefound that weak interactions do not destroy the quan-tum boomerang effect, in the sense that the center ofmass of a wave packet launched with a finite velocity isstill retro-reflected after a few scattering times, slowlyreturning towards its initial position. The boomerang ef-fect is only partial though, as the quantum particle doesnot return to its initial position but stops on the wayback. We have interpreted this phenomenon as a conse-quence of the destruction of interference between multiplescattering paths induced by a time-dependent nonlinearphase acquired along a path. To characterize this phe-nomenon, we have introduced a break time – and thecorresponding break energy – the characteristic time be-yond which the destruction of interference prevents thewave packet to further move back to its initial position.We have finally shown that different initial states andinteraction strengths can all be described by means of asingle parameter, the nonlinear energy estimated at thebreak time.Our analysis is limited to the regime of weak disor-der and weak interactions. When the disorder strengthis increased, the quantitative description becomes morecomplicated, but the overall conclusions are expected tobe qualitatively identical, provided interactions remainweak. Indeed, a wave-packet in a relatively strong disor-der contains many energy components, each energy beingcharacterized by a scattering time and a mean free path.In the non-interacting limit, each energy component willdisplay a boomerang effect described by Eq. (7), but thesuperposition of various energy components will lead toa complicated (cid:104) x ( t ) (cid:105) function. In the presence of weakinteractions, each energy component will display a par-tial boomerang effect, to that (cid:104) x (cid:105) ∞ is again likely to benonzero. Whether an effective break time can be definedin such a case is an open question. For weak disorderand stronger interactions, the break time is likely to de-crease until it becomes comparable to the scattering time.Whether a single parameter also controls this regime isan interesting question left for future studies. Anotherimportant question – especially if one envisions experi-ments with ultra-cold atoms – is to know whether theobserved softening of the boomerang effect due to inter-actions remains valid beyond the mean-field description.Studying the full many-body boomerang effect is a chal-lenging task. ACKNOWLEDGMENTS
We kindly thank PL-Grid Infrastructure for provid-ing computational resources. JJ and JZ acknowledgesupport of National Science Centre (Poland) under theproject OPUS11 2016/21/B/ST2/01086. NC acknowl- t . . . . . . . R e Σ ( g ) τ g = 0 . g = 15 . g = 30 . g = 45 . g = 60 . g = 75 . FIG. 10. Numerically calculated real parts of the self energyΣ ( g ) (solid lines) for plane waves vs. time for several valuesof the interaction strength. In the plot we additionally shownonlinear energies E NL ( t ) (dashed lines), and indicate as dot-ted lines the value 2 gρ . Energies are expressed in units of (cid:126) k /m and time in units of m/ (cid:126) k . The legend indicatescurves from bottom to top. The value of the mean scatteringtime τ is indicated by the small black arrow near the timeaxis. edges financial support from the Agence Nationale de laRecherche (grant ANR-19-CE30-0028-01 CONFOCAL).We also acknowledge support of Polish-French bilateralproject Polonium 40490ZE. Appendix A: Short time behavior
We have seen that interactions modify the long-timedynamics of the center of mass position, and that achange of either parameters σ or g can be encompassedin the use of the nonlinear energy. In this appendix, weshow that the very same concept – the nonlinear time-dependent energy – also accurately describes the short-time dynamics of the system, through a change of thereal part of the self-energy and of the scattering time.
1. Self energy in interacting systems
The self energy is a key concept in the description ofdisordered quantum systems. It is a complex valued func-tion of the state energy E and its momentum (cid:126) k describ-ing the disorder-induced energy shift and the exponentialdecay of average Green’s functions in configuration space.It is noted Σ E ( k ) and defined by: G RE ( k ) = 1 E − E − Σ E ( k ) , (A1)where G R denotes the averaged retarded Green’s functionand E = (cid:126) k / m the disorder-free energy. The self- energy vanishes in the absence of disorder and is muchsmaller than E in the weak-disorder limit (by a factor1 /k (cid:96) ) [34].The evolution operator is the temporal Fourier trans-form of the Green’s function. If the self-energy is asmooth function of E , one obtains for the evolution ofa plane wave: (cid:104) ψ | ψ ( t ) (cid:105) = e − i ( E k +Σ E ( k ) ) t/ (cid:126) = e − i ( E k +Re Σ E ( k ) ) t/ (cid:126) e Im Σ E t/ (cid:126) , (A2)so that Re Σ is an energy shift and − Im Σ the decay rateinduced by the disorder. At the Born approximation, wehave in 1D: Σ (0) E ( k ) = − i (cid:126) τ , (A3)where the superscript (0) refers to zero interactions.In the presence of interactions, the situation is in gen-eral much more complicated. Because of the nonlinearityof the GPE, the notion of evolution operator no longerexists and the overlap Eq. (A2) has no reason to be anexponential function of time. However, it is possible todefine an effective self-energy using Eq. (A2), the left-hand-side of the equation being computed numericallyfrom the solution | ψ ( t ) (cid:105) of the GPE. The obtained self-energy Σ ( g ) E ( k ) depends on time.To analyze the impact of interactions on the self en-ergy, we introduce its nonlinear part Σ ( g ) , defined as:Σ ( g ) = Σ ( g ) E ( k ) − Σ (0) E ( k ) , (A4)where both Σ ( g ) E ( k ) and Σ (0) E ( k ) are calculated numer-ically. The real part of this quantity is plotted in Fig.10 in the plane-wave limit σ = ∞ . Σ ( g ) increases over afew mean scattering times and then saturates at roughlytwice its initial value.It is easy to compute Σ ( g ) at t = 0 from the GPE. Theresult is: Σ ( g ) ( t = 0) = E NL ( t = 0) = gρ , (A5)where the first equality is valid for any initial state, whilethe second holds only for a plane wave. At time longerthan the scattering time, the randomization of the wavefunction phenomenon described in Sec. IV C is respon-sible for a doubling of the nonlinear energy. It is thusvery natural, and fully confirmed by the numerical re-sults in Fig. (10) as well as by a theoretical approach [40]to have: Σ ( g ) ( t (cid:29) τ ) = E NL ( t (cid:29) τ ) = 2 gρ . (A6)This close connection between the nonlinear energy andthe nonlinear part of the self-energy also exists at in-termediate times, as shown in Fig. 10, where we alsoplot numerically computed nonlinear energies. After aninitial growth, both Σ ( g ) ( t ) and E NL ( t ) saturate around .
00 0 .
02 0 .
04 0 .
06 0 . E shortNL . . . . . . τ σ = ∞ σ = 40 σ = 20 σ = 10 Eq. (A7)
FIG. 11. Fitted values of the mean scattering time τ for wavepackets of different initial widths vs. the nonlinear energyaveraged over the fit time window, E shortNL . Time is in unitsof m/ (cid:126) k and energy in units of (cid:126) k /m . The black dashedline shows the prediction of Eq. (A7), where the wave packetwith σ = 40 /k is used for τ . gρ and follow a close evolution (even though the growthrate of Re Σ ( g ) is slightly lower than for the nonlinear en-ergy). Altogether, it suggests that Re Σ ( g ) and E NL mayhave a similar status for the problem of interacting disor-dered systems. This corroborates the conclusion of Sec.IV C, since Re Σ is typically involved in the calculationof any observable, in particular of the CMP.
2. Modification of the scattering time
We now show that the nonlinear energy, Eq. (14),which governs the long-time behavior of the quantumboomerang effect, also controls the change in the meanscattering time.As described in Sec. III, during the first part of thetime evolution, precisely in the range t < − τ , seeFig. 1, the CMP is essentially not modified by interac-tions. The only difference between the interacting andnon-interacting cases is that the mean scattering timeand mean free path are increased. We have used the the- oretical prediction for (cid:104) x ( t ) (cid:105) [25] in the non-interactinglimit, which depends on τ and (cid:96) , to fit the data in theinteracting case in the short-time region, and thus accessthe dependence of τ and (cid:96) on the interaction strength g . Fits have been performed including weights inverselyproportional to the square of the statistical errors in thetime window t ∈ [0 , t fit ]. Our choice is t fit = 20 τ . Thisvalue is chosen such that the fits enclose the whole ballis-tic motion and the beginning of the reflection. The valueof t fit slightly influences the fitted parameters, but thechanges are smaller than the error originating from thefitting procedure.The fits return both τ and (cid:96) for a given interactionstrength g . This allows us to calculate the average veloc-ity v = (cid:96)/τ . From our data, we observe that the averagevelocity remains almost unaffected across all studied in-teraction values, although it is a little higher than thepredicted value (cid:126) k /m at the Born approximation. Thisapparent discrepancy is caused by higher order correc-tions to the Born approximation and is of the order of1 /k (cid:96) (cid:28)
1. We can thus restrict the analysis to themean scattering time τ only.The fitted values of τ are shown in Fig. 11 as a functionof the nonlinear energy averaged over the fit time window, E shortNL . To explain these curves, we expand the meanscattering time τ ( E (cid:39) E + E NL ) to leading order in E NL (cid:28) E , using the Born approximation, Eq. (3). Thisyields: τ (cid:39) τ + (cid:126) k γ E NL . (A7)The linear increase of τ is well visible in Fig. 11. Thefact that curves at different σ are slightly shifted upwardsis due to the (small) dependence of the g − independentpart, τ , on σ , which shows up beyond the Born approxi-mation. Eq. (A7) is shown in Fig. 11, where for τ we usethe numerical value of the scattering time for σ = 40 /k and g = 0. The agreement between Eq. (A7) and thedata is very good. At larger values of E shortNL , the curvesstart to deviate from the linear behavior, bending up-wards. This effect is smaller for decreasing wave-packetwidths σ . The observed change corresponds to a rela-tively strong boomerang breakdown, so that the use ofthe non-interacting theoretical prediction for the fit be-comes less reliable. [1] P. W. Anderson, “Absence of diffusion in certain randomlattices,” Phys. Rev. , 1492–1505 (1958).[2] Ferdinand Evers and Alexander D. Mirlin, “Andersontransitions,” Rev. Mod. Phys. , 1355–1417 (2008).[3] AA Chabanov, M Stoytchev, and AZ Genack, “Statis-tical signatures of photon localization,” Nature , 850(2000).[4] Tal Schwartz, Guy Bartal, Shmuel Fishman, andMordechai Segev, “Transport and Anderson localization in disordered two-dimensional photonic lattices,” Nature , 52 (2007).[5] Hefei Hu, A Strybulevych, JH Page, Sergey E Skipetrov,and Bart A van Tiggelen, “Localization of ultrasound ina three-dimensional elastic network,” Nature Physics ,945 (2008).[6] Juliette Billy, Vincent Josse, Zhanchun Zuo, AlainBernard, Ben Hambrecht, Pierre Lugan, David Cl´ement,Laurent Sanchez-Palencia, Philippe Bouyer, and Alain Aspect, “Direct observation of Anderson localization ofmatter waves in a controlled disorder,” Nature , 891(2008).[7] Fred Jendrzejewski, Alain Bernard, Killian Mueller,Patrick Cheinet, Vincent Josse, Marie Piraud, LucaPezz´e, Laurent Sanchez-Palencia, Alain Aspect, andPhilippe Bouyer, “Three-dimensional localization of ul-tracold atoms in an optical disordered potential,” NaturePhysics , 398 (2012).[8] Isam Manai, Jean-Fran¸cois Cl´ement, Radu Chicireanu,Cl´ement Hainaut, Jean Claude Garreau, Pascal Szrift-giser, and Dominique Delande, “Experimental obser-vation of two-dimensional Anderson localization withthe atomic kicked rotor,” Phys. Rev. Lett. , 240603(2015).[9] Giulia Semeghini, Manuele Landini, Patricia Castilho,Sanjukta Roy, Giacomo Spagnolli, Andreas Trenkwalder,Marco Fattori, Massimo Inguscio, and Giovanni Mod-ugno, “Measurement of the mobility edge for 3d Ander-son localization,” Nature Physics , 554 (2015).[10] Donald H White, Thomas A Haase, Dylan J Brown,Maarten D Hoogerland, Mojdeh S Najafabadi, John LHelm, Christopher Gies, Daniel Schumayer, andDavid AW Hutchinson, “Observation of two-dimensionalAnderson localisation of ultracold atoms,” arXiv preprintarXiv:1911.04858 (2019).[11] Serge Aubry and Gilles Andr´e, “Analyticity breaking andAnderson localization in incommensurate lattices,” Ann.Israel Phys. Soc , 18 (1980).[12] Giacomo Roati, Chiara D’Errico, Leonardo Fallani,Marco Fattori, Chiara Fort, Matteo Zaccanti, GiovanniModugno, and Massimo Inguscio, “Anderson localiza-tion of a non-interacting Bose–Einstein condensate,” Na-ture , 895 (2008).[13] I. V. Gornyi, A. D. Mirlin, and D. G. Polyakov, “Inter-acting electrons in disordered wires: Anderson localiza-tion and low- t transport,” Phys. Rev. Lett. , 206603(2005).[14] Denis M Basko, Igor L Aleiner, and Boris L Alt-shuler, “Metal–insulator transition in a weakly interact-ing many-electron system with localized single-particlestates,” Ann. Phys. (N.Y.) , 1126–1205 (2006).[15] David A. Huse, Rahul Nandkishore, and VadimOganesyan, “Phenomenology of fully many-body-localized systems,” Phys. Rev. B , 174202 (2014).[16] Fabien Alet and Nicolas Laflorencie, “Many-body local-ization: An introduction and selected topics,” ComptesRendus Physique , 498 – 525 (2018), Quantum simu-lation / Simulation quantique.[17] Michael Schreiber, Sean S. Hodgman, Pranjal Bordia,Henrik P. L¨uschen, Mark H. Fischer, Ronen Vosk, EhudAltman, Ulrich Schneider, and Immanuel Bloch, “Obser-vation of many-body localization of interacting fermionsin a quasirandom optical lattice,” Science , 842–845(2015).[18] Jae-Yoon Choi, Sebastian Hild, Johannes Zeiher, Pe-ter Schauß, Antonio Rubio-Abadal, Tarik Yefsah, VedikaKhemani, David A Huse, Immanuel Bloch, and Chris-tian Gross, “Exploring the many-body localization tran-sition in two dimensions,” Science , 1547–1552(2016).[19] J. ˇSuntajs, J. Bonˇca, T. Prosen, and L. Vidmar, “Quan-tum chaos challenges many-body localization,” arXivpreprint , arXiv:1905.06345 (2019). [20] D. A. Abanin, J. H. Bardarson, G. De Tomasi,S. Gopalakrishnan, V. Khemani, S. A. Parameswaran,F. Pollmann, A. C. Potter, M. Serbyn, and R. Vasseur,“Distinguishing localization from chaos: challenges infinite-size systems,” arXiv preprint , arXiv:1911.04501(2019).[21] Piotr Sierant, Dominique Delande, and Jakub Za-krzewski, “Thouless time analysis of Anderson and many-body localization transitions,” Phys. Rev. Lett. ,186601 (2020).[22] R. K. Panda, A. Scardicchio, M. Schulz, S. R. Taylor,and M. ˇZnidariˇc, “Can we study the many-body localisa-tion transition?” EPL (Europhysics Letters) , 67003(2020).[23] Henrik P. L¨uschen, Sebastian Scherg, ThomasKohlert, Michael Schreiber, Pranjal Bordia, XiaoLi, S. Das Sarma, and Immanuel Bloch, “Single-particlemobility edge in a one-dimensional quasiperiodic opticallattice,” Phys. Rev. Lett. , 160404 (2018).[24] Jan Major, Giovanna Morigi, and Jakub Zakrzewski,“Single-particle localization in dynamical potentials,”Phys. Rev. A , 053633 (2018).[25] Tony Prat, Dominique Delande, and Nicolas Cherroret,“Quantum boomeranglike effect of wave packets in ran-dom media,” Phys. Rev. A , 023629 (2019).[26] A. S. Pikovsky and D. L. Shepelyansky, “Destructionof Anderson localization by a weak nonlinearity,” Phys.Rev. Lett. , 094101 (2008).[27] Ch. Skokos, D. O. Krimer, S. Komineas, and S. Flach,“Delocalization of wave packets in disordered nonlinearchains,” Phys. Rev. E , 056211 (2009).[28] S. Flach, D. O. Krimer, and Ch. Skokos, “Universalspreading of wave packets in disordered nonlinear sys-tems,” Phys. Rev. Lett. , 024101 (2009).[29] Nicolas Cherroret, Benoˆıt Vermersch, Jean Claude Gar-reau, and Dominique Delande, “How nonlinear inter-actions challenge the three-dimensional Anderson transi-tion,” Phys. Rev. Lett. , 170603 (2014).[30] Ihor Vakulchyk, Mikhail V. Fistul, and Sergej Flach,“Wave packet spreading with disordered nonlineardiscrete-time quantum walks,” Phys. Rev. Lett. ,040501 (2019).[31] Krzysztof Sacha, Cord A. M¨uller, Dominique Delande,and Jakub Zakrzewski, “Anderson localization of soli-tons,” Phys. Rev. Lett. , 210402 (2009).[32] Dominique Delande, Krzysztof Sacha, Marcin P(cid:32)lodzie´n,Sanat K Avazbaev, and Jakub Zakrzewski, “Many-bodyAnderson localization in one-dimensional systems,” NewJournal of Physics , 045021 (2013).[33] Lev Pitaevskii and Sandro Stringari, Bose-Einstein con-densation and superfluidity , Vol. 164 (Oxford UniversityPress, 2016).[34] Eric Akkermans and Gilles Montambaux,
Mesoscopicphysics of electrons and photons (Cambridge universitypress, 2007).[35] Lifshits Il’ya Mikhailovich, Sergeyi Andreevich Gre-deskul, and Leonid Andreevich Pastur,
Introductionto the theory of disordered systems (Wiley-Interscience,1988).[36] VL Berezinskii, “Kinetics of a quantum particle in a one-dimensional random potential,” Soviet Journal of Exper-imental and Theoretical Physics , 620 (1974).[37] Note that there is no constraint on the size of the wavepacket σ compared to the mean free path (cid:96) . The expres- sion (6) does not depend on the ratio of the these twolength scales.[38] Nicolas Cherroret and Thomas Wellens, “Fokker-Planckequation for transport of wave packets in nonlinear dis-ordered media,” Physical Review E , 021114 (2011).[39] Tobias Geiger, Andreas Buchleitner, and ThomasWellens, “Microscopic scattering theory for interactingbosons in weak random potentials,” New Journal ofPhysics , 115015 (2013).[40] Thomas Wellens, private communication (2017-2019).[41] H. Tal-Ezer and R. Kosloff, J. Chem. Phys. , 3967(1984).[42] C. Leforestier, R.H. Bisseling, C. Cerjan, M.D. Feit,R. Friesner, A. Guldberg, A. Hammerich, G. Jolicard,W. Karrlein, H.-D. Meyer, N. Lipkin, O. Roncero,and R. Kosloff, “A comparison of different propagationschemes for the time dependent schr¨odinger equation,”J. Comput. Phys. , 59 (1991).[43] S. Roche and D. Mayou, “Conductivity of quasiperiodicsystems: A numerical study,” Phys. Rev. Lett. , 2518–2521 (1997).[44] Holger Fehske, Jens Schleede, Gerald Schubert, GerhardWellein, Vladimir S. Filinov, and Alan R. Bishop, “Nu-merical approaches to time evolution of complex quan-tum systems,” Physics Letters A , 2182 – 2188 (2009).[45] A preliminary study of this partial destruction of theboomerang effect can be found on the arXiv preprint [51].However, these preliminary results are not part of the published paper [25].[46] Marco Tavora, Achim Rosch, and Aditi Mitra, “Quenchdynamics of one-dimensional interacting bosons in adisordered potential: Elastic dephasing and criticalspeeding-up of thermalization,” Phys. Rev. Lett. ,010601 (2014).[47] Nicolas Cherroret, Tomasz Karpiuk, Benoˆıt Gr´emaud,and Christian Miniatura, “Thermalization of matterwaves in speckle potentials,” Phys. Rev. A , 063614(2015).[48] E. Abrahams, P. W. Anderson, D. C. Licciardello, andT. V. Ramakrishnan, “Scaling theory of localization: Ab-sence of quantum diffusion in two dimensions,” Phys.Rev. Lett. , 673–676 (1979).[49] Joseph W Goodman, Speckle phenomena in optics: the-ory and applications (Roberts and Company Publishers,Greenwood Village (Co), 2007).[50] The averaging must in principle be performed over manydisorder realizations, while the interaction energy maybe different for each disorder realization. In practice, aspatial averaging over the wave packet size is equivalent,provided the wave packet contains many speckle grains,which is the case if σ (cid:29) /k0