Thermalization in a quasi-1D ultracold bosonic gas

We study the collisional processes that can lead to thermalization in one-dimensional systems. For two body collisions excitations of transverse modes are the prerequisite for energy exchange and thermalzation. At very low temperatures excitations of transverse modes are exponentially suppressed, thermalization by two body collisions stops and the system should become integrable. In quantum mechanics virtual excitations of higher radial modes are possible. These virtually excited radial modes give rise to effective three-body velocity-changing collisions which lead to thermalization. We show that these three-body elastic interactions are suppressed by pairwise quantum correlations when approaching the strongly correlated regime. If the relative momentum $k$ is small compared to the two-body coupling constant $c$ the three-particle scattering state is suppressed by a factor of $(k/c)^{12}$, which is proportional to $\gamma ^{12}$, that is to the square of the three-body correlation function at zero distance in the limit of the Lieb-Liniger parameter $\gamma \gg 1$. This demonstrates that in one dimensional quantum systems it is not the freeze-out of two body collisions but the strong quantum correlations which ensures absence of thermalization on experimentally relevant time scales.


I. INTRODUCTION
One-dimensional (1D) systems [1,2] are a model to study the fundamental processes of dynamics and (de)coherence in interacting many-body quantum systems. Ultracold atoms in strongly elongated traps with ω r ≫ ω z (ω r , ω z being the frequencies of the radial and longitudinal confinement, respectively) offer the possibility to implement 1D quantum physics if both the temperature T and chemical potential µ are small compared to the energy scale given by the transverse confinement: 1D systems of ultra-cold atoms were implemented in both optical lattices [3] and atom chips [4]. In the limit of zero temperature they are a realization of the Lieb-Liniger model [5] of spinless bosons with contact (point-like) interaction, a prime example of an integrable system. An important parameter characterizing an 1D system of bosons with point-like interactions described by the (three dimensional) s-wave scattering length α s is the Lieb-Liniger parameter [5] where m is the mass ot the bosonic atom, n 1D the linear density of the atoms in the 1D trap characterized by the transverse confinement frequency ω r , and l r is the fundamental length scale of the localization of an atom in the transversal direction given by l r = h/(mω r ).
The limit γ ≪ 1 corresponds to a weakly-interacting regime, whereas γ ≫ 1 signifies strongly-interacting, strongly correlated (Tonks-Girardeau) regime. In an integrable system [6,7] the number of their integrals of motion equals exactly the number of their degrees of freedom. Thus such a system always "remembers" its initial state in the course of its dynamical evolution, and thermalization does not occur. In an integrable system the finite spread of initial energy may lead only to relaxation towards the generalized Gibbs (or fully constrained thermodynamic) ensemble [8]. Strictly speaking, there is no thermalization in any closed system, but for non-integrable systems the eigenstate thermalization hypothesis [9] holds, enabling dephasing to mimic the relaxation to the thermal equilibrium.
Recently strong inhibition of thermalization, was observed on an optical lattice experiment with bosons deep in the 1D regime [10]. On the other hand interference experiments on atom chips with pairs of weakly interacting Bose gases fulfilling the conditions of Eq. (1) are in a good agreement with the thermal-equilibrium description of the 1D atomic ensembles [11][12][13].
In the present paper we investigate the collisional properties of Bose gases in a 1D geometry and how they contribute to thermalization. We follow thereby a procedure outlined in our two previous works [14,15], and give a more indepth description of the underlying theoretical considerations. We start with the calcualtion of the freeze-out of thermalization providing two body collisions. We then proceed to show that virtual excitations of excited states lead to effective three-body collisions, which lead to a term in the Hamiltonian that breaks integrability and enables thermalization. We then proced to estimte the effects of quantum correlations in 1D and show how they suppress the three-body term for strongly correlated 1D systems, thereby extending the time scale, on which a quasi-1D system can be approximately described as integrable.
In our theoretical considerations we consider identical bosons in a tight 1D wave guide with radial confinement given by a 2D harmonic oscillator with frequency ω r (we set ω z = 0). The contact interaction is described by the pseudopotential 4πh 2 m −1 α s δ(r − r ′ ). The Hamiltonian of the 1D system is: Thereby the field operatorsψ(r) are assumed to vanish for x 2 + y 2 → ∞ and to be periodic along z with the period L. For the solutions to Eq. 5 we make the Ansatz : where L is the quantization length and the atomic annihilation and creation operatorsâ {n,ℓ} k andâ † {n,ℓ} k obey the standard bosonic commutation rules. ϕ n,ℓ (x, y) is the normalized eigenfunction of both the radial confinement Hamiltonian,Ĥ (r) ϕ n,ℓ (x, y) = (n + 1)hω r ϕ n,ℓ (x, y) and the z-projection of the orbital momentum, −i[x(∂/∂y) − y(∂/∂x)]ϕ n,ℓ (x, y) = ℓϕ n,ℓ (x, y).

II. TWO-BODY COLLISIONS
We first look at collisions of two identical bosonic atoms that are initially in the transverse ground state of the radial confinement. If the collision is restricted to 1D, that is both atoms remain after the collision in the transverse ground state, then there can be no energy exchange and cinsequently no thermalization. For such two-body collisions to contribute to energy exchange and thermalization, they have to lead to a change in transverse excitation. By symmetry ∆n 1 + ∆n 2 must be even. For atoms in the transverse ground state (n 1 = n 2 = 0) ∆n 1 = ∆n 2 band following the above considerations their orbital-momentum quantum numbers after collision are restricted to −ℓ and +ℓ. The rate of populating the radially excited modes by pairwise atomic collisions Γ 2b , can then be estimated for a non-degenerate Bose gas, using Fermi's golden rule. For k B T <hω r this rate is The dimensionless quantity ζ = n 1D α 2 s /l r combines two dimensionless parameters characterizing a 1D system. n 1D α s ∝ μ hωr is a measure how much the interaction energy (the chemical potential µ) is smaller than the energy scale given by the transverse confinement. α s /l r characterizes the relation between the transverse confinement and the strength of the contact interaction. Its importance can be seen when looking at how the effective 1D coupling constant g 1D of pairwise interacting bosonic atoms in a waveguide changes with confinement due to virtual excitation of the radial modes. Following Olshanii [16] g 1D = 2hω r α s /[1 − C ′ αs √ 2lr ], C ′ ≈ 1.46 and increases as α s /l r grows. This points to the ratio α s /l r ≪ 1 as the measure how much the 1D approximation is violated. In a general case we have to change in all the following expressions α s to α s /[1 − C ′ αs √ 2lr ]. Eq. (7) has a also transparent physical interpretation: The rate Γ 2b is related to the 3D atomic density (∼ n 1D /l 2 r ), times the s-wave scattering cross-section (∼ α 2 s ), times the exponential Boltzmann factor for the fraction of atoms fast enough to scatter into higher radial modes, times the corresponding velocity of the collision (∼h/(ml r )).
Looking at the scaling of Eq. (7) one immedeately sees that the collision rate leading to thermalization (Γ 2b ) rapidly diminishes when the temperature approaches T ∼hω r and is suppressed by more then a factor of 50 (e −4 ) for T = 1 2h ω r . Estimating the numbers for recent experients [12] in 87 Rb: α s = 5.3 nm, n 1D = 50 µm −1 , ω r /(2π) = 3 kHz, T = 30 nK (ζ ≈ 0.007) one obtains a collision rate of Γ 2b ∼ 0.02 s −1 . Ths is at least one order of magnitude too small for two body collisions to be responsible for the thermalization required in the evaporative cooling process leading to these low temperatures.

III. THREE-BODY COLLISIONS
If the kinetic energy of the collision is less than 2hω r , then the radial modes can be excited only virtually. Such processes contribute to the system dynamics in the second and higher orders of perturbation theory.
The simplest case is when after the collision the radial motion state is |{n ′ . Then only one more collision is enough to de-excite the radial mode and bring the system back on the energy shell [see Fig.  1(a)]. Such a process yields an effective three-body collision already in the second order of perturbation theory.

A. Perturbative approach
We will now calculate the matrix elements for the process shown in figure 1(a) which leads to effective three-body collisions. In our perturbation calculation the small parameter is n 1D α s , that is the mean field interaction µ is much smaller then the energy scalehω r connected to the transverse confinement. In addition, to avoid complications related to the confinement-induced resonance in 1D scattering [16] we assume α s ≪ l r . We can then rewrite the Hamiltonian (4) asĤ where all the terms irrelevant to the process under discussion (Fig. 1a) are gathered inR, and the radial matrix element f 2p,0;0,0 0,0;0,0 is given by: It connects to two atoms in the ground state of the incoming channel, to one atom remaining in the same state, and the other being excited to a state with zero orbital-momentum quantum number and even main quantum number n = 2p, p = 0, 1, 2, ... (remember for Boson n and ℓ are required to have the same parity). To evaluate Eq. (9), we recall that the normalized radial wave functions ϕ n,ℓ (x, y) of interest are real and can be expressed through Laguerre polynomials L p ϕ n=2p,ℓ=0 (x, y) = (πl 2 r ) −1 exp − Since L p (0) = 1, we obtain [17] dx dy dx ′ dy ′ ϕ n=2p,ℓ=0 x independently of p. Then we easily obtain the necessary matrix element as f 2p,0;0,0 0,0;0,0 = C 2p,0;0,0 2p,0;0,0 /2πl 2 r , where the coefficient C 2p,0;0,0 2p,0;0,0 is defined by the expansion Comparing the coefficients in front of (x 2 + y 2 ) p in the left-and right-hand sides of Eq. (12), we obtain C 2p,0;0,0 2p,0;0,0 = 2 −p and f 2p,0;0,0 0,0;0,0 = In our consideration we are only interested in the case where the collision energy of the two atoms is always much smaller thanhω r . Then, using the matrix element (13) and adiabatically eliminating the radially excited mode operators, we obtain from the original Hamiltonian Eq. (4) an effective 1D Hamiltonian: where we write for simplicityâ k ≡â {0,0} k and the numerical constant ξ is given by Note that the relative contribution of the virtual states with the excitation energy higher than 2hω r given by (ξ − 1)/ξ is remarkably small. The summation in the last term of Eq. (14) is taken over all the kinetic momenta obeying the conservation law Introducing the field operatorψ(z) = L −1/2 kâ k exp(ikz), we can rewrite Eq. (14) aŝ The first and second terms in Eqs. (14) or (16) correspond to the Lieb-Liniger model. The third (cubic) term stems from the effective three-body collisions mediated by virtually excited states 1 . This third (cubic) term in Eq. 16 violates the integrability in the 1D system. The fact that the effective three-body interactions are dominated by virtual excitations of the lowest even-parity excited state may seem surprising, since the correct calculation of the effective two-body 1D coupling constant requires taking into account of the infinite number of states [16]. However, in the latter case one deals with the calculation of the two-body wave function, which has in 1D a 1/z singularity (stemming from the 3D boundary condition at r → 0 that provides the correct asymptotic form of the scattered s-wave). The removal of this singlularity yields the regular part of the two-body wave function and, through this regular part, the scattering amplitude and, hence, the effective coupling in 1D. On the other hand, if one tries to obtain the effective 1D interaction by adiabatic elimination of all the excited states, one gets a divergent series ∞ n=0 n −1/2 in the expression for the effective 1D coupling constant. The aforementioned regularization of the wave functions formally corresponds to the renormalization of this divergent series via substituting it by a finite expression lim s→∞ [16]. In our case, processes related to three-body collisions do not give rise to additional singularities in the many-body wave function, and no additional regularization is needed. The sum over all excited states thus converges. The convergence is rapid enough to ensure fair estimation of the whole sum by its first term.
Accurate calculation of the effective three-body interaction potential U 3b yields where Obviously, the sum in Eq. (18) converges and gives a sharp-peaked function rapidly (exponentially) decreasing at distances much larger than l r . Since Eq. (1) holds, we consider scattering events with transferred momenta much less thanh/l r . In this case we can use approximation . Then, by taking the matrix element of the effective interaction U 3b and dividing it by 3! (the number of permutation of three identical particles) we obtain the last term in the second-quantized Hamiltonian (14). Before continuing we want to point out similarities with other recent works: (1) The mechanism discussed here is to a certain extent similar to the virtual association of atoms to a molecular dimer [18]. In our discussion here, virtual excitation of radial modes during a two-atom collision temporarily localize the interatomic distance on the length scale ∼ l r . Scattering a third atom on such a transient structure of finite size and mass 2m leads to an effective threebody collision. In Ref. [18], collisions of a third atom bring "virtual" dimers, enhanced in size by a Feshbach resonance down to the energy shell, thus bringing about "quantum chemistry" in 1D. (2) In a similar way effective three-body interactions between polar molecules emerge, due to virtual transitions to an off-resonant internal state [19].

B. Variational approach
The cubic term in Eq. (16) is negative and thus supports no bound ground state. Therefore we have to consider it as a first correction term to the purely pairwise interaction energy in the effective 1D Hamiltonian. Considering the mean-field limit of Eq. (16),ψ ≈ √ n 1D exp(iθ), we get the energy density (per unit length) this expansion is correct in the limit of the small linear density n 1D a s ≪ 1. If we neglect the contribution of the radial levels with the main quantum number larger than 2 by setting ξ ≈ 1, we see that Eq. (19) is the expansion up to the cubic term of the energy density obtained by the variational method by Salasnich, Parola and Reatto [20] E var =h Here one assumes the wave function of the transversal atomic motion to be ∝ exp[−(

IV. CALCULATIONS OF THE COLLISION RATES
We now turn to the collision rates for these effective three body collisions. We start with calculating the rate Γ k1k2k3 for the decay of a specific state |k 1 , k 2 , k 3 ≡ a † k1 a † k2 a † k3 |vac , |vac being the vacuum state of the atomic field wherein atoms are absent (should not be confused with the vacuum of elementary excitations) due to three-body collisions. The final states of the decay are |k ′ |vac . To make the calcualtion simple we assume the 1D bosonic gas being strongly non-degenerate (k B T much higher than the chemical potential) and weakly-interacting. This enables us to neglect the probability double (and higher) occupation of any k-mode, therefore assuming all involved atomic momenta to be different, and assume the elementary excitations coinciding with the atomic plane waves with the free-particle (quadratic) dispersion law. Then Fermi's golden rule yields To account for the condition k ′ we add an additional integration over k ′ 1 with the necessary delta-function: Since the atoms are indistinguishable, we need to integrate over the volume W ′ in the k ′ -space that corresponds to a unique ordering of the variables, e.g., k ′ 1 > k ′ 2 > k ′ 3 . The integral over W ′ of a function fully symmetric over permutations of k ′ 1 , k ′ 2 , k ′ 3 amounts to 1/3! of the integral over the whole k ′ -space. If all the involved momenta are different (as is the case for a system far from degeneracy) all possible ways of ordering three bosonic creation operators with lower indices k ′ 1 , k ′ 2 , k ′ 3 and three bosonic annihilation operators with lower indices k 1 , k 2 , k 3 finally yield To evaluate the integral in Eq. (23), we perform an orthogonal transformation from k ′ 1 , k ′ 2 , k ′ 3 to Jacobi co-ordinates (here we deal with the 1D analog of the hyperspherical co-ordinates, which are used in the three-body problem [21]) in the wavenumber space: Further we introduce the hyperangle χ ′ via and express Eq. (23) as and accurate evaluation of the integral yields (cf. [14]) To calculate the three-body collision rate Γ 3b per atom, we need to multiply Γ k1k2k3 by the product of populations N f kj of the states |k j = a † kj |vac , j = 1, 2, 3 (the occupation probabilities are normalized to unity, ∞ −∞ dk f k = 1) and integrate the whole k 1 , k 2 , k 3 -space, divide by 3! to take into account the indistinguishability of the bosons. Then we obtain the number of three-body collisions per unit time in the whole system. This number should be divided by N to obtain the rate per atom: The result of Eq. (28) may seem counterintuitive at first: the collision rate is independent of temperature, and it is proportional to the dimensional parameter ζ 2 and the radial confinement ω r . The physics behind the first observation is related to the fact that the collision kinetic energy is small compared to the virtual excitation energy. This was one of our assumptions in deriving the effective three body collisions and is required by the condition (k B T <hω r ) to be fulfilled when building a 1D system (Eq. (1)). Consequently the composite matrix element of the second-order process should not depend in leading order on the velocities (energies) of the colliding particles and hence on temperature (see Eq. (9)). In addition the phase space volume for the scattered particles is independent on the incoming momenta k 1 , k 2 , and k 3 .
The other terms can be motivated the following basic physics considerations: Since effective three-body elastic scattering is the dominant process the scattering rate must be proportional to the 3D density squared: (n 1D /l 2 r ) 2 . Furthermore, the scattering rate contains the square of the matrix element corresponding to the diagram in Fig. 1(a), where each vertex is proportional to α s , therefore this rate is proportional to α 4 s . The factorh/m provides the correct dimensionality (s −1 ).
We can now compare the scattering rates for the two-body collisions Γ 2b (Eq. (7)) or effective three-body collisions Γ 3b (Eq. (28)) and evaluate their contributions to thermalization and the breakdown of integrability in 1D systems. For k B T <hω r we find a simple scaling: The relative importance of two-body collisions and the effective three-body collions mediated by virtual excitations is determined by the dimensionless quantity ζ e Ratio between the scattering rates for the two routes to thermalization in quasi-1D systems: Γ 3b for the effctive threebody collisions and Γ 2b for two-body collisions leading to excited transverse states. The points represent the ratios evaluated for various sets of experimental parameters from [11] (points), [12] (crosses), [13] (triangle), and [23] (diamonds). In these experiments the parameter ζ was often close to 0.007 (the ratio Γ 3b /Γ 2b for ζ = 0.007 exactly is shown by the solid curve). For comparison, we plot also Γ 3b /Γ 2b for ζ = 0.002 (dashed curve) and 0.02 (dot-dashed curve). Units on the axes are dimensionless.
scattering rate due to virtual excitations dominates, and can lead to thermalization even when the thermalization processes due to two-body collisions are frozen out. For example, in typical atom chip experiment [11][12][13] the three body reate Γ 3b dominates above the two-body collisions at k B T ≤ 1 2h ω r . A detailed comparison of the two rates Γ 2b and Γ 3b and their relation to typical experimental parameters is given in Fig. 2. The scattering rate due to virtual excitations of the radial modes can dominate over real excitations for typical parameters of the recent experiment [12].
The above calculation was for a non-degenerate ultracold gas. In a degenerate gas we need to consider the Bogoliubov-type spectrum of elementary excitations [22], that is phononic in the long-wavelength regime, as well as the relation between atoms and elementary excitations via the Bogoliubov transformations and bosonic amplification of scattering to modes, which are initially occupied.
Taking into account all these factors, we find the rate of damping of a fast particle in a quasicondensate (see also Ref. [24]): Here the momenta of the scattered elementary excitations are defined by the expressions reciprocal to Eq. (24,25): The energy of a mode with the momentumhk is ε k =h 2 η k /(2m) with η k = k 2 (k 2 + 8n 1D α s /l 2 r ).
The static structure factor of a quasicondensate is In equilibrium the population of the elementary mode with the momentumhk is given by the Bose-Einstein statistics with the mean occupation number for the mode with the momentumhk The initial kinetic energyh 2 k 2 0 /(2m) of the fast atom is assumed here to be large compared to both the mean-field interaction energy per particle in the quasicondensate and the temperature: Under condition (32) one of the scattered particles is always fast, and one of the three structure factors appearing in Eq. (30) is always very close to 1 (and the corresponding occupation number is close to 0). In the most scattering events the other two particles are also fast, and, hence, S k ′ j ≈ 1 and n k ′ j ≈ 0 for all three particles, j = 1, 2, 3. Only for the scattering events with small transferred momentum two of the structure factors are significantly less than 1 and/or the corresponding populations approach the high-temperature limit k B T /ε k ′ j . However, in the practically interesting case where k B T ∼hω r n 1D α s the contribution of scattering events with small transferred momentum is relatively small, and The result of Eq. (33) for a fast atom in a quasicondensate is also obtained by Tan, Pustilnik and Glazman [24].

V. CALCULATIONS OF THE THERMALIZATION RATES
We now turn to quantify thermalization in tightly confined Bosons in a quasi-1D geometry by both two body collisions and the effective interaction (14). We again for simplicity consider a non-degenerate, weakly-interacting (the Lieb-Liniger parameter [5] γ = 2α s /(n 1D l 2 r ) being much less than 1) gas of bosonic atoms. The assumptions of weak interaction and non-degeneracy enable us to express the three-particle distribution function through the product of single-particle distribution functions f k . In contrast, calculation of relaxation via three-body collisions of low-energy excited states in the stronger interacting regime and especially for γ ≥ 1, requires to consider the (strong) quantum correlations in the quasi-1D bosonic system. We will discuss the effects of correlations on scattering rate Γ 3b and on thermalization in section VI We start be writing the Boltzmann equation with a three-body collision integral [25], taking into account the indistinguishability of the particles: with Eq. (34) can be easily understood: After integration over k ′ and k ′′ the loss term in Eq. (34) is simply −Γ 3b f k , which is the elastic three-body collision rate per atom. On the other hand, the three-atom state |k, k ′ , k ′′ is populated by elastic three-body collisions from those states |K 0 , K −1 , K +1 which have the same center-of-mass momentum: Since the kinetic energy of the relative motion,h 2k2 /(2m) is conserved, the states |K 0 , K −1 , K +1 (from where the state |k, k ′ , k ′′ can be populated from) can be fully parametrized by the hyperangle γ.
We now use the following Ansatz for the perturbed momentum distribution to solve Eq. (34). Therby k th = √ 2mk B T /h and H 4 is the Hermite polynomial of the 4th order. This Ansatz is the simplest nontrivial perturbation that retains dk kf k = 0. We then proceed to linearizing Eq. (34) with respect to the perturbation amplitude ε 4 (t) and obtain an exponential solution ε 4 (t) = ε 4 (0) exp(−Γ 3b [4] t) with Γ 3b [4] = C [4]h m n 1D l 2 r 2 α 4 s = C [4] ω r ζ 2 (38) with the numerical constant To estimate the validity of our Ansatz we note that using a higher-order Hermite polynomial H n in Eq. (37) leaving the functional dependence on the parameters of the system unchanged and leads only to a minor modification of the numerical prefactor. For example, for n = 5 and 6 the thermalization rates are given by 10 27 C 3b ω r ζ 2 and 34 81 C 3b ω r ζ 2 , respectively. It is interesting to note that in these 1D systems the thermalization rate due to the effective three-body collisions (Γ 3b [4] ) is about a factor 3 smaller then the collision rate (Γ 3b ). This suggests that in 1D systems thermalization requires also about 3 collisions, similar to 3D [26]. Fig. 3 shows numerical values of Γ 3b [4] as a function of the 1D density of 87 Rb atoms and the radial trapping frequency.
For comparison we calculate numerically the thermalization rate Γ 2b [4] for two-body collisions involving the real transitions between the ground and excited radial states. We follow hereby the same Ansatz and perturb the velocity distribution of atoms in the ground and excited state as given by Eq. (37), the Boltzmannian distribution of overall populations between the levels being kept intact. In the parameter range of interest we find numerically Γ 2b The ratio of the thermalization rates for the two-body and three-body processes is therefore very close to the respective ratio of the collision rates, shown in Fig. 2.
It is interesting to note that we find for both processes that thermalization in 1D needs about 3 collisions capable to distribute energy. This is very close to the 2.7 collisions required for thermalization in 3D [26].
For the typical parameters of an ultracold 87 Rb gas on an atom chip [12] (ω r ≈ 2π × 3 kHz, n 1D ≈ 50 µm −1 ) we obtain Γ 3b [4] ≈ 2 s −1 . This thermalization rate is temperature-independent and much larger than the one calculated from the simple two-body collisions with the energy sufficient to excite radial modes Γ 2b [4] ≈ 3 × 10 −3 s −1 at the lowest temperatures measured (30 nK). The estimated Γ 3b [4] is consistent with the time needed for evaporative cooling of a 87 Rb gas on an atom chip well belowhω r [11,12].

VI. SUPPRESSION OF THERMALIZATION BY ATOMIC CORRELATIONS
The thermalization rate Γ 3b [4] given by Eq. (38) was calculated for a weakly-interacting, non-degenerate gas. Calculation the thermalization rate Γ 3b G [4] in a general case requires to take into account additional physics. First we need to consider the effects of quantum degeneracy and second the fact that the disperison relations for the elementary excitations in a 1D quantum (degenerate) system may differ significantly compared to a free particle, especially for phonon-like excitations. These effects, together with the bosonic amplification of the scattering into thermally populated modes, tend to accelerate thermalization. A third observation is that the three body rates require three particles to be close to the same location. Quantum mechanically this is characterized by the third order correlation function g 3 (0). A full consideration of the above competing effects will require extended numerical analysis of many particular cases and transcends beyond the scope of the present manuscript. We will give here physical arguments of what to expect.
We start by pointing out that the form of the secondary-quantized Hamiltonian Eq. (16) allows us to give a simple estimate of the ratio of these two rates Here ̺ is a phase-space factor accounting for the dispersion law of elementary excitation, which changes from freeparticle-like to phonon-like. But it changes the thermalization rate less dramatically than the second factor associated with the local three-body correlation function and will become the dominating factor when approaching the strongly correlated regime (γ > 1). For a non-degenerate weakly-interacting Bose gas g 3 (0) = 3! = 6. For a degenerate 1D Bose gas g 3 (0) has been recently calculated for the whole range of atomic repulsion strength (0 < γ < ∞) by Cheianov et al. [27]. In the zero-temperature limit g 3 (0) rapidly decreases from 1 to 16π 6 /(15γ 6 ) as γ grows from 0 to values γ ≫ 1. We now turn to the above conjecture on suppression of the thermalization by atomic correlations. A detailed calculation can be found in Ref. [15], here we sketch the basic physics argument. To look at the correlations we start by considering N identical bosons in 1D configuration with the Hamiltonian 16, which, after rescaling of units, takes the formĤ (41) c = 2α s /l 2 r is the strength of interaction of two atoms in the tight waveguide with ground state size l r . U 3b is obtained by by adiabatic elimination of transverse modes virtually excited by the 3D short-range pairwise atomic interaction [14]. The explicit form of U 3b is given by Eq. (17) within a numerical prefactor, We follow now our detailed calcultions in Ref. [15] and estimate the three-body scattering amplitude in the presence of the delta-functional pairwise repulsive interactions. The stronger the pairwise interparticle repulsion, the smaller is the probability of a close encounter of three particles which will result in a suppression of the three-body scattering amplitude. The simplest case is to analyze the Hamiltonian (41) for N = 3 particles. For that purpose we express it in hyperspherical coordinates R, χ defined as [21] and obtain for the Hamiltonian The corresponding Schrödinger equation for the three-particle wave function iŝ The wavenumbers k j are defined from the set of transcendental equations [5], provided that the periodic boundary conditions are set on the interval of the length L. By setting L → ∞, we obtain a continuous spectrum, where k j 's are real for repulsive interaction (c > 0). We can now separate the center-of-mass motion and describe the relative motion in hyperspherical coordinates. This leads to the Ansatz : where the kinetic energy of the relative motion is given by . In the adiabatic hyperspherical approximation [30], which holds in the long-wavelength limit and neglects coupling of different scattering channels as well as accumulation of phase shifts of the scattered wave due to non-adiabatic effects we get The hyperangular part B 0 (R, χ) of this wave function is the eigenfunction of the Hamiltonian (43) with fixed R, corresponding to the lowest eigenvalue λ 0 (R), which is the smallest positive root of the transcendental equation Using the regular hexagon symmetry group [31] of the Hamiltonian (43), we write the hyperangular part of Eq. (45) as whereB 0 (R, χ) = cos[λ 0 (R)(π/6 − |χ|)], |χ| ≤ π/3 0, otherwise After integrating out the hyperangular variable, the Schrödinger equation in the adiabatic hyperspherical approximation reduces to with the boundary conditions requiring F 0 to be finite for both R = 0 and R → ∞, for the "partial wave" corresponding to the lowest eigenvalue λ 0 (R), whose asymptotic expressions are We can solving now Eq. (49) analytically in two regions, cR ≪ 1 and cR ≫ 1, with λ 0 (R) approximated by Eq.
(51), and tailoring the solutions by quasiclassical expressions for F 0 (R) in the intermediate hyperradius range. The scattering amplitudef 0 (for its definition in planar geometry see [32][33][34][35]) can then be obtained from the asymptotic form of the wave function at R → ∞ where H 3 (z) = J 3 (z) + iY 3 (z) is the Hankel function of the first kind and J 3 (z), Y 3 (z) are the third-order Bessel functions. The behavior of F 0 at R → 0 is defined by the details of the potential U 00 (R), but the result can be finally expressed via the effective vertex of the three-body elastic collisions, thus yielding where Ω ≈ 1 is a numerical constant and ξ ≈ 1.15 is defined by Eq. (15). We use thereby the fact that U 00 (R) differs significantly from zero on the length scale l r , over which a virtually excited particle can propagate, and cl r ≪ 1. From Eq. (53) we conclude that the three-body scattering amplitude decreases in proportion to (k/c) 6 as k/c → 0, i.e., when the interaction is strong enough to induce significant atomic correlations. The three-body scattering rate in a 1D system of bosons in the case of strong pairwise interaction is suppressed by a factor ∼ (k/c) 12 . By averaging over collision momenta in a moderately-excited strongly-interacting state we obtain the scattering rate suppression factor ∼ (k/c) 12 ∼ γ −12 .
We can now compare our result with the zero-distance three-particle correlation function g 3 (0). In the strong interaction limit γ ≫ 1 g 3 (0) ∝ γ −6 [27,36] and gives a direct physical motivation of our original conjecture [14] that the pairwise interactions and the quantum correlations induced by them in a strongly-interacting 1D bosonic system suppress the three-body elastic scattering rate, and, hence, thermalization, by a factor ∝ g 2 3 (0). In other words, strong quantum correlations extend the time scale, on which a quasi-1D system approaching the Tonks-Girardeau regime can be considered as approximately integrable.
Moreover, we can corroborate this conjecture by observing that the thermalization rate is proportional to the square of the matrix element of the transition operator that is proportional to dzψ †ψ †ψ †ψψψ . Since the energies of the products of the elastic three-body process are low (of about k B T ), we may assume that the correlations in the initial and the final states are the same, and the transition matrix element can be regarded as proportional to g 3 (0) that yields again the g 2 3 (0) scaling of the rate. In contrast to this, the inelastic three-body processes are accompanied by a large energy release, and after an inelastic collision the newly formed dimer molecule and the fast atom leave the system almost immediately. Therefore the inelastic three-body relaxation rate in a 1D ultracold Bose gas is proportional to the first power of g 3 (0) [37].
We now compare the calculated thermalization rates to the quantum Newton's cradle experiment [10] where a lower boundary for the damping time of the 1D motion towards a Gaussian profile was estimated. For three different Lieb-Liniger parameters γ = 1.4, 3.2, 18 Kinoshita et al estimate lower bounds to the thermalization time from the consistency of the the experimental momentum distributions with the experimentally observed heating during the time interval of 0.5, 0.5, 1.0 seconds probed. They find one sigma lower limits of 2.6 s, 25 s, and 13 s respectively. In their experiment the motion of two groups of 87 Rb atoms was excited at the relative velocity equal to 4 recoil velocities, which is far above the width of the ground-state velocity distribution. We therefore can not expect the collision rate to be suppressed in proportion to g 2 3 (0) that, as described above for slow collisions. Instead, we have to apply the estimates for damping of a fast particle discussed in the end of section IV. In a strongly-correlated system one has to take the particle correlations into account and eq. 33 has to be multiplied by a factor g 2 (0) denoting the two-particle correlation function at zero distance (Tan, Pustilnik and Glazman [24]): Substituting the experimental parameters of Ref. [10] and taking the values for g 2 (0) from Ref. [38], we obtain Γ damp k0 ≈ 15 s −1 , 1.7 s −1 , and 7 × 10 −3 s −1 for γ = 1.4, 3.2, and 18, respectively. To compare this calculated damping rates to the experiment in [10] one has to consider that (1) the two colliding clouds overlap only for a very short time during each oscillation and (2) that the thermalization rate is a factor 3 longer (section V). Taking this into account we estimate the respective thermalization times of 2.6 s, 35 s and > 1000 s. These rates are consistent with the experimental findings of Ref. [10]. For a more detailed comparison one would need longer time scale experiments with lower intrinsic heating and more detailed calculations of the dynamics including the damping due to three body collisions discussed here and in [24].

VII. CONCLUSION
A radially confined atomic gas is never perfectly 1D, and radial motion can be excited, either in reality or virtually even if Eq. (1) holds. This possibility leads to effective three-body collisions, which arise in the second order of perturbation theory and can be associated with virtual excitation of radial modes. These processes lead to thermalization even when two body collisions are frozen out at k B T ≪hω r and provide a mechanism to break integrability in 1D systems. In other words, the freeze-out of the radial modes is only a necessary, but not sufficient condition for integrability in 1D systems. Our estimations of the relaxation rates for weakly interacting quasi-1D Bose gases are consistent with recent experimental observations for weakly-interacting quasicondensates [11,12].
These effective three-body collisions can be suppressed by quantum correlations caused by strong pairwise repulsions. If they dominate, as in a strongly correlated 1D Tonks-Girardeau gas, they suppress the influence of the integrabilitybreaking interaction term. The thermalization rate decreases in proportion to g 2 3 (0) as the system enters the regime of strong correlations (γ ≫ 1), and the system (remaining non-integrable in the strict sense) behaves like (almost) integrable on time scales short compared to the inverse thermalization rate.
The effective three-body collisions, and their suppression by quantum correlations should be accessible in experiments looking at the damping of fast, particle-like excitations in systems with γ varying in a broad range of values from less than 1 to ∼ 10. This work is supported by the EC (STREP MIDAS) and the FWF.