Prethermalization to thermalization crossover in a dilute Bose gas following an interaction ramp
Mathias Van Regemortel, Hadrien Kurkjian, Iacopo Carusotto, Michiel Wouters
PPrethermalization to thermalization crossover in a dilute Bose gas following aninteraction ramp
Mathias Van Regemortel, ∗ Hadrien Kurkjian, and Michiel Wouters
TQC, Universiteit Antwerpen, Universiteitsplein 1, B-2610 Antwerpen, Belgium
Iacopo Carusotto
INO-CNR BEC Center and Dipartimento di Fisica,Universit`a di Trento, via Sommarive 14, 38123 Povo, Italy (Dated: January 16, 2019)The dynamics of a weakly interacting Bose gas at low temperatures is close to integrable due to theapproximate quadratic nature of the many-body Hamiltonian. While the short-time physics afteran abrupt ramp of the interaction constant is dominated by the integrable dynamics, integrabilityis broken at longer times by higher-order interaction terms in the Bogoliubov Hamiltonian, inparticular by Beliaev-Landau scatterings involving three quasiparticles. The two-stage relaxationprocess is highlighted in the evolution of local observables such as the density-density correlationfunction: an integrable dephasing mechanism leads the system to a prethermal stage, followedby true thermalization conveyed by quasiparticle collisions. Our results bring the crossover fromprethermalization to thermalization within reach of current experiments with ultracold atomic gases.
I. INTRODUCTION
Ever since the development of quantum mechanics inthe early days, it has been a central question how the uni-tary time evolution of a quantum wavefunction of manyparticles may generate a seemingly thermal ensemble inthe long-time limit – at least in the eyes of an exper-imenter with limited tools to probe the system. The eigenstate thermalization hypothesis (ETH) [1, 2] aimsto address this question by stating that expectation val-ues of macroscopic observables computed with respect toa single generic eigenstate of energy E are the same asthe microcanonical average around the corresponding en-ergy. The hypothesis has been verified numerically for awide series of chaotic quantum systems [3, 4].Since it relies on the hypothesis of ergodicity, ETHis not expected to hold for integrable quantum systems.There, an extensive number of conserved quantities re-stricts the full quantum dynamics to a small subspace ofthe total phase space, thereby preventing thermalization.The long-time states of integrable systems can still bestatistically described by a stationary generalized Gibbsensemble (GGE) [5], that incorporates all the conservedcharges, as recently seen in a cold atom experiment [6].Another seminal experimental example is the quantumNewton cradle [7]. In the same spirit as ETH, a repre-sentative eigenstate of the integrable Hamiltonian can beidentified based on these conserved charges, which cor-rectly reproduces expectation values of local observables[8].Similarly, approximate integrable systems can gothrough a dephasing stage, after which they are left ina prethermal state [9], also described by a GGE withall the approximately conserved quantities. Nevertheless, ∗ [email protected] prethermalized thermalized
Figure 1. (a) A pictorial image of the separation of time scalesconsidered in this work. The gas is brought out of equilib-rium by abruptly ramping up the interaction constant g ina time τ s . Then, the approximate integrability of Hamilto-nian (1) leads the system through a prethermalization stage(blue shades) on a time scale τ preth set by the chemical poten-tial µ . Finally, a thermal equilibrium is reached on a vastlylonger time scale τ th through Beliaev-Landau collisions (redshades). Alternatively, the red shades can be seen as repre-senting the growth of thermodynamic entropy. (b) Diagramsof the predominant non-integrable collisions that drive thesystem toward full thermalization; Beliaev decay (up) andLandau scattering (down). at longer times true thermalization sets in, conveyed byhigher-order relaxation processes, such as illustrated inFig. 1(a).As of now, the literature on the crossover from aprethermalized to a thermalized state after a globalquench has been mostly restricted to toy models. It hasbeen studied how a 1D chain [10] or liquid [11] of spinlessfermions with weak integrability breaking first relaxes toa prethermal state, after which a kinetic picture allowsus to understand the full thermalization dynamics of themodel. In this article, we aim to bring the study of thetwo distinct relaxation mechanisms within the realm ofcurrent experiments. a r X i v : . [ c ond - m a t . qu a n t - g a s ] J a n II. THE HAMILTONIAN
The weakly interacting Bose gas provides a naturalexample of an experimentally relevant system that isclose to being integrable. It is described by the standardHamiltonian in 3D (we use units of (cid:126) = 1)ˆ H = (cid:88) k k m ˆ a † k ˆ a k + g V (cid:88) p , k , q ˆ a † p + q ˆ a † k − q ˆ a k ˆ a p . (1)Here, V is the volume of the gas, m is the particle massand g is the effective interaction constant, found from the s -wave scattering length a s as g = 4 πa s /m .Starting from the ground state of an ideal gas withdensity n , we perform an abrupt ramp of the interactionconstant g i → g f (with g i = 0 < g f ) within a nonzerotime window τ s (see Fig. 1(a)), and study the subsequentdynamics under Hamiltonian (1). Experimentally, thiscan be done with a Feshbach resonance [12] by suddenlyramping up an external magnetic field. Recently, thismechanism was utilized to probe the analog of cosmicSakharov oscillations in a 2D bosonic gas [13]. In lowdimensions, the interaction constant can also be modifiedby varying the transverse confinement [14].When interactions are weak (small na s ) and the tem-perature well below the critical temperature, almost allparticles are found in the k = 0 mode, justifying thereplacement ˆ a → (cid:104) ˆ a (cid:105) ≡ √ n V , where n ≈ n is thecondensate density . The dynamics of the bosonic gasafter an interaction quench was studied on the level of aquadratic approximation in the fluctuation operators ˆ a k ( k (cid:54) = 0) [17, 18], and later the departure from the prether-malized state was considered [19] and the damping of theoscillations was added by hand [20].We, however, seek to explicitly retain terms contain-ing three fluctuation operators as well, so as to describehigher-order (non-integrable) scatterings that eventuallylead the system toward thermalization. In the literatureon superfluidity, these are commonly studied in the con-text of Beliaev decay and Landau damping [21], wherethey are responsible for the damping of a phonon [22–25].By truncating (1) to third order in fluctuation opera-tors, we find the approximate Hamiltonian,ˆ H ≈ E + ˆ H + ˆ H . (2)The quadratic part can be diagonalized with the standardBogoliubov transformation ˆ a k = u k ˆ χ k + v k ˆ χ †− k , with u k , v k = ± (cid:115) k / m + gn ω k ± , (3) We perform a simplified approach, where the number of particlesis not conserved, but number-conserving approaches [15] wouldresult in exactly the same Hamiltonians ˆ H and ˆ H [16] k ( ξ − ) -1 n k -6 -4 -2 τ − s ( µ ) T ( µ ) Figure 2. Effect of an interaction quench g i = 0 → g f .(a) The momentum distribution n ( χ ) k of quasiparticles afterthe interaction ramp for decreasing switching times τ s = { , . , . , . } × µ − from bottom to top; the thickblack line is the limiting case τ s →
0, for which we have re-sult (13). n ( χ ) k constitutes a conserved charge of ˆ H and onlyevolves under ˆ H . (b) Upon a decrease in the switching time τ s more energy is injected into the system, which then trans-lates into a higher equilibrium temperature T . For τ s → T = O ( τ − / s ). and the quasiparticle frequency ω k = (cid:114) k m (cid:16) k m + 2 gn (cid:17) . (4)In terms of the Bogoliubov operators, Hamiltonian (2) isthen expressed as [16]ˆ H = (cid:88) k ω k ˆ χ † k ˆ χ k , (5)ˆ H = g (cid:114) n V (cid:88) k , q (cid:16) A k , q ˆ χ † k ˆ χ † q ˆ χ k + q + B k , q ˆ χ k ˆ χ q ˆ χ − k − q + h.c (cid:17) , (6)with the matrix elements of ˆ H A k , q = u k u q u k + q + v k v q v k + q + (cid:0) u k + q + v k + q (cid:1)(cid:0) u k v q + u q v k (cid:1) , (7) B k , q = 13 ( u k u − k − q v q + ( u q + v q )( u k v − k − q + u − k − q v k )+ v k v − k − q u q ) . (8)Upon taking the thermodynamic limit and rescaling thewave numbers with the healing length after the quench ξ = (cid:112) /mµ , k → ˜ k = kξ one notes that the densityof states times the matrix elements of ˆ H squared scalesas 1 / ( nξ ) = (cid:112) (4 π ) na s , exactly like the condensatedepletion n − n . Therefore, if the number of depletedparticles is sufficiently small, the dynamics under the in-tegrable Hamiltonian ˆ H occurs on a substantially fastertime scale than the ergodic dynamics of ˆ H . A similarreasoning to compare third and fourth order terms of theHamiltonian justifies the omission of ˆ H in (2). III. THE EQUATIONS OF MOTION
We start by looking at the short-time dynamics, gener-ated by ˆ H , such as studied in [17]. In particular, the in-teraction ramp takes place within a nonzero time window,short enough so that we can safely neglect any effects ofˆ H during the quench. We return to the basis of particleoperators and find that the dynamics of the quadraticcorrelation functions n ( a ) k = (cid:104) ˆ a † k ˆ a k (cid:105) and c ( a ) k = (cid:104) ˆ a k ˆ a − k (cid:105) is governed by [17] ∂ t n ( a ) k = − (cid:2) g ( t ) n c ( a ) k (cid:3) , (9) i∂ t c ( a ) k = 2 (cid:16) k m + g ( t ) n (cid:17) c ( a ) k + g ( t ) n (2 n ( a ) k + 1) . (10)This system of equations is readily integrated numeri-cally for a given temporal profile g ( t ) and with appro-priate initial conditions. It has been intensely studied inthe context of the dynamical Casimir effect [26], wherea modulation of the interaction constant or condensatedensity causes a change of vacuum for the quasiparti-cle operators ˆ χ k [14, 18, 27, 28]. The correlations ofthe quasiparticles, in turn, are evaluated with the lineartransform n ( χ ) k = (cid:0) u k + v k (cid:1) n ( a ) k − u k v k Re (cid:2) c ( a ) k (cid:3) + v k , (11) c ( χ ) k = u k c ( a ) k + v k c ( a ) ∗ k − u k v k n ( a ) k − u k v k . (12)In the limit of instantaneous switching time τ s → n ( a ) k and c ( a ) k are zero just after the quench and we find thecorrelation functions after the quench as n ( χ ) k = v k , c ( χ ) k = − u k v k . (13)In Fig. 2(a), we show the quasiparticle momentum dis-tribution for different τ s and see that it converges to (13)for shorter τ s .We now stick to the basis of Bogoliubov operators ˆ χ k .Their quadratic correlation functions evolve trivially un-der ˆ H as n ( χ ) k ( t ) = n ( χ ) k and c ( χ ) k ( t ) = ˜ c ( χ ) k e − iω k t , mak-ing n ( χ ) k and ˜ c ( χ ) k conserved quantities of ˆ H related tothe integrable dynamics. However, they do experience avariation under the full Hamiltonian (2), which breaksthe integrability. Via Heisenberg’s equation of motion,we derive their dynamics under ˆ H : ∂ t n ( χ ) k = 2 g (cid:114) n V Im (cid:26) (cid:88) q B k , − q R ∗ k , q +2 A k , q − k M q , k + A q , k − q M ∗ k , q (cid:27) (14) i∂ t ˜ c ( χ ) k = 2 g (cid:114) n V (cid:88) q (cid:26) B − k , q M k , q +2 A k , − q M ∗ q , k + A q , k − q R k , q (cid:27) e iω k t , (15) where we have introduced the correlation functions ofthree quasiparticles M k , q = (cid:68) ˆ χ † k − q ˆ χ † q ˆ χ k (cid:69) , R k , q = (cid:68) ˆ χ q − k ˆ χ − q ˆ χ k (cid:69) . (16)We next evaluate the equation of motion for these third-order correlators i∂ t M k , q = ( ω k − ω q − ω k − q ) M k , q + g (cid:114) n V F ( M ) k , q , (17) i∂ t R k , q = ( ω k + ω q + ω q − k ) R k , q + g (cid:114) n V F ( R ) k , q . (18)Here, the matrices F ( M,R ) k , q contain correlators of four op-erators. More generally, a connected correlator of p op-erators couples to correlators of p + 1 operators on theright-hand side, making this an ever-growing hierarchy[29]. However, as explained in [30], fourth-order corre-lators in F ( M,R ) k , q can be approximately factorized intoproducts of second-order correlators using Wick’s theo-rem , thus establishing a truncated hierarchy of correla-tions functions. After this factorization, we find the driveterm in Eqs. (17)–(18) as (we drop the superscript · ( χ ) for ease of notation) F ( M ) k , q = 2 A k , − q (cid:16) c ∗ q ( n k − q − n k ) − c ∗ k − q c k (cid:17) +2 A k , q − k (cid:16) c ∗ k − q ( n q − n k ) − c k c ∗ q (cid:17) +2 A q , k − q (cid:16) n k − q ( n q − n k ) − n k ( n q + 1) (cid:17) +3 B k , − q (cid:16) c ∗ k − q c ∗ q − c k n k − q (cid:17) +3 B k , q − k (cid:16) c ∗ k − q c ∗ q − c k ( n q + 1) (cid:17) − B q , k − q c k (cid:16) n q + n k − q + 1 (cid:17) , (19)and F ( R ) k , q = 2 A k , − q (cid:16) c k − q ( n k + n q + 1) + c k c q (cid:17) +2 A k , q − k (cid:16) c q ( n k + n k − q + 1) + c k − q c k (cid:17) +2 A q , k − q (cid:16) c k ( n q + n k − q + 1) + c q c k − q (cid:17) +3 B k , − q (cid:16) ( n k − q + 1)( n k + n q + 1) (cid:17) +3 B k , q − k (cid:16) n q ( n k + n k − q + 1) + n k − q + 1 (cid:17) +3 B q , k − q (cid:16) n k ( n q + n k − q + 1) (cid:17) . (20)As such, we establish a closed set of differential equationsfor correlators up to order three, which approximatelydescribes the dynamics of the bosonic gas after the in-teraction ramp, provided n (cid:29) n − n , ensuring thatconnected correlators of higher order have a decreasingmagnitude. Since the average value of linear operators (cid:104) ˆ χ k (cid:105) is zero to zerothorder in ˆ H , the products of a first and third order correlatorsare negligible compared to the terms in (19) and (20). IV. THE KINETIC EQUATIONS
In the long-time limit, the coupled system of equa-tions (14)–(15) and (17)–(18) reproduces the well-knownkinetic equations. This can be seen by formally solving(17) as M k , q ( t ) = − ig (cid:114) n V (cid:90) t ds F ( M ) k , q ( s ) e i ( ω k − ω q − ω k − q )( s − t ) , (21)and similar for R k , q ( t ) in (18). These expressions cannow be plugged into (14)–(15), after which we obtain ef-fective dynamics by ( i ) sending the integration boundary t → ∞ in (21), thus singling out non-oscillating terms inthe integral over s , and ( ii ) time averaging Eq. (14) toremove the contributions that oscillate rapidly with time t . The result is that the evolution of quasiparticle occu-pation numbers is governed by the kinetic (or quantumBoltzmann) equations ∂ t n k = 4 π g n V (cid:26) (cid:88) q A q , k − q δ (cid:0) ω k − ω q − ω k − q (cid:1) × (cid:16) n k − q n q − n k ( n q + n k − q + 1) (cid:17) +2 (cid:88) q A k , q − k δ (cid:0) ω q − ω k − ω q − k (cid:1) × (cid:16) n q ( n k + n q − k + 1) − n k n q − k (cid:17)(cid:27) . (22)Within the kinetic approximation, the oscillation fre-quencies from the evolution of M k , q have been trans- lated into δ -functions imposing energy conservation forthe redistribution of quasiparticle occupation numbers.In our method, the kinetic equations come as a limitingbehavior, so that deviations from them can be studiedquantitatively, as we do in Fig. 3; this to our knowledgehas not been done previously in 3D.With (22), we rederive the kinetic equation that isknown from the literature on Beliaev-Landau scattering,where it is commonly established with Fermi’s goldenrule [16]. The first term represents the redistribution ofquasiparticles through Beliaev decay, where a quasipar-ticle with high momentum k decays into (or is formedfrom) two with q and k − q . The second term, in turn,describes the Landau process of absorption (or emission)of the quasiparticle with momentum k by a quasiparti-cle q − k (or q ). Notice that the Landau term comeswith an additional factor 2 from the two possible scat-tering channels [23]. See Fig. 1(b) for the correspondingdiagrams.Through the same analysis, we obtain the evolution ofpair correlations, i∂ t c k = (2 ω k + 2 δω k − iγ k ) c k + I k ( { c q } ) (23)The first term ur this equation contains the evolution of c k under the Bogoliubov frequency z k = ω k + δω k − iγ k / H in second order perturbation the-ory in the instantaneous Fock state | n k , { n q } q (cid:54) = k (cid:105) (seethe Appendix for the explicit derivation). This containsthe Landau-Beliaev damping rate γ k = 4 π g n V (cid:88) q (cid:2) A q , k − q ( n q + n k − q + 1) δ (cid:0) ω k − ω q − ω k − q (cid:1) + 2 A k , q − k ( n q − k − n q ) δ (cid:0) ω q − ω k − ω q − k (cid:1) , (cid:3) (24)(where the Beliaev and Landau parts are respectively the first and second summation) and the frequency shift δω k = g n V P (cid:88) q (cid:34) A q , k − q ( n q + n k − q + 1) ω k − ω q − ω k − q + 4 A k , q − k ( n q − k − n q ) ω k + ω q − k − ω q − B k , q ( n q + n − k − q + 1) ω k + ω q + ω − k − q (cid:35) , (25)where P denotes the Cauchy principal value and we have used the fact that our momentum distribution remainssymmetric n ( χ ) k = n ( χ ) − k . The collisional integral I k ( { c q } ) accounts for the fact that the distribution { c q } in modes q (cid:54) = k changes dynamically with c k ; it is given by I k ( { c q } ) = 2 g n V (cid:34) (cid:88) q A q , k − q c q c k − q ω q + ω k − q − ω k + i + + 4 A k , q − k c q c ∗ q − k ω q − ω k − ω q − k + i + + 18 B k , q ˜ c ∗− k − q ˜ c ∗ q − ω k − ω q − ω − k − q + i + (cid:35) , (26)Eq. (23) describes how the coherence between modes k and − k evolves under three-body scatterings, and in particular I k describes how it is affected by the coherence in other modes q . Remark that c k may show a temporal evolution evenif the populations n k are prepared at thermal equilibrium, reflecting the underlying Landau-Beliaev scatterings thatmaintain equilibrium. To our knowledge, this equation was not found in the literature, unlike the kinetic equation(22) on n k [31].In Fig. 3, we perform a quantitative comparison be- tween the full integration of the truncated hierarchy t ( µ − ) n k k=0.313 ξ -1 t ( µ − ) n k × -4 k=3.98 ξ -1 t ( µ − ) n k × -3 k=2.5 ξ -1 t ( µ − ) n k × -6 k=5 ξ -1 Figure 3. A comparison of the quasiparticle occupation num-bers n k as produced by the integration of the hierarchyof correlation functions (blue lines) and the derived kineticequation in the adiabatic limit (red lines) for the quench g i = 0 → g f = 0 . µξ (corresponding to nξ = 20) in τ s = 0 . µ − for different momenta (see the initial momen-tum distribution of quasi-particles right after the quench onFig. 2). and the approximate kinetic description (22). We showthe evolution of the quasiparticle occupation numbers atshort times t ∼ /µ . We observe that the curve of n k ( t )predicted by the kinetic equation differs in two distinctways from that of the hierarchy: ( i ) the evolution at veryshort times is not well captured by the kinetic equation,which results in a small offset (controlled by the inter-action strength na s ) between the two curves, an offsetthen conserved all along the evolution and ( ii ) contraryto the kinetic description, the hierarchy of correlationsretains high-frequency components in n k ( t ). Those twodifferences are directly related to the approximations ( i )and ( ii ) detailed in the main text below (21), on whichkinetic equations are based.In the long-time limit, we find that (22) and (23) con-verge to the values in a thermal ensemble. The momen-tum distribution of quasiparticles approaches the Bose-Einstein distribution n th k = 1 e βω k − , (27)with β = 1 /k B T the inverse temperature set by the totalinjected energy, while the anomalous correlations vanish.The energy after the quench on the level of the quadraticHamiltonian, E = E + (cid:80) k ω k n ( χ ) k , is conserved un-der the kinetic equations. However, using the value ofthe mode occupation number n ( χ ) k for an infinitely fastquench [see Eq. (13)] leads to an ultraviolet divergenceof this injected energy. This divergence is regularizedby a finite switching time: this sets an effective cutoff in energy ω k max ∝ /τ s , corresponding to a momentumcutoff k max ∝ / √ τ s in the limit of fast quench τ s → E − E ∝ / √ τ s .This enables us to fix the total injected energy with theswitching time τ s and, consequently, the final equilibriumtemperature of the gas by matching this energy with theenergy of a thermal ensemble. When k B T > µ , we havethat E − E ∝ T / , such that we derive the asymp-totic scaling T ∝ τ − / s . In Fig. 2(b), we show the fullvariation of equilibrium temperature with switching time τ s . V. THE DENSITY-DENSITY CORRELATIONFUNCTION
Finally, we investigate the behavior of macroscopic ob-servables in real space, which are expected to exhibit thetwo distinct relaxation stages. We concentrate on dis-tances of the order of the (equilibrium) thermal wave-length 2 π/k th , with ω k th = k B T . We choose τ − s = 0 . µ ,such that k B T = 0 . µ and therefore the thermal wave-length is of the same order as the healing length ξ . For k ∼ /ξ , the kinetic equations are accurate for times t > /µ , and introduce an offset of the order of (cid:112) na ,our small parameter. Therefore, they correctly describethe dynamics of spatial correlations at length scales ∼ ξ in the weakly interacting limit, as we have also checkednumerically.The first relaxation stage of local observables to theirprethermal value is caused by a dephasing mechanismwhere all k -modes interfere desructively. We there-fore define the annihilation operator in position spaceˆ a ( r ) = 1 / √ V (cid:80) k e i k · r ˆ a k . Our analysis is now focused onthe evolution of the density-density correlation function,defined as g (2) ( r − r (cid:48) ; t ) = (cid:104) : ˆ n ( r ) ˆ n ( r (cid:48) ) : (cid:105) t / (cid:104) ˆ n ( r ) (cid:105) t for ahomogeneously distributed gas, where ‘:’ denotes normalordering and ˆ n ( r ) = ˆ a † ( r )ˆ a ( r ) is the local density opera-tor. The density-density correlation function has provenits importance previously in the context of analog gravity[32, 33], where the correlation pattern shows a fingerprintof the analog of Hawking radiation at an acoustic blackhole’s horizon [34].On the Gaussian level, the density correlation functioncan be simplified to g (2) ( r − r (cid:48) ; t ) = 1 + 2 n (cid:0) n ( r − r (cid:48) ; t ) + Re (cid:8) m ( r − r (cid:48) ; t ) (cid:9)(cid:1) , (28)where we defined n ( r − r (cid:48) ; t ) = (cid:104) ˆ a † ( r )ˆ a ( r (cid:48) ) (cid:105) t = 1 V (cid:88) k (cid:54) =0 e i k · ( r − r (cid:48) ) n ( a ) k ( t ) , (29)and analogous for m ( r − r (cid:48) ) = (cid:104) ˆ a ( r )ˆ a ( r (cid:48) ) (cid:105) . The quadraticcorrelations of fluctuations, n ( a ) k and c ( a ) k , can be obtainedfrom the quasiparticle correlations n ( χ ) k and c ( χ ) k throughthe inverse of the transformation (11)–(12). -1 g ( ) ( x ; t ) − × -3 -2-1.5-1-0.5 x=0.94 ξ -1 × -4 -50510 x=1.9 ξ t ( µ − ) -1 g ( ) ( x ; t ) − × -4 -20246 x=2.8 ξ t ( µ − ) -1 × -4 -10123 x=3.8 ξ Figure 4. The evolution of the density correlation functionafter the quench g i = 0 → g f = 0 . µξ in τ s = 0 . µ − for varying distances x = | r − r (cid:48) | according to the kineticpicture Eqs. (22) and (23). The horizontal solid lines indi-cate the asymptotic values for the prethermal and thermal(quasi)stationary ensemble; the temperature T = 0 . µ isfound from the initial state. At short distances, we clearlynotice a relaxation to a prethermal plateau on a time scale ofthe order of τ preth = µ − (for ramps τ s ∼ µ − and x ∼ ξ ),this is due to a dephasing mechanism in ˆ H . In g (2) ( x ; t ), thisis manifested as a fast oscillation at short times, which thendiminishes due to a destructive interference between all k -modes once the light-cone correlation peak has moved awayfrom the considered distance x [18]. Then, at much latertimes, τ therm ∼ µ − , a new equilibrium value is found thatcorresponds to the value in the thermal ensemble throughthe much slower dynamics of ˆ H . The difference between theprethermal and thermal value vanishes for increasing separa-tion x as the correlation function drops to zero in this limit. In Fig. 4 we show the evolution of the density corre-lation function after a ramp g i = 0 → g f = 0 . µξ (sothat na s = 1 . · − ) at different distances x = | r − r (cid:48) | .We observe a clear first relaxation, approximately to theprethermal value on a time scale τ preth ∼ µ − after aninitial oscillation due to the light-cone peak [18, 35] thatdies out due to dephasing once this has traveled away;this is governed by ˆ H . At much longer times, the scat-terings contained in ˆ H cause a new relaxation, this timeto the thermal value. We find that the thermalizationtime τ therm ∼ µ − is in qualitative agreement withthe Beliaev-Landau lifetime of the thermal wavenumber1 /γ BLk th ∼ /µ (cid:112) na s for k B T ∼ µ [16]. VI. CONCLUSIONS
We have illustrated that the crossover from a prether-malized to a thermalized state can be witnessed in a coldatomic gas by probing the density correlations after a sudden interaction ramp. The switching time of the rampdetermines the final temperature in the equilibrium en-semble. While a simple dephasing mechanism, treatedon the level of the (integrable) quadratic Hamiltonian,causes local observables to relax to a prethermal value,a more sophisticated approach is needed to describe thethermalization stage. Here, third-order interaction pro-cesses, known as Beliaev-Landau collisions, are the pre-dominant mechanism to lead the system away from inte-grability and, eventually, to thermal equilibrium. Whenfocusing on most relevant length scales of the order of theequilibrium thermal wavelength, a kinetic description issufficient to describe the final relaxation. In principle,our predictions are within reach of current experimentswith ultracold atomic gases.
ACKNOWLEDGMENTS
MVR gratefully acknowledges support in the form of aPh. D. fellowship of the Research Foundation - Flanders(FWO) and hospitality at the BEC Center in Trento.HK is supported by the FWO and the European UnionH2020 program under the MSC Grant Agreement No.665501. MW acknowledge financial support from theFWO-Odysseus program. IC was funded by the EU-FET Proactive grant AQuS, Project No. 640800, andby Provincia Autonoma di Trento, partially through theproject “On silicon chip quantum optics for quantumcomputing and secure communications (SiQuro)”.
APPENDIX: PERTURBED BOGOLIUBOVENERGY IN AN ARBITRARY EXCITED STATE
To recover Eqs. (24–25), we treat ˆ H as a perturbationof ˆ H and we recall [36, 37] that the complex poles E ofthe resolvent (or equivalently of the Green function) in agiven state | ψ (cid:105) are given to second order in the pertur-bation by E = E + (cid:104) ψ | ˆ H | ψ (cid:105) + (cid:104) ψ | ˆ H ˆ QE − ˆ H ˆ H | ψ (cid:105) (30)= E + (cid:104) ψ | ˆ H | ψ (cid:105) + (cid:88) λ |(cid:104) ψ | ˆ H | λ (cid:105)| E − (cid:104) λ | ˆ H | λ (cid:105) (31)where E = (cid:104) ψ | ˆ H | ψ (cid:105) is the unperturbed energy, ˆ Q =1 − | ψ (cid:105)(cid:104) ψ | projects orthogonally to | ψ (cid:105) and the states | λ (cid:105) are therefore orthogonal to | ψ (cid:105) . We apply Eq. (31) to theFock states | n k , { n q } q (cid:54) = k (cid:105) and | n k − , { n q } q (cid:54) = k (cid:105) whoseperturbed energies, respectively E ( n k ) and E ( n k − z k ≡ E ( n k ) − E ( n k −
1) (32)Changing the sum over the intermediate Fock states λ into a sum over the scattered momentum q (taking careto avoid double countings) and replacing in the denom-inator E ( n k ) by its zeroth-order approximation E ( n k ) we get z k = z k , + g n V (cid:88) q (cid:34) A q , k − q [(1 + n k − q )(1 + n q ) − n k − q n q ] z k , − ω q − ω k − q + 4 A k , q − k [(1 + n q ) n q − k − n q (1 + n q − k )] z k , + ω q − k − ω q + 18 B k , q [ n − k − q n q − (1 + n − k − q )(1 + n q )] z k , + ω q + ω − k − q (cid:35) (33)where z k , = E ( n k ) − E ( n k −
1) = ω k + i + is theunperturbed Bogoliubov frequency. The counting fac-tors originate in our symmetric writing (6) of ˆ H wherethe same operator appears more than once. Using thePlemelj formula 1 / ( x + i + ) = P (1 /x ) − iπδ ( x ) to splitthe real and imaginary parts finally yields Eqs. (24–25). [1] J. M. Deutsch, Physical Review A , 2046 (1991).[2] M. Srednicki, Physical Review E , 888 (1994).[3] M. Rigol, V. Dunjko, and M. Olshanii, Nature , 854(2008).[4] M. Rigol, Physical Review A , 053607 (2009).[5] M. Rigol, V. Dunjko, V. Yurovsky, and M. Olshanii,Physical review letters , 050405 (2007).[6] T. Langen, S. Erne, R. Geiger, B. Rauer, T. Schweigler,M. Kuhnert, W. Rohringer, I. E. Mazets, T. Gasenzer,and J. Schmiedmayer, Science , 207 (2015).[7] T. Kinoshita, T. Wenger, and D. S. Weiss, Nature ,900 (2006).[8] J.-S. Caux and F. H. Essler, Physical review letters ,257203 (2013).[9] J. Berges, S. Bors´anyi, and C. Wetterich, Phys. Rev.Lett. , 142002 (2004).[10] B. Bertini, F. H. Essler, S. Groha, and N. J. Robinson,Physical review letters , 180601 (2015).[11] M. Buchhold, M. Heyl, and S. Diehl, Phys. Rev. A ,013601 (2016).[12] C. Chin, R. Grimm, P. Julienne, and E. Tiesinga, Re-views of Modern Physics , 1225 (2010).[13] C.-L. Hung, V. Gurarie, and C. Chin, Science , 1237557(2013).[14] J.-C. Jaskula, G. B. Partridge, M. Bonneau, R. Lopes,J. Ruaudel, D. Boiron, and C. I. Westbrook, Phys. Rev.Lett. , 220401 (2012).[15] Y. Castin and R. Dum, Physical Review A , 3008(1998).[16] A. Sinatra, Y. Castin, and E. Witkowska, Physical Re-view A , 033614 (2009).[17] S. S. Natu and E. J. Mueller, Physical Review A ,053607 (2013).[18] I. Carusotto, R. Balbinot, A. Fabbri, and A. Recati, TheEuropean Physical Journal D , 391 (2010).[19] G. Menegoz and A. Silva, Journal of Statistical Mechan-ics: Theory and Experiment , P05035 (2015). [20] A. Ran¸con and K. Levin, Physical Review A , 021602(2014).[21] S. Beliaev, Sov. Phys. JETP , 299 (1958).[22] L. Pitaevskii and S. Stringari, Physics Letters A , 398(1997).[23] S. Giorgini, Physical Review A , 2949 (1998).[24] F. Chevy, V. Bretin, P. Rosenbusch, K. W. Madison, andJ. Dalibard, Phys. Rev. Lett. , 250402 (2002).[25] N. Katz, J. Steinhauer, R. Ozeri, and N. Davidson, Phys-ical review letters , 220401 (2002).[26] V. Dodonov, Physica Scripta , 038105 (2010).[27] S. Koghee and M. Wouters, Physical review letters ,036406 (2014).[28] C. Wilson, G. Johansson, A. Pourkabirian, M. Simoen,J. Johansson, T. Duty, F. Nori, and P. Delsing, Nature , 376 (2011).[29] N. P. Proukakis, Journal of Physics B: Atomic, Molecularand Optical Physics , 4737 (2001).[30] M. Van Regemortel, W. Casteels, I. Carusotto, andM. Wouters, Phys. Rev. A , 053854 (2017).[31] E. Lipschitz and L. Pitaevskii, in Landau and LifshitzCourse of Theoretical Physics , Vol. 10 (Pergamon Press,New York, 1981) Chap. VII.[32] R. Balbinot, A. Fabbri, S. Fagnocchi, A. Recati, andI. Carusotto, Physical Review A , 021603 (2008).[33] I. Carusotto, S. Fagnocchi, A. Recati, R. Balbinot, andA. Fabbri, New Journal of Physics , 103001 (2008).[34] J. Steinhauer, Nature Physics , 959 (2016).[35] P. Calabrese and J. Cardy, Journal of Statistical Mechan-ics: Theory and Experiment , P04010 (2005).[36] C. Cohen-Tannoudji, J. Dupont-Roc, and G. Gryn-berg, Processus d’interaction entre photons et atomes (InterEditions et ´Editions du CNRS, Paris, 1988).[37] H. Kurkjian, Y. Castin, and A. Sinatra, Annalen derPhysik529