Thermalization of a quantum Newton's cradle in a one-dimensional quasicondensate

We study the nonequilibrium dynamics of the quantum Newton's cradle in a one-dimensional (1D) Bose gas in the weakly interacting quasicondensate regime. This is the opposite regime to the original quantum Newton's cradle experiment of [Kinoshita et al., Nature 440, 900 (2006)], which was realized in the strongly interacting 1D Bose gas. Using finite temperature $c$-field methods, we calculate the characteristic relaxation rates to the final equilibrium state. Hence, we identify the different dynamical regimes of the system in the parameter space that characterizes the strength of interatomic interactions, the initial temperature, and the magnitude of the Bragg momentum used to initiate the collisional oscillations of the cradle. In all parameter regimes, we find that the system relaxes to a final equilibrium state for which the momentum distribution is consistent with a thermal distribution. For sufficiently large initial Bragg momentum, the system can undergo hundreds of repeated collisional oscillations before reaching the final thermal equilibrium. The corresponding thermalization timescales can reach tens of seconds, which is an order of magnitude smaller than in the strongly interacting regime.

Figure 1

Evolution of the position-space (left) and momentum-space (right) density distributions in the quantum Newton's cradle setup of a harmonically trapped 1D quasicondensate, generated from $c$-field simulations. The three examples displayed all have parameter values of $\overline{\omega }=0.0696$ and ${\chi }_0=0.0562$ for the initial thermal equilibrium state, whereas the wave number for the Bragg pulse varies as follows: (a) $q_0=1.16$, (b) $q_0=2.06$, and (c) $q_0=3.67$.

Figure 2

Momentum distributions of a harmonically trapped 1D quasicondensate in the quantum Newton's cradle setup. In each example, we show the initial (blue, dotted line), relaxed (red, full), and thermal fitted distributions (black, dashed). The initial and relaxed distributions displayed in (a)–(c) are, respectively, for the same parameters as in Fig. 1, i.e., for $({\chi }_0,\overline{\omega )=(0.0562,0.0696)}$ initially and three different values of $q_0$ as shown, whereas the thermal fitted distributions, plotted for comparison with the relaxed ones, are for the following fitting parameters: (a) $({\chi }_0,\overline{\omega })=(0.250,0.045)$, (b) $({\chi }_0,\overline{\omega })=(0.238,0.0462)$, and (c) $({\chi }_0,\overline{\omega })=(0.500,0.150)$.

Figure 3

Comparison of the total internal energy difference per atom (circles) between the initial and relaxed states with the kinetic energy per atom (dashed line) of ${\hslash }^2{\left(2k_0\right)}^2/2m$ [or $4q_0^2$, in units of ${\hslash }^2/\left(2m{x}_0^2\right)]$, transferred by the symmetric Bragg momentum kick of $\pm 2\hslash k_0$, as functions of the dimensionless momentum $q_0$, for $({\chi }_0,\overline{\omega })=(0.0562,0.0696)$. The data points represent the numerically calculated values of the average internal energy, averaged over 128 stochastic realizations, with the error bars indicating one standard deviation from the average; the points labeled (a)–(c) correspond to the examples of Fig. 1. All energies have been scaled by ${\hslash }^2/\left(2m{x}_0^2\right)$; when multiplied by the total number of atoms $N$, this translates to $N{\hslash }^2/\left(2m{x}_0^2\right)=\frac{1}{2}{k}_{B}T\overline{N}$, where $\overline{N}=\int {\left|\phi (\xi ,\tau )\right|}^{2}d\xi =4\sqrt{2}/\left(3{\chi }_0\overline{\omega }\right)$ and $T$ is the initial temperature.

Figure 4

Power spectra of the population of the lowest energy harmonic oscillator mode for a single realization of the simulations of Figs. 1–1, respectively. We have used ${\mathrm{\Omega }}_{\tau }$ to denote the dimensionless Fourier frequency component of the dimensionless time $\tau$ and expressed it in units of the dimensionless trap frequency $\overline{\omega }$. It can be seen that all spectra are continuous.

Figure 5

Diagram of the dynamical regimes for the quasicondensate Newton's cradle in the ${\chi }_0$ versus $q_0$ parameter space, with $\overline{\omega }=0.0696$. The crossover boundaries (solid and dashed lines) between the different regimes, corresponding to fast (regions I and II) and slow (region III) relaxation, equate to the right-hand sides of Eqs. (14) and (15). The circles labeled (a)–(c) represent, respectively, the parameter combinations used in the examples displayed in Figs. 1, 2, and 8 (see text). The colored rectangle corresponds to the scanned parameter space of ${\chi }_0$ and $q_0$, for which we extracted the dimensionless decay rate $\overline{\mathrm{\Gamma }}$ from the numerical simulations (see Sec. 3d). As we see, the crossover boundary between the regions II and III predicted from simple arguments is in good agreement with the quantitative picture emerging from the simulations of the decay rate.

Figure 6

Momentum distributions averaged over one ($m\mathrm{th}$) oscillation period, $n_m\left(k\right)$, for parameter values corresponding to $({\chi }_0,q_0,\overline{\omega })=(0.0562,3.67,0.0696)$, i.e., example (c) in Fig. 1. Distributions for $m=1,30,100,200,$ and 300 are shown. The average distributions for the earlier oscillations, $m=1$ and 30, still have two separate peaks at the initial momentum kicks $\pm 2q_0$ and a smaller third peak at $q=0$, indicating there is still some oscillatory dynamics present in the system, while the large-$m$ averages approach stationary thermal distributions as those in Fig. 2.

Figure 7

Momentum distributions averaged over one ($m\mathrm{th}$) oscillation period, $n_m\left(k\right)$, for parameter values corresponding to $({\chi }_0,q_0,\overline{\omega })=(0.0562,3.67,0.0696)$ as in Fig. 6. However, here we initialize the system with 10% of particles remaining unperturbed by the Bragg pulse (see text). We see that in this case the distinct three-peaked distribution, observed in Refs. [9, 32], does occur initially ($m=1$), in contrast to Fig. 6, however, the three peaked structure eventually dampens out as can be seen from the example for $m=30$.

Figure 8

Relaxation of the rms width $W_m$ of the average momentum distribution over the $m\mathrm{th}$ oscillation period as a function of oscillation number $m$. The symbols correspond to the numerically evaluated values of the width $W_m$, for the same parameters as in Fig. 1, whereas the solid, dashed, and dash-dotted lines are based on exponential fits of Eq. (18); the widths are rescaled by their respective values of $W_1$ to be easily visible on the same graph. Note that the exponential fits begin after the respective prethermalization periods which are are not well approximated by a single exponential. The long-time relaxation rates $\overline{\mathrm{\Gamma }}$ extracted from these fits are as follows: (a) $\overline{\mathrm{\Gamma }}=0.06$, (b) $\overline{\mathrm{\Gamma }}=0.025$, and (c) $\overline{\mathrm{\Gamma }}=0.0065$.

Figure 9

Comparison of the power spectra of $|a_0{\left(t\right)|}^2$ for GPE simulations with $M=16$ for increasing normalization, corresponding to the infinite temperature simulations as described by Bland et al. [38]. Note we have used ${\mathrm{\Omega }}_{\tilde{t}}$ to denote the corresponding Fourier dimension of the dimensionless time $\tilde{t}$. The normalizations are (a) $N=15$, (b) $N=30$, (c) $N=50$, and (d) $N=150$. For $N=15$, the spectra is quasidiscrete, but for larger values of $N$ the spectra are continuous.

Figure 10

Comparison of different randomized initial conditions for the same normalization of $N=15$, with $M=16$. Left: Comparison of the power spectra of $|a_0{\left(t\right)|}^2$ for two different randomized initial conditions (black and grey lines) as described by Bland et al., with normalization $N=15$. Right: Comparison of the respective initial densities (top) and phases (bottom).