Thermal Fluctuations of Vortex Matter in Trapped Bose–Einstein Condensates
arXiv:cond-mat/0604416v3 [cond-mat.stat-mech] 26 Oct 2006
S. Kragset1, E. Babaev2,1, and A. Sudb?1
Department of Physics, Norwegian University of Science and Technology, N-7491 Trondheim, Norway 2 Laboratory of Atomic and Solid State Physics, Cornell University, Ithaca, NY 14853-2501, USA We perform Monte Carlo studies of vortices in three dimensions in a cylindrical con?nement, with uniform and nonuniform density. The former is relevant to rotating 4 He, the latter is relevant to a rotating trapped Bose–Einstein condensate. In the former case we ?nd dominant angular thermal vortex ?uctuations close to the cylinder wall. For the latter case, a novel e?ect is that at low temperatures the vortex solid close to the center of the trap crosses directly over to a tension-less vortex tangle near the edge of the trap. At higher temperatures an intermediate tensionful vortex liquid located between the vortex solid and the vortex tangle, may exist.
Topological vortex excitations are hallmarks of quantum ordered states such as superconductors and super?uids. Vortices are important at vastly di?erent length scales ranging from the dynamics of neutron stars to transport properties in superconductors and rotational response of ultracold atomic gases in optical traps [1, 2]. Particularly rich physics is associated with various orderings of quantum vortices, and the transitions between them. In certain systems, even numerous aggregate states of vortex matter are possible . The simplest, yet extremely important, temperature-induced phase transition in vortex matter in a one-component superconductor or super?uid is the one from an ordered vortex line lattice (VLL) to a disordered vortex liquid (VL) state. This has been a subject of intense research in the context of high-Tc superconductivity. Remarkable progress has recently been made in the physics of vortices in Bose-Einstein condensates (BEC) of ultra-cold atoms [1, 2]. These systems are extremely clean and have physical parameters controllable in a wide range. This has led to many suggestions for testing a number of general physical concepts in BECs. Since individual vortices, as well as their ordering patterns, can be resolved in ultra-cold gases [1, 2], a natural question arises if in these systems a thermally induced phase transition from VLL to VL can be observed, and if it can shed more light on fundamental properties of the VL states in quantum ?uids in general. Detailed theoretical studies of a model uniform system, with density corresponding to the values attainable in the center of a trap, indicate that in present experiments rotation is too slow and/or particle numbers are too high to obtain a VLL melting if the approximation of a uniform system holds . We also mention studies of quantum VLL melting in a one dimensional optical lattice . In a harmonic trap, the condensate density will gradually be depleted from the center of the trap towards its edge. This suggests that thermal ?uctuation e?ects will be enhanced close to the edges. We thus focus on investigating possible crossover states of vortex matter in a spatially inhomogeneous system such as a trapped BEC. Melting of vortex lines in a bosonic condensate can
be modeled by the uniformly frustrated lattice 3D XY model (see e.g. [6, 7, 8]). We study a system which takes into account the presence of a harmonic trap by using the Hamiltonian H=?
Pij cos(θj ? θi ? Aij ),
where θi is the phase of the condensate wave function at position i. The factor Pij ≡ P (rij ) = 1?(rij /R)2 if rij ≤ R, where R is the cloud size and rij is the radial distance from the cloud center. For rij > R we set Pij = 0. We account for circulation by introducing the potential j Aij = i dl · A. Here, ? × A = (0, 0, 2πf ), and f is the number of rotation-induced vortices per plaquette in the xy-plane. We perform Monte Carlo (MC) simulations on Eq. (1) with the Metropolis algorithm. We have considered cubic numerical grids of size L3 , with L = 36, 72, and along the z-direction we impose periodic boundary conditions to model an elongated system. The ?lling fraction is f = 1/36, and temperatures T ∈ (0.30, 3.0) . The temperature at which the Bose condensation takes place at zero rotation in a bulk system, is given by T = 2.2 in our units. We note that a rapidly rotating system should be better described by the Lowest Landau Level (LLL) approximation , characterized for instance by the fact that the intervortex separation is comparable to the vortex-core size. However, several experiments on elongated systems are clearly outside the LLL regime. For instance, in Ref. , the inter-vortex separation is given as 5.0?m, while the healing length (vortex-core radius ξ) is 0.2?m. A parameter estimating the validity of the LLL approximation is the ratio λ of the interaction energy scale to the level spacing of the transverse harmonic con?nement , λ = 4π? 2 as n/(M ? ω⊥ ), where as h h is the s-wave scattering length, M is the particle mass, and ω⊥ is the trap frequency. The LLL approximation requires λ ? 1, while the parameters of Ref.  correspond to at least λ ∈ 1 ? 100. Under such conditions, we expect Eq. (1) to be adequate. (See also caption of Fig. 3).
In the problem of vortex matter in a trap we encounter two speci?c circumstances, namely a ?nite-size situation and an inhomogeneous density pro?le. We begin by examining the consequences of ?nite size, by studying the system in a cylindrical container with a uniform Pij . Such a situation is indeed relevant for the physics of liquid 4 He. Vortex orderings in such a geometry at zero temperature were studied in , however the VLL melting for this case was addressed only for a planar geometry (see e.g. ). Since the melting process in three dimensional vortex matter is very di?erent from that in two dimensions due to the importance of vortex line bending, this problem warrants careful consideration. Fig. 1 shows the results of simulations of vortex matter in a cylindrical container with Pij = 1, rij ≤ R. At low temperatures, the simulations reproduce orderings with circular distortions of VLL near the container wall, as predicted for 4 He in a zero-temperature treatment of the problem . For a large number of vortices the system reacquires the hexagonal lattice symmetry away from the wall, see Fig. 1 (bottom row). Increasing the temperature in the case of small number of vortices (top row of Fig. 1) the dominant vortex ?uctuations are associated with angular displacements, while radially the vortex density remains ordered. For a larger number of vortices (bottom row of Fig. 1) we ?nd dominance of angular ?uctuations only for the vortices situated close to container wall, while the center of the system does not display this phenomenon. The crossover to a uniformly molten vortex system occurs in both cases only at a higher temperature. The two-step thermal crossover in the vortex pattern we ?nd is analogous to that in two dimensions where vortices are point-like objects (see e.g. ). There is, however, a principal di?erence in our case, since in three dimensions the VLL melting is accompanied by signi?cant vortex bending ?uctuations.
T = 0.50 T = 1.00 T = 1.25 T = 1.67
1). In a uniform and in?nite system, the ?uctuations can cause either VLL melting via a ?rst-order phase transition or a second order transition associated with a thermally induced proliferation of closed vortex loops near the critical temperature where the vortices loose their line-tension [8, 13]. In a ?nite-size inhomogeneous system, the situation is di?erent. A density gradient in a trap may e?ectively be viewed as a temperature gradient in a uniform system. It is clear that for low, but ?nite, temperatures there will be a ?nite area near the edge of the cloud which e?ectively would be at a high enough temperature to feature an annulus of tension-less tangle of vortices. This is a phase where the vortex-line tension (free energy per unit length) has vanished through the proliferation of vortex loops . The boundary where this tangle sets in, marks the true boundary of the BEC. An issue to be adressed is whether we encounter an appreciable VL (i.e. tensionful but disordered vortex state) layer between the tension-less vortex tangle and the VLL. Fig. 2 shows snapshots of vortex con?gurations in the model Eq. 1 generated by MC simulations at T = 0.5 and T = 1.0. For better visualization we choose the vortex radius to be 0.4 of the grid spacing so it should not be associated with the core size. The sharp bends at short length scales result from the presence of a numerical grid. Nonetheless, this model has proved to be accurate for describing vortex ?uctuations at scales larger than the grid spacing [6, 7, 8].
T = 0.50
T = 1.00
T = 1.25
T = 1.67
FIG. 2: (Color online) Snapshots of vortex con?gurations in a rotating trapped BEC at T = 0.5 (left ?gures) and at T = 1.0 (right ?gures). The top row shows a selection of 16 out of 72 layers in z direction. The bottom row shows smaller selections in the xy plane, but 32 out of 72 layers in z direction. Fluctuations are minimal in the trap center, and increase towards the edge of the trap. A distinct front separating regions of ordered and disordered vortices is easily identi?ed.
FIG. 1: xy positions of vortices in a cylindrical container integrated over z-direction, and averaged over every tenth of a total of 5 · 105 MC sweeps. Top and bottom rows have L = 36 and L = 72, respectively. At T = 0.5, 1.0, we discern circular ordering close to the cylinder wall combined with a hexagonally ordered state closer to the center. At T = 1.25, 1.67 we observe dominance of angular ?uctuations closest to the edge.
Let us now turn to the case of a harmonic trap (Eq.
We next locate the vortex liquid layer between the ordered vortex state and the tension-less vortex tangle closest to the edge of the trap. To this end, we take snapshots of vortex matter at di?erent temperatures, in each case integrating over the z-direction. In a resulting image, straight vortex lines will be seen as bright spots while bent vortices will be seen as smeared out spots. This may be related to experiments, where at least for
non-equilibrated vortex systems the z-integration renders vortices essentially indistinguishable [14, 15]. The results are given in the upper row of Fig. 3 for di?erent temperatures. There we can identify regions of rather straight and ordered vortex lines and a smeared region. To obT = 0.50 T = 1.00 T = 1.25 T = 1.67
in z direction, de?ned in a selected region between two cylinders of radii R1 and R2 . We do so by applying a twist ?(rij ) ≡ ?ij = ?? if R1 ≤ rij < R2 , z 0 otherwise, (2)
to the model Eq. 1 and de?ning the modi?ed helicity modulus as follows, 1 ? Υz (R1 , R2 ) ≡ ′ N 1 ? T N′
Pij cos(θj ? θi ? Aij ) [Pij sin(θj ? θi ? Aij )]
T = 0.50
T = 1.00
T = 1.25
T = 1.67
T = 0.50
T = 1.00
T = 1.25
T = 1.67
FIG. 3: xy-positions of vortices in a trapped BEC integrated over z-direction. Top row are snapshots, while middle and bottom rows are averages of 105 and 5 · 105 MC sweeps, respectively. Every tenth con?guration has been sampled. This provides information on the stability of the ordered region and the evolution of the disordered region as T varies. Intervortex distance 2r0 = 2 ? /M ? and healing h √ length ξ = ? / 2M g2D n ≤ ca/2, where a is lattice constant h and c is a constant of order, but less than, unity. Hence, 2? ?/g2D n ≤ c2 /9, which according to Coddington et al.  h validates our approach.
tain further insight into the vortex matter in this case, we also perform a thermal averaging as in Fig. 1. This is shown in the second and third rows in Fig. 3. By averaging over di?erent number of snapshots we identify a well-de?ned boundary between the ordered and disordered regions. Indeed, in a ?nite system, averaging will eventually produce a complete smearing even in the center of the trap. We observe signatures of this e?ect in the clear di?erences between the third picture in second and third row of Fig. 3, where averaging was made over 105 and 5 · 105 MC sweeps, respectively. Thus, the time scale of the ?uctuations in the ordered regions are dramatically larger than those related to the ?uctuations in the disordered regions. We next investigate the character of the vortex state in the disordered region. It is known that in the VL the helicity modulus, or equivalently the super?uid density, is zero in any direction [7, 8]. Monitoring of the helicity modulus could be employed to identify a region of a possible tensionful VL in the above pictures, as explained below. In a trapped system the global helicity modulus Υz  has no rigorous meaning, due to the non-uniform Pij . However, we introduce a modi?ed helicity modulus
Here, ′ is over all links where ?ij is nonzero (depending on R1 and R2 ) and N ′ is the number of these links. In a uniform extended system the proliferation of vortex loops happens via a second order phase transition . When loops proliferate, the condensate is destroyed. The temperature of vortex-loop proliferation decreases with increasing rotation (see Fig. 12 in ). Alternatively, the destruction of the phase-coherence along the z-axis is caused by destruction of the lattice order via a ?rst order transition. This scenario, if it is realized in a trapped system, should be manifest in the shape of the helicity modulus as the transition is approached, in that it should be signi?cantly di?erent from the case without rotation. Namely, one should see a remnant of a ?rst order phase transition with a near-discontinuity in the helicity modulus, rather than the continuous variation characterizing a second order transition driven by a proliferation of vortex loops. If one were to observe no appreciable di?erence in the temperature dependence of the helicity with and without rotation, one would conclude that the demarcation line seen in the images separates an ordered region from tension-less vortex tangle, with no discernible VL region. This scenario would imply a well de?ned and regular structure of the boundary of the VLL. ? The results for Υz (R1 , R2 ) are shown in Fig. 4. These measurements indeed show that the presence of a rota? tion signi?cantly reduces the temperature at which Υz vanishes. This reduction relative to the case of no rotation decreases with increasing R1 , R2 . That is, panel (d) is similar to panel (c) (no trap), whereas in panel ? (g) there is little di?erence between Υz (R1 , R2 ) with and without rotation. Thus, for the latter case the presence ? of vortices essentially does not in?uence Υz , and the destruction of super?uid density is driven by the proliferation of vortex loops. Panel (g) is connected to the three leftmost panels in Fig. 3 in the following sense. The distinct demarcation line in these leftmost images separates a vortex solid from a tension-less vortex tangle with no visible tensionful VL region. This is consistent with the experiments showing a very regular edge structure for systems with large number of vortices [1, 2]. Note also
1 0.8 0.6 0.4 0.2 0 0
No rotation (a)
With rotation (b)
the absence of circular distortion for the vortices at the edge of the system. On the other hand, by comparing ? panels (c) and (d), we see that Υz and Υz are quite similar. This points toward a possibility of the existence of a tensionful VL phase close to the center of the trap.
? Υz (r)
0.4 0.6 r/R
0.4 0.6 r/R
0.8 0.6 0.4 0.2 0 0.8 0.6 0.4 0.2 0 ? Υ((R1 , R2 ), T ) 0.6 0.4 0.2 0 0.4 0.2 0
We have considered vortex matter in the model Eq. 1 with uniform and nonuniform density. The uniform case features dominant angular vortex ?uctuations near the wall of the cylinder. Vortex matter in a trapped BEC is more complicated, due to density gradients. We have identi?ed a number of inhomogenous vortex states. Notabley, we ?nd a direct crossover from a vortex solid to tension-less vortex tangle with no discernible intermediate tensionful vortex liquid at low temperatures as the edge of the trap is approached. This explains very regular edges of VLL at ?nite temperatures in many experiments. At higher temperatures, a possible tensionful vortex liquid state located between a vortex solid at the center and a tension-less vortex tangle closer to the edge is identi?ed. Our simulations indicate strong bending ?uctuations of vortices in this region which may obscure its visibility in experiments. An experimental observation of the discrepancy between visible vortex numbers and rotation frequency (as we predict for disordered vortex states) has in fact been observed, albeit in an anharmonic trap . This work was supported by the Research Council of Norway, Grant Nos. 158518/431, 158547/431 (NANOMAT), 167498/V30 (STORFORSK), the National Science Foundation, Grant No. DMR-0302347, Nordforsk Network on Low-Dimensional Physics. Discussions and correspondence with J. Dalibard, A. K. Nguyen, V. Schweikhard, E. Sm?rgrav, M. Zwierlein and especially with Erich J. Mueller, are gratefully acknowledged.
0 0.5 1
T = 0.50
? FIG. 4: Results for Υz (R1 , R2 ). The two top panels show thermal depletion of the super?uid density in the model Eq. 1 (r is the distance from the center of the trap) at the temperatures T = 2.50 (the lowermost curve), T = 2.00, T = 1.67, T = 1.25, T = 1.00, T = 0.50. In panels (c)–(g), the upper curve (+) is the helicity modulus without rotation, while the lower curve (×) the helicity modulus with rotationinduced vortices with ?lling fraction f = 1/36 as functions of temperature. Panel (c) shows Υz for a cubic uniform system with periodic boundary conditions. The upper curve (+) has the properties of a second order transition (helicity modulus vanishes because of vortex-loop proliferation), whereas the lower curve (×) has the ?nite-size appearance of a ?rst order transition (suggesting that the helicity modulus vanishes because of vortex lattice melting) [7, 8]. The ? ? remaining panels (d)-(g) show Υz (0, R/4), Υz (R/4, 2R/4), ? z (2R/4, 3R/4), and Υz (3R/4, R), respectively. The vortex ? Υ plots on the left (obtained at T = 0.50) de?nes the radii R1 and R2 as white circles. Taking parameters from Ref. 2, and using ? = (h/M )Nv /2πR2  where Nv is the number of vortices in the trap, we ?nd ? ? 100Hz. Since ω⊥ ? 500Hz , this puts us well outside the LLL regime.
 K. W. Madison et. al., Phys. Rev. Lett. 84, 806 (2000); I. Coddington et. al., Phys. Rev. Lett. 91, 100402 (2003); V. Schweikhard et. al., ibid 92, 040404 (2004); N. L. Smith et al., ibid, 93, 080406 (2004); I. Coddington et al., Phys. Rev. A73, 063607 (2004); S. R. Muniz et al., Phys. Rev. A73, 041605 (2006).  J. R. Abo-Shaeer, C. Raman, J. M. Vogels, and W. Ketterle, Science 292, 476 (2001).  E. Babaev, A. Sudb?, and N. W. Ashcroft, Nature, 431, 666 (2004).  S. A. Gi?ord and G. Baym, Phys. Rev. A 70, 033602 (2004).  M. Snoek and H. T. C. Stoof, cond-mat/0601695.  R. E. Hetzel, A. Sudb?, and D. A. Huse, Phys. Rev. Lett. 69, 518 (1992).  X. Hu, S. Miyashita, and M. Tachiki, Phys. Rev. Lett. 79, 3498 (1997).  A. K. Nguyen and A. Sudb?, Phys. Rev. B 60, 15307 (1999).
 At least 105 MC sweeps were discarded for equilibration at each temperature and up to 5 · 105 sweeps were used for calculating averages.  G. Watanabe, G. Baym, and C. J. Pethick, Phys. Rev. Lett. 93, 190401 (2004); N. R. Cooper, S. Komineas, and N. Read, Phys. Rev. A 70, 033604 (2004).  L. J. Campbell and R. M. Zi?, Phys. Rev. B 20, 1886 (1979).
 Yu. Lozovik and E. Rakoch, Phys. Rev. B 57, 1214 (1998).  Z. Tesanovic, Phys. Rev. B 59, 6449 (1999).  J. R. Abo-Shaeer, C. Raman, and W. Ketterle, Phys. Rev. Lett. 88, 070409 (2002).  V. Bretin, S. Stock, Y. Seurin, and J. Dalibard, Phys. Rev. Lett. 92, 050403 (2004).