Two-dimensional dynamics of expansion of a degenerate Bose gas Igor E. Mazets1,2 1Vienna Center for Quantum Science and Technology, Atominstitut, TU Wien, 1020 Vienna, Austria 2Ioffe Physico-Technical Institute of the Russian Academy of Sciences, 194021 St.Peterburg, Russia Expansion of a degenerate Bose gas released from a pancake-like trap is numerically simulated underassumption of separation of motion in theplane of theloose initial trapping and the motion in the direction of the initial tight trapping. The initial conditions for the phase fluctuations are 2 generated using the extension to the two-dimensional case of the description of the phase noise by 1 the Ornstein-Uhlenbeck stochastic process. The numerical simulations, taking into account both 0 thefinitesizeofthetwo-dimensionalsystemandtheatomicinteractions,whichcannotbeneglected 2 on the early stage of expansion, did not reproduce the scaling law for the peaks in the density fluctuation spectra experimentally observed by Choi, Seo, Kwon, and Shin [Phys. Rev. Lett. 109, t c 125301(2012)]. Thelatterexperimentalresultsresultsmaythusrequireanexplanationbeyondour O current assumptions. 5 PACSnumbers: 67.85.-d,03.75.Lm,03.75.Hh ] s a Correlations in degenerate Bose gases with repul- being equal to (ωx, ωy, ωz) = 2π (3.0, 3.9, 370) Hz. × g sive interactions manifest themselves most apparently Thechemicalpotentialwasaboutµ =2π¯h 260Hz,i.e., 0 - via phase coherence. Even low-dimensional degenerate less than the spacing between the discrete×levels of the t n bosonic systems (quasicondensates) demonstrate phase trapping Hamiltonian in the tighgtly confined direction. a coherence over finite distances. After release from the Thecoherentfractionofatomicensembledecreasedfrom u trap, phase fluctuations are converted, in the course of 0.78to 0.12 as the temperaure T grew from 20 nK up to q free expansion, into density-density correlations(density the BKT transition temperature T 70 nK. . BKT t ≈ a ripples). Theory for expansion of both one-dimensional The experimental results [5] turned out to be quite m (1D) and two-dimensional (2D) trapping geometries has surprising. Instead of recovering the scaling (2), they been developed [1] and experimentally proven for 1D brought about another dependence, - d trapped gases [2] (similar behavior has been previously n found in expanding clouds, which were in the three- ¯hqn2te/(2πm)=αnγ, (3) o dimensional Bose-Einstein condensate rather in the 1D with0.2<α<0.45and0.7<γ <1for expansiontimes c [ regime inside a very elongated trap [3]). in the range 10 ms < te < 25 ms. In other words, the Incontrastto1DBosegases,wherecorrelationsdecay spectral peaks of the density fluctuations were not only 1 exponentially and their characteristic length λT deter- shifted with respect to Eq. (2), but the spacing between v mines the typical wavelengthscale of the density ripples the ajacent peaks was significantly smaller than Eq. (2) 3 1 and the timescale mλ2T/¯h of their emergence (m being predicts. Choi et al. ruled out the effects of the cloud fi- 8 themassoftheatom),correlationsindegenerate2DBose nitesizeandsuggestedthatsuchanunexpectedbehavior 1 gases decay with the distance according to a power law, might be intrinsic to the expansion dynamics [5]. . if the temperature is below the point of the Beresinskii- To check this assumption, we performed numerical 0 1 Kosterlitz-Thouless (BKT) transition [4]. Such a power- simulations of 2D expansion of a degenerate gas. Since 2 lawdecayrendersnospecificcorrelationlength,thespec- thetemperatureand,hence,thecoherentfractiondonot 1 trumofthedensityfluctuationsofa2Dgasreleasedfrom affect significantly the peak positions q [5], we find it n : atrapevolvesinaself-similarwayatasymptoticallylarge safe to neglect thermal population of the excited lev- v i expansion times te. In particular, after averaging the els of motion in the tightly trapped (z) direction and X power spectrum of density fluctuations model the system’s dynamics by the Gross-Pitaevskii r equation (GPE), taking into account the interaction ef- a 2 P = d2ρδn(ρ)e−iqρ , (1) fects, which are important at the initial stage of expan- q sion. Furthemore, we separate the motion in z-direction (cid:12)Z (cid:12) (cid:12) (cid:12) and in the (x,y)-plane. This would be impossible if the where δn(ρ) is the lo(cid:12)cal density fluctua(cid:12)tion, over the di- (cid:12) (cid:12) Beliaev and Landau damping of elementary excitations rections of the 2D wave vector q, the nth peak position in expanding but still dense enough atomic cloud were q , n=1,2,3, ..., is expected to satisfy the relation [1] n giving rise to significant population of modes with non- ¯hq2t /(2πm)=n 1/2. (2) zero z-components of kinetic momentum. However, the n e − estimation of the ratio of the damping rate of an ele- Recently, Choi et al. [5] reported their experimental mentary excitation to its frequency yields [6] the value results on probing fluctuations in a 2D gas of bosonic a1/2/(n1/2λ2 ) 10−2, where a is the atomic s-wave (23Na) atoms in free expansion after a sudden release ∼scatstering3Dlengthth,∼n is the threse-dimensional atomic 3D fromapancake-likeopticaltrap,thetrappingfrequencies number density, and λ is the de Broglie wave length of th 2 atomsattemperatureT. Evenforelementaryexcitations ofour numericalsimulations. Ona (2M+1) (2M+1) × with the energy close to the chemical potential the typi- square grid Eq. (8) takes the form calrelaxationtime(and,hence,thecharacteristictimeof M M scatteringintoz-direction)inathree-dimensionaldegen- ε Z = const dϕ exp w einraotuerg2aDswsyitshtetmheisbouflkthdeenorsditeyreoqfu1a0l0tomtsh,ewpheicahkidsemnsuicthy lx=Y−Mly=Y−MZ lxly (− 2 lxly × longerthanthetypicalexpansiontime. Thereforewecan separatethemotionindifferentdirections. Inzdirection, [(ϕlxly −ϕlxly−1)2+(ϕlxly −ϕlx−1ly)2]). (9) we have a free expansion of a wave packet, which is ini- tially, at t = 0, the wave function of the ground state Here the pair of subscripts l and l denotes a vari- x y of the harmonic trapping potential with the fundamen- able at the point with coordinates x = l ∆s and x tal frequency ω . Then the motion in the (x,y) plane is z y = l ∆s [note that the grid step ∆s does not ap- y described by the GPE pear in Eq. (9) in our 2D problem]. The density profile and its peak value are characterized by dimen- ∂ ¯h2 ∂2 ∂2 i¯h∂tΨ=−2m(cid:18)∂x2 + ∂y2(cid:19)Ψ+µ(t)|Ψ|2Ψ. (4) asinodnleεss=qu¯han2nti2tDie(s0,w0l)x/l(ym=kBTn2),D(rlexs∆pesc,tliyve∆lys.)/nP2eDr(i0o,d0ic) The fast (but still not infinitely fast) evolution of boundary conditions are set, i.e., ϕlxM+1 = ϕlx−M and the quasicondensate profile in z-direction affects the ϕM+1ly =ϕ−Mly. slow motion in the perpendicular plane) only through Assumeweknowthe2M+1phasesϕlx−1ly. Thenwe cancalculatethephasesϕ inthenextrowofthegrid the time-dependent nonlinearity in Eq. (4). The lxly using the transfer matrix approach[9]. time-dependent effective 2D coupling constant µ(t) Weassumethatthegridstepissmallenoughandthus takes into account the growth of the transverse pro- file width 1+(ωzt)2. The normalization condition we can write wlxly/wlx0 ≈ wlx−1ly/wlx−10. The trans- ∝ formation [10] max Ψ(x,y,t=0) = Ψ(x = 0,y = 0,t=0) 1 fixes | p | | | ≡ the prefactor, and we obtain M w µ(t)= 1+µ(0ωzt)2. (5) Alx−1j =lyX=−MVlx;jlyrwllxxl0yϕlx−1ly (10) diagonalizes the non-negative quadratic form Initialshape of Ψ2 is tphe inverted-parabolicThomas- | | Fermi profile corresponding to the chemical potential µ0 M w 2M of the trapped gas of sodium atoms. The phase ϕ = wlxly(ϕlx−1ly −ϕlx−1ly−1)2 = JlxjA2lx−1j. argΨ at t=0 can be generated using the generalization lyX=−M lx0 Xj=0 of the stochastic Ornstein-Uhlenbeck process, previously (11) used to model random phase distributions of phases in BothV andJ canfoundnumericallybystandard lx;jly lxj 1D quasicondensates [7, 8]. methods of linear algebra. Then Eq. (9) factorizes: Z = AssumethatthetemperatureT of2DBosegasisbelow const 2M Z , where each factor j=0 j the BKT transition, phase and density fluctuations are small and, therefore, the Hamiltonian can be linearized: Q M ε H = ¯h2 d2ρ n2D(ρ)( ⊥ϕ)2+g˜δn2+ (∇⊥δn)2 . Zj = lx=Y−MZ dAlxj exp"− 2lx(Alx−1j −Alxj)2− 2m ∇ 4n (ρ) Z (cid:20) 2D (cid:21)(6) εlxJlxjA2 (12) 2 lxj Here n = dzn is the 2D number density, g˜ is the # 2D 3D dimensionlesscouplingstrength(g˜ 0.013intheexperi- R ≈ corresponds to a 1D evolution (along the x-coordinate) ment[5]), ⊥ isthe2Dgradientinthe(x,y)-plane. Then ∇ of a A describable by the Ornstein-Uhlenbeck stochas- the partition function of a trapped gas is j tic process (which can be easily seen from the anal- ogy with the relative phase between two tunnel-coupled Z = ϕ δn exp[ H/(k T)]. (7) D D − B 1D quasicondensates at thermal equilibrium [7]). Here Z Z ε =εw . The formula [11] lx lx0 The density fluctuation variable can be integrated out, aZnd=wceonosbttainDϕ exp −2m¯hk2 T d2ρn2D(ρ)(∇⊥ϕ)2 . Alxj =Alx−1je−√Jlxj +vuu12−εlex−2√JlJxljxjς, (13) Z (cid:20) B Z (cid:21) t p (8) where ς is a (pseudo)random number obeying the Gaus- We can use the thermodynamic expression(8) to gen- siandistribution with zero mean andunity variance,up- erateinitialphasedistributionforindividualrealizations dates the variable A on the next step. Then we obtain j 3 7 instead of recovering Eq. (3), we found that the peak positions are well described by the formula 6 ¯hq2t /(2πm)=β(n ℓ), (15) Ls 5 n e − unit 4 where b. Har 3 β ≈0.9, ℓ≈0.7 (16) q P 2 forthe wholerangeofexperimentallyrelevantexpansion times (see Fig. 2). 1 Note, that the transversal trapping frequency ω = z 0 2π 370 Hz of Ref. [5] is approximately equal to the 0.0 0.2 0.4 0.6 0.8 1.0 1.2 × radial trapping frequency in the experiment of Ref. [3], q HΜm-1L but the cigar-shaped configuration of the atomic cloud resulted in the latter case in a faster ( t−2 instead of FIG.1: Solidline: thepowerspectrumofdensityfluctuations t−1 at expansion times >1 ms) decre∝ase of the three- fort =11msobtainedfromEq. (4),i.e.,takingintoaccount ∝dimensional density of ato∼ms than in the 2D case. Ra- e interaction of atoms in an expanding atomic cloud. Parame- dial trapping frequencies in recent experiments with 1D ters of the trapped degenerate cloud are taken from Ref. [5] ultracold gases on atom chips [2, 8] were by an order of (see details in the text). Dashed line: spectrum of density magnitudehigher,thusfurtherdecreasingthetimescale, fluctuations obtained for a purely ballistic expansion during on which atomic interactions still play a role. Therefore thesamet . Thespectraareaveragedover5realisationseach. e expansion of atomic clouds in the experiments [2, 3, 8] wasveryclosetoballistic,incontrasttothecasemodeled in our present work. the initial phase values in the next grid row: Our numerical results do not reproduce the observa- M tions by Choi et al. [5]. For example, our calculations w ϕlxly = Vlx;jly wlx0 Alxj. (14) predict the values of the 5th peak position q5, which are lyX=−M r lxly at relevant expansion times larger by a factor 2–3 than the observed ones. In general, our calculations predict Subsequent steps are done using the new localvalues for the linear depenfdence of q2 on the peak number n, in n the density profile. contrasttothenonlinearfittingformula(3)foundinRef. Eq. (13), together with chosing the (pseudo)random [5]. initial values of the Aj variables according to a Gaus- Such a conspicuous disagreement between the theory sian distribution with zero mean and variance hA2lxji = that assumes separation of motion in the former trap 1/(2ε J ), enables us to prepare the initial condi- lx lxj tions for eachrealizationof our numericalmodeling [12]. p æ We solved Eq. (4) using the split-step pseudospectral 0.8 method for the typical parameters of Ref. [5] and ob- tained the power spectra of density fluctuations. As the æ à density fluctuation δn in Eq. (1) we took the difference 0.6 between the local densities given at the same expansion L ´ timet bytwosolutions,onecontaininginitialphasefluc- -1 à tuatioens and other having constant phase everywhere at HΜm 0.4 æ ´ Dƒ t=0. AnexampleofsuchaspectrumPq,thatisPqaver- 2qn à ƒ agedoverthe directionofthe wavevectorq, isdisplayed ´ D in Fig. 1. The variation of the temperature from 60 nK æ ƒ 0.2 à D to 20 nK (corresponding to the increase of ε from 15 to ´ Dƒ 45) changes the height of the peaks of P only, leaving the peak positions q unmodified, aswe caqnexpect from æà´Dƒ n 0.0 theexperimentalresults[5]. ThesolutionoftheGPE(4) 1 2 3 4 5 is also juxtaposed with the density fluctuations spectra n for purely ballistic expansion [i.e., by setting µ(t) 0]. ≡ In the case of purely ballistic expansion, we found that FIG.2: Thepeakpositionsforthepowerspectrumofdensity Eq. (2) is satisfied with 1% accuracy for n 2; for the fluctuations (the peak number n is dimensionless) of a 23Na first peak position we obtained ¯hq12te/(2πm≥) ≈ 0.4 for ultracold gasinteractingduringtheexpansion forte=11ms 10 ms t 23 ms. (filled circles), 14 ms (filled squares), 17 ms (crosses), 20 ms e As w≤e can≤see fromFig.1, interactionofatoms during (open squares), and 23 ms (triangles). Straight lines (solid, expansion is not negligible and shifts the peak positions long-dashed, dotted, short-dashed, and dot-dashed, respec- tively) display thefittingby Eqs. (15, 16). with respect to the ballistic expansion case. However, 4 plane and in the direction of the tight confinement be- times as strongly as at longer ones. fore the release of the atom from the trap and the ex- The factthatthe experimentalspectra ofdensity fluc- perimental results [5] calls for explanantion. We do not tuations are describable by Eq. (3) when no vortices expect the scattering of phonons into z-direction in the are present in the 2D system [5] allow us to rule out expandingquasicondensatetobe areasonfor sucha dis- anyvortex-relatedexplanationofthisdiscrepancy,which crepancy, since the estimations based on the theory of therefore still remains puzzling. Ref. [6] show that this process is too slow compared to the typical expansion times. Moreover, the influence of This work was supported by the FWF (project sucha scatteringwouldgrowwith t , but our results de- P 22590-N16). The author thanks J. Schmiedmayer and e viate strongly from the experimental findins at shorter Yong-il Shin for helpful discussions. [1] A. Imambekov, I. E. Mazets, D. S. Petrov, V. Gritsev, [7] H.-P. Stimming, N. J. Mauser, J. Schmiedmayer, and I. S. Manz, S. Hofferberth, T. Schumm, E. Demler, and J. E. Mazets, Phys.Rev.Lett. 105, 015301 (2010). Schmiedmayer,Phys.Rev. A 80, 033604 (2009). [8] T. Betz, S. Manz, R. Bu¨cker, T. Berrada, Ch. Koller, [2] S. Manz, R. Bu¨cker, T. Betz, Ch. Koller, S. Hoffer- G.Kazakov,I.E.Mazets, H.-P.Stimming,A.Perrin,T. berth, I. E. Mazets, A. Imambekov, E. Demler, A. Per- Schumm, and J. Schmiedmayer, Phys. Rev. Lett. 106, rin,J.Schmiedmayer,andT.Schumm,Phys.Rev.A81, 020407 (2011). 031610(R) (2010). [9] J. Zinn-Justin,Quantum Field Theory and Critical Phe- [3] S.Dettmer,D.Hellweg,P.Ryytty,J.J.Arlt,W.Ertmer, nomena, 4th edition (Clarendon, Oxford, 2002), Ch. 2 and K. Sengstock, D. S. Petrov, G. V. Shlyapnikov, H. and 23. [4] KLVe.rteLtu..t8Bzm7e,rae1nz6nin0,4sLk0i6.i,S(Sa2on0vt0o.1sP),.haynsd.JME.TLPew34en,s6t1e0in(,1P97h2y)s;.JR.eMv.. [10] Th2je=M0rVellxa;tjiloynVslx;jPl′yMly==−δlMyl′yVlhxo;jldly.Vlx;j′ly = δjj′ and Kosterlitzand D.J. Thouless, J.Phys.C61181 (1973). [11] PD. T. Gillespie, Phys. Rev.E 54, 2084 (1996). [12] As a consistency check of our method, we can demon- [5] Jae-yoonChoi,SangWonSeo,WooJinKwon,andYong- il Shin,Phys. Rev.Lett. 109, 125301 (2012). strate that the description via the Ornstein-Uhlenbeck [6] L.P. Pitaevskii and S. Stringari, Phys. Lett. A 235, 398 stochastic process yields the correct correlation function −1/(2πε) (1997); C. J. Pethick and H. Smith, Bose-Einstein Con- hexp[iϕ(0,0)−iϕ(x,y)]i∝ x2+y2 . densation in Dilute Gases, 2nd edition (Cambridge Uni- versity Press, Cambridge, 2008), p.313. (cid:16)p (cid:17)