Interacting ultracold atomic kicked rotors: loss of dynamical localization

We study the fate of dynamical localization of two quantum kicked rotors with contact interaction, which relates to experimental realizations of the rotors with ultra-cold atomic gases. A single kicked rotor is known to exhibit dynamical localization, which takes place in momentum space. The contact interaction affects the evolution of the relative momentum k of a pair of interacting rotors in a non-analytic way. Consequently the evolution operator U is exciting large relative momenta with amplitudes which decay only as a power law 1/k4. This is in contrast to the center-of-mass momentum K for which the amplitudes excited by U decay superexponentially fast with K. Therefore dynamical localization is preserved for the center-of-mass momentum, but destroyed for the relative momentum for any nonzero strength of interaction.

The quantum kicked rotor (QKR) model is a canonical model to explore quantum chaos 1,2 . It describes a quantum rotor degree of freedom which is periodically kicked by a force periodic in the angle. The QKR enjoys dynamical localization (DL) -i.e. the arresting of the growth of the momentum despite the absence of a cutoff in the frequency of the kick drive. DL was first discovered numerically by Casati, Chirikov, Ford, and Izrailev 3 and later confirmed experimentally for Rydberg atoms in a microwave field 4,5 and ultracold atomic gases in a modulated standing wave of a near-resonant laser 6 . A recent work reports on the experimental observation of DL with laser-kicked molecular rotors 7 . Coupled kicked rotors can also be realized by placing atoms in pulsed, incommensurate optical lattices 8 . If the driving period is an irrational multiple of 2π, the rotor is localized in the momentum space, even though the classical counterpart shows diffusive momentum growth. This happens because classical chaotic diffusion is suppressed by quantum interference effects. The mechanism of DL was described in a seminal paper by Fishman, Grempel and Prange 9 . These authors demonstrated that the kicked rotor model maps directly to an Anderson-like model with a quasi-periodic potential, which originates from the irrational driving periods. Therefore DL is closely related to Anderson localization of waves in truly random (uncorrelated) potentials.
The original quantum kicked rotor corresponds to a single quantum particle problem. The effect of interactions on Anderson localization has been attracting a lot of interest recently and several theoretical studies considered various versions of interacting kicked rotors. In ref. 10 a similar problem was studied for a simpler, integrable model of linear rotors 9 , where localization can survive in the presence of interactions due to integrability. The authors of ref. 11 analyzed coupled relativistic rotors which might be applicable to fermions in pulsed magnetic fields, and reported that DL can be destroyed by suitable parameter tuning. In ref. 12, two kicked rotors with product sinusoidal interaction at the kick were studied with respect to temporal fluctuations in the reduced density matrix. In ref. 13, the coupling was sinusoidal depending on the two rotors relative coordinates: recovering of the chaotic behavior was found above some kicking threshold in the semi-classical approximation. In ref. 14, the interaction at the kick of the kicked rotors contained both product and relative coordinate dependent sinusoidal terms. Localization was found for weak coupling and quasi-diffusive regime was found for stronger interaction with a complex intermediate regime.
From the experimental perspective, interaction between rotors is negligible for Rydberg atoms and laser-kicked molecular rotors. However the interaction between ultracold atoms in a Bose-Einstein condensate (BEC) can be substantial, and even tunable using Feshbach resonances 15 , which is particularly true for sodium atoms used in ref. 6. The atom-atom interaction in this case is typically of a contact type, i.e. the atoms interact through a δ(x 1 − x 2 ) potential 15 . For the experimental realization in ref. 6, this interaction of BEC atoms persists at all times -in contrast to the kick potential, and in contrast to the theoretical studies discussed above, which Scientific RepoRts | 7:41139 | DOI: 10.1038/srep41139 consider a kicked (time-dependent) interaction. A δ(x) interaction is long-ranged in momentum space, and can therefore have a qualitatively strong impact on DL for interacting ultracold atoms. Will DL survive, or not?
In this work, we address this question. We consider two bosons interacting via a δ-function potential that are driven by a periodic kicking potential. The wave function for two δ-function interacting bosons is computed. We use the center of mass and relative coordinates with appropriate periodic boundary conditions. A repulsive δ-interaction is considered, that does not lead to the appearance of a bound state. In the chosen basis the matrix elements of the time evolution operator show different decay rates along the center of mass (superexponential) and the relative momentum (algebraic) directions. Due to this qualitative difference in the decay properties of the matrix elements, dynamical localization (i.e. exponential localization with a finite localization length) is destroyed for the relative momentum, while being preserved for the center-of-mass momentum.

Model and Methods
We consider two bosons moving on a ring [0, 2π) with δ-function interaction and periodic kicking potential. The Hamiltonian of the model is given by: 2 is the two-body Hamiltonian of the Lieb-Liniger model 16 2 , λ is the interaction strength, ξ is the kicking strength. In the experimental setup of ref. 6, the kicking potential, a standing wave of a near resonant light, is cigar-shape and not on a ring. That will in general change the boundary conditions, and momentum selection rules, as compared to a ring structure with periodic boundary conditions. These changes are however not of qualitative nature for non-interacting particles, while multiple collisions for interacting particles could yield additional sensitivity to boundary conditions. However, our main results will be obtained from analyzing the impact of only one kick. Furthermore, recent toroidal Bose-Einstein condensates have been experimentally realized which emulate periodic boundary conditions 17 . If DL is destroyed in just a two-body interacting model, we expect it also be destroyed in the experimental many-body interaction case. Therefore we can view our model as a first step towards the consideration of a many-body interacting system and a simple paradigmatic case of just two interacting atoms which is the building block of reaching out to many-body interactions.
We start by computing the wave function of the two bosons system with δ-function interaction. It can be represented in a center of mass and relative coordinates frame -(y 1 , y 2 ), where In this frame, the first part of (1) reads It splits in two parts, describes a free moving particle and H y 2 describes a single particle with δ-function potential. The wave function of the complete system is φ φ φ = y y y y ( , ) ( ) ( ) 1 2 1 1 2 2 -the product of two single particle wave functions φ 1 and φ 2 , that satisfy . The total eigenenergy of the system is = + E E E y y 1 2 . Because of the periodicity, the complete wave function satisfies φ φ π = + = y y y y ( , ) ( 4 , ) . This can be simplified into three identities: , which serve as the periodic boundary conditions for φ 1 and φ 2 .
The wavefunction for the free moving particle is: select the quantized values of K: e i4Kπ = 1, giving , and the eigenenergy describes a massive particle on a one-dimensional ring of circumference 4π and a δ-function singularity at y 2 = 0. The eigenstates of this problem can be either symmetric or antisymmetric around y 2 = 0. Since rotors are bosons, the wave function should be invariant under permutation x 1 ↔ x 2 and only symmetric functions φ 2 (y 2 ) = φ 2 (− y 2 ) are allowed. Then the derivative is an antisymmetric function The wave function is continuous at 0: φ 2 (+0) = φ 2 (−0), but its derivative has a jump: , it follows that φ 2 (y 2 ) = e i2Kπ φ 2 (y 2 + 2π) and φ 2 (y 2 ) = φ 2 (y 2 + 4π). It is worth noting that the center of mass and relative momenta do not decouple completely, due to the boundary conditions.
With these boundary conditions it is easy to compute the wave function φ 2 (see Supplementary material for more details): Scientific RepoRts | 7:41139 | DOI: 10.1038/srep41139 with A λ = ħ 2 /Mλ being a dimensionless inverse interaction strength and the eigenenergy We use A λ to measure the strength of the interaction, to which it is inversely proportional. We attach the momenta K and k to the full wavefunction -φ y y ( , ) K k 1 2 . Now the eigenstates of two bosons with δ-function interaction can be written as φ K k , the corresponding wave function is The eigenenergy is given by . The wavefunctions are symmetric with respect to k: if k is a solution of Eq. (6) then − k is also a solution. As follows from Eq. (4), the wave function with k and − k are exactly the same for integer K or differ by a global phase π for half integer K. Consequently φ − y y ( , ) K k 1 2 is equivalent to φ y y ( , ) K k 1 2 , reflecting the bosonic nature of the rotors. In the following discussion, we only consider the k > 0 case.
In the presence of a periodic driving potential in Eq. (1), the dynamics is described by the time evolution operator over one period U (Floquet propagator) 18 . The matrix of U is computed below in the basis of φ K k (see also the Supplementary material). Starting from some initial state |ψ(0)〉 of the rotors, the final state |ψ(NT)〉 after N driving periods is obtained the initial state by matrix multiplication:

Results and Discussion
For periodically kicked interacting rotors, the Floquet operator reads where H δ is given by Eq. (2) and H k = 2ξ cos(y 1 /2) cos(y 2 /2). In the basis of φ K k , the matrix elements of U become . For a fixed A λ and large enough k and r (such that λ kA , λ  r A 1), the matrix elements of U decay as with the asymptotics along the center-of-mass momentum direction following from the asymptotic properties of Bessel functions 19 (Notice that the dependence on |K − R| is through the index of the Bessel function). Therefore, for large k and r, the matrix element U KR kr decays super-exponentially fast with the center of mass K momentum (i.e. faster than exponential and asymptotically as a factorial), due to the scaling of 2  which is controlled by the first order derivative of the Bessel function (13) (the same conclusion holds for any values of r and k, see Supplementary material). The decay along the relative momentum k direction however is a power law k −4 , reflecting the presence of a singular δ-function interaction. This is our key result: super-exponential decay of the matrix element ensures the survival of dynamical localization for the center-of-mass momentum, while the power-law decay destroys it for the relative momentum (in the sense of exponential localization with a finite localization length). Smearing the δ-function of the interaction into a smooth (analytic) one will introduce an energy (wave number) cutoff k max and corresponding exponential matrix element decay for energies beyond that cutoff, ensuring survival of the DL for the relative momenta in this case. This follows from the observation that the matrix element is essentially a Fourier transform of a function related to the interaction potential with respect to k (see Supplementary material). It is non-analytic in the case of the δ-interaction (since its first derivative has a jump), while it is analytic for smooth (analytic) interaction potentials. The asymptotic convergence of a Fourier transform is exponential for analytic functions, but only a power law for the non-analytic ones.
The comparison of the numerical results to the asymptotic behavior is presented in Fig. 1. The top figure shows the decay of matrix element U KR kr with relative momentum k for several values of the coupling A λ indicated by colors. The other momenta r, K, R are fixed. The power law fit = − − . U k log 4 log 0 2 KR kr 10 10 (the black line) agrees well with the numerical values of the matrix elements U KR kr with a given K, R, r and  k r. The small-k dependence is not sensitive to the inverse interaction strength A λ . The power-law decay of the matrix elements is not monotonic with A λ : initially (blue, red and green curves) the prefactor is decreasing, however upon further decrease of A λ (magenta curve) it starts to increase and there is a non-monotonic intermediate part. This non-monotonicity can be explained from Eqs (11 and 12): for very small and very large A λ , f kr is small, behaving as 2rA λ and 1/(2kA λ ) respectively. Therefore, with decreasing A λ , and for given k, K, r, R, the matrix element of U (11) will first increase from a small value and then decrease back, which is precisely the non-monotonicity observed. For small A λ the power-laws have the same prefactors: this follows from (12) -for a fixed r and This is what we observe for the smallest values of A λ in Fig. 1. The bottom plot in Fig. 1 shows the decay of the matrix elements as a function of K: a faster than exponential decay is observed, agreeing with the asymptotic behavior (13). The prefactor also shows non-monotonicity -the matrix elements initially increase with increasing A λ and later decrease -and has the same origin as above.
The impact of the decay properties of the matrix elements of U is observed in the evolution of an initial state ψ φ = (0) K k 0 0 with fixed momenta K 0 and k 0 according to Eq. (8). This choice of initial state is well-suited to detect delocalization in the momentum space. The final state after N driving periods is are the expansion coefficients over the eigenstates of the two-body Lieb-Liniger problem. Figure 2 shows the final state after N = 100 -left column -and N = 5000 -right column -driving periods for two different values A λ , corresponding to moderately strong interaction between the rotors. The final state gets extended along the relative momentum k direction. Also the extension is more pronounced with decreasing A λ . Figure 3 shows the amplitude distribution of the final state after N = 5000 driving periods for several values of A λ . Figure 3(a) shows the final state for the case of two essentially non-interacting rotors (A λ = 10 14 ). The final state is localized in both K and k momenta directions, displaying the dynamical localization of the non-interacting kicked rotor model. As the strength of the interaction is increasing, Fig. 3(b,c), the final state starts to extend along the k direction. The interaction between the two rotors delocalizes the state in the relative momentum k direction. The localization length along the center of mass momentum K direction also increases, as seen in Fig. 3(b,c).   final state have values similar to the (c) case, since the matrix elements of U become independent of A λ as we have discussed earlier.
In order to quantify the spreading of the initial state with N, we compute the evolution of the variance of the momenta: Figure 4 shows the evolution of the variance with N for several strengths of the interaction A λ . For all the values of the interaction except A λ = 10 the variance has a clearly increasing trend spanning several orders of magnitude in N, therefore signaling delocalization along the relative momentum k direction. The non-monotonic dependence of the variance on A λ has the same origin as the non-monotonic dependence of the matrix elements of U that we discussed above.

Conclusions
In conclusion, we studied dynamical properties of two interacting kicked quantum rotors. Due to the non-analyticity of the δ-function interaction, the matrix elements of the time evolution operator exhibit different decay behaviors in center-of-mass and relative momentum directions. Along the center of mass momentum direction, matrix elements decay super-exponentially. Along the relative momentum direction, matrix elements decay as a power-law with the exponent 4. As a result the center-of-mass motion remains localized like in the non-interacting case, while the relative motion becomes extended. The destruction of dynamical localization should be observable in toroidal traps 17 and setups similar to the one used in ref. 6 using Feshbach resonances. Furthermore, while contact interaction is a good approximations, our conclusions hold for other regularized (analytic) longer ranged interactions as well, at least up to the corresponding energy cutoff which will separate the convergence criteria of the Fourier series of the contact interaction potential from its regularized version in momentum space. The extension to many interacting particles is also an interesting path. To analyze this, we need to consider many interacting atoms and study the highly complex case of many-body interactions for quantum kicked rotors. This is a challenging task, as can be observed from the complexity of results on mean field treatments in refs 20-22 which range from complete destruction of dynamical localization to modifications of Anderson transitions. However, we strongly believe, that in the many-body version of our model, this mechanism still holds, leading to the destruction of DL, and is even more efficient due to the increased number of relaxation channels.