Roton in a few-body dipolar system

We solve numerically exactly the many-body 1D model of bosons interacting via short-range and dipolar forces and moving in the box with periodic boundary conditions. We show that the lowest energy states with fixed total momentum can be smoothly transformed from the typical states of collective character to states resembling single particle excitations. In particular, we identify the celebrated roton state. The smooth transition is realized by simultaneous tuning short-range interactions and adjusting a trap geometry. With our methods we study the weakly interacting regime as well as the regime beyond the range of validity of the Bogoliubov approximation.


I. INTRODUCTION
In the 30s of the last century unusual properties of the Helium-II were discovered. The subsequent results of Allen and Misener [1], Kapitza [2] were simulating the development of theoretical models [3][4][5][6][7]. The qualitative theory of superfluidity is due to Landau [5][6][7]. He deduced from the measurement of the specific heat [8] and the second sound velocity [9] that the excitations in the Helium-II must have a peculiar spectrum, with the local minimum [7]. The excitation at the local minimum has been called a "roton". Later Feynman alone [10] and with Cohen [11] formulated the very first, yet semiquantitative microscopic model explaining the origin of this local minimum. Finally, in Helium the roton was observed experimentally [12], but rather unsatisfactory agreement between theory and measurement suggested that the exact nature of the rotonic excitation was still missing. It was finally understood many years later by means of subtle ansatzes for the roton's wave function [13,14]. It should be emphasized that liquid Helium-II is a strongly correlated (with a small condensate fraction) system, where roton's characteristic momentum scales as the interatomic distance. There are still active studies of the roton state in this regime [15].
At the beginning of XXI century the roton-maxon spectrum was predicted in completely different physical systems -dipolar gas of ultracold atoms in constrained geometries [16,17]. Unlike in Helium-II, in this case the interactions are weak and a condensate fraction is dominant. Therefore, one can use the mean field or Bogoliubov description and find the roton state as a Bogoliubov quasi particle [16][17][18][19][20][21][22][23][24][25][26][27][28][29][30][31][32][33][34]. The dispersion curve of such systems is related to a specific k-dependence of an effective interaction potential rather than to strong correlations. Possibility of changing the particles polarization as well as almost free tuning of the short range interactions combined with the trap geometry modifications enables unprecedented flexibility in the study of the roton spectrum in dipolar gases [16][17][18][19][20][21][22][23][24][25][26][27][28][29][30][31][32][33][34] ending with a recent experimental confirmation of the phenomenon [35]. Usually the dipolar system is studied within the Bogoliubov approximation, so that there is no access to the detailed structure of the low lying excitations. Only a few many-body investigations were performed for the roton state using different techniques [36][37][38][39]. A good attempt can be made by a numerically exact solution to a many-body problem with a rotonic characteristic. Even if found for relatively small number of particles, modern experiments with a precise control over only a few atoms in optical lattices or single traps (see for instance [40][41][42][43][44]) allow to test its physical predictions.
In this work we present numerically exact results for a quasi-1D model, which admits the roton excitation. In a number of recent papers low excitation of 1D interacting bosons were already investigated (see for instance [45][46][47][48][49][50] and references therein). The historically earliest example is the famous Lieb-Liniger model [51,52] comprising of N contact interacting bosons moving on the circle. Their seminal analytical solution predicts two branches of elementary excitations, which was also observed experimentally [42]. The upper type-I excitation branch was immediately recognized as the Bogoliubov excitation spectrum [53]. The states of the lower type-II excitation branch, were identified later with grey and dark solitons arising in the mean field theory of ultracold gases [54][55][56][57]. Little is known about the classification of exact many-body elementary excitations in the dipolar gas. For (quasi)-1D model with bosons interacting only by repulsive dipolar interactions the lowest energy states resemble rather type II excitations from the Lieb-Liniger model [48,58] and for at least weakly interacting regime the picture with two branches of elementary excitations is also expected in this case [59]. On the other hand, in the dipolar systems, well understood Bogoliubov spectrum may exhibit a local minimum identified as a roton [16,17]. When the interatomic forces are of the attractive character on the short-scale, whereas the long-range part of potential is repulsive, the interplay of these two interactions may lower the energy of the roton mode even to the ground state level. It opens a significant question: is it possible in a dipolar analogue of the Lieb-Liniger model that the two branches cross, such that it is a type-I Bogoliubov excitation, in particular the roton, which would appear in the lower branch?
It is a purpose of this work to show that by tuning short-range interactions and adjusting a ring geometry one can continuously change [60] the character of the lowest energy state for a given total momentum of the system from a type-II excitation to the roton mode. We also analyze a numerically exact roton's wave function in the weakly interacting regime and its position and momentum properties.

II. MODEL
We consider N dipolar bosons confined in both transverse directionsŷ andẑ with a tight harmonic trap of a frequency ω ⊥ . Multi-particle wave-function is approximately the Gaussian in tight directions for all variables. It requires the chemical potential µ much smaller than energy of the first excited state in the transverse direction, µ ≪ ω ⊥ . In the longitudinal directionx the space is assumed to be finite, with the length L and with the periodic boundary conditions imposing quantisation of momenta in that direction. All atoms are polarized along theẑ axis. The above system corresponds to atoms moving on the circumference of a circle, having the dipole moments perpendicular to the circle-plane. Hence, in analogy with nuclear physics [61,62] and following [54,55] we call the lowest energy states of a given total momentum of the system, the yrast states. Our quasi-1D system is governed by Hamiltonian withâ k (â † k ) anihilating (creating) a boson with momentum k. The effective potential consists of the longrange dipolar part and the short-range part, namely The quasi-1D dipolar potential reads V dd (k) = where d is an atomic dipole moment and µ 0 is the vacuum permeability. This effective quasi-1D potential comes from integration of the full 3D dipolar interaction over both transverse variables. The singular part coming from this integration is incorporated with the short range interaction. The function f which appears in Eq. 1 is equal to f (u) = u e u Ei(u), where Ei is the exponential integral [63]. Stability of our calculations requires smoothing of a usual short range interaction model used in the ultracold physics, the delta function. We choose a Gaussian model [64][65][66][67][68][69][70][71][72], namely V sr (k) = V 0 e − 1 2 k 2 r 2 with r standing for the potential range and |V 0 | for its depth. This step makes our model more realistic, imitating the attractive van der Waals interaction. For convenience we set V 0 = 2 a ml 2 ⊥ with a mimicking an usual scattering length. The relation between Gaussian model and the real scattering length can be found in [73] and references therein. Below we use box units where L/2π, 2π /L and 4π 2 2 /mL 2 are the units of length, momentum and energy respectively.
Our effective potential V eff (k) corresponds to calculating the interactions along the circumference positions with the periodicity of the system accounted for. However, we checked that it would not be changed significantly if a bit elegant but more realistic geometric distance over the chord was used (see Appendix A).
We access the many-body eigenstates of Hamiltonian (1) by exact diagonalization using the Lanczos algorithm [74]. Our calculations are performed in the Fock space spanned by the plane-wave basis with a maximum total kinetic energy of the system E max = k 2 max /2 -with single-particle momentum k max ≫ 1/r -sufficiently high to assure convergence. Here we employ the fact that the total momentum of the systemK = k kâ † kâ k is conserved, Ĥ ,K = 0, so its eigenvalues K are good quantum numbers, used here together with the total number of atoms N to label different eigenstates |N, K, i enumerated by i with i = 0 corresponding to an yrast state. We remind the Reader that for finite systems on the ring it suffices to consider the eigenstates only up to K/N = 1/2 [52,55]. This comes from the presence of the so called umklapp process [52]. Any eigenstate with a total momentum ) may be understood as the state with a total momentum K with a shifted center-of-mass momentum (see Fig. 5 in [55]). Note that such shifting does not change the internal structure of the state.

III. RESULTS
In the following paragraphs of this work we consider two different situations, namely with weak interactions where the depletion (given by P (K = 0) in Fig. 1d and 2d) of a ground state is less than 5% and stronger interactions where its value is around 20%. Note, that the latter case is still far from the Helium-II regime. We present in Fig. 1 our analysis of yrast states for the first situation. We consider N = 16 dysprosium atoms with a dd = 132 a 0 and the potential range r = 182 a 0 , where a 0 is the Bohr radius. We initially set a and ω ⊥ corresponding to the usual situation where the yrast states energies clearly do not follow the Bogoliubov spectrum (black dashed line) given by: and rather resemble the lowest excitation branch from the Lieb-Liniger model [58,59](black squares in Fig. 1a).
Then we continuously change a and ω ⊥ (a simillar effect would be observed if one changed the length of the box) keeping V eff (0) = const. We finally end with the profoundly different spectrum (red points in Fig. 1a) closer to a corresponding Bogoliubov dispersion relation (red dashed line), in particular with the characteristic inflection point for K = 2. Our result suggest that at least some of yrast states may change their character from collective type-II excitations [58] to type-I ones. Moreover, we argue that the inflection point can be identified with the roton-like state.
To test our hypothesis about the change of the character of the yrast state for K = 2 we compare it together with the first excited state with the number conserving Bogoliubov approximation [75] sketched here briefly. The spectrum in the Bogoliubov approximation is explained by the concept of quasiparticles that, in our case, has to be rewritten in terms of Fock states in particle basis. We use the following Ansatz [76] for the Bogoliubov vacuum (K = 0): where |vac is the particle vacuum and u k , v k = To trace a continuous transformation of the yrast state from type-II to type-I excitation we evaluate the fidelities | N, K, i|N, K B | 2 , which is depicted in Fig. 1b. For the initial values of the parameters a and ω ⊥ the first excited state is a Bogoliubov excitation, whereas the yrast state remains a type-II excitation [58]-a fact observed in the Lieb-Liniger result as well. Then we observe a gradual exchange of the states' character as we modify the effective potential ending with a complete role reversal of the two first states. Note, that the sum of the fidelities (black dotted line in Fig. 1 b) is almost equal to 1 at any stage of the transition. It means that Bogoliubov excitation, to a good approximation, remains in a plane spanned by the two lowest eigenstates.
To show the qualitative change of the yrast state for K = 2 we calculate the normalized second order correlation function g 2 (z) := Ψ † (z)Ψ † (0)Ψ(0)Ψ(z) / Ψ † (z)Ψ(z) Ψ † (0)Ψ(0) (Fig. 1c), which can be measured in experiments with ultracold atoms, see for instance [77][78][79][80]. We observe a dramatic difference between two yrast states for border cases from Fig. 1b (marked as black and red points). Namely that almost flat distribution typical for type-II excitation in weakly interacting regime is replaced by a function exhibiting an enhanced regular modulation with the number of maxima given by K rot , which is the roton momentum.
Note, that for a small number of particles we are able to find stable solutions corresponding to realistic, physical gas parameters (N V sr (0) = −23.98, N V dd (0) = 26.98) only for the Bogoliubov spectrum with the inflection, not to the one with the characteristic local minimum. Using the gas parameters for which the Bogoliubov spectrum has the local minimum implies much stronger interactions for our few-body system. Our result would approach Bogoliubov's predictions in the limit of N → ∞ (see Appendix B). Instead of going to larger systems, we turn to the strong interactions scenario with N = 10, which is beyond the Bogoliubov approximation. In Fig. 2 we summarize our findings, where the characteristic local minimum for K = 3 is present. In this case the spectrum is calculated with an accuracy of several percent only. Our previous conclusions hold also for this situation. However, the role reversal of the lowest states is more subtle because of higher momentum of the roton. At the end of the transition we stay with the yrast state, which still has the overlap with Bogoliubov (50%) and at the same time exhibits the enhanced regular modulation in second order correlation function and the local minimum in the spectrum. It is the roton-like state in a regime between weak interaction and Helium-II regime.
Note that the presence, position, and depth of the roton minimum for both cases considered in this work are tunable by varying the number of atoms N , trapping frequency ω ⊥ and the short-range coupling strength as it was predicted for the roton state in the meanfield studies of ultracold dipolar gases [16]. We choose K rot /N < 1/2 to minimize the impact of the umklapp process [52], discussed earlier in this work, on the eigenstates.
To fully comprehend the difference between the two types of low-energy excitations we study the probability P (k) = 1 N â † kâ k of finding a single-particle moving with momentum k for yrast states with various K. For both weak and strong interactions the type-II yrast states (black markers in Fig. 1d and 2d) for K > 1 beside k = 0 mainly consists of k = 1 states, which is more visible as we increase K. It corresponds to a dominant role of one of the Dicke states (exactly K atoms with k = 1 and N − K with k = 0) in their many-body wave function [58,59], especially for weak interactions. On the other hand, in the rotonic cases (red markers in Fig. 1d and 2d) we observe a local maximum of P (k) for k = K rot for the yrast states with K rot , which clearly resembles recently published result by F. Ferlaino's group [35]. It means that the yrast state for K rot has a single particle excitation character rather than a collective one, so that within our, experimentally achievable, procedure one can completely change the character of the low-energy excitations.

IV. DISCUSSION
We find with our numerically exact treatment that all the properties of the roton state discussed earlier can be understood by analysing contributions of different Fock states to its wave function. In both cases of interactions studied in this work, we find that the dominant contribution to the roton states comes from the so called W state |0 −kmax , ... (N −1) 0 , 0 1 , 1 2 , ..., 0 kmax , as one would expect for the Bogoliubov excitations [59]. The latter state is important from the fundamental point of view, as representative of an entanglement class [81], and applied side -it can be used to beat the standard quantum limit for the metrological tasks [82]. The state was recently produced via non-demolition measurement [83]. According to our earlier findings [59], which holds also for purely dipolar repulsion [58], the low-lying excitations of weakly interacting bosons are highly-entangled states dominated by the Dicke state, a result of the bosonic statistics mainly. However, the interplay between the short-range and longrange interactions of the opposite sign can promote the excitation with the dominant W state as a low-lying excitation for K > 1 in the system.

V. CONCLUSION
To summarize, we showed that manipulating physical parameters in our model one can continuously alter the character of a given yrast state from type-II excitation to the roton mode. We emphasise the fact that the effect is already present in relatively small systems enabling use of the simplest exact diagonalization of the whole Hamiltonian. All interesting properties of the roton-like mode both in the momentum and the position representations come from the fact, that the W state plays the dominant role in the roton state in the plane wave basis. It is in the stark contrast to the weakly repulsive bosons, where the dominant role of the Dicke states is observed [58,59]. We show that the normalized second order correlation function, accessible in experiments, displays characteristic enhanced regular modulation for the roton state. Within our many-body model we access stronger regimes between the weakly interacting one and the Helium-II scenario, finding the roton mode also in this case. Our results open new questions concerning quasi-1D systems with both long-range and short-range interactions. Is it possible to fully replace type-II branch with type-I branch as low-lying excitations? Would solitonic branch still exist in the spectrum? The thermodynamic properties of dipolar bosons were investigated only approximately, in the weakly interacting regime [84][85][86]. The results presented in this paper can motivate further research in this direction, but using full many-body approach accounting for the lower branch and the transitions discussed here. In the main text we use the effective potential (in the momentum representation) V eff (k), that originates as follow. For quasi-1D model on the infinite line, the effective potential in the space representation U 1D (x) takes the form: Erfc |q| √ 2 . (A1) As we consider a finite system with the periodic boundary conditions, we introduce U periodic (x) = n∈Z U 1D (x + nL). From Poisson summation formula it satisfies . However, if one wants to deal with a real ring with atoms moving on its circumference the effective potential should rather depend on the geometric distance over the chord U ring (x) = U 1D L π sin( πx L ) . In Fig. 3 we compare both approaches. As we see both curves are almost indistinguishable in the regions where the value of the effective potential is meaningful. A very small difference in all cases from Fig. 3 is observed only on the potential tail.  In the Bogoliubov approximation one operates with the gas parameters N V sr (0), N V dd (0) (in the box units defined in the main text) with the assumption of weak interactions and large number of atoms N . Obviously, in the many-body approach, where N is finite, the energy of the pairwise interactions is significantly higher. Then, one can ask how many atoms (how weak interactions) one needs to converge with the many-body solution towards N → ∞ limit. To answer it, we study energies of a series of yrast states (left panel of Fig. 4) and their overlaps with the corresponding Bogoliubov excitations given by fidelities (right panel of Fig. 4) defined in the main text. We obtain both the spectrum and the fidelities for different number of atoms N ranging from 7 to 16. The parameters for different N are chosen to always produce the same Bogoliubov excitation spectrum with the inflection point as for red dashed line in Fig. 1a in the main text. We see that even for small number of atoms N = 16 we obtain very good overlap with the Bogoliubov approximation, especially for K ≤ 2. However, our numerically exact solution includes all the possible correlations between atoms, hence it cannot be fully reproduced by single Bogoliubov excitation.