Quasi-one-dimensional Bose gases with large scattering length

合集下载

Bose-Einstein condensation of scalar fields on hyperbolic manifolds

a rXiv:h ep-th/9213v11O ct1992Bose-Einstein condensation of scalar ﬁelds on hyperbolic manifolds Guido Cognola and Luciano Vanzo Dipartimento di Fisica -Universit`a di Trento ∗,Italia and Istituto Nazionale di Fisica Nucleare,Gruppo Collegato di Trento september 1992PACS numbers:03.70Theory of quantized ﬁelds05.90Other topics in statistical physics and thermodynamics1IntroductionBose-Einstein condensation for a non relativistic ideal gas has a long history[1].The physical phenomenon is well described in many text books (see for example ref.[2])and a rigorous mathematical discussion of it was given by many authors[3,4].The generalization to a relativistic idel Bose gas is non trivial and only recently has been discussed in a series of papers[5,6,7].It is well known that in the thermodynamic limit(inﬁnite volume andﬁxed density)there is a phase transition of theﬁrst kind in correspondence of the critical temperature at which the condensation manifests itself.At that temperature,the ﬁrst derivative of some continuous thermodynamic quantities has a jump.If the volume is keepedﬁnite there is no phase transition,nevertheless the phenomenon of condensation still occurs,but the critical temperature in this case is not well deﬁned.For manifolds with compact hyperbolic spatial part of the kind H N/Γ,Γbe-ing a discrete group of isometries for the N-dimensional Lobachevsky space H N, zero temperature eﬀects as well asﬁnite temperature eﬀects induced by non-trivial topology,have been recently studied in some detail[8,9,10,11,12,?,14,15,16]. To our knowledge,a similar analysis has not yet been carried out for non compact hyperbolic manifolds.Hyperbolic spaces have remarkable properties.For example,the continuous spec-trum of the Laplace-Beltrami operator has a gap determined by the curvature ra-dius of H N,implying that masslessﬁelds have correlation functions exponentially decreasing at inﬁnity(such a gap is not present for the Dirac operator).For that reason,H4was recently proposed as an excellent infrared regulator for massless quantumﬁeld theory and QCD[17].Critical behaviour is even more striking.In two ﬂat dimension vortex conﬁgurations of a complex scalarﬁeld,the XY model for He4ﬁlms have energy logarithmically divergent with distance,while on H2it isﬁnite. This implies that the XY model is disordered at anyﬁnite temperature on H2.Even quantum mechanics on H2has been the subject of extensive investigations[18].The manifold H4is also of interest as it is the Euclidean section of anti-de Sitter space which emerges as the ground state of extended supergravity theories.The stress tensor on this manifold has been recently computed for both boson and fermion ﬁelds using zeta-function methods[19].In the present paper we shall discussﬁnite temperature eﬀects and in particular the Bose-Einstein condensation for a relativistic ideal gas in a3+1dimensional ultrastatic space-time M=R×H3.We focus our attention just on H3,because such a manifold could be really relevant for cosmological and astrophysical applications.To this aim we shall derive the thermodynamic potential for a charged scalarﬁeld of mass m on M,using zeta function,which on H3is exactly known.We shall see that the thermodynamic potential has two branch points when the chemical potentialµriches±ωo,ω2o=κ+m2being the lower bound of the spectrum of the operator L m=−△+m2and−κthe negative constant curvature of H3.The values±ωo will be riched byµ=µ(T)of course for T=0,but also for T=T c>0.This is the critical temperature at which the Bose gas condensates.The paper is organized as follows.In section2we study the elementary properties of the Laplace-Beltrami operator on H3;in particular we derive its spectrum and build up from it the related zeta-function.In section3we brieﬂy recall how zeta-function can be used in order to regularize the partition function and we derive the regularized expression for the thermodynamic potential.In section4we discuss the Bose-Einstein condensation and derive the critical temperatures in both the cases of low and high temperatures.In section5we consider in detail the low and high temperature limits and derive the jump of theﬁrst derivative of the speciﬁc heat.The paper end with some considerations on the results obtained and some suggestions for further developments.2The spectrum and the zeta function of Laplace-Beltrami operator on H3For the aims of the present paper,the3-dimensional Lobachevsky space H3can be seen as a Riemannian manifold of constant negative curvature−κ,with hyperbolic metric dl2=d̺2+sinh2̺(dϑ2+sin2ϑdϑ2)and measure dΩ=sinh2̺d̺dΣ,dΣbeing the measure on S2.For convenience,here we normalize the curvature−κto−1.In these coordinates,the Laplace-Beltrami operator△reads△=∂2∂̺+1so,in order to derive it,it is suﬃcient to study radial wave functions of−△,that is solution of equationd2ud̺+λu=0(2) which reduces tod2vνsinh̺(4) Now,the L2(dΩ)scalar product for uν(̺)is(uν,uν′)=4πν2δ(ν−ν′)(5)from which the density of states̺(ν)=Vν2/2π2directly follows.As usual,we have introduced the large,ﬁnite volume V to avoid divergences.When possible,the limit V→∞shall be understood.At this point the computation of zeta function is straightforward.As we shall see in the following,what we are really interested in,is the zeta function related to the operators Q±=L1/2m±µ.The eigenvalues of L m areω2(ν)=ν2+a2=ν2+κ+m2,then we getζ(s;Q±)=V(4π)3/2Γ(s−1/2)(2a)s−3F(s+1,s−3;s−12a)(6)where F(α,β;γ;z)is the hypergeometric function.For its properties and its integral representations see for example ref.[20].It has to be noted that eq.(6)is the very same one has on aﬂat space for a massiveﬁeld with mass equal to a.Here in fact, the curvature plays the role of an eﬀective mass.As we see from eq.(6),the zeta function related to the pseudo-diﬀerential oper-ators Q±has simple poles at the points s n=3,2,1,−1,−2,−3,...with residues b n(±µ)=Res(ζ(s;Q±),s n)given byb3(±µ)=Vπ2;b1(±µ)=V2a);(8)c−n=(−1)n nV(2a)n+3dz F(α,β;γ;z)=αβgdV(11)where Lµ=−(∂τ−µ)2+L m and the Wick rotationτ=ix0has to be understood. In eq.(11)the integration has to be taken over allﬁeldsφ(τ,x a)withβ-periodicity with respect toτ.The eigenvalues of the whole operator Lµ,sayµn,νreadµn,ν= 2πnlog det(ℓ−2Lµ)=−1β[logℓ2ζ(0;Lµ)+ζ′(0;Lµ)](13)βℓbeing an arbitrary normalization parameter coming from the scalar path-integral measure.Note thatℓ,which has the dimensions of a mass,is necessary in order to keep the zeta-function dimensionless for all s.Theﬁnite temperature andµdependent part of the thermodynamic potential does not suﬀer of the presence of such an arbitrary parameter.On the contrary,ℓenters in the regularized expression of vacuum energy and this creates an ambiguity[12],which is proportional to the heat kernel expansion coeﬃcient K N(L m)related to L m(in general,K N(L m)=0). When the theory has a natural scale parameter,like the mass of the particle or the constant curvature of the manifold,the ambiguity can be removed by an”ad hoc”choice ofℓ[23].Here we would like to study the behaviour of thermodynamic quantities,then we are only interested in theµand T dependent part of the thermodynamic potential; that is a well deﬁned quantity,which does not need regularization.To compute it, it is not necessary to use all the analytic properties of zeta function(for a careful derivation of vacuum energy see for example refs.[12,14]).Then we can proceed in a formal way and directly compute log det Lµdisregarding the vacuum energy divergent term.First of all we observe that∞ n=−∞log(ω2+(2πn/β+iµ)2)=∞ n=−∞ dω24ω cothβ2(ω−µ) dω2(14)=−log 1−e−β(ω+µ) −log 1−e−β(ω−µ) −βωUsing eq.(13),recalling thatω2=ν2+a2and by integrating overνwith the state density that we have derived in the previous section,we get the standard result1Ω(β,µ)=−2π2β log 1−e−β(ω(ν)+µ) +log 1−e−β(ω(ν)−µ) ν2dν(15)V+∞ n=1cosh nβµK2(anβ)π2∞ n=0µ2n2;L m)β−s ds(17)πi1E(β,µ)=−where K2is the modiﬁed Bessel function,c is a suﬃciently large real number and ζR(s)is the usual Riemann zeta-function.The integral representations(17)and (18),which are valid for|µ|<a,are useful for high temperature expansion.On the contrary,the representation(16)in terms of modiﬁed Bessel functions is more useful for the low temperature expansion,since the asymptotics of Kνis well known.4Bose-Einstein condensationIn order to discuss Bose-Einstein condensation we have to analyze the behaviour of the charge density∂Ω(V,β,z)ρ=zV(expβωj−z)(20) and the activity z=expβµhas been introduced.Theωj in the sum are meant to be the Dirichlet eigenvalues for any normal domain V⊂H3.That is,V is a smooth connected submanifold of H3with non empty piecewise C∞boundary.By the inﬁnite volume limit we shall mean that a nested sequence of normal domains V k has been choosen together with Dirichlet boundary conditions and such that V k≡H3.The reason for this choice is the following theorem due to Mac Kean (see for example[24]):—ifωok denotes the smallest Dirichlet eigenvalue for any sequence of normal domains V kﬁlling all of H3thenωok≥a and lim k→∞ωok=a.(Although the above inequality is also true for Neumann boundary conditions,the existence of the limit in not assured to the authors knowledge).Now we can show the convergence of theﬁnite volume activity z k to a limit point ¯z as k→∞.Toﬁx ideas,let us supposeρ≥0:then z k∈(1,expβωok).Since ρ(V,β,z)is an increasing function of z such thatρ(V,β,1)=0andρ(V,β,∞)=∞, for eachﬁxed V k there is a unique z k(¯ρ,β)∈(1,βexpωok)such that¯ρ=ρ(V k,β,z k).By compactness,the sequence z k must have at least one ﬁxed point ¯z and as ωok →a 2as k goes to inﬁnity,by Mc Kean theorem,¯z ∈[1,exp βa ].From this point on,the mathematical analysis of the inﬁnite volume limit exactly parallels the one in ﬂat space for non relativistic systems,as it is done in various references [25,26,4].Inparticular,there is a critical temperature T c over which there are no particles in the ground state.T c is the unique solution of the equation̺=sinh βa cosh βaV ∂µ= ∞0 1e β(ω(ν)−µ)−1 ν22π2 ∞0ν2dνκ+m 2.Thevery diﬀerence between ﬂat and hyperbolic spaces occurs for massless particles.We shall return on this important point in a moment.Solutions of eq.(21)can be easily obtained in the two cases βa ≫1and βa ≪1(in the case of massive bosons these correspond to non relativistic and ultrarela-tivistic limits respectively).We have in fact̺≃T 3e x 2/2a −1= aT 2π2 ∞0x 2dx3;βa ≪1(24)from which we get the corresponding critical temperaturesT c =2πζR (3/2)2/3;βa ≫1(25)T c = 3̺∂µ̺(T,µ)(29)and since ∂µ̺diverges for µ=a we obtain µ′(T +c )=0.This is not the case of µ′′.In fact we shall see that µ′′(T +c )is diﬀerent from zero and therefore µ′′(T )isadiscontinuous function of temperature.This implay that theﬁrst derivative of the speciﬁc heat C V has a jump for T=T c given bydC VdT T−c=µ′′(T+c)∂U(T,µ)∂TT=T+c(30)U(T,µ)being the internal energy,which can be derived by means of equation U(β,µ)=−µ̺V+∂πΓ(k+1)Γ(−k+5/2)(2s)−k(32) Then,for small T we haveE(β,µ)≃−a4Vanβ5/2∞ k=0Γ(k+5/2) 2π 3/2∞ n=1e−nβ(a−|µ|)T A2;T∂̺2̺(35)where A=2.363and C=−2.612are two coeﬃcients of the expansion ∞ n=1e−nxNow,using eq.(30),we have the standard resultdC VdT T−c=3̺C2T c(37) The high temperature expansion could be obtained by using eq.(17),like in ref.[15].Here we shall use eq.(18),because for the aim of the present paper it is more ing the properties ofζ(s;Q±),which we have discussed in section2,we see that the integrand function in eq.(18)ζR(s+1)Γ(s)[ζ(s;Q+)+ζ(s;Q−)]β−(s+1)(38) has simple poles at s=3,1,0,−3,−5,−7,...and a double pole at s=−1.In-tegrating this function on a closed path containing all the poles,we get the high temperature expansion,valid for T>T c(hereγis the Euler-Mascheroni constant)E(β,µ)≃−V π212β2(a2−2µ2)+(a2−µ2)3/224π2(3a2−µ2)+a44π+γ−3(−2π)n(2n+1)where we have used the formulaζ′R(−2n)=Γ(2n+1)ζR(2n+1)3+µT(a2−µ2)1/212π2(41)−2(−2π)nForµ=a,the leading term of this expression gives again the result(24).From eq.(41),by a strightforward computation and taking only the leading terms into account,one getsµ′′(T+c)≃−12π2andﬁnally,from eqs.(30)and(31) dC VdT T−c≃−32̺π2pands adiabatically,can represent a manifold of the form we have considered.The problem we have studied then canﬁnd physical applications in the standard model of the universe.References[1]A.Einstein.Berl.Ber.,22,261,(1924).[2]K.Huang.Statistical Mechanics.J.Wiley and Sons,Inc.,New York,(1963).[3]H.Araki and E.J.Woods.J.Math.Phys.,4,637,(1963).[4]ndau and m.Math.Phys.,70,43,(1979).[5]H.E.Haber and H.A.Weldom.Phys.Rev.Lett.,46,1497,(1981).[6]H.E.Haber and H.A.Weldom.J.Math.Phys.,23,1852,(1982).[7]H.E.Haber and H.A.Weldom.Phys.Rev.D,25,502,(1982).[8]A.A.Bytsenko and Yu.P.Goncharov.Mod.Phys.Lett.A,6,669,(1991).[9]Yu.P.Goncharov and A.A.Bytsenko.Class.Quantum Grav.,8,L211,(1991).[10]A.A.Bytsenko and Y.P.Goncharov.Class.Quantum Grav.,8,2269,(1991).[11]A.A.Bytsenko and S.Zerbini.Class.Quant.Gravity,9,1365,(1992).[12]G.Cognola,L.Vanzo and S.Zerbini.J.Math.Phys.,33,222,(1992).[13]A.A.Bytsenko,L.Vanzo and S.Zerbini.Mod.Phys.Lett.A,7,397,(1992).[14]A.A.Bytsenko,G.Cognola and L.Vanzo.Vacuum energy for3+1dimesionalspace-time with compact hyperbolic spatial part.Technical Report,Universit´a di Trento,UTF255,(1992).to appear in JMP.[15]G.Cognola and L.Vanzo.Thermodynamic potential for scalarﬁelds in space-time with hyperbolic spatial part.Technical Report,Universit´a di Trento,UTF 258,(199).to be published.[16]A,A.Bytsenko,L.Vanzo and S.Zerbini.Zeta-function regularization approachtoﬁnite temperature eﬀects in Kaluza-Klein space-time.Technical Report,Uni-versit´a di Trento,UTF259,(1992).to be published.[17]C.G.Callan and F.Wilczek.Nucl.Phys.B,340,366,(1990).[18]N.Balasz and C.Voros.Phys.Rep.,143,109,(1986).[19]R.Camporesi and A.Higuchi.Phys.Rev.D,45,3591,(1992).[20]I.S.Gradshteyn and I.M.Ryzhik.Table of integrals,series and products.Aca-demic press,Inc.,New York,(1980).[21]A.Actor.Phys.Lett.B,157,53,(1985).[22]m.Math.Phys.,55,133,(1977).[23]J.S.Dowker and J.P.Schoﬁeld.Nucl.Phys.B,327,267,(1989).[24]I.Chavel.Eingenvalues in Riemannian Geometry.Accademic Press,(1984).[25]R.Ziﬀ,G.E.Uhlenbeck and M.Kac.Phys.Rep.,32,169,(1977).[26]J.T.Lewis and J.V.Pul´m.Math.Phys.,36,1,(1974).[27]J.M.Blatt and S.T.Butler.Phys.Rev.Lett.,100,476,(1955).。

Bose-Einstein condensation in a one-dimensional interacting system due to power-law trappin

a r X i v :c o n d -m a t /9810197v 1 [c o n d -m a t .s t a t -m e c h ] 16 O c t 1998Accepted to PHYSICAL REVIEW A for publicationBose-Einstein condensation in a one-dimensional interacting system due to power-lawtrapping potentialsM.Bayindir,B.Tanatar,and Z.GedikDepartment of Physics,Bilkent University,Bilkent,06533Ankara,TurkeyWe examine the possibility of Bose-Einstein condensation in one-dimensional interacting Bose gas subjected to conﬁning potentials of the form V ext (x )=V 0(|x |/a )γ,in which γ<2,by solving the Gross-Pitaevskii equation within the semi-classical two-ﬂuid model.The condensate fraction,chemical potential,ground state energy,and speciﬁc heat of the system are calculated for various values of interaction strengths.Our results show that a signiﬁcant fraction of the particles is in the lowest energy state for ﬁnite number of particles at low temperature indicating a phase transition for weakly interacting systems.PACS numbers:03.75.Fi,05.30.Jp,67.40.Kh,64.60.-i,32.80.PjI.INTRODUCTIONThe recent observations of Bose-Einstein condensation (BEC)in trapped atomic gases [1–5]have renewed inter-est in bosonic systems [6,7].BEC is characterized by a macroscopic occupation of the ground state for T <T 0,where T 0depends on the system parameters.The success of experimental manipulation of externally applied trap potentials bring about the possibility of examining two or even one-dimensional Bose-Einstein condensates.Since the transition temperature T 0increases with decreasing system dimension,it was suggested that BEC may be achieved more favorably in low-dimensional systems [8].The possibility of BEC in one -(1D)and two-dimensional (2D)homogeneous Bose gases is ruled out by the Hohen-berg theorem [9].However,due to spatially varying po-tentials which break the translational invariance,BEC can occur in low-dimensional inhomogeneous systems.The existence of BEC is shown in a 1D noninteracting Bose gas in the presence of a gravitational ﬁeld [10],an attractive-δimpurity [11],and power-law trapping po-tentials [12].Recently,many authors have discussed the possibility of BEC in 1D trapped Bose gases relevant to the magnetically trapped ultracold alkali-metal atoms [13–18].Pearson and his co-workers [19]studied the in-teracting Bose gas in 1D power-law potentials employing the path-integral Monte Carlo (PIMC)method.They have found that a macroscopically large number of atoms occupy the lowest single-particle state in a ﬁnite system of hard-core bosons at some critical temperature.It is important to note that the recent BEC experiments are carried out with ﬁnite number of atoms (ranging from several thousands to several millions),therefore the ther-modynamic limit argument in some theoretical studies [15]does not apply here [8].The aim of this paper is to study the two-body interac-tion eﬀects on the BEC in 1D systems under power-law trap potentials.For ideal bosons in harmonic oscillator traps transition to a condensed state is prohibited.It is anticipated that the external potentials more conﬁn-ing than the harmonic oscillator type would be possible experimentally.It was also argued [15]that in the ther-modynamic limit there can be no BEC phase transition for nonideal bosons in 1D.Since the realistic systems are weakly interacting and contain ﬁnite number of particles,we employ the mean-ﬁeld theory [20,21]as applied to a two-ﬂuid model.Such an approach has been shown to capture the essential physics in 3D systems [21].The 2D version [22]is also in qualitative agreement with the results of PIMC simulations on hard-core bosons [23].In the remaining sections we outline the two-ﬂuid model and present our results for an interacting 1D Bose gas in power-law potentials.II.THEORYIn this paper we shall investigate the Bose-Einstein condensation phenomenon for 1D interacting Bose gas conﬁned in a power-law potential:V ext (x )=V 0|x |κF (γ)G (γ)2γ/(2+γ),(2)andN 0/N =1−TF (γ)=1x 1/γ−1dx1−x,(4)and G (γ)=∞x 1/γ−1/2dxNk B T 0=Γ(1/γ+3/2)ζ(1/γ+3/2)T 01/γ+3/2.(6)Figure 1shows the variation of the critical temperature T 0as a function of the exponent γin the trapping po-tential.It should be noted that T 0vanishes for harmonic potential due to the divergence of the function G (γ=2).It appears that the maximum T 0is attained for γ≈0.5,and for a constant trap potential (i.e.V ext (x )=V 0)the BEC disappears consistent with the Hohenberg theorem.0.00.5 1.0 1.5 2.0γ0.00.20.40.6k B T 0 (A r . U n .)FIG.1.The variation of the critical temperature T 0withthe external potential exponent γ.We are interested in how the short-range interactioneﬀects modify the picture presented above.To this end,we employ the mean-ﬁeld formalism and describe the col-lective dynamics of a Bose condensate by its macroscopictime-dependent wave function Υ(x,t )=Ψ(x )exp (−iµt ),where µis the chemical potential.The condensate wavefunction Ψ(x )satisﬁes the Gross-Pitaevskii (GP)equa-tion [24,25]−¯h 2dx 2+V ext (x )+2gn 1(x )+g Ψ2(x )Ψ(x )=µΨ(x ),(7)where g is the repulsive,short-range interaction strength,and n 1(x )is the average noncondensed particle distribu-tion function.We treat the interaction strength g as a phenomenological parameter without going into the de-tails of actually relating it to any microscopic descrip-tion [26].In the semi-classical two-ﬂuid model [27,28]the noncondensed particles can be treated as bosons in an eﬀective potential [21,29]V eﬀ(x )=V ext (x )+2gn 1(x )+2g Ψ2(x ).(8)The density distribution function is given byn 1(x )=dpexp {[p 2/2m +V eﬀ(x )−µ]/k B T }−1,(9)and the total number of particles N ﬁxes the chemical potential through the relationN =N 0+ρ(E )dE2mgθ[µ−V ext (x )−2gn 1(x )],(12)where θ[x ]is the unit step function.More precisely,the Thomas-Fermi approximation [7,20,30]would be valid when the interaction energy ∼gN 0/Λ,far exceeds the kinetic energy ¯h 2/2m Λ2,where Λis the spatial extent of the condensate cloud.For a linear trap potential (i.e.γ=1),a variational estimate for Λis given by Λ= ¯h 2/2m (π/2)1/22a/V 0 1/3.We note that the Thomas-Fermi approximation would breakdown for tem-peratures close to T 0where N 0is expected to become very small.The above set of equations [Eqs.(9)-(12)]need to be solved self-consistently to obtain the various physical quantities such as the chemical potential µ(N,T ),the condensate fraction N 0/N ,and the eﬀective potential V eﬀ.In a 3D system,Minguzzi et al .[21]solved a simi-lar system of equations numerically and also introduced an approximate semi-analytical solution by treating the interaction eﬀects perturbatively.Motivated by the suc-cess [21,22]of the perturbative approach we consider aweakly interacting system in1D.To zero-order in gn1(r), the eﬀective potential becomesV eﬀ(x)= V ext(x)ifµ<V ext(x)2µ−V ext(x)ifµ>V ext(x).(13) Figure2displays the typical form of the eﬀective po-tential within our semi-analytic approximation scheme. The most noteworthy aspect is that the eﬀective poten-tial as seen by the bosons acquire a double-well shape because of the interactions.We can explain this result by a simple argument.Let the number of particles in the left and right wells be N L and N R,respectively,so that N=N L+N R.The nonlinear or interaction term in the GP equation may be approximately regarded as V=N2L+N2R.Therefore,the problem reduces to the minimization of the interaction potential V,which is achieved for N L=N R.FIG.2.Eﬀective potential V eﬀ(x)in the presence of in-teraction(x0=(µ/V0)1/γa).Thick dotted line represents external potential V ext(x).The number of condensed atoms is calculated to beN0=2γa√ze x−1+ 2µ/k B Tµ/k B TH(γ,µ,xk B T)(2µ/k B T−x)1/γ−1/2dxexp[(E−µ)/k B T]−1=κ(k B T)1/γ+1/2J(γ,µ,T),(18) whereJ(γ,µ,T)= ∞2µ/k B T x1/γ+1/2dxze x−1.and Ecis the energy of the particles in the condensateE c=g(1+γ)(2γ+1)gV1/γ.(19)The kinetic energy of the condensed particles is neglected within our Thomas-Fermi approximation to the GP equa-tion.III.RESULTS AND DISCUSSIONUp to now we have based our formulation for arbitrary γ,but in the rest of this work we shall present our re-sults forγ=1.Our calculations show that the results for other values ofγare qualitatively similar.In Figs. 3and4we calculate the condensate fraction as a func-tion of temperature for various values of the interaction strengthη=g/V0a(at constant N=105)and diﬀerent number of particles(at constantη=0.001),respectively. We observe that as the interaction strengthηis increased, the depletion of the condensate becomes more apprecia-ble(Fig.3).As shown in the correspondingﬁgures,a signiﬁcant fraction of the particles occupies the ground state of the system for T<T0.The temperature depen-dence of the chemical potential is plotted in Figs.5and 6for various interaction strengths(constant N=105) and diﬀerent number of particles(constantη=0.001) respectively.0.00.20.40.60.8 1.0T/T 00.00.20.40.60.81.0N 0/NN 0/N=1−(T/T 0)3/2η=10−5η=10−3η=10−1η=10FIG.3.The condensate fraction N 0/N versus temperature T /T 0for N =105and for various interaction strengths η.Eﬀects of interactions on µ(N,T )are seen as large de-viations from the noninteracting behavior for T <T 0.In Fig.7we show the ground state energy of an interacting 1D system of bosons as a function of temperature for dif-ferent interaction strengths.For small η,and T <T 0, E is similar to that in a noninteracting system.As ηincreases,some diﬀerences start to become noticeable,and for η≈1we observe a small bump developing in E .This may indicate the breakdown of our approxi-mate scheme for large enough interaction strengths,as we can ﬁnd no fundamental reason for such behavior.It is also possible that the Thomas-Fermi approximation em-ployed is violated as the transition to a condensed state is approached.0.00.20.40.60.8 1.0T/T 00.00.20.40.60.81.0N 0/NN 0/N=1−(T/T 0)3/2N=108N=105N=103N=101FIG.4.The condensed fraction N 0/N versus temperature T /T 0for η=0.001and for diﬀerent number of particles N .0.00.20.40.60.8 1.0 1.2T/T 0−100100200300400µ/V 0η=1η=0.1η=0.001η=0.00001FIG.5.The temperature dependence of the chemical potential µ(N,T )for various interaction strength and for N =105particles.Although it is conceivable to imagine the full solution of the mean-ﬁeld equations [Eq.(9)-(12)]may remedy the situation for larger values of η,the PIMC simulations [19]also seem to indicate that the condensation is inhibited for strongly interacting systems.The results for the spe-ciﬁc heat calculated from the total energy curves,i.e.C V =d E /dT ,are depicted in Fig.8.The sharp peak at T =T 0tends to be smoothed out with increasing in-teraction strength.It is known that the eﬀects of ﬁnite number of particles are also responsible for such a be-havior [20].In our treatment these two eﬀects are not disentangled.It was pointed out by Ingold and Lam-brecht [14]that the identiﬁcation of the BEC should also be based on the behavior of C V around T ≈T 0.0.00.20.40.60.8 1.0 1.2T/T 0−5050100µ/V 0N=107N=105N=103N=101FIG.6.The temperature dependence of the chemical po-tential µ(N,T )for diﬀerent number of particles N and for η=0.001.0.00.20.40.60.8 1.0 1.2T/T 00.00.20.40.60.8<E >/N k B T 0η=0η=0.001η=0.1η=1Maxwell−BoltzmannFIG.7.The temperature dependence of the total energy of 1D Bose gas for various interaction strengths ηand N =105particles.Our calculations indicate that the peak structure of C V remains even in the presence of weak interactions,thus we are led to conclude that a true transition to a Bose-Einstein condensed state is predicted within the present approach.0.00.20.40.60.81.01.2T/T 00.00.20.40.60.81.0C V /N k Bη=0η=0.001η=0.1Maxwell−BoltzmannFIG.8.The temperature dependence of the speciﬁc heat C V for various interaction strengths ηand N =105particles.IV.CONCLUDING REMARKSIn this work we have applied the mean-ﬁeld,semi-classical two-ﬂuid model to interacting bosons in 1D power-law trap potentials.We have found that for a range of interaction strengths the behavior of the thermo-dynamic quantities resembles to that of non-interactingbosons.Thus,BEC in the sense of macroscopic occu-pation of the ground state,occurs when the short-range interparticle interactions are not too strong.Our results are in qualitative agreement with the recent PIMC sim-ulations [19]of similar systems.Both 2D and 1D sim-ulation results [19,23]indicate a phase transition for a ﬁnite number system,in contrast to the situation in the thermodynamic limit.Since systems of much larger size can be studied within the present approach,our work complements the PIMC calculations.The possibility of studying the tunneling phenomenon of condensed bosons in spatially diﬀerent regions sepa-rated by a barrier has recently attracted some attention [31–34].In particular,Dalfovo et al .[32]have shown that a Josephson-type tunneling current may exist for bosons under the inﬂuence of a double-well trap potential.Za-pata et al .[34]have estimated the Josephson coupling energy in terms of the condensate density.It is inter-esting to speculate on such a possibility in the present case,since the eﬀective potential in our description is of the form of a double-well potential (cf.Fig.2).In our treatment,the interaction eﬀects modify the single-well trap potential into one which exhibits two minima.Thus if we think of this eﬀective potential as the one seen by the condensed bosons and according to the general ar-guments [31–34]based on two weakly connected systems we should have an oscillating ﬂux of particles when the chemical potential in the two wells is diﬀerent.Any con-ﬁguration with N L =N R which is always the case for odd number of bosons will result in an oscillatory mo-tion.It would be interesting to explore these ideas in future work.ACKNOWLEDGMENTSThis work was supported by the Scientiﬁc and Techni-cal Research Council of Turkey (TUBITAK)under Grant No.TBAG-1736and TBAG-1662.We gratefully ac-knowledge useful discussions with Prof.C.Yalabık and E.Demirel.[5]D.J.Han,R.H.Wynar,Ph.Courteille,and D.J.Heinzen,Phys.Rev.A57,R4114(1998).[6]I.F.Silvera,in Bose-Einstein Condensation,Ed.by A.Griﬃn,D.W.Snoke,and S.Stringari(Cambridge Uni-versity Press,Cambridge,1995).[7]F.Dalfovo,S.Giorgini,L.P.Pitaevskii,and S.Stringari,preprint,cond-mat/9806038(to be published in Reviews of Modern Physics);A.S.Parkins and D.F.Walls,Phys.Rep.303,1(1998).[8]W.Ketterle and N.J.van Druten,Phys.Rev.A54,656(1996).[9]P.C.Hohenberg,Phys.Rev.158,383(1967).[10]A.Widom,Phys.Rev.176,254(1968).[11]L.C.Ioriatti,Jr.,S.G.Rosa,Jr.,and O.Hipolito,Am.J.Phys.44,744(1976).[12]V.Bagnato and D.Kleppner,Phys.Rev.A44,7439(1991).[13]T.Haugset and H.Haugerud,Phys.Rev.A57,3809(1998).[14]G.-L.Ingold and mbrecht,Eur.Phys.J.D1,29(1998).[15]W.J.Mullin,J.Low Temp.Phys.110,167(1998)[16]W.Deng and P.M.Hui,Solid State Commun.104,729(1997).[17]H.Monien,M.Linn,and N.Elstner,preprint,cond-mat/9711178.[18]F.Brosens,J.T.Devreese,and L.F.Lemmens,SolidState Commun.100,123(1996).[19]S.Pearson,T.Pang,and C.Chen,Phys.Rev.A58,1485(1998).[20]S.Giorgini,L.P.Pitaevskii,and S.Stringari,Phys.Rev.A,54,4633(1996);J.Low Temp.Phys.109,309(1997).[21]A.Minguzzi,S.Conti,and M.P.Tosi,J.Phys.Cond.Matter9,L33(1997).[22]M.Bayindir and B.Tanatar,Phys.Rev.A58,3134(1998).[23]S.Heinrichs and W.J.Mullin,preprint,cond-mat/9807331.[24]E.P.Gross,Nuovo Cimento20,454(1961).[25]L.P.Pitaevskii,Zh.Eksp.Teor.Fiz.40,646(1961)[Sov.Phys.JETP13,451(1961)].[26]M.Olshanii,Phys.Rev.Lett.81,938(1998).[27]R.J.Dodd,K.Burnett,M.Edwards,and C.W.Clark,Acta Phys.Pol.A93,45(1998).[28]M.Naraschewski and D.M.Stamper-Kurn,Phys.Rev.A58,2423(1998).[29]V.Bagnato,D.E.Pritchard,and D.Kleppner,Phys.Rev.A35,4354(1987).[30]G.Baym and C.J.Pethick,Phys.Rev.Lett.76,6(1996).[31]S.Grossmann and M.Holthaus,Z.Naturforsch.50a,323(1995).[32]F.Dalfovo,L.P.Pitaevskii,and S.Stringari,Phys.Rev.A54,4213(1996).[33]A.Smerzi,S.Fantoni,S.Giovanazzi,and S.R.Shenoy,Phys.Rev.Lett.79,4950(1997).[34]I.Zapata,F.Sols,and A.J.Leggett,Phys.Rev.A57,28(1998).。

Bethe Ansatz Solutions and Excitation Gap of the Attractive Bose-Hubbard Model

a r X i v :c o n d -m a t /0108314v 1 [c o n d -m a t .s t a t -m e c h ] 20 A u g 2001Bethe Ansatz Solutions and Excitation Gap of the Attractive Bose-Hubbard ModelDeok-Sun Lee and Doochul KimSchool of Physics,Seoul National University,Seoul 151-747,KoreaThe energy gap between the ground state and the ﬁrst excited state of the one-dimensional attractive Bose-Hubbard Hamiltonian is investigated in connection with directed polymers in random media.The excitation gap ∆is obtained by exact diagonalization of the Hamiltonian in the two-and three-particle sectors and also by an exact Bethe Ansatz solution in the two-particle sector.The dynamic exponent z is found to be 2.However,in the intermediate range of the size L where UL ∼O (1),U being the attractive interaction,the eﬀective dynamic exponent shows an anomalous peak reaching high values of 2.4and 2.7for the two-and the three-particle sectors,respectively.The anomalous behavior is related to a change in the sign of the ﬁrst excited-state energy.In the two-particle sector,we use the Bethe Ansatz solution to obtain the eﬀective dynamic exponent as a function of the scaling variable UL/π.The continuum version,the attractive delta-function Bose-gas Hamiltonian,is integrable by the Bethe Ansatz with suitable quantum numbers,the distributions of which are not known in general.Quantum numbers are proposed for the ﬁrst excited state and are conﬁrmed numerically for an arbitrary number of particles.I.INTRODUCTIONThe dynamics of many simple non-equilibrium sys-tems are often studied through corresponding quantum Hamiltonians.Examples are the asymmetric XXZ chain Hamiltonian and the attractive Bose-Hubbard Hamilto-nian for the single-step growth model [1]and the directed polymers in random media (DPRM)[2],respectively.The single-step growth model is a Kardar-Parisi-Zhang (KPZ)universality class growth model where the inter-face height h (x,t )grows in a stochastic manner under the condition that h (x ±1,t )−h (x,t )=±1.The process is also called the asymmetric exclusion process (ASEP)in a diﬀerent context.The evolution of the probability distri-bution for h (x,t )is generated by the asymmetric XXZ chain Hamiltonian [3].The entire information about the dynamics is coded in the generating function e αh (x,t ) .Its time evolution,in turn,is given by the modiﬁed asym-metric XXZ chain Hamiltonian [4–6],H XXZ (α)=−L i =1e 2α/L σ−i σ+i +1+12L i =1(b i b †i +1+b †i b i +1−2)−UL i =1b †ib i (b †ib i −1)4Lρ(1−ρ)and −n√4Lρ(1−ρ)≫1and the density of particles is ﬁ-nite in the limit L →∞,∆(α)behaves as ∆(α)∼L −1.However,when α∆(α)∼L−3/2[3,11].The dynamic exponent z=3/2is a characteristic of the dynamic universality class of the KPZ-type surface growth.When the number of par-ticles isﬁnite and the density of particles is very low, it is known that∆(α)∼L−2[12].However,whenα<0,which corresponds to the ferromagnetic phase, most Bethe Ansatz solutions are not available althoughthe Bethe Ansatz equations continue to hold.Asαbe-comes negative,the quasi-particle momenta appearing inthe Bethe Ansatz equations become complex,so solutions are diﬃcult to obtain analytically.The attractive Bose-Hubbard Hamiltonian is expected to have some resemblance to the ferromagnetic phaseof the asymmetric XXZ chain Hamiltonian consider-ing the equivalence ofαand−n.The equivalence isidentiﬁed indirectly by comparing the two scaling vari-ablesα LU under the relation U= 4ρ(1−ρ)or the two generating functions exp(αh(x,t) and Z(x,t)n under the relation Z(x,t)=e−h(x,t).In contrast to the asymmetric XXZ chain Hamiltonian,theBose-Hubbard Hamiltonian does not satisfy the Bethe Ansatz except in the two-particle sector[13].Instead, the attractive delta-function Bose-gas Hamiltonian,H D(n)=−1∂x2i−Ui<jδ(x i−x j),(4)which is the continuum version of the attractive Bose-Hubbard Hamiltonian,is known to be integrable by the Bethe Ansatz.The attractive delta-function Bose gas has been studied in Refs.[14]and[15].The ground-state energy is obtained from the Bethe Ansatz solution by us-ing the symmetric distribution of the purely imaginary quasi-particle momenta.However,the structure of the energy spectra is not well known for the same reason as in the asymmetric XXZ chain Hamiltonian withα<0. The unknown energy spectra itself prevents one from un-derstanding the dynamics of DPRM near the stationary state.In this paper,we discuss in Section II the distribu-tion of the quantum numbers appearing in the Bethe Ansatz equation for theﬁrst excited state of the attrac-tive delta-function Bose-gas Hamiltonian,the knowledge of which is essential for solving the Bethe Ansatz equa-tion.In Section III,the excitation gap of the attractive Bose-Hubbard Hamiltonian with a small number of par-ticles is investigated through the exact diagonalization method.We show that the gap decays as∆∼L−2,i.e., z=2,but that the exponent becomes anomalous when U∼L−1.The emergence of the anomalous exponent is explained in connection with the transition of theﬁrst excited state from a positive energy state to a negative energy state.The Bethe Ansatz solutions in the two-particle sector show how the behavior of the gap varies with the interaction.We give a summary and discussion in Section IV.II.QUANTUM NUMBER DISTRIBUTION FOR THE FIRST EXCITED STATEIn this section,we study the Bethe Ansatz solutions for the ground state and theﬁrst excited state of the attrac-tive delta-function Bose-gas Hamiltonian.The eigenstate of H D(n),Eq.(4),is of the formφ(x1,x2,...,x n)= P A(P)exp(ik P1x1+ik P2x2+···+ik P n x n),(5)where P is a permutation of1,2,...,n and x1≤x2≤...≤x n with no three x’s being equal.The quasi-particle momenta k j’s are determined by solving the Bethe Ansatz equations,k j L=2πI j+ l=jθ(k j−k l2+j,(j=1,2,...,n),(7)and the quasi-particle momenta are distributed symmet-rically on the imaginary axis in the complex-k plane. Care should be taken when dealing with theﬁrst excited state.For the repulsive delta-function Bose-gas Hamilto-nian,where U is replaced by−U in Eq.(4),the quantum numbers for one of theﬁrst excited states areI j=−n+12.(8)However,for the attractive case,by following the move-ment of the momenta as U changes sign,weﬁnd that the quantum numbers for theﬁrst excited state should be given byI j=−n−12(=I1).(9)That is,the two quantum numbers I1and I n become the same.Such a peculiar distribution of I j’s does not ap-pear in other Bethe Ansatz solutions such as those for the XXZ chain Hamiltonian or the repulsive delta-function Bose-gas Hamiltonian.We remark that even though the two I j’s are the same,all k j’s are distinct;otherwise,the wavefunction vanishes.Such a distribution of quan-tum numbers is conﬁrmed by the consistency between the energies obtained by diagonalizing the Bose-HubbardHamiltonian exactly and those obtained by solving the Bethe Ansatz equations with the above quantum num-bers for very weak interactions,for which the two Hamil-tonians possess almost the same energy spectra.When there is no interaction(U=0),all quasi-particlemomenta,k j’s,are zero for the ground state while for theﬁrst excited state,all the k j’s are zero except thelast one,k n=2π/L.In the complex-k plane,as the very weak repulsive interaction is turned on,the n−1momenta are shifted inﬁnitesimally from k=0withk1<k2<···<k n−1,and the n th momentum is shifted inﬁnitesimally to the left from k=2π/L.All the mo-menta remain on the real axis.When the interaction is weakly attractive,the n−1momenta become complexwith Im k1<Im k2<···<Im k n−1and Re k j≃0for j=1,2,...,n−1,and the n th momentum remains on thereal axis,but is shifted to the left.Figure1shows the dis-tribution of the quantum numbers and the quasi-particlemomenta in the presence of a very weak attractive in-teraction.The quasi-particle momenta are obtained by solving Eq.(6).Knowledge of the distribution of the quantum num-bers is essential for solving the Bethe Ansatz equations of the attractive delta-function Bose-gas Hamiltonian.For the original attractive Bose-Hubbard Hamiltonian,the Bethe Ansatz solutions are the exact solutions for the two-particle sector only,but are good approximate so-lutions in other sectors provided the density is very low and the interaction is very weak.This is because the Bethe Ansatz for the Bose-Hubbard Hamiltonian fails once states with sites occupied by more than three parti-cles are included.Thus,for the sectors with three or more particles,the Bethe Ansatz solutions may be regarded as approximate eigenstates provided states with more than three particles at a site do not play an important role in the eigenfunctions.In Ref.[13],it is shown that the error in the Bethe Ansatz due to multiply-occupied sites (occupied by more than three particles)is proportional to U2,where U(>0)in Ref.[13]corresponds to−U in Eq.(2).This applies to the attractive interaction case also.For the repulsive Bose-Hubbard Hamiltonian,the Bethe Ansatz is a good approximation when the density is low and the interaction is strong because the strong re-pulsion prevents many particles from occupying the same site[16].For the attractive Bose-Hubbard Hamiltonian, the Bethe Ansatz is good when the density is low and the interaction is weak because a weak attraction is better for preventing many particles from occupying the same site and because the error is proportional to U2.III.POWER-LA W DEPENDENCE ANDANOMALOUS EXPONENTWe are interested in the scaling limit L→∞with the scaling variable n√byE 0=−4sinh 2 κ2Lcosh2q 2L sinh2q2sinh κ−U,(12)andqL =logU +2cos(πU −2cos(π2s κ−s U,(14)which gives s κ≃1.151.When the size of the system L is increased by δL with U =U∗,the changes of κand q ,δκand δq ,are,from Eqs.(12)and (13),δκ=−πs κ(4s κ2−s U 2)L 2≡−πΓδL(4/π)s U −s U 2+4δLL 2.(15)The perturbative expansion ∆(L +δL )≃∆(L )(1−z (δL/L )),under the assumption that ∆(L )∼L −z ,gives the value of z eﬀat U ∗:z eﬀ=21+s κΓ+Σlog((L −1)/(L +1))(17)by using the solutions of Eqs.(12)and (13)for suﬃciently large L .As discussed above,the exponent z eﬀshows an anomalous peak near U =U ∗or UL/π=s U and ap-proaches 2.0as UL/π→0or ∞.Figure 6shows a plot of z eﬀversus the scaling variable UL/πat L =10000.IV.SUMMARY AND DISCUSSIONAs the asymmetric XXZ chain generates the dy-namics of the single-step growth model,the attractive Bose-Hubbard Hamiltonian governs the dynamics of the DPRM.We studied the attractive Bose-Hubbard Hamil-tonian and its continuum version,the attractive delta-function Bose-gas Hamiltonian concentrating on the be-havior of the excitation gap,which is related to the char-acteristics of DPRM relaxing into the stationary state.For the attractive delta-function Bose gas Hamiltonian,The quantum numbers for the ﬁrst excited state in the Bethe Ansatz equation are found for the attractive delta-function Bose gas Hamiltonian,and the distribution of the quasi-particle momenta is discussed in the presence of a very weak attractive interaction.Our result is the start-ing point for a further elucidation of the Bethe Ansatz solutions.We show that the excitation gap depends on the size of the system as a power law,∆∼L −z ,and that the exponent z can be calculated by using an exact diag-onalization of the attractive Bose-Hubbard Hamiltonian in the two-and the three-particle sectors and by using the Bethe Ansatz solution in the two-particle sector.The exponent z is 2.0.However,for the intermediate region where UL ∼O (1),the eﬀective exponent z eﬀshows a peak.The equivalence of the diﬀerential equations govern-ing the single-step growth model and DPRM implies some inherent equivalence in the corresponding Hamil-tonians.The power-law behavior of the excitation gap,∆∼L −2,for the attractive Bose-Hubbard Hamiltonian with a very weak interaction is the same as that for the asymmetric XXZ chain Hamiltonian with a small num-ber of particles,which is expected considering the rela-tion U =4ρ(1−ρ).The fact that the excitation gap behaves anomalously for U ∼L −1implies the possibility of an anomalous dynamic exponent z for a ﬁnite scaling variable n√[1]M.Plischke,Z.Racz,and D.Liu,Phys.Rev.B 35,3485(1987).[2]M.Kardar,Nucl.Phys.B 290[FS20],582(1987).[3]L.H.Gwa and H.Spohn,Phys.Rev.A 46,844(1992).[4]B.Derrida and J.L.Lebowitz,Phys.Rev.Lett.80,209(1998).[5]D.-S.Lee and D.Kim,Phys.Rev.E 59,6476(1999).[6]B.Derrida and C.Appert,J.Stat.Phys.94,1(1999).[7]J.Krug and H.Spohn,in Solids Far from Equilibrium ,edited by C.Godr´e che (Cambridge University Press,Cambridge,1991),p.412.[8]B.Derrida and K.Mallick,J.Phys.A 30,1031(1997).[9]S.-C.Park,J.-M.Park,and D.Kim,unpublished.[10]E.Brunet and B.Derrida,Phys.Rev.E 61,6789(2000).[11]D.Kim,Phys.Rev.E 52,3512(1995).[12]M.Henkel and G.Sch¨u tz,Physica A 206,187(1994).[13]T.C.Choy and F.D.M.Haldane,Phys.Lett.90A ,83(1982).[14]E.H.Lieb and W.Liniger,Phys.Rev.130,1605(1963).[15]J.G.Muga and R.F.Snider,Phys.Rev.A 57,3317(1998).[16]W.Krauth,Phys.Rev.B 44,9772(1991).(a)ω-0.10.1-0.50.5(b)Re k Im k 0FIG.1.For the ﬁrst excited state,(a)the quantum num-bers I j ’s are depicted in the complex-ωplane with ω=e 2πiI/L and (b)the quasi-particle momenta k j ’s are shown in the complex-k plane.Here,the size of the system L is 20,the number of particles n is 10,and the attractive interaction U is 0.0025.The ﬁlled circle in (a)is where the two quantum numbers overlap.0.5 1102030EL n =2 U =0.05ground state first excited state0.5 1102030EL n =2 U =0.5ground state first excited state-3.4-3.3-3.2 102030EL n =2 U =5ground state first excited state0.5 11020 30EL n =3 U =0.05ground state first excited state-0.5 0 0.5 1020 30ELn =3 U =0.5ground state first excited state-12.17-12.16-12.15 1020 30ELn =3 U =5ground state first excited stateFIG.2.Ground-state energies and ﬁrst excited-state ener-gies are plotted versus the size of the system L (4≤L ≤30)for U =0.05,0.5,and 5in the two-and the three-particle sectors.The dotted line represents E =0.For all values of U and L ,the ground-state energy is negative.On the other hand,when U =0.5,the excited-state energy becomes nega-tive near L ≃14in the two-particle sector and L ≃6in the three-particle sector.The signs of the excited-state energies for U =0.05and 5do not change in the range of L shown here.0.0010.010.1110102030∆LU=0.05U=0.5U=5FIG.3.Log-log plot of the excitation gaps (∆)versus the size of the system (L )in the two-particle sector.Data for U =0.05and 5approach straight lines with slope z =2.0,but those for U =0.5show a strong crossover before approach-ing the asymptotic behavior.The solid line for U =0.5is that ﬁtted in the range 14≤L ≤18,and shows an eﬀective z ≃2.4.0.00010.001 0.010.1110102030∆LU=0.05U=0.5U=5FIG.4.Same as in Fig.3,but for the three-particle sec-tor.The ﬁtted solid line used the data for 8≤L ≤12,and has a slope of approximately 2.7.(a)k (b)k FIG.5.Distributions of the quasi-particle momenta,k j ’s,for the ground state (ﬁlled circles)and the ﬁrst excited state (open circles)are shown in the complex-k plane for n =2.The size of the system L is 100and the interaction U is (a)0.001and (b)0.1.22.22.4s U510z e f fU L/πFIG.6.Eﬀective exponent z eﬀin the two-particle sector versus the scaling variable UL/πat L =10000.The interac-tion U varies from 0.0001to 0.001.At UL/π=s U ≃2.181,z eﬀ≃2.401.。

Comments on ``Observation of Strong Quantum Depletion in a Gaseous Bose-Einstein Condensate

Comments on “Observation of Strong Quantum Depletion in a Gaseous Bose-Einstein Condensate” [cond-mat/0601184]
Sang-Hoon Kim
Division of Liberal Arts, Mokpo National Maritime University, Mokpo 530-729, Korea (Dated: February 6, 2008)
int N0 (a, T ) ∼ N0 (T )e−βU(a,N0 ) .
g3 (z) +
∞ l=1 l n
2 3
kB T h ¯ω
2
g2 (z) + ..., (1)
(5)
where gn (z) = z /l is the Bose functions. z = eβ(µ−ε0 ) is the fugacity between 0 and 1, ε0 = (3/2)¯ ω h is the zero-point energy of the harmonic oscillator, and µ(T ) is the chemical potential. N0 is the ground state occupation number given by N0 = 1/(z −1 − 1). The others are common notations. For interacting system there is a shift of the zero point energy by the interaction. Then, the new ground state occupation number is written as 1 int , (2) N0 = β(ε +U−µ) e 0 −1 where U = U (a, N0 ) is the interaction energy between atoms and may obtained from the solution of the GPE. We consider the positive scattering length only, so does

19 Bose-Einstein Condensation

Chapter19Bose-Einstein CondensationAbstract Bose-Einstein condensation(BEC)refers to a prediction of quantum sta-tistical mechanics(Bose[1],Einstein[2])where an ideal gas of identical bosons undergoes a phase transition when the thermal de Broglie wavelength exceeds the mean spacing between the particles.Under these conditions,bosons are stimulated by the presence of other bosons in the lowest energy state to occupy that state as well,resulting in a macroscopic occupation of a single quantum state.The con-densate that forms constitutes a macroscopic quantum-mechanical object.BEC was ﬁrst observed in1995,seventy years after the initial predictions,and resulted in the award of2001Nobel Prize in Physics to Cornell,Ketterle and Weiman.The exper-imental observation of BEC was achieved in a dilute gas of alkali atoms in a mag-netic trap.Theﬁrst experiments used87Rb atoms[3],23Na[4],7Li[5],and H[6] more recently metastable He has been condensed[7].The list of BEC atoms now includes molecular systems such as Rb2[8],Li2[9]and Cs2[10].In order to cool the atoms to the required temperature(∼200nK)and densities(1013–1014cm−3) for the observation of BEC a combination of optical cooling and evaporative cooling were employed.Early experiments used magnetic traps but now optical dipole traps are also common.Condensates containing up to5×109atoms have been achieved for atoms with a positive scattering length(repulsive interaction),but small con-densates have also been achieved with only a few hundred atoms.In recent years Fermi degenerate gases have been produced[11],but we will not discuss these in this chapter.BECs are now routinely produced in dozens of laboratories around the world. They have provided a wonderful test bed for condensed matter physics with stunning experimental demonstrations of,among other things,interference between conden-sates,superﬂuidity and vortices.More recently they have been used to create opti-cally nonlinear media to demonstrate electromagnetically induced transparency and neutral atom arrays in an optical lattice via a Mott insulator transition.Many experiments on BECs are well described by a semiclassical theory dis-cussed below.Typically these involve condensates with a large number of atoms, and in some ways are analogous to describing a laser in terms of a semiclassi-cal meanﬁeld.More recent experiments however have begun to probe quantum39739819Bose-Einstein Condensation properties of the condensate,and are related to the fundamental discreteness of the ﬁeld and nonlinear quantum dynamics.In this chapter,we discuss some of these quantum properties of the condensate.We shall make use of“few mode”approxi-mations which treat only essential condensate modes and ignore all noncondensate modes.This enables us to use techniques developed for treating quantum optical systems described in earlier chapters of this book.19.1Hamiltonian:Binary Collision ModelThe effects of interparticle interactions are of fundamental importance in the study of dilute–gas Bose–Einstein condensates.Although the actual interaction potential between atoms is typically very complex,the regime of operation of current exper-iments is such that interactions can in fact be treated very accurately with a much–simpliﬁed model.In particular,at very low temperature the de Broglie wavelengths of the atoms are very large compared to the range of the interatomic potential.This, together with the fact that the density and energy of the atoms are so low that they rarely approach each other very closely,means that atom–atom interactions are ef-fectively weak and dominated by(elastic)s–wave scattering.It follows also that to a good approximation one need only consider binary collisions(i.e.,three–body processes can be neglected)in the theoretical model.The s–wave scattering is characterised by the s–wave scattering length,a,the sign of which depends sensitively on the precise details of the interatomic potential [a>0(a<0)for repulsive(attractive)interactions].Given the conditions described above,the interaction potential can be approximated byU(r−r )=U0δ(r−r ),(19.1) (i.e.,a hard sphere potential)with U0the interaction“strength,”given byU0=4π¯h2am,(19.2)and the Hamiltonian for the system of weakly interacting bosons in an external potential,V trap(r),can be written in the second quantised form asˆH=d3rˆΨ†(r)−¯h22m∇2+V trap(r)ˆΨ(r)+12d3rd3r ˆΨ†(r)ˆΨ†(r )U(r−r )ˆΨ(r )ˆΨ(r)(19.3)whereˆΨ(r)andˆΨ†(r)are the bosonﬁeld operators that annihilate or create a par-ticle at the position r,respectively.19.2Mean–Field Theory —Gross-Pitaevskii Equation 399To put a quantitative estimate on the applicability of the model,if ρis the density of bosons,then a necessary condition is that a 3ρ 1(for a >0).This condition is indeed satisﬁed in the alkali gas BEC experiments [3,4],where achieved densities of the order of 1012−1013cm −3correspond to a 3ρ 10−5−10−6.19.2Mean–Field Theory —Gross-Pitaevskii EquationThe Heisenberg equation of motion for ˆΨ(r )is derived as i¯h ∂ˆΨ(r ,t )∂t = −¯h 22m ∇2+V trap (r ) ˆΨ(r ,t )+U 0ˆΨ†(r ,t )ˆΨ(r ,t )ˆΨ(r ,t ),(19.4)which cannot in general be solved.In the mean–ﬁeld approach,however,the expec-tation value of (19.4)is taken and the ﬁeld operator decomposed asˆΨ(r ,t )=Ψ(r ,t )+˜Ψ(r ,t ),(19.5)where Ψ(r ,t )= ˆΨ(r ,t ) is the “condensate wave function”and ˜Ψ(r )describes quantum and thermal ﬂuctuations around this mean value.The quantity Ψ(r ,t )is in fact a classical ﬁeld possessing a well–deﬁned phase,reﬂecting a broken gauge sym-metry associated with the condensation process.The expectation value of ˜Ψ(r ,t )is zero and,in the mean–ﬁeld theory,its effects are assumed to be small,amounting to the assumption of the thermodynamic limit,where the number of particles tends to inﬁnity while the density is held ﬁxed.For the effects of ˜Ψ(r )to be negligibly small in the equation for Ψ(r )also amounts to an assumption of zero temperature (i.e.,pure condensate).Given that this is so,and using the normalisationd 3r |Ψ(r ,t )|2=1,(19.6)one is lead to the nonlinear Schr¨o dinger equation,or “Gross–Pitaevskii equation”(GP equation),for the condensate wave function Ψ(r ,t )[13],i¯h ∂Ψ(r ,t )∂t = −¯h 22m ∇2+V trap (r )+NU 0|Ψ(r ,t )|2 Ψ(r ,t ),(19.7)where N is the mean number of particles in the condensate.The nonlinear interaction term (or mean–ﬁeld pseudo–potential)is proportional to the number of atoms in the condensate and to the s –wave scattering length through the parameter U 0.A stationary solution forthe condensate wavefunction may be found by substi-tuting ψ(r ,t )=exp −i μt ¯h ψ(r )into (19.7)(where μis the chemical potential of the condensate).This yields the time independent equation,40019Bose-Einstein Condensation−¯h2 2m ∇2+V trap(r)+NU0|ψ(r)|2ψ(r)=μψ(r).(19.8)The GP equation has proved most successful in describing many of the meanﬁeld properties of the condensate.The reader is referred to the review articles listed in further reading for a comprehensive list of references.In this chapter we shall focus on the quantum properties of the condensate and to facilitate our investigations we shall go to a single mode model.19.3Single Mode ApproximationThe study of the quantum statistical properties of the condensate(at T=0)can be reduced to a relatively simple model by using a mode expansion and subsequent truncation to just a single mode(the“condensate mode”).In particular,one writes the Heisenberg atomicﬁeld annihilation operator as a mode expansion over single–particle states,ˆΨ(r,t)=∑αaα(t)ψα(r)exp−iμαt/¯h=a0(t)ψ0(r)exp−iμ0t/¯h+˜Ψ(r,t),(19.9) where[aα(t),a†β(t)]=δαβand{ψα(r)}are a complete orthonormal basis set and {μα}the corresponding eigenvalues.Theﬁrst term in the second line of(19.9)acts only on the condensate state vector,withψ0(r)chosen as a solution of the station-ary GP equation(19.8)(with chemical potentialμ0and mean number of condensate atoms N).The second term,˜Ψ(r,t),accounts for non–condensate atoms.Substitut-ing this mode expansion into the HamiltonianˆH=d3rˆΨ†(r)−¯h22m∇2+V trap(r)ˆΨ(r)+(U0/2)d3rˆΨ†(r)ˆΨ†(r)ˆΨ(r)ˆΨ(r),(19.10)and retaining only condensate terms,one arrives at the single–mode effective Hamil-tonianˆH=¯h˜ω0a †a0+¯hκa†0a†0a0a0,(19.11)where¯h˜ω0=d3rψ∗0(r)−¯h22m∇2+V trap(r)ψ0(r),(19.12)and¯hκ=U02d3r|ψ0(r)|4.(19.13)19.5Quantum Phase Diffusion:Collapses and Revivals of the Condensate Phase401 We have assumed that the state is prepared slowly,with damping and pumping rates vanishingly small compared to the trap frequencies and collision rates.This means that the condensate remains in thermodynamic equilibrium throughout its prepara-tion.Finally,the atom number distribution is assumed to be sufﬁciently narrow that the parameters˜ω0andκ,which of course depend on the atom number,can be re-garded as constants(evaluated at the mean atom number).In practice,this proves to be a very good approximation.19.4Quantum State of the CondensateA Bose-Einstein condensate(BEC)is often viewed as a coherent state of the atomic ﬁeld with a deﬁnite phase.The Hamiltonian for the atomicﬁeld is independent of the condensate phase(see Exercise19.1)so it is often convenient to invoke a symmetry breaking Bogoliubovﬁeld to select a particular phase.In addition,a coherent state implies a superposition of number states,whereas in a single trap experiment there is aﬁxed number of atoms in the trap(even if we are ignorant of that number)and the state of a simple trapped condensate must be a number state(or,more precisely, a mixture of number states as we do not know the number in the trap from one preparation to the next).These problems may be bypassed by considering a system of two condensates for which the total number of atoms N isﬁxed.Then,a general state of the system is a superposition of number difference states of the form,|ψ =N∑k=0c k|k,N−k (19.14)As we have a well deﬁned superposition state,we can legitimately consider the relative phase of the two condensates which is a Hermitian observable.We describe in Sect.19.6how a particular relative phase is established due to the measurement process.The identiﬁcation of the condensate state as a coherent state must be modiﬁed in the presence of collisions except in the case of very strong damping.19.5Quantum Phase Diffusion:Collapsesand Revivals of the Condensate PhaseThe macroscopic wavefunction for the condensate for a relatively strong number of atoms will exhibit collapses and revivals arising from the quantum evolution of an initial state with a spread in atom number[21].The initial collapse has been described as quantum phase diffusion[20].The origins of the collapses and revivals may be seen straightforwardly from the single–mode model.From the Hamiltonian40219Bose-Einstein CondensationˆH =¯h ˜ω0a †0a 0+¯h κa †0a †0a 0a 0,(19.15)the Heisenberg equation of motion for the condensate mode operator follows as˙a 0(t )=−i ¯h [a 0,H ]=−i ˜ω0a 0+2κa †0a 0a 0 ,(19.16)for which a solution can be written in the form a 0(t )=exp −i ˜ω0+2κa †0a 0 t a 0(0).(19.17)Writing the initial state of the condensate,|i ,as a superposition of number states,|i =∑n c n |n ,(19.18)the expectation value i |a 0(t )|i is given byi |a 0(t )|i =∑n c ∗n −1c n √n exp {−i [˜ω0+2κ(n −1)]t }=∑nc ∗n −1c n √n exp −i μt ¯h exp {−2i κ(n −N )t },(19.19)where the relationship μ=¯h ˜ω0+2¯h κ(N −1),(19.20)has been used [this expression for μuses the approximation n 2 =N 2+(Δn )2≈N 2].The factor exp (−i μt /¯h )describes the deterministic motion of the condensate mode in phase space and can be removed by transforming to a rotating frame of reference,allowing one to writei |a 0(t )|i =∑nc ∗n −1c n √n {cos [2κ(n −N )t ]−isin [2κ(n −N )t ]}.(19.21)This expression consists of a weighted sum of trigonometric functions with different frequencies.With time,these functions alternately “dephase”and “rephase,”giving rise to collapses and revivals,respectively,in analogy with the behaviour of the Jaynes–Cummings Model of the interaction of a two–level atom with a single elec-tromagnetic ﬁeld mode described in Sect.10.2.The period of the revivals follows di-rectly from (19.21)as T =π/κ.The collapse time can be derived by considering the spread of frequencies for particle numbers between n =N +(Δn )and n =N −(Δn ),which yields (ΔΩ)=2κ(Δn );from this one estimates t coll 2π/(ΔΩ)=T /(Δn ),as before.From the expression t coll T /(Δn ),it follows that the time taken for collapse depends on the statistics of the condensate;in particular,on the “width”of the initial distribution.This dependence is illustrated in Fig.19.1,where the real part of a 0(t )19.5Quantum Phase Diffusion:Collapses and Revivals of the Condensate Phase403Fig.19.1The real part ofthe condensate amplitudeversus time,Re { a 0(t ) }for an amplitude–squeezed state,(a )and a coherent state (b )with the same mean numberof atoms,N =250.20.40.60.81-11234560b a is plotted as a function of time for two different initial states:(a)an amplitude–squeezed state,(b)a coherent state.The mean number of atoms is chosen in each case to be N =25.The timescales of the collapses show clear differences;the more strongly number–squeezed the state is,the longer its collapse time.The revival times,how-ever,are independent of the degree of number squeezing and depend only on the interaction parameter,κ.For example,a condensate of Rb 2,000atoms with the ω/2π=60Hz,has revival time of approximately 8s,which lies within the typical lifetime of the experimental condensate (10–20s).One can examine this phenomenon in the context of the interference between a pair of condensates and indeed one ﬁnds that the visibility of the interference pat-tern also exhibits collapses and revivals,offering an alternative means of detecting this effect.To see this,consider,as above,that atoms are released from two conden-sates with momenta k 1and k 2respectively.Collisions within each condensate are described by the Hamiltonian (neglecting cross–collisions)ˆH =¯h κ a †1a 1 2+ a †2a 22 ,(19.22)from which the intensity at the detector follows asI (x ,t )=I 0 [a †1(t )exp i k 1x +a †2(t )expi k 2x ][a 1(t )exp −i k 1x +a 2(t )exp −i k 2x ] =I 0 a †1a 1 + a †2a 2+ a †1exp 2i a †1a 1−a †2a 2 κt a 2 exp −i φ(x )+h .c . ,(19.23)where φ(x )=(k 2−k 1)x .If one assumes that each condensate is initially in a coherent state of amplitude |α|,with a relative phase φbetween the two condensates,i.e.,assuming that|ϕ(t =0) =|α |αe −i φ ,(19.24)40419Bose-Einstein Condensation then one obtains for the intensityI(x,t)=I0|α|221+exp2|α|2(cos(2κt)−1)cos[φ(x)−φ].(19.25)From this expression,it is clear that the visibility of the interference pattern under-goes collapses and revivals with a period equal toπ/κ.For short times t 1/2κ, this can be written asI(x,t)=I0|α|221+exp−|α|2κ2t2,(19.26)from which the collapse time can be identiﬁed as t coll=1/κ|α|.An experimental demonstration of the collapse and revival of a condensate was done by the group of Bloch in2002[12].In the experiment coherent states of87Rb atoms were prepared in a three dimensional optical lattice where the tunneling is larger than the on-site repulsion.The condensates in each well were phase coherent with constant relative phases between the sites,and the number distribution in each well is close to Poisonnian.As the optical dipole potential is increased the depth of the potential wells increases and the inter-well tunneling decreases producing a sub-Poisson number distribution in each well due to the repulsive interaction between the atoms.After preparing the states in each well,the well depth is rapidly increased to create isolated potential wells.The nonlinear interaction of(19.15)then determines the dynamics in each well.After some time interval,the hold time,the condensate is released from the trap and the resulting interference pattern is imaged.As the meanﬁeld amplitude in each well undergoes a collapse the resulting interference pattern visibility decreases.However as the meanﬁeld revives,the visibility of the interference pattern also revives.The experimental results are shown in Fig.19.2.Fig.19.2The interference pattern imaged from the released condensate after different hold times. In(d)the interference fringes have entirely vanished indicating a complete collapse of the am-plitude of the condensate.In(g),the wait time is now close to the complete revival time for the coherent amplitude and the fringe pattern is restored.From Fig.2of[12]19.6Interference of Two Bose–Einstein Condensates and Measurement–Induced Phase405 19.6Interference of Two Bose–Einstein Condensatesand Measurement–Induced PhaseThe standard approach to a Bose–Einstein condensate assumes that it exhibits a well–deﬁned amplitude,which unavoidably introduces the condensate phase.Is this phase just a formal construct,not relevant to any real measurement,or can one ac-tually observe something in an experiment?Since one needs a phase reference to observe a phase,two options are available for investigation of the above question. One could compare the condensate phase to itself at a different time,thereby ex-amining the condensate phase dynamics,or one could compare the phases of two distinct condensates.This second option has been studied by a number of groups, pioneered by the work of Javanainen and Yoo[23]who consider a pair of statisti-cally independent,physically–separated condensates allowed to drop and,by virtue of their horizontal motion,overlap as they reach the surface of an atomic detec-tor.The essential result of the analysis is that,even though no phase information is initially present(the initial condensates may,for example,be in number states),an interference pattern may be formed and a relative phase established as a result of the measurement.This result may be regarded as a constructive example of sponta-neous symmetry breaking.Every particular measurement produces a certain relative phase between the condensates;however,this phase is random,so that the symme-try of the system,being broken in a single measurement,is restored if an ensemble of measurements is considered.The physical conﬁguration we have just described and the predicted interference between two overlapping condensates was realised in a beautiful experiment per-formed by Andrews et al.[18]at MIT.The observed fringe pattern is shown in Fig.19.8.19.6.1Interference of Two Condensates Initially in Number States To outline this effect,we follow the working of Javanainen and Yoo[23]and consider two condensates made to overlap at the surface of an atom detector.The condensates each contain N/2(noninteracting)atoms of momenta k1and k2,respec-tively,and in the detection region the appropriateﬁeld operator isˆψ(x)=1√2a1+a2exp iφ(x),(19.27)whereφ(x)=(k2−k1)x and a1and a2are the atom annihilation operators for the ﬁrst and second condensate,respectively.For simplicity,the momenta are set to±π, so thatφ(x)=2πx.The initial state vector is represented simply by|ϕ(0) =|N/2,N/2 .(19.28)40619Bose-Einstein Condensation Assuming destructive measurement of atomic position,whereby none of the atoms interacts with the detector twice,a direct analogy can be drawn with the theory of absorptive photodetection and the joint counting rate R m for m atomic detections at positions {x 1,···,x m }and times {t 1,···,t m }can be deﬁned as the normally–ordered averageR m (x 1,t 1,...,x m ,t m )=K m ˆψ†(x 1,t 1)···ˆψ†(x m ,t m )ˆψ(x m ,t m )···ˆψ(x 1,t 1) .(19.29)Here,K m is a constant that incorporates the sensitivity of the detectors,and R m =0if m >N ,i.e.,no more than N detections can occur.Further assuming that all atoms are in fact detected,the joint probability density for detecting m atoms at positions {x 1,···,x m }follows asp m (x 1,···,x m )=(N −m )!N ! ˆψ†(x 1)···ˆψ†(x m )ˆψ(x m )···ˆψ(x 1) (19.30)The conditional probability density ,which gives the probability of detecting an atom at the position x m given m −1previous detections at positions {x 1,···,x m −1},is deﬁned as p (x m |x 1,···,x m −1)=p m (x 1,···,x m )p m −1(x 1,···,x m −1),(19.31)and offers a straightforward means of directly simulating a sequence of atom detections [23,24].This follows from the fact that,by virtue of the form for p m (x 1,···,x m ),the conditional probabilities can all be expressed in the simple formp (x m |x 1,···,x m −1)=1+βcos (2πx m +ϕ),(19.32)where βand ϕare parameters that depend on {x 1,···,x m −1}.The origin of this form can be seen from the action of each measurement on the previous result,ϕm |ˆψ†(x )ˆψ(x )|ϕm =(N −m )+2A cos [θ−φ(x )],(19.33)with A exp −i θ= ϕm |a †1a 2|ϕm .So,to simulate an experiment,one begins with the distribution p 1(x )=1,i.e.,one chooses the ﬁrst random number (the position of the ﬁrst atom detection),x 1,from a uniform distribution in the interval [0,1](obviously,before any measurements are made,there is no information about the phase or visibility of the interference).After this “measurement,”the state of the system is|ϕ1 =ˆψ(x 1)|ϕ0 = N /2 |(N /2)−1,N /2 +|N /2,(N /2)−1 expi φ(x 1) .(19.34)That is,one now has an entangled state containing phase information due to the fact that one does not know from which condensate the detected atom came.The corre-sponding conditional probability density for the second detection can be derived as19.6Interference of Two Bose–Einstein Condensates and Measurement–Induced Phase 407n u m b e r o f a t o m s n u m b e r o f a t o m s 8position Fig.19.3(a )Numerical simulation of 5,000atomic detections for N =10,000(circles).The solid curve is a least-squares ﬁt using the function 1+βcos (2πx +ϕ).The free parameters are the visibility βand the phase ϕ.The detection positions are sorted into 50equally spaced bins.(b )Collisions included (κ=2γgiving a visibility of about one-half of the no collision case.From Wong et al.[24]40819Bose-Einstein Condensationp (x |x 1)=p 2(x 1,x )p 1(x 1)=1N −1 ˆψ†(x 1)ˆψ†(x )ˆψ(x )ˆψ(x 1) ˆψ†(x 1)ˆψ(x 1) (19.35)=12 1+N 2(N −1)cos [φ(x )−φ(x 1)] .(19.36)Hence,after just one measurement the visibility (for large N )is already close to 1/2,with the phase of the interference pattern dependent on the ﬁrst measurement x 1.The second position,x 2,is chosen from the distribution (19.36).The conditional proba-bility p (x |x 1)has,of course,the form (19.32),with βand ϕtaking simple analytic forms.However,expressions for βand ϕbecome more complicated with increasing m ,and in practice the approach one takes is to simply calculate p (x |x 1,···,x m −1)numerically for two values of x [using the form (19.30)for p m (x 1,...,x m −1,x ),and noting that p m −1(x 1,...,x m −1)is simply a number already determined by the simu-lation]and then,using these values,solve for βand ϕ.This then deﬁnes exactly the distribution from which to choose x m .The results of simulations making use of the above procedure are shown in Figs 19.3–19.4.Figure 19.3shows a histogram of 5,000atom detections from condensates initially containing N /2=5,000atoms each with and without colli-sions.From a ﬁt of the data to a function of the form 1+βcos (2πx +ϕ),the visibil-ity of the interference pattern,β,is calculated to be 1.The conditional probability distributions calculated before each detection contain what one can deﬁne as a con-000.10.20.30.40.50.60.70.80.91102030405060number of atoms decided 708090100x=0x=1x=2x=4x=6Fig.19.4Averaged conditional visibility as a function of the number of detected atoms.From Wong et al.[13]19.7Quantum Tunneling of a Two Component Condensate40900.51 1.520.500.5Θz ο00.51 1.520.500.5Θx ο(b)1,234elliptic saddle Fig.19.5Fixed point bifurcation diagram of the two mode semiclassical BEC dynamics.(a )z ∗,(b )x ∗.Solid line is stable while dashed line is unstable.ditional visibility .Following the value of this conditional visibility gives a quantita-tive measure of the buildup of the interference pattern as a function of the number of detections.The conditional visibility,averaged over many simulations,is shown as a function of the number of detections in Fig.19.4for N =200.One clearly sees the sudden increase to a value of approximately 0.5after the ﬁrst detection,followed by a steady rise towards the value 1.0(in the absence of collisions)as each further detection provides more information about the phase of the interference pattern.One can also follow the evolution of the conditional phase contained within the conditional probability distribution.The ﬁnal phase produced by each individual simulation is,of course,random but the trajectories are seen to stabilise about a particular value after approximately 50detections (for N =200).19.7Quantum Tunneling of a Two Component CondensateA two component condensate in a double well potential is a non trivial nonlinear dynamical model.Suppose the trapping potential in (19.3)is given byV (r )=b (x 2−q 20)2+12m ω2t (y 2+z 2)(19.37)where ωt is the trap frequency in the y –z plane.The potential has elliptic ﬁxed points at r 1=+q 0x ,r 2=−q 0x near which the linearised motion is harmonic withfrequency ω0=q o (8b /m )1/2.For simplicity we set ωt =ω0and scale the length in units of r 0= ¯h /2m ω0,which is the position uncertainty in the harmonic oscillatorground state.The barrier height is B =(¯h ω/8)(q 0/r 0)2.We can justify a two mode expansion of the condensate ﬁeld by assuming the potential parameters are chosen so that the two lowest single particle energy eigenstates are below the barrier,with41019Bose-Einstein Condensation the next highest energy eigenstate separated from the ground state doublet by a large gap.We will further assume that the interaction term is sufﬁciently weak that, near zero temperature,the condensate wave functions are well approximated by the single particle wave functions.The potential may be expanded around the two stableﬁxed points to quadratic orderV(r)=˜V(2)(r−r j)+...(19.38) where j=1,2and˜V(2)(r)=4bq2|r|2(19.39) We can now use as the local mode functions the single particle wave functions for harmonic oscillators ground states,with energy E0,localised in each well,u j(r)=−(−1)j(2πr20)3/4exp−14((x−q0)2+y2+z2)/r20(19.40)These states are almost orthogonal,with the deviation from orthogonality given by the overlap under the barrier,d3r u∗j(r)u k(r)=δj,k+(1−δj,k)ε(19.41) withε=e−12q20/r20.The localised states in(19.40)may be used to approximate the single particle energy(and parity)eigenstates asu±≈1√2[u1(r)±u2(r)](19.42)corresponding to the energy eigenvalues E±=E0±R withR=d3r u∗1(r)[V(r)−˜V(r−r1)]u2(r)(19.43)A localised state is thus an even or odd superposition of the two lowest energy eigenstates.Under time evolution the relative phase of the superposition can change sign after a time T=2π/Ω,the tunneling time,where the tunneling frequency is given byΩ=2R¯h=38ω0q20r2e−q20/2r20(19.44)We now make the two-mode approximation by expanding theﬁeld operator asˆψ(r,t)=c1(t)u1(r)+c2(t)u2(r)(19.45) where。

玻色爱因斯坦凝聚体中暗孤子动力学研究_刘超飞

收稿日期：2016-01-03基金项目：国家自然科学基金资助项目（11304130，11365010，61565007）；江西省教育厅资助项目（GJJ150685）；江西省科技厅资助项目（20151BAB212002）；江西理工大学清江青年英才支持计划资助作者简介：刘超飞（1981-），男，博士，副教授，主要从事玻色爱因斯坦凝聚等方面的研究，E-mail:liuchaofei0809@.江西理工大学学报ＪｏｕｒｎａｌｏｆＪｉａｎｇｘｉＵｎｉｖｅｒｓｉｔｙｏｆＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙ第37卷第5期2016年10月Vol.37,No.5Oct.20160引言在原子凝聚体中，原子间的排斥相互作用导致了许多有趣的非线性现象，一个最典型的例子就是暗孤子的形成[1-8].实验上，通过相印记[1-2]和扰动原子密度[3]方法，暗孤波已在稀释玻色-爱因斯坦凝聚体中产生.在相印记方法中，凝聚体的一个部分被远失谐激光束短时间照射，使它获得了相移，但没有产生重大的密度扰动.根据相映射实验，将出现一个暗孤立波，以及作为副产品的声波.在随后阶段的相映射实验中，由于其固有的不稳定性和横向激发[4]，人们观察到暗孤波退化为涡旋环.另一项实验涉及到通过慢光技术从凝聚体中突然清除一个盘状区域，它将产生暗孤波（暗孤波也将蜕化成涡旋环）和相反传播的声波[3].暗孤子的大小接文章编号：2095-3046（2016）05-0102-10DOI:10.13265/ki.jxlgdxxb.2016.05.016玻色爱因斯坦凝聚体中暗孤子动力学研究刘超飞a ，潘小青a ，张赣源b（江西理工大学，a.理学院；b.应用科学院，江西赣州341000）摘要：文章从侦测暗孤子能量的角度，全面介绍了玻色爱因斯坦凝聚体中暗孤子发生声波辐射以及与声波相互作用的动力学行为.虽然暗孤子-声波相互作用导致暗孤子能量变化，暗孤子在简谐势阱中，以及简谐势阱受到扰动时，都具有类似粒子运动的动力学行为.考虑Rabi 耦合时，暗孤子还可以转化为稳定传播的矢量暗孤子.文章对玻色爱因斯坦凝聚体中暗孤子动力学研究的进行全面总结，将加深人们对暗孤子现象的认识.此外，类似暗孤子-声波相互作用的行为，也将出现在其他孤子动力学研究中.关键词：玻色爱因斯坦凝聚；暗孤子；声波；孤子声波相互作用；简谐势阱中图分类号：O469文献标志码：AKinetic study of dark soliton in Bose-Einstein condensateLIU Chaofei a ,PAN Xiaoqing a ,ZHANG Ganyuan b(a.Faculty of Science;b.Faculty of Applied Science,Jiangxi University of Science and Technology,Ganzhou 341000,China)Abstract ：By exploring the energy of dark soliton,this paper systematically introduces the dynamical behavior of sound-emission of dark soliton and dark soliton-sound interaction in Bose-Einstein condensate.Although the dark soliton-sound interaction leads to the change of the dark soliton ′s energy,dark soliton displays the particle-like behavior very well in the harmonic potential and even in the periodic perturbed harmonic trap.Under the Rabi coupling,dark soliton can transfer into the vector dark soliton,which propagates stably in the condensates.This paper provides a full overview of the kinetic study of dark soliton,and it will greatly increase people ′s knowledge about the dark soliton.Furthermore,similar behaviors of the dark soliton-sound interaction will occur in the dynamical investigation of other soliton.Key words ：Bose-Einstein condensate;dark soliton;sound waves;soliton-sound interaction;harmonic trap刘超飞，等：玻色爱因斯坦凝聚体中暗孤子动力学研究近当前成像技术的极限.在这些实验中，通过一再释放势阱中的凝聚体，让其膨胀，然后采取一个光学吸收图像来实现成像.产生暗孤子的更先进的方法也已经提出[9-10]，包括合并相映射和密度工程方法[5-6].暗孤子也可用两个凝聚碰撞产生[11]，以及用凝聚体的布拉格光学晶格反映[12-13]产生.迄今为止，在几何形状上，暗孤子的实验可以从球对称[4]到高度拉长的情况（长宽比大于30[2]）.这些系统在性质上仍然是三维，致使暗孤子容易由于横向不稳定性而被破坏，从而迅速衰变为涡旋.这是个关键因素，它限制了观察到的孤子寿命，其值约为数十毫秒.然而，最近的实验发现，在准一维凝聚体中，暗孤波将是亚稳的，其寿命可大大延长到直至数秒[14].理论上，根据零温平均场理论，稀释原子形成的玻色爱因斯坦凝聚体可由Gross-Pitaevskii方程描述，它是一个很好的非线性系统.在外势场为零时，Gross-Pitaevskii方程支持暗孤子解.与光学系统的一个重要的区别是，势阱能导致凝聚体有不均匀的背景密度.对于凝聚体中的暗孤波，其性质与光学系统中的暗孤波类似.三维暗孤波由于横向激发，是不稳定的，因此，暗孤子容易衰变弯曲形成涡旋.在实验中，可将凝聚体在横向高度压缩，从而构成所谓的一维体系，而暗孤子的运动则由不均匀的纵向密度控制.总的说来，暗孤子是一个局部的密度缺陷，就像一个凹槽，其周围充满凝聚体.并且，暗孤子的两边存在位相差，它是散焦的色散效应与聚焦的非线性相互作用之间达到平衡的结果.因此，暗孤子的主要特性之一，是能在传播中保持其局域化的形状不变[11,15-17].通常，研究孤子的文章主要是给出新的孤子解，或者探索孤子的不稳定性等.在这篇文章中，系统性的介绍玻色爱因斯坦凝聚体中暗孤子动力学研究.与大多数论文不同，本文的研究重点是暗孤子能量的计算和动力学演化中的能量侦测.暗孤子受到外界环境的干扰后，发生声波辐射，暗孤子能量降低，速度变快，同时，外界的声波又能反作用于暗孤子，使其能量增加.这种研究暗孤子的方式，最初在文献[18-19]上介绍.实际上，这种研究思路完全可以推广到其他孤子解的动力学研究上.1暗孤子的数值解凝聚体背景密度均匀，且为.则含有速度为和位置在的暗孤子的一维的形式的凝聚体的波函数为：ψs（z，t）=n姨exp（-iμ攸t）·i v+1-v2姨tanh[1-v2姨（z-vt）姨姨]（1）这里ξ＝攸n gm姨为凝聚体的愈合长度，它可用来刻画暗孤子的尺度.暗孤子的速度依赖于密度n d和通过其中心的相移S.并且，v/c=1-（n d/n）姨= cos（S/2），孤子速度的最大值由Bogoliubov声速度/ C=ng/m姨决定.这里有两个极限情况：①固定暗孤子完全是黑色的，即有一个零密度节点，以及π相滑移；②孤子速度为c时，无相滑移，也与背景的密度无差异，因此难以分辨.图1展示了各种速度暗孤子的密度和相位.一个重要特点是，固定孤子的能量最高，在v=c时，孤子能量基本上为零，这导致人们认为孤子具有负有效质量的设想.由于不会耗散，孤子往往类似于粒子.事实上，对于一阶弱作用力，暗孤子就像一个有效质量为负的经典粒子[11,20-22].这意味着，例如，在一个谐势阱的凝聚体中，暗孤子将趋于来回震荡.S图1各种速度暗弧子密度和相位（b）暗孤子的相位Sπ/20.0-π/2-505z/ξv/c=0v/c=0.25v/c=0.50v/c=0.751.00.50.0-505z/ξ（a）暗孤子的密度nn/nv/c=0v/c=0.25v/c=0.50v/c=0.75第３７卷第５期103与纵向均匀系统（如非线性光学纤维）不同，原子凝聚体沿孤子的运动方向有束缚.例如，暗孤子在谐势阱束缚下的凝聚体中，将趋于在势阱中来回震荡.由于这一空间束缚，孤子一般会与其他激发共同存在，例如声波.因此，文中说的孤子并非纯数学意义中的孤子，而是指一个空间区域（即“孤子地区”），该区域存在密度凹陷和相位差，以及其他可能的激发.这可以粗略地视为一个受到扰动的孤子[23].2暗孤子的能量在无限大体系中，重整化的一维暗孤子能量（即除去背景流体的贡献部分）由式（2）给出，E s ol ＝4攸n 3/21-vc23/2（2）然而，在非均匀凝聚体中，只有当暗孤子在密度局部均匀区域，该方程式才有效.这将在后面的文章中进一步加以说明.还有一种得出一维暗孤子能量方法，它基于对Gross-Pitaevskii 方程的数值积分，即使在密度不均匀时仍然有效，ε（ψ）=攸22m荦ψ2+V ψ2+12g ψ4（3）暗孤子能量E s 是通过对孤子位置Z s 积分一个距离z ints ，然后减去对应的时间独立性背景密度n TI 的贡献，即E s =Z s +z intsZ s -z ints乙ε（ψ）d z -Z s +z intsZ s -z ints乙ε（n TI 姨）d z （4）“孤子区域”必须足够大，以包含暗孤子能量的绝大部分.实际上，暗孤子的速度和背景密度都影响暗孤子密度凹陷处的宽度.图2说明了积分后的各种速度的孤子能量（实线），分别与该区域的大小和背景密度的函数关系.当z ints >5ξ时，积分得到的能量值几乎与从公式（4）（虚线）渐近预测的数值完全相等，所以我们选择“孤子区域”为（Z s ±5ξ）.在时间依赖性模拟中，“孤子区域”能含有声波.通常，很难能区分孤子能和声波能，但后者的数值，在孤子的速度不是很大时，是非常小的.3简谐势阱中的暗孤子假设凝聚体在一维谐势阱中，现在考虑暗孤子在凝聚体中的动力学行为，简谐势阱为：V （z ）=12ω2z z 2（5）这种势阱通常是由磁场形成.系统在空间上是有限的，这是体系的一个重要特点.因此，凝聚体的大小可用托马斯-费米半径刻画.对于一个束缚频率为ωz 的势阱，托马斯-费米密度分布是一个倒抛物线型，其中n TF =（1-ω2z z 2/2），托马斯-费米半径为R TF =2/ω2z 姨.图3（a ）展示了一个速度为v =0.5c 的暗孤子在势阱中，其初始位置为z=0，凝聚体的密度峰值为n 0，纵向束缚频率为ωz =2姨×10-2（μ/攸）.实际上，除了在边界附近由于小动能贡献导致‘尾巴’状热云外，托马斯-费米密度分布与真实实验状况吻合得很好（图3（a ））.因此，凝聚体的实际大小刚好大于托马斯-费米半径100ξ.在该系统中原子密度的时间演化由图3（b ）显示，其纵坐标为位置，横坐标为时间.暗孤子是一个局部的密度极小.它在向势阱壁移动过程中减速，当其密度极小处触及零密度时，孤子的运动方向改变.众所周知，当孤子远离势阱中心时，其相滑移增加，并在最大振幅处达到π.随后，孤子改变其运动方向，孤子的相滑移变到-π.图2通过积分得到的暗孤子能量E s与积分区间宽度z int 的函数关系1.51.00.50.0012345（a ）固定背景密度n 0时z int /ξv =0.5cv =0.75cv =0.25cv =0c E s /μ（b ）固定孤子速度v =0.5c 时1.00.80.60.40.20.0012345n =n 0n =0.8n 0n =0.6n 0n =0.4n 0n =0.2n 0z int/ξE s /μ江西理工大学学报2016年10月104在谐势阱中的暗孤子，其振荡频率近似为ωz /2姨[14,15,24-28].这是由分析托马斯-费米密度分布得到的，并且人们已用数值模拟证实了这一结论.在图3（b ）中，孤子振荡周期约为T s =630（ξ/c ），而势阱的周期约为T z =444（ξ/c ），这与理论预测结果相同.孤子的运动扰动背景流体，致使流体振荡，其幅度约为2%n 0.可以用下式定义背景凝聚体的偶极振荡，D =乙z ψ（z ）2d z（6）凝聚体的偶极运动D 和孤子路径Z s 由图3（c ）给出.暗孤子的振荡频率为ωz /2姨（实线），它诱发了势阱中背景流体的偶极振荡（虚线），其频率ωz 为[27].在一定程度上，孤子行为就像搅拌器，搅拌着流体.暗孤子将在势阱中加速，孤子由于辐射声波而衰减.在这种情况下，暗孤子的深度将变浅，而其速度将变快，从而更进一步逼近势阱壁，并导致了与反阻尼类似的现象.这与阻尼谐振子相比，结果正好相反.对于阻尼谐振子，其振荡幅度随时间减小.但是，简谐势阱中看不到任何的孤子震荡幅度的净变化，因此，不能推断出孤子的衰变.在图4（a ）中仅仅观察到的孤子振幅的小周期调制，例如，围绕其平均幅度，大约有1％的最大调制幅度变化.类似的结果也出现在了孤子能量的进化中，如图4（b ）所示.孤子能量的平均值保持不变，但有振荡调制.3.1速度对暗孤子运动的影响图5（a ）显示了不同初始速度的暗孤子在一个固定的简谐势阱中的演化路径.增加孤子的速度，其主要结果是孤子振幅增加.但即使孤子速度高达0.7c ，孤子的振荡频率仍保持在预期值ωs =ωz /2姨的周围.但对于运动速度非常快的孤子，如v =0.9c ，如图5（c ）中（点虚线），我们看到其值略有偏差，它趋向于数值更高的振荡频率.这很可能是因为这个快速运动的浅孤子进入了凝聚体的边界造成.在边界处，凝聚体的密度偏离于托马斯-费米密度分布.相比之下，较慢的孤子的振荡束缚在势阱中心，而势阱中心的密度分布几乎与托马斯-费米密度分布相同.对于托马斯-费米密度分布，在简谐势阱中的暗孤子的振荡幅度与孤子的初始速度成正比.为比较各种速度的中心孤子，我们可以使用这一关系，重整化孤子的位置.图5（b ）显示了孤子位置的重整化图.各速度下的重整化位置随孤子速度的增加而幅度增大.对于不同速度的孤子的能量的振荡演化，这种效应也被观察到，如图5（c ）所示.对于快速运动的孤子（例如，v =0.9c ）（点虚线），其能量调制延伸到了最初能量的0.4倍.1.00.50.0-1001001.00.50.0z /ξn /n 0V /μ（a ）凝聚体在简谐势阱中的密度（左轴，实线），其中ωz =2姨×10-2（μ/攸）（右轴，虚线）.速度为v =0.5c 的暗孤子在势阱中心.图4简谐势阱中暗孤子振幅与能量关系E s /μ0.860.850.845000100001500020000t /（ξ·c -1）（b ）暗孤子能量E s 的演化100-100z /ξ（b ）凝聚体随时间演化的重整化图500-50500-505001000150020002500t /（ξ·c -1）（c ）暗孤子位置Z s （实线，左轴）和凝聚体的偶极运动D （虚线，右轴）图3简谐势阱中的暗孤子运动z s /ξD /ξ2n 05001000150020002500t /（ξ·c -1）5150z S/ξ（a ）暗孤子离开势阱中心的距离Z s5000100001500020000t /（ξ·c -1）刘超飞，等：玻色爱因斯坦凝聚体中暗孤子动力学研究第３７卷第５期105刚才我们已经看到，这些小的位置和能量调制，由振荡孤子对背景流体的干扰产生，随后反馈到孤子上.孤子速度越快，有效质量越小，所以背景流体的振荡将会对它们有更大的反馈作用，从而引起更大的调制.3.2束缚势阱强度对暗孤子行为的影响增加势阱的纵向强度，同时保持凝聚体密度峰值固定，这将减少凝聚体的空间范围.对于固定速度的孤子，振荡振幅的绝对值下降，但仍然与托马斯-费米半径形成一个近似的常数比值.图6（a ）显示了各种强度的简谐势阱中，暗孤子位置的变化，其位置已根据托马斯-费米半径重新标度，而时间单位也调整为ω-1z .对于较低的势阱频率，例如，ωz =2姨×10-2（μ/攸）（黑色实线），暗孤子的振荡频率即为其分析预测值ωz /2姨.而对于较高频率的势阱，例如，ωz =62姨×10-2（μ/攸）（点虚线），暗孤子的振荡频率比预测值大.图6（c ）给出了各种初始速度的孤子的振荡频率与势阱频率的函数.对于弱势阱有ωs /ωz ≈2姨，分析值与预测值吻合（虚线）.然而，增加势阱的强度，ωs /ωz 的比值偏离预测值，并随势阱强度的增加而单调增加.这种偏差是不可忽略的，在这里的势阱频率范围里，这个值可高达10％.造成此偏差的原因，是因为对孤子频率的分析预测值，假定了凝聚体为托马斯-费米密度分布.在我们的数值方法里，其密度峰值保持固定，这一假定仅对弱的简谐势阱有效（大量的粒子）.图6（d ）显示了整个凝聚体在各种势频率中的轴向密度分布.随着势阱频率的增加，密度越来越背离倒抛物线型的托马斯-费米密度分布.这偏差在凝聚体的边界处最为明显.甚至可以看到"尾巴"状的低密度伸展通过托马斯-费米半径.为了突出这种偏差，文章还在同一图中绘制了实际密度与托马斯-费米密度分布的差值.文章认为，这一偏差可以解释暗孤子的振荡频率与势阱频率的变化有关.当势阱强度增加时，暗孤子位置（图6（a ））和能量（图6（b ））的调制，由于暗孤子与偶极振荡相互作用的增大而增加.在这里，凝聚体的尺度减小，从而其有效质量降低，而孤子基本保持相同的大小.因此，振荡暗孤子诱发背景凝聚体相对更大的扰动，从而导致孤子的动力学调制更大.4在周期性扰动势阱中的暗孤子通常，人们假设凝聚体在静态简谐势阱中，势阱为V har （x ）=m ω2x 2/2，ω是势阱的频率.在这里，我（c ）孤子能量的演化，该结果经过了由最初的孤子能量E inits 的重整10002000300040005000t /（ξ·c -1）图5不同速度暗孤子在简谐势阱中位置与能量演化1.00.80.60.4E s /E s i n t10610410210098z S /v s /ξ（b ）暗孤子在势阱中的距离（用孤子速度重整后的结果）10002000300040005000t /（ξ·c -1）100500-50-100z s /ξ（a ）暗孤子在无限深简谐势阱（ωz =2姨×10-2（μ/攸））中位置的演化.其中初始速度v /c =0.1（实线），0.5（虚线），0.7（点线）和0.9（点虚线）10002000300040005000t /（ξ·c -1）江西理工大学学报2016年10月106们考虑整个势阱存在扰动，这种势阱可写为V Ext （x ，t ）=m ω2[x +h sin （ωd t ）]2（7）h 和ωd 分别是扰动的幅度和频率[29].我们用数值模拟对以上模型进行研究.图7显示了在各种振幅的扰动下，孤子能量的演变.显然，孤子的运动方向与扰动的运动方向的耦合决定了孤子的演化.在没有扰动的情况下（h =0），该模型将退化为孤子在简谐势阱中的振荡试验.在这种情况下，孤子在背景密度不均匀的凝聚体中传播，其外形变得不对称，同时它还会向相反方向辐射声波[18].众所周知，在简谐势阱中的孤子的振动频率为ωsol =ω/2姨，孤子会发射和重新接收声波.总体而言，孤子不断受到孤子自身带来的流体的扰动，但并不会衰退.当势阱的运动方向与孤子的运动方向相反时（h >0），孤子往往首先获得能量直至达到峰值，然后孤子能量减小到原来的值.这种能量的增益损失周期性的重复.能量变化的周期为1516ξ/c .因此在大量的时间里孤子能量比其初始值大.增加扰动幅度，可以提高孤子在其能量循环中获得和失去能量的能力.相反，当势阱与孤子的移动方向相同时（h <0），孤子往往首先失去能量直至最低值，然后它重新获得能量，恢复其原始值.这个过程构成一个损失获得循环，其周期为1516ξ/c .同样，增加扰动的幅度，孤子失去更多的能量，然后恢复至初始值.因此在大量的时间里孤子能量比其初始值小.实际上，图7比较了在受周期性扰动和不受周期性扰动的简谐势阱中，暗孤子的演变.图8显示了在各种振幅的扰动下，相应的暗孤子的位置的演变.一般来说，凝聚体会伴随势阱运动.由于势阱的振幅和振荡频率都非常小，势阱振荡导致孤子的轨道与没有受扰动的情况（h =0）发生偏离.孤子振荡周期出现波动.整体而言，孤子的振荡频率仍然围绕着ωsol =ω/2姨.因此，这一特征也确保了势阱移动的方向与孤子移动方向的耦合.当势阱与孤子有相同的运动方向时，凝聚体伴随势阱运动.因此，孤子被携带着运动，它偏离振荡中心更远.从而使孤子能量往往比原来的值小（见图8）.但是，当势阱运动与孤子运动方向相反时，凝聚体的运动相对地减小了孤子的震荡幅度，孤子能量出现增益-损失循环.一般来说，如果势阱扰动频率等于天然粒子的振荡频率，很可能引发共鸣.虽然很多研究显示孤子具有粒子状特性，但是，孤子毕竟不是正常的粒子，因此即便势阱振荡与孤子震荡能很好的匹配，也无法造成孤子的共振行为.图6暗孤子在不同频率简谐势阱中位置与能量演化0.90.80.70204060t /ω-1（b ）对应的孤子能量的进化E s /μ（d ）图（a ）和（b ）凝聚体的密度（左轴），以及其与托马斯-费米值的密度偏离（n -n TF ）1.00.750.500.250.000.5 1.0Z /R TFn /n 00.040.020.00-0.02n -n T F /n 00.500.250.00-0.25-0.50z s /R T F（a ）暗孤子位置与托马斯-费米半径R TF 的比值.其中势阱强度ωz =ω0=2姨×10-2（μ/攸）（实线），2ω0（虚线），4ω0（点线）和6ω0（点虚线）204060t /ω-1刘超飞，等：玻色爱因斯坦凝聚体中暗孤子动力学研究第３７卷第５期（c ）孤子振荡频率ωs 与势阱频率ωz 的函数关系，其中.虚线为分析值ωs =ωz /2姨0.800.750.700.10.2ωs /ωzωz /（c ·ξ-1）v /c =0.25v /c =0.50v /c =0.751075暗孤子动力学研究展望实验中，暗孤子的确能展示良好的粒子状特性.例如，暗孤子在简谐势阱中，通常会来回振荡[11,15].如果一个静态的暗孤子最初位置不在势阱中心，它将受到势阱的外力作用，使之加速向势阱中心运动.最近，Parker 和他的同事考虑了对简谐势阱做一些修正，充分展示了可能出现的暗孤子行为[18,19,30].将一个紧束缚的势阱嵌入一个弱束缚的简谐势阱中，这样就可控制声波的逃逸[18].将光晶格势阱加入简谐势阱中，就可用于干扰暗孤子[30].此外，还可以考虑了参数驱动以及阻尼机制[19].而参考文献[31]中，有限温度效应对暗孤子的影响受到了系统性的探讨.类似的，Bilas 和Pavloff 研究了准一维玻色爱因斯坦凝聚体中，随机势对运动暗孤子的影响[32]，还研究了暗孤子在传播途中遇到障碍的情况[33].除了上述单成分凝聚体中暗孤子的工作，随着对BECs 中孤子的深入研究，人们在多元凝聚体混合物中还发现了矢量孤子、如亮-暗矢量孤子[34-37]、亮-亮矢量孤子、暗-暗矢量孤子[38-40].和单分量凝聚体相比，玻色爱因斯坦凝聚体的二元混合物已经被显示具有迷人的宏观量子现象，如复杂的空间结构[41-43]、亚稳态[44-46]、对称破缺不稳定性[47-49].迄今为止，我们已经知道种间相互作用系数对凝聚体混合物的基态结构起决定作用.当不等式g 12≤g 1g 2姨满足时，两种凝聚体是易融合的；当g 12>g 1g 2姨时，由于很强的种间排斥相互作用，凝聚体为不可融合的.就考察孤子来说，多个分量这个自由度的引入给系统带来了丰富的非线性现象，比如：孤子链、孤子对、多模激发等.除此之外，一种新型的孤子，即共生孤子，在两分量87Rb 和85Rb 的凝聚体中被发现.此时，只要分量间原子吸引力足够强，便能够克服各自分量原子内的排斥力而起到一个有效吸引的作用，从而在两分量玻色―爱因斯坦凝聚中形成亮孤子.其实，早在1993年，Kivshar 等[50]通过求图8扰动中的暗孤子震荡（b ）图（a ）的放大图.长度单位为ξ=攸/m μ姨（a ）在各种扰动幅度下，孤子轨迹随时间的变化403020100-10-20-30-40200040006000X /qt /（ξ·c -1）403020100-10-20-30-40300600900h =3ξh =2ξh =1ξh =0h =-1ξh =-2ξh =-3ξt /（ξ·c -1）X /q1.281.241.201.151.121.083000600090001200015000t /（ξ·c -1）（a ）在各种扰动幅度下，孤子能量随时间的变化.孤子初始速度为0.3c 和初始位置为x =0E /μ 1.281.241.201.151.121.0850010001500t /（ξ·c -1）（b ）图（a ）的放大图.体系参数为ω=2姨/100（c /ξ），ωd =ω/2姨h =3ξh =2ξh =1ξh =0h =-1ξh =-2ξh =-3ξ图7孤子能量受扰动幅度的影响E /μ江西理工大学学报2016年10月108解两个耦合非线性薛定谔方程，显示了矢量暗孤子存在的可能性.近来，在耦合的一维非线性薛定谔方程的框架内，双组分凝聚体的矢量暗孤子得到了相应的研究[51].然而，这些研究基于稳定的媒质.无论是凝聚体的种类，还是凝聚体两成分的比率，都是固定的.最近，Rabi 耦合[52]被用于将凝聚体从某一成分向其他成分的凝聚体转化.这就暗示着暗孤子在动态凝聚体媒质中运行是可能的.对于理想情况，即凝聚体种间相互作用与种内相互作用强度相同时，我们可以看到由一种成分构成的暗孤子可以转化为另一成分的暗孤子[53]（如图9所示）.并且，暗孤子转变为动态的矢量暗孤子后，其运动轨迹不受Rabi 耦合强度的影响.而种间相互作用与种内相互作用强度不相同时，矢量暗孤子的运动轨迹受到Rabi 耦合的影响比较明显.但是，这一长时间模拟所说明的最主要的结论为：暗孤子可以在具有Rabi 耦合的凝聚体中存在.这一结果在特定凝聚体比率的矢量暗孤子的设计上，具有非常重要的意义.将来，还可以通过控制Rabi 耦合，如在特定时间终结Rabi 耦合，从而得到特定比率的凝聚体混合物，以及矢量暗孤子.当然，矢量暗孤子稳定存在的内在机制等还有待进一步的探索.相信该研究成果将给实验上认识凝聚体中暗孤子等激发行为提供理论支持.6结论文章根据详细的能量计算分析，对玻色爱因斯坦凝聚体中暗孤子动力学行为的数值研究进行了系统的介绍.从分析暗孤子的数值解、暗孤子的能量计算开始，重点探讨了暗孤子在简谐势阱中的动力学行为，以及简谐势阱出现扰动时的情况.玻色爱因斯坦凝聚体中的暗孤子具有类似粒子的运动行为，在非均匀密度的凝聚体中，暗孤子发生声波图9玻色爱因斯坦凝聚体二元混合物中的矢量暗孤子的震荡行为.Rabi 耦合强度为0.025，g 1=g 2=g 12=1（a ）凝聚体成分1的演化和暗孤子震荡行为100500-50-1006001200180024003000t /（ξ·c -1）X0.13750.27500.41250.55000.58750.82500.96251.100ψ12（b ）凝聚体成分2的演化与暗孤子的震荡行为6001200180024003000X0.13750.27500.41250.55000.58750.82500.96251.100ψ22100500-50-100（c ）两凝聚体的密度和6001200180024003000X0.13750.27500.41250.55000.58750.82500.96251.100ψ12+ψ22100500-50-100t /（ξ·c -1）t /（ξ·c -1）刘超飞，等：玻色爱因斯坦凝聚体中暗孤子动力学研究第３７卷第５期109。

Quantum group invariant, nonextensive quantum statistical mechanics

∗
Electronic address:ubriaco@
1
II. NONEXTENSIVE STATISTICAL MECHANICS AND QUANTUM GROUP INVARIANCE
In this section we study the consequences of introducing quantum group invariance in nonextensive statistical mechanics. The density matrix ρ ˆ for nonextensive thermodynamics was reformulated in [12,13] in terms of normalized q -expectation values. The recent formalism deﬁnes the expectation values of an operator as Aq =
Quantum group invariant, nonextensive quantum statistical mechanics
Marcelo R. Ubriaco∗ Laboratory of Theoretical Physics Department of Physics University of Puerto Rico P. O. Box 23343, R´ ıo Piedras PR 00931-3343, USA
R
,
(1)
where pR is the probability for the ensemble to be in the state R and q is parameter which in principle can be any real number. The Tsallis entropy Sq possesses the properties of the Shanon’s entropy except that of additivity. Given two independent systems Σ and Σ′ , the entropy Sq satisﬁes the pseudoadditivity property

Exciton-phonon droplets with Bose-Einstein condensate transport and optical properties

2
D. Roubtsov et al.
where ncoh (x, t) ≈ ncore (x − vt) is the ballistic (superﬂuid) part of the packet, ncore (x − vt) ≃ |Ψ0 |2 (x − vt), (4) and ∆n(x, t) is the noncondensed part of it. The following decomposition can be written for the out-of-condensate part: ˆ† δ ψ ˆ(x, t) ≈ δncloud (x, t) + δntail (x, t). ∆n(x, t) = δ ψ The challenging problem is how to describe the spatially inhomogeneous state of the droplet with the exciˆ(x, t), where tonic BEC inside in terms of Ψ0 (x, t) and δ ψ ˆ δ ψ is the “ﬂuctuating” part of the exciton Bose-ﬁeld. For example, within the quasistationary approximation, one has to calculate ncore (x/L0 ) and ∆n(x) ≃ δno, cloud (x/Lch ) + δntail (5) and understand how the diﬀerent characteristic lengthes and coherence properties appear in the theory with the Bose-condensate. Note that if the excitonic packet moves in a crystal (or another semiconductor structure), it interacts with thermal phonons, noncondensed excitons, impurities and other imperfections of the lattice, etc.. Then, the coherent core of the packet can be found in a quasi-stable state, and the ﬂuctuations of φo (x − vt) and, especially, ϕc (x, t) can be of a great importance for possible experimental veriﬁcations of their existence. 2 Exciton-Phonon Condensate To obtain the necessary density of excitons nx in the excitonic cloud and, thus, meet the BEC conditions, the crystals are irradiated by laser pulses with h ¯ ωL ≫ Egap , and the temperature of the crystal is T ≃ 1 ∼ 5 K. If the cross-section area S of an excitation spot on the surface of the crystal can be made large enough, such as S ≈ Ssurf , the hot droplet of paraexcitons can acquire an average momentum during its thermolization process (T ∗ (t) → T ). Indeed, the phonon wind, or the ﬂow of nonequilibrium phonons, blows unidirectionally from the surface into the bulk [5] and transfers the nonzero momentum to the excitonic cloud, Px ≃ Nx ¯ h k0 = 0 and Px ⊥ Ssurf , see Fig. 1. As a result, the packet of moving excitons and nonequilibrium phonons of the phonon wind (Nph ≃ Nx ) is actually the system that undergoes the transition toward developing the Bose-Einstein correlations at T ∗ < Tc . Let us assume that the condensate has been already formed inside the moving excitonic droplet, and the following representation of the exciton Bose-ﬁeld holds: ˆ = Ψ0 + δ ψ ˆ. For the displacement ﬁeld of the crystal ψ ˆ , we introduce a nontrivial coherent part too, i.e., u ˆ= u ˆ , u0 = 0. The important property of the excitonu0 + δ u phonon condensate is a kind of the self-consistency condition, roughly, ∂ u0 (x − vt) ∝ |Ψ0 (x − vt)|2 . (6) In these terms, the moving packet contains both the macroscopically occupied exciton-phonon condensate, or

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

a r X i v :c o n d -m a t /0308585v 1 [c o n d -m a t .s o f t ] 27 A u g 2003Quasi-one-dimensional Bose gases with large scattering lengthG.E.Astrakharchik (1,2),D.Blume (3),S.Giorgini (1),and B.E.Granger (4)(1)Dipartimento di Fisica,Universit`a di Trento and BEC-INFM,I-38050Povo,Italy(2)Institute of Spectroscopy,142190Troitsk,Moscow region,Russia(3)Department of Physics,Washington State University,Pullman,WA 99164-2814,USA(4)Institute for Theoretical Atomic,Molecular and Optical Physics,Harvard-Smithsonian CFA,Cambridge,MA 02138,USA(Dated:February 2,2008)Bose gases conﬁned in highly-elongated harmonic traps are investigated over a wide range of interaction strengths using quantum Monte Carlo techniques.We ﬁnd that the properties of a Bose gas under tight transverse conﬁnement are well reproduced by a 1d model Hamiltonian with contact interactions.We point out the existence of a unitary regime,where the properties of the quasi-1d Bose gas become independent of the actual value of the 3d scattering length a 3d .In this unitary regime,the energy of the system is well described by a hard rod equation of state.We investigate the stability of quasi-1d Bose gases with positive and negative a 3d .PACS numbers:03.75.FiIn recent years the study of quasi-1d quantum Bose gases has attracted a great deal of interest.Intriguing properties of quasi-1d gases,such as the exact mapping between interacting bosons and non-interacting fermions,have been predicted [1,2,3].A bosonic gas that behaves as if it consisted of spinless fermions,a so-called Tonks-Girardeau (TG)gas,cannot be described within mean-ﬁeld theory since it exhibits strong correlations;instead,a many-body framework is called for.While experimen-tal evidence of quasi-1d behavior has been reported for bosonic atomic gases under highly-elongated harmonic conﬁnement [4],TG gases have not been observed yet.It has been suggested,however,that TG gases can be realized experimentally for either low atomic densities or strong atom-atom interaction strengths.The 3d s -wave scattering length a 3d ,and hence the strength of atom-atom interactions,can be tuned to essentially any value,including zero and ±∞,by utilizing a magnetic atom-atom Feshbach resonance [5,6].Utilizing a two-body Feshbach resonance,3d degener-ate gases with large scattering length a 3d have been stud-ied experimentally and theoretically.For a 3d →±∞,it is predicted that the behavior of the strongly-correlated gas is independent of a 3d [7].For homogeneous 3d Bose gases,this unitary regime can most likely not be reached ex-perimentally since three-body recombination is expected to set in when a 3d becomes comparable to the average interparticle distance.Three-body recombination leads to cluster formation,and hence “destroys”the gas-like state.The situation is diﬀerent for Fermi gases,for which the unitary regime has already been achieved experimen-tally [6].In this case,the Fermi pressure stabilizes the system even for large |a 3d |.It has been predicted that three-body recombination processes are suppressed for strongly interacting 1d Bose gases [8].These studies raise the question whether a highly-elongated inhomoge-neous Bose gas,that is,an inhomogeneous quasi-1d Bose gas,is stable as a 3d →±∞.This Letter investigates the properties of a quasi-1d Bose gas at zero temperature over a wide range of in-teraction strengths within a microscopic,highly accurate many-body framework.We ﬁnd that the system i)is well described by a 1d model Hamiltonian with contact inter-actions and renormalized coupling constant [2]for any value of the 3d scattering length a 3d ;ii)behaves like a TG gas for a critical positive value of a 3d ;iii)reaches a unitary regime for large values of |a 3d |,where the proper-ties of the quasi-1d Bose gas become independent of the actual value of a 3d and are similar to those of a hard-rod gas;and iv)becomes unstable against cluster formation for a critical negative value of a 3d .Our study is based on the 3d Hamiltonian H 3d ,H 3d =N i =1−¯h 22ω2ρρ2i +ω2z z 2i +i<jV (r ij ),(1)which describes N spin-polarized mass m bosons un-der highly-elongated conﬁnement with ωz =λωρ,whereλ≪1.The coordinates ρi =¯h/mωρ.To simulate the behavior of a 3d near a ﬁeld-dependent Feshbach resonance,we vary the well depth V 0and consequently the scattering length a 3d .FIG.1:g 1d [dashed line,Eq.(3)]and a 1d [solid line,Eq.(4)]as a function of a 3d .A vertical arrow indicates the value of a 3d where g 1d diverges,a c 3d /a ρ=0.9684.Horizontal arrows indicate the asymptotic values of g 1d and a 1d ,respectively,as a 3d →±∞(g 1d =−1.9368a ρ¯h ωρand a 1d =1.0326a ρ).Inset:a 3d as a function of the well depth V 0for V SR (r ).Importantly,a 3d diverges for particular values of V 0(the inset of Fig.1shows one such divergence).At each di-vergence,a new two-body s -wave bound state is pulled in.Our numerical calculations (see below)are performed for the well depths V 0shown in the inset of Fig.1.We consider situations where the bosonic gas described by Eq.(1)is in the 1d regime for any value of the 3d scattering length a 3d ,which implies Nλ≪1.If the range of V (r )is much smaller than a ρ,it is predicted [2]that the properties of the 3d gas are well described by the 1d contact-interaction Hamiltonian H 1d ,H 1d =N i =1−¯h 2∂z 2i+mma 2ρ12=1.0326.Alternatively,g 1d canbe expressed through the eﬀective 1d scattering length a 1d ,g 1d =−2¯h 2/(ma 1d )[2],wherea 1d =−a ρa ρFIG.2:3d DMC energy per particle calculated using V HS(diamonds)and V SR(asterisks),respectively,together with1d DMC energy per particle calculated using g1d[squares,Eq.(3)]and g01d(pluses),respectively,as a function of a3d[(a)linear scale;(b)logarithmic scale]for N=5andλ=0.01.The statistical uncertainty of the DMC energies is smallerthan the symbol size.Dotted and solid lines show the1denergy per particle calculated within the LDA for g01d(usingthe LL equation of state)and for g1d,Eq.(3)(using the LLequation of state for g1d>0,and the HR equation of statefor g1d<0),respectively.A dashed horizontal line indicatesthe TG energy,and a vertical arrow the position where g1d,Eq.(3),diverges.V SR(asterisks),respectively.For small a3d/aρ,the en-ergies for these two two-body potentials agree within thestatistical uncertainty.For a3d>∼aρ,however,clear dis-crepancies are visible.The DMC energies for V SR crossthe TG energy per particle(indicated by a dashed hor-izontal line),E/N−¯hωρ=¯hωρλN/2,very close to thevalue a c3d=0.9684aρ[indicated by a vertical arrow inFig.2(b)],while the energies for V HS cross the TG en-ergy per particle at a notably smaller value of a3d.To compare our results obtained for the3d Hamilto-nian,H3d,with those for the1d Hamiltonian,H1d,wealso solve the Schr¨o dinger equation for H1d,Eq.(2).Forpositive values of the coupling constant g1d we calculatethe many-body ground state energy by the exact DMCmethod.For g1d<0,however,the1d Hamiltonian sup-ports many-body bound states,and,as in the3d case,we use the FN-DMC method to describe the lowest-lyinggas-like state.For N=2and g1d<0,theﬁrst excitedstate of the Schr¨o dinger equation for H1d,Eq.(2),has anode at z12=a1d.To solve the many-body Schr¨o dingerequation for H1d with negative g1d for the lowest gas-likestate by the FN-DMC technique,we parameterize ourmany-body nodal surface in terms of a1d.This many-body nodal surface is expected to be good if the densityof the gas is low.In addition to the3d energy per particle,Fig.2showsthe resulting1d energy per particle obtained by solvingthe Schr¨o dinger equation for H1d,Eq.(2),for the renor-malized coupling constant g1d[squares,Eq.(3)],and theunrenormalized coupling constant g01d(pluses),respec-tively.The1d energies calculated using the two diﬀerentcoupling constants agree well for small a3d,while cleardiscrepancies become apparent for larger a3d.Impor-tantly,the1d energies calculated using the renormalized1d coupling constant g1d agree well with the3d energiescalculated using the short-range potential V SR(aster-isks)up to very large values of the3d scattering lengtha3d,and also for negative a3d.In contrast,at large a3dthe1d energies deviate clearly from the3d energies calcu-lated using the hard-sphere potential V HS(diamonds).We conclude that the renormalization of the eﬀective1dcoupling constant g1d and the1d scattering length a1dare crucial to reproduce the results of the3d Hamilto-nian H3d when a3d>∼aρand when a3d negative.Smalldeviations between the1d energies calculated using therenormalized1d coupling constant g1d and the3d ener-gies calculated using the short-range potential V SR re-main;we attribute these to theﬁnite range of V SR.If the size of the cloud is much larger than the har-monic oscillator length a z,where a z=2 1+128√45π2√aρ+··· .(5)Theﬁrst term corresponds to the energy per particle inthe TG regime.Similarly,the linear density in the centerof the cloud is to lowest order given by the TG result,n1d=√4 ticle for the LL equation of state(g1d>0)as well asthe HR equation of state(g1d<0)calculated withinthe LDA.Remarkably,the LDA energies nearly coin-cide with the1d many-body DMC energies(pluses andsquares,respectively);ﬁnite-size eﬀects play a role onlyfor a3d≪aρ.Our calculations establish for theﬁrsttime that a simple treatment,i.e.,a HR equation ofstate treated within the LDA,describes trapped quasi-1dgases with negative coupling constant g1d well over a widerange of the3d scattering length a3d.For a3d→−0,that is,for large a1d,the HR equation of state using the LDA, cannot properly describe trapped quasi-1d gases,which are expected to become unstable against formation of cluster-like many-body bound states for a1d≈1/n1d.We hence investigate the regime with negative a3d in more detail within a many-body framework.By comparing with results for the3d Hamiltonian, Eq.(1),we have shown above that the1d HamiltonianH1d,Eq.(2),provides an excellent description of quasi-1d gases.We hence base our stability analysis of quasi-1d gases with large eﬀective1d scattering length a1d onH1d.We solve the many-body Schr¨o dinger equation for H1d using the variational quantum Monte Carlo(VMC)method.Our variational N-particle Bijl-Jastrow-typewave function consists of one-and two-body terms.The one-body terms are written as a function of a single vari-ational parameterα,which determines the size of the atomic gas.The two-body term is parameterized by a1d,and explicitly accounts for correlations.Figure3shows the resulting VMC energy per parti-cle for N=5andλ=0.01as a function of the vari-ational parameterα(Gaussian width)for four diﬀerenta1d.For a1d/aρ=1and2,Fig.3shows a local minimumatαmin≈a z.The minimum VMC energy nearly coin-cides with the essentially exact DMC energy,which sug-gests that our variational wave function provides a highlyaccurate description of the quasi-1d many-body system. The energy barrier atα≈0.2a z decreases with increas-ing a1d,and disappears for a1d/aρ≈3.We interpret this vanishing of the energy barrier as an indication of insta-bility[15].For small a1d,the energy barrier separates the gas-like state from cluster-like bound states.For larger a1d,this energy barrier disappears and the gas-like state becomes unstable against cluster formation.We additionally performed variational calculations forlarger N and diﬀerentλ.Weﬁnd that the onset of in-stability of quasi-1d Bose gases can be described by the product of the1d scattering length a1d and the linear density at the trap center n1d.To be speciﬁc,our many-body calculations suggest that a quasi-1d gas is stable for a1d n1d<∼0.35,and becomes unstable for a1d n1d>∼0.35. Our analysis suggests that the quasi-1d unitary regime can be reached experimentally.By tuning the3d scat-tering length,it is further possible to investigate the on-set of instability.By reducingλone should be able to stabilize relatively large quasi-1d systems.FIG.3:VMC energy per particle as a function of the varia-tional parameterαfor N=5,λ=0.01,and a1d/aρ=1.0326 (pluses),2(asterisks),3(diamonds)and4(triangles).In conclusion,we investigated the energetics of a Bose gas under highly-elongated harmonic conﬁnement over a wide range of the3d scattering length.Weﬁnd that thequasi-1d gas can be described by a many-body1d model Hamiltonian with zero-range interactions and renormal-ized coupling constant.For a3d→±∞,the quasi-1d gas enters a unitary regime,where all properties of the sys-tem are independent of a3d.In the vicinity of the unitaryregime,the quasi-1d system behaves like a gas of HRs. For negative a3d,quasi-1d gases become unstable against cluster formation for a critical value of a1d n1d.GEA and SG acknowledge support by the Ministero dell’Istruzione,dell’Universit`a e della Ricerca(MIUR). DB acknowledges support by the NSF(grant0331529) and by the BEC Center at the University of Trento,and BEG by the NSF through a grant to ITAMP.[1]M.Girardeau,J.Math.Phys.1,516(1960).[2]M.Olshanii,Phys.Rev.Lett.81,938(1998).[3]D.S.Petrov,G.V.Shlyapnikov and J.T.M.Walraven,Phys.Rev.Lett.85,3745(2000).[4]A.G¨o rlitz et al.,Phys.Rev.Lett.87,130402(2001);F.Schreck et al.,ibid.87,080403(2001);M.Greiner et al., ibid.87,160405(2001).[5]S.Inouye et al.,Nature392,151(1998);S.L.Cornishet al.,Phys.Rev.Lett.85,1795(2000);T Loftus et al., Phys.Rev.Lett.89,203201(2002).[6]K.M.O’Hara et al.,Science298,2179(2002);T.Bourdelet al.,Phys.Rev.Lett.91,020402(2003).[7]S.Cowell et al.,Phys.Rev.Lett.88,210403(2002);H.Heiselberg,Phys.Rev.A63,043606(2001).[8]D.M.Gangardt and G.V.Shlyapnikov,Phys.Rev.Lett.90,010401(2003).[9]E.H.Lieb and W.Liniger,Phys.Rev.130,1605(1963);E.H.Lieb,ibid.130,1616(1963).[10]J.B.McGuire,J.Math.Phys.5,622(1964).[11]D.Blume,Phys.Rev.A66,053613(2002);G.E.As-5trakharchik and S.Giorgini ibid.66,053614(2002). [12]J.Boronat and J.Casulleras,Phys.Rev.B49,8920(1994).[13]P.J.Reynolds et al.,J.Chem.Phys.77,5593(1982).[14]V.Dunjko et al.,Phys.Rev.Lett.86,5413(2001).[15]J.L.Bohn,B.D.Esry and C.H.Greene,Phys.Rev.A58,584(1998);H.T.C.Stoof,J.Stat.Phys.87,1353(1997).。