Critical temperature of interacting Bose gases in two and three dimensions

0.46 MB·English

by  S. Pilati

Critical temperature of interacting Bose gases in two and three dimensions.
S. Pilati and S. Giorgini
Dipartimento di Fisica, Universit`a di Trento and CNR-INFM BEC Center, I-38050 Povo, Trento, Italy
N. Prokof’ev
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
Theoretische Physik, ETH Zu¨rich, CH-8093 Zu¨rich, Switzerland and
Russian Research Center “Kurchatov Institute”, 123182 Moscow, Russia
8 We calculate the superfluid transition temperature of homogeneous interacting Bose gases in
0 threeand two spatial dimensions using large-scale Path Integral Monte Carlo simulations (with up
0 to N = 105 particles). In 3D we investigate the limits of the universal critical behavior in terms
2 of the scattering length alone by using different models for the interatomic potential. We find that
this type of universality sets in at small values of the gas parameter na3 . 10−4. This value is
n different from the estimate na3 . 10−6 for the validity of the asymptotic expansion in the limit of
Ja vanishingna3. In2DwestudytheBerezinskii-Kosterlitz-Thoulesstransitionofagaswithhard-core
interactions. Forthissystemwefindgoodagreementwiththeclassicallattice|ψ|4 modeluptovery
2
large densities. We also explain the origin of the existing discrepancy between previous studies of
2
thesame problem.

]
r PACSnumbers:
e
h
ot The theoreticaldeterminationofthe superfluidtransi- calculate the higher-order logarithmic correction to (1).
. tion temperature in homogeneous, interacting Bose sys- Continuous-space studies of a gas of hard spheres,
t
a tems is a fine example of a many-body problem that can based on the conventional PIMC algorithm [14], were
m
be quantitatively addressed only by “exact” numerical carried out in Refs. [7, 13]. Both calculations suffered
- techniques. This fact is well understood in the case of from two shortcomings: first, the number of particles in
d 4
strongly interacting quantum fluids, such as liquid He thesimulationwasonlyfewhundredsmakingtheextrap-
n
o where the critical temperature in bulk [1] as well as in olation to the thermodynamic limit difficult; second, the
c two-dimensional configurations [2, 3] was calculated us- algorithmisknowntobeinefficientforsimulationsofthe
[ ing path integral Monte Carlo (PIMC) simulations, but superfluid density. The results of Ref. [13] for T agree
c
1 at first glance is surprising in the case of dilute gases. with the asymptotic law (1) at small densities, in con-
v However, in three dimensions (3D) the presence of any trasttothe significantlylowervaluesobtainedinRef.[7]
4 finite interaction changes the universality class of the (see Fig. 1). Moreover, there is strong (about ten stan-
2
transition from the Gaussian complex-field model, cor- dard deviations!) discrepancy between Refs. [7] and [13]
4
responding to the ideal gas Bose-Einstein condensation athigher densities calling for further investigationof the
3
1. t(hBeEcCr)itticeamlpteemraptuerraetTucr0e, tTo tchaantnooftthbeeXobYtaminoeddeflr.oTmhuTs0, proInble2mD. the BKT transition temperature of a weakly
0 c c
perturbatively[4]. Intwodimensions(2D)thesuperfluid interacting gas is written in the form [15]
8
0 transition, which belongs to the Berezinskii-Kosterlitz- 2π~2n2D 1
v: Tinhteorualcetsison(BeKffTec)tsunainvdertshaelirteyicslansos [u5n,p6e]r,tuisrbinedduccreidticbayl TBKT = mkB ln(ξ/4π)+lnln(1/n2Da22D) , (2)
i
X temperature to start with. where n2Da22D is the corresponding gas parameter, and
ar In a 3D weakly repulsive gas the critical temperature a2Disthe2Dscatteringlength. Thenumericalcoefficient
shift is fixed by the s-wave scattering length a (a > 0), ξ was calculated in Ref. [16] from lattice Monte Carlo
which characterizes interatomic interactions at low tem- simulations of the classical ψ 4 model, similarly to the
| |
perature [4, 7, 8, 9, 10, 11, 12, 13], 3D case, yielding the value ξ =380 3.
±
An important question concerns the range of applica-
T =T0 1+c(an1/3) . (1) bilityofEqs.(1)-(2)sincetheywerederivedforthelimit
c c h i
of vanishingly small interaction strength. More broadly,
Here n is the gas density and T0 = one is interested in knowing up to what value of the gas
c
(2π~2/mk )[n/ζ(3/2)]2/3 with m the particle mass parameter it is possible to express T as a function of
B c
and ζ(3/2) 2.612. The numerical coefficient c in na3 alone and ignore the dependence on the interaction
≃
Eq. (1) was calculated in Refs. [9, 10] by solving the potential details. These questions are particularly rele-
effective 3D classical ψ 4 model using lattice Monte vant for the 2D case where experimental determinations
| |
Carlosimulations. The reportedvalue is c=1.29 0.05. ofthecriticaltemperatureintrappedconfigurationshave
±
The same classical model was employed in Ref. [11] to become available [17, 18].
2
InthisLetterwereportPIMCresultsforthesuperfluid 1.08
PRA04
transition temperature in interacting Bose gases both in PRL97
1.07 Our results: HS
3D and in 2D. We carry out large-scale simulations of Our results: SS
homogeneoussystemswithuptoN =105 particles. The 1.06
1.05
simulations are based on the worm algorithm, recently
extended to continous-space systems, which allows for a 0C 1.04
T
reliable and efficient calculation of the superfluid den- /C 1.03
T
sity [3]. With new data the extrapolation of the critical 1.02
temperature to the thermodynamic limit can be done
1.01
with the level of accuracy that was unreachable in pre-
1
vious attempts. In 3D we determine the dependence of
T on the gas parameter na3, reporting good agreement 0.99
c
with the expansion (1) in the dilute regime and signif- 10-7 10-6 10-5 10-4 10-3 10-2 10-1
na3
icant deviations from previous studies [7, 13] at higher
densities. Furthermore, we carry out simulations with
both a hard- and a soft-sphere interatomic potential to FIG. 1: (color online). Critical temperature in 3D as a func-
investigate the universal behavior in terms of the scat- tionofthegasparameterna3. ThesymbolslabeledbyPRA04
tering length. In 2D we calculate TBKT for a hard-disk correspond to the results of Ref. [13], the ones labeled by
PRL97correspond toRef.[7]. Thedashed line(green) isthe
gas as a function of the interaction strength, finding re-
expansion (1) of Ref. [9] and the dotted line (black) is the
sultsinexcellentagreementwiththeprediction(2)upto
expansion of Ref. [11] including logarithmic corrections. The
regimes of surprisingly high density.
solid line (red) is a guide to theeye.
We consider a 3D system of N particles with peri-
odic boundary conditions described by the Hamiltonian
1.4
H = (~2/2m) N 2 + V(r r ), where r
denote−s the coorPdiin=a1te∇siof tPhei<i-jth p|airt−iclej.| Two-bodyi 1.3 NN==2150624
N=4096
interactionsaremodeledby the followingpotentials: the N=16386
1.2
hard-sphere (HS) potential, VHS(r) = + if r < a
and zero otherwise, and the soft-sphere (S∞S) potential, ρ/ 1.1
S
VSS(r) = V0 (V0 > 0) if r < R0 and zero otherwise. In 3ρ 1
theHScasethes-wavescatteringlengthacoincideswith 1/N
therangeofthepotential,whileintheSScaseitisgiven 0.9
by a = R0[1 tanh(K0R0)/K0R0] with K0 = √V0m/~. 0.8
−
We always use the range R0 = 5a and adjust V0 to ob-
tain the desired value of a. We notice that the HS and 0.7
SS model represent two extreme cases of repulsive po- 0.99 0.995 1 1.005 1.01
tentials: in the HS case the energy is entirely kinetic, T/T0C
while for the SS case the Born approximation result
aB = (m/~2)R0∞V(r)r2dr accounts for more than 80% FIG. 2: (color online). Results for a 3D non-interacting gas.
of the value of the scattering length. Scaled superfluid fraction as a function of temperature for
different system sizes. The lines are linear fits from Eq. (3).
InaPIMCsimulationoneobtainsaveragesofphysical
quantities over a set of stochastically generated configu-
rations R = (r1,...,rN) sampled from a probability dis- ticle number N using the scaling ansatz[21]
tributionproportionaltothedensitymatrixρ(R,R,β)=
Re−βH R , where β = 1/k T is the inverse temper- N1/3ρ (t,N)/ρ=f(tN1/3ν)=f(0)+f′(0)tN1/3ν +....
B s
h | | i
ature. The superfluid density ρ is obtained from the (3)
s
winding number estimator [1], which accounts for long Here, t = (T T )/T is the reduced temperature, ν is
c c
permutation cycles of identicalparticles occurring in the thecriticalexp−onentofthecorrelationlengthξ(t) t−ν,
system. The calculation of ρ(R,R,β) is based on the andf(x)isauniversalanalyticfunctionwhichallo∼wsfor
pair-productdecomposition,wherethetwo-bodydensity a linear expansion around x=0.
matrix associated with the relative motion of the pair is We first consider the non-interacting case (see Fig. 2).
determined exactly both for the HS and the SS poten- The scaling curves all intersect at the same value of the
tial [19, 20]. The critical temperature T is determined reduced temperature according to Eq. (3) with T /T0 =
c c c
fromcalculationsofthesuperfluidfractionρ /ρ(ρ=mn 1.0005(4) and ν = 0.96(3). The value of ν is consistent
s
isthetotalmassdensity)forsystemswithincreasingpar- withthepredictionν =1ofthecomplexGaussianmodel.
3
2 1.072 TABLEI:TransitiontemperatureTc ofa3Dhard-spheregas
1.071 fordifferentvaluesofthegasparameter. Theresultswiththe
1.07
ρ/ 11..68 0T/TCC 1111....000066666789 lab1el×SnSa103c−o7rrespond1.0tTo0c6/a9T(sc10o0ft)-sphere2g×ansa1.03−3 1.T06c/2T4(c04)
S 1.065
3ρ 1.4 1.0x104 3.0x104 5.0x104 7.0x104 9.0x104 5×10−7 1.0091(7) (SS)2×10−3 1.0277(4)
1/N N 1×10−6 1.0127(7) 5×10−3 1.0652(4)
1.2 N=64 1×10−5 1.0214(7) 1×10−2 1.0627(4)
1 NN==2150624 1×10−4 1.0351(7) (SS)1×10−2 0.9880(5)
NN==420199652 (SS) 1×10−4 1.0359(12) 5×10−2 1.0060(5)
N=100000
0.8
1.058 1.06 1.062 1.064 1.066 1.068 1.07 1.072 1.074
T/T0
C
1
FIG.3: (coloronline). Resultsfora3Dhard-spheregaswith
N=64
na3 = 5×10−3. Scaled superfluid fraction as a function of N=256
N=1024
temperature for different system sizes. The lines are linear 0.8 N=4096
fits from Eq. (3). The inset shows the N dependence of the
intersectionpointbetweenlinescorrespondingtopairsofcon-
0.6
secutive system sizes. 1
ρ
/S 0.9
ρ
0.8
0.4
ρ 0.7
The correspondingresults for the hard-spheregas at the ρ/S 0.6 T/T*=0.167 TBKT
density na3 = 5 10−3 are reported in Fig. 3. In this 0.2 0.5 TT//TT**==00..2214
× 0.4 T/T*=0.231
casetheintersectionpointbetweencurvesforconsecutive 1×103 2×103 3×103 4×103
N
systemsizesclearlydriftstowardslowertemperaturesfor 0
0.1 0.15 0.2 0.25
largerN (see insetofFig.3). Thiseffectarisesfromcor- T/T*
rectionsto the scalinglaw[the righthand side ofEq.(3)
hastobe multiplied by(1+CN−ω/3+...)withω 0.8 FIG. 4: (color online). Results for a 2D hard-disk gas with
and C of order unity]. A reliable extrapolation t≈o the na2 = 0.01. Superfluid density as a function of tempera-
thermodynamic limit using Eq. (3) requires large-scale ture for different system sizes. The solid line (red) shows
systems on order of N & 104. The temperature T and the extrapolation to the thermodynamic limit. The dashed
c
line (black) corresponds to the BKT universal jump. In the
the exponentν aredeterminedfromEq.(3)by consider-
inset we show the dependence of the superfluid density on
ing series of data corresponding to systems with N suf-
the system size for different temperatures. The dotted lines
ficiently large. We find T /T0 = 1.0652(4) (see Table I)
c c are fits using the Kosterlitz-Thouless renormalization group
and ν = 0.70(4), in agreement with the value ν = 0.67 equations.
corresponding to the universality class of the XY model
in three dimensions [22]. By using the same procedure
3
we calculate T as a function of the gas parameter na
c
both for the HS and the SS potential. The results are and the flow of intersection points with the system size
shown in Fig. 1 and are reported in Table I. (which is substantial between N =256 and N 20000).
∼
Except for the largestdensities, our results for the HS ThediscrepancywithRef.[13]reducesatverysmallden-
gasaresystematicallyhigherthantheonesofRef.[7]and sities na3 . 10−6, where finite-size corrections to T ap-
c
lower compared to the ones of Ref. [13] [23]. With much pear to be less relevant and there is agreement with the
better accuracy for ρ and larger system sizes we are in expansion(1). ThecriticaltemperatureT firstincreases
s c
the position to explain the origin of discrepancy. It is with na3, then goes through a maximum and for larger
two-fold. First,thenumberofimaginarytimeslicesused valuesofthegasparameterdecreasesbelowtheT0value.
c
inRef.[13]was15,afactorof3largerthaninRef.[7]. We For example, T /T0 = 0.70 in liquid 4He corresponding
c c
find,bydoingsimulationswithupto96slices,thatabout to the effective gas parameter na3 0.21 [24]. The HS
≃
25 slices have to be used to ensure that the correspond- gas is finally expected to become a solid at the freezing
3
ing systematic errorsarenegligible. Second, Refs. [7, 13] density na 0.25 [24]. Concerning the comparison be-
underestimatederrorbarsbyafactoroftenbecausemul- tween differe≃nt model potentials, at na3 = 10−4 we find
tiple intersections between scaling curves (i) render the very good agreement between the HS and SS gas, while
procedure of locating the intersection point ambiguous, for higher densities deviations start to become evident
and(ii)preventonefromdetectingcorrectionstoscaling and the SS results are significantly smaller.
4
0.26 Inconclusion,wehavecarriedoutanumericalstudyof
thesuperfluidtransitiontemperatureofinteractingBose
0.24 gasesin 3D and2D andestablishedthe limits ofvalidity
of the asymptotic expansions (1) and (2), and the uni-
versal description in terms of the scattering length. We
0.22
* havealsoexplainedandresolvedthediscrepancybetween
T
/C previous studies.
T
0.2
0.18
[1] E.L. Pollock and D.M. Ceperley, Phys. Rev. B 36, 8343
0.16 10-4 10-3 10-2 10-1 (1987).
n a2 [2] D.M. Ceperley and E.L. Pollock, Phys. Rev. B 39, 2084
2D 2D
(1989).
[3] M. Boninsegni, N. Prokof’ev and B. Svistunov, Phys.
FIG. 5: (color online). Critical temperature of a 2D gas of Rev. Lett.96, 070601 (2006).
hard disks as a function of thegas parameter na2. The solid [4] G. Baym et al., Phys. Rev. Lett. 83, 1703 (1999); Eur.
line is the result (2). Phys. J. B 24, 104 (2001).
[5] V.L. Berezinskii, Sov.Phys. JETP 34, 610 (1972).
[6] J.M. Kosterlitz and D.J. Thouless, J. Phys. C 6, 1181
(1973); J.M. Kosterlitz, J. Phys. C 7, 1046 (1974).
The simulations of the 2D Bose gas are carried out [7] P. Gru¨ter, D. Ceperley and F. Lalo¨e, Phys. Rev. Lett.
in a similar way, using a hard-core potential with range 79, 3549 (1997).
[8] M.HolzmannandW.Krauth,Phys.Rev.Lett.83,2687
a2D (hard-diskpotential)tomodelthe interactions. The
(1999).
resultsforthesuperfluidfractionasafunctionoftemper-
[9] V.A. Kashurnikov, N.V. Prokof’ev and B.V. Svistunov,
ature are reported in Fig. 4 for the value n2Da22D =0.01 Phys. Rev.Lett. 87, 120402 (2001).
of the gas parameter. The extrapolation to the thermo- [10] P. Arnold and G. Moore, Phys. Rev. Lett. 87, 120401
dynamic limit is carriedout by fitting numerical data at (2001).
finite N to the solutionof the Kosterlitz-Thoulessrecur- [11] P. Arnold, G. Moore and B. Tom´aˇsik, Phys. Rev. A 65,
sion relations[6] 013606 (2001).
[12] B. Kastening, Phys. Rev.A 68, 061601(R) (2003).
df(N)/dlnN = y2(N)f2(N)/2 [13] K. Nho and D.P. Landau, Phys. Rev. A 70, 053614
− (2004).
dy(N)/dlnN = [1 f(N)]y(N), (4) [14] D.M. Ceperley, Rev.Mod. Phys. 67, 1601 (1995).
−
[15] V.N. Popov, Functional Integrals and Collective Excita-
wheref(N)=(~2πn2D/2mkBT)ρs(N)/ρisadimension- tions (Cambridge University Press, Cambridge, 1987);
less function proportional to the superfluid fraction and D.S. Fisher and P.C. Hohenberg, Phys.Rev. B 37, 4936
y(N) is proportional to the vortex fugacity. The start- (1988).
ing values f0 and y0 of the recursion relations (4) are [16] N. Prokof’ev, O. Ruebenacker and B. Svistunov, Phys.
determined from a best fit to the results correspond- Rev. Lett.87, 270402 (2001).
[17] Z.Hadzibabicet al.,Nature(London)441, 1118(2006).
ing to different values of N and system temperature
[18] P. Kru¨ger, Z. Hadzibabic and J. Dalibard, Phys. Rev.
in close proximity of the transition (see Fig. 4). From Lett. 99, 040402 (2007).
these initial values one determines the critical temper-
[19] S. Pilati, K. Sakkos, J. Boronat, J. Casulleras and S.
ature TBKT in the thermodynamic limit. Here TBKT is Giorgini, Phys. Rev.A 74, 043621 (2006).
written in units of T∗ = 2π~2n2D/(2mkB), which pro- [20] Formoredetailsontheimplementationofthewormalgo-
vides the natural temperature scale for quantum degen- rithmusingthepairproductdecompositionseeS.Pilati,
eracy in 2D. In Fig. 4 we also show the prediction of PhD Thesis (Trento), unpublished.
[21] E.L. Pollock and K.J. Runge, Phys. Rev. B 46, 3535
the universaljump ofthe superfluid fractionat the tran-
sition ρs/ρ = 2mkBTBKT/(π~2n2D) [25], in nice agree- [22] J(1.C99.2L)e. Guillou and J.Zinn-Justin,Phys.Rev.Lett.39,
ment with the temperature dependence of our extrapo-
95 (1977).
latedcurve. TheresultsforTBKTasafunctionofn2Da22D [23] We notice that by using small system sizes similar to
are reported in Fig. 5, where they are compared with the ones of Ref. [13], the value of Tc increases becoming
the prediction obtained from Eq. (2). The agreement is consistent with the findingsof that paper.
surprisingly good up to very large values of the gas pa- [24] M.H. Kalos, D.LevesqueandL. Verlet,Phys.Rev.A 9,
risamcleotseer.toAthtethfreeehziignhgesptoidnetnnsi2tDyan222DD≃a22D0.3=2 o0f.1a, gwahsicohf [25] D2112.70R81.((N1199el77s47o))n.. and J.M. Kosterlitz, Phys. Rev. Lett. 39,
harddisks[26],the dilutegasresult(2)isonlyabout7% [26] L. Xing, Phys. Rev.B 42, 8426 (1990).
larger than the PIMC value.