Probing an interacting Bose gas in a quasi-two-dimensional
trapbyShih-Kuang Tung
B.A., National Tsing Hua University, 1999
M.S., National Tsing Hua University, 2001

A thesis submitted to the
Faculty of the Graduate School of the
University of Colorado in partial fulfillmentof the requirements for the degree ofDoctor of Philosophy
Department of Physics
2010This thesis entitled:
Probing an interacting Bose gas in a quasi-two-dimensional trap
written by Shih-Kuang Tung
has been approved for the Department of Physics

E. A. Cornell

D. S. Jin

Date

The final copy of this thesis has been examined by the signatories, and we find that boththe content and the form meet acceptable presentation standards of scholarly work in theabove mentioned discipline.

iii

Tung, Shih-Kuang (Ph.D., Physics )Probing an interacting Bose gas in a quasi-two-dimensional trap
Thesis directed by Prof. E. A. Cornell

A fluctuating regime associated with the Berenzinskii-Kosterlitz-Thouless (BKT) transition is studied in a two-dimensional system of near-degenerate bosons (87 Rb atoms). Two
imaging procedures, the in-situ imaging method and the focusing-imaging method, are implemented to image the coordinate-space and momentum-space density distributions. The
scaled compressibility is extracted from the coordinate-space density distributions obtained
from in-situ images. Comparing the measured compressibility to the prediction from oursemiclassical mean-field model, an onset of a regime which is beyond the mean-field prediction is identified and first resolved at phase-space density nλ2 & 3. Information about
the coherence of the system is extracted from the momentum-space distributions obtained
from the focusing images. The spatial extent of the coherence at a size of the high density region in the system only appears at a distinguishably higher phase-space density ofnλ2 ≈ 8. Therefore, a very interesting regime that is beyond mean-field prediction, but not
yet a quasicondensate, is identified. This regime, which does not exist in three dimensions,is a product of the enhanced interactions associated with reduced dimensionality.

Dedication

To Claire, Tse-Chi, Mei-Shia, and Shih-Hong.

v

Acknowledgements

Until this day, I still feel lucky that I can be one of those guys working at JILA.
Therefore, the first person I should give my gratitude to is my boss, Eric. Not only he gaveme this opportunity to work with all these wonderful people, but also reminded me howmuch I love baseball. As a graduate student, one of the worst things you can do is probably
rooting for a baseball team that has a long rivalry history with your boss’ team. Throughthose going-back-and-forth emails and conversations, I learned how to talk baseball in theproper way. Here, I also have to give credit to Aaron and Russ. They also helped equallywith my ability of formal baseball talking. All I can say is that this was all beyond myimagination, and it was fun. So thank you all for the great experience.
In the lab, I would like to give my special thanks to three people, Peter Engels, Volker
Schweikhard, and Giacomo Lamporesi. Peter was my first teacher in the lab and helped mesurvive through my first year. Volker told me everything I need to know to run the machineand showed me what the “German engineering” means. I had been working with Volker
for almost five years, and with Volker, we had done some works that I am very proud of.
Giacomo joined the lab, when I just took over the lab from Volker. Giacomo is always very
patient with what I have to say about physics whether it makes sense or not. Giacomo andI had done some very interesting works. Some of those works are still very advance andvery interesting. Also, thanks to Giacomo, I finally knew what the “Italian dinning” reallymeans. In my last year at JILA, Dan and Lin joined the lab. I wish them good luck in the
lab.

vi

I would like also thank Russ, Laura, and Huanqian at the EDM labs, Juan and Rob
at the 85 Rb lab, Jason and John at the 40 K lab, John and Ty in the 87 Rb-40 K mixer lab forloaning me equipments and discuss physics with me. Also, I would like to thank Dan, Matt,
Phil, and Tyler for bringing my English conversation skill up to another level.
At JILA, I learned from all the best, Eric, Debbie, Jun, and many others. I learned
from them how to be a successful physicist, and what it takes to be successful. Eric always
has the magic to make complicated things easy to understand. Debbie always has the magic
to asks the right and tough questions to make me think deeper and push myself forward.
Although I didn’t have a very direct connection with Jun in physics, I still learned from him
by watching when lights in his office goes off. Jun has this magic to stay awake and stayenthusiastic for physics. He works so hard to make me want to work even harder. So, I
thank these great teachers.
People at JILA shops is always happy to help. Here I would like to thank several
people whom I often asked for help. They are Carl at electronic shop. Hans and Tom at themachine shop. J. R. at the JILA computing team. Those are the guys who makes my work
at JILA much easier.
Very special thanks to Russ and Juan. They tried very hard to make this place likemy second home, and I surely felt that way. In the end, I would like to give my gratitude to
my lovely wife, Claire, for bringing me happiness outside the physics world.

vii

Contents

Chapter
1 Introduction

1

1.1

My work at JILA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1

1.2

This thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

32 Interacting bosons in two dimensions

5

2.1

Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

2.2First-order correlation function . . . . . . . . . . . . . . . . . . . . . . . . .

7

2.2.1

The ideal Bose gas . . . . . . . . . . . . . . . . . . . . . . . . . . . .

9

2.2.2

The weakly interacting Bose gas . . . . . . . . . . . . . . . . . . . . .

10

2.3

Berezinskii-Kosterlitz-Thouless (BKT) transition . . . . . . . . . . . . . . . .

13

2.4

Experiments on 2D ultracold bosonic atoms . . . . . . . . . . . . . . . . . .

15

3 Experimental details

18

3.1

Optical lattice . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

18

3.2

Time-averaged orbiting magnetic trap . . . . . . . . . . . . . . . . . . . . . .

23

3.3

Imaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

26

3.3.1

In-situ imaging method . . . . . . . . . . . . . . . . . . . . . . . . .

30

3.3.2

Focusing-imaging method . . . . . . . . . . . . . . . . . . . . . . . .

31

viii

4 Analysis of in-situ images
4.1

4.2

The semiclassical mean-field model . . . . . . . . . . . . . . . . . . . . . . .

38

4.1.1

Local-density approximation . . . . . . . . . . . . . . . . . . . . . . .

39

4.1.2

Semiclassical approximation . . . . . . . . . . . . . . . . . . . . . . .

40

4.1.3

Mean-field approximation . . . . . . . . . . . . . . . . . . . . . . . .

41

4.1.4

The “no-condensate” model . . . . . . . . . . . . . . . . . . . . . . .

42Compressiblity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

43

4.2.1

48

Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5 Analysis of focusing images
5.1

38

Mean-field model for momentum-space distribution . . . . . . . . . . . . . .

6 Conclusion and outlook

53
55
64

6.1

Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

64

6.2

Outlook . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

65

Bibliography

67

Appendix
A Observation of vortex pinning in Bose-Einstein condensates [1, 2]

75

B Vortex proliferation in the Berezinskii-Kosterlitz-Thouless regime on a two-dimensionallattice of Bose-Einstein condensates [1, 3]

87

Chapter 1
Introduction

1.1

My work at JILA
In my early work at JILA (2003–2006), I participated in experiments studying rotating

Bose-Einstein condensates (BEC). My first experiment at JILA studied the dynamics oftwo-component vortex lattices [4]. In this experiment, we first prepared our rotating BEC in
one hyperfine state, then used a microwave pulse to transfer some of the atoms into anotherstate. Because of interactions between two different vortex lattices, interlaced vortex lattices
were observed. In the process of achieving equilibrium, the structure of each vortex latticechanged from a triangular one to a square one. This interaction-induced structural transition
was studied in the experiment. In the meantime, I also helped in the experiment studying
the equilibrium properties of rotating BECs [5].
In another experiment, we used a two-dimensional (2D) optical lattice to manipulate
a single-component vortex lattice [2]. The optical lattice was made to co-rotate with avortex lattice and to provide pinning forces to pin the vortices. After carefully stabilizing
and aligning our experiment, we were able to make the rotating optical lattice co-rotate
such that a vortex lattice saw a stationary pinning potential in its rotating frame. By
using a square optical lattice to pin a triangular vortex lattice, we observed the structureof vortices changes from a triangular lattice to a square lattice. The structural transitionin this experiment was caused by the pinning force provided by the optical lattice. In this2

experiment, we studied the dynamics of a pinned vortex lattice and interactions between avortex lattice and an optical lattice.
Following the pinning experiment, we changed our focus from the physics happening
in the rotating frame to the physics in the stationary frame. We loaded a three-dimensional(3D) stationary BEC into a 2D trapping array that was provided by the optical lattice.
By changing the height of the optical potential in each trapping site, we were able to gofrom the condition of EJ À ET to the condition of EJ ¿ ET , where EJ is the tunnelingenergy between adjacent trapping sites, and ET is the thermal energy. In this experiment[3], we observed how the phase of the small BECs in each trapping site became random bywatching how the small BECs reconnected once the optical lattice was turned off. In thecase of EJ À ET , the atoms in different condensates communicated with each other via
tunneling. Therefore, they maintained their phase coherence. When EJ ¿ ET , the thermal
fluctuations in the system were large enough to destroy phase coherence in the system. In
this case, when we reconnected the small BECs, vortices resulted. In contrast to previous
experiments, vortices occurred in the stationary frame. This phenomenon is related to theBerezinskii-Kosterlitz-Thouless (BKT) transition in the X-Y model [6, 7]. My work in thisexperiment was more than changing the mirrors and tweaking the knobs. I developed a 3D
Gross-Pitaveskii equation (GPE) solver to calculate the tunneling energy between adjacent
trapping sites.
In 2007–2008, Giacomo and I tried a couple of experiments. First, we tried to rotate astationary condensate using a rotating optical lattice. We watched its rotation decaying asa function of temperature. Second, we tried to rotate a 2D condensate in a one-dimensional(1D) optical lattice. In this experiment, we successfully rotated a BEC, observed tens ofvortices, and kept a similar number of vortices in the optical lattice “more than 10 seconds.”
In the end, we didn’t submit the work for publication. We decided to switch the focus of
the experiment to studying the BKT transition in a continuous system.

3

For the work described in the main body of this thesis, we used a similar setup to theParis group [8, 9] to create many isolated 2D systems. We then used a microwave-pumping
scheme to image in-situ one layer (or one 2D system) in the trap. We then compared
the measured result with the one predicted from the semiclassical Hartree-Fock mean-fieldmodel. This model is similar to the one reported in [10]. The amazing thing about thesemiclassical Hartree-Fock mean-field model is that there is really nothing amazing about it.
This model discards all the interaction-related magic such as phonons and vortices. As muchas one needs the phonons and vortices to account for BKT physics, which are beyond themean-field description, we were able to identify the onset of a deviation in our measurementfrom the semiclassical Hartree-Fock mean-field model.
We also used a complementary imaging technique known as the focusing-imaging
method [11, 12] to measure the momentum distribution of the system more accurately thanprevious measurements [13, 14]. From the momentum distribution of the system, we observed
the onset of a long spatial extent of coherence in the system. By comparing both the in-situ
and focusing-imaging methods, we were able to identify rather accurately a very interesting
regime that was beyond the mean-field description, but not quite yet a quasicondensate.
Using our data, we also tested the critical phase-space density proposed in [15].
1.2This thesis
The organization of my thesis is as follows: In chapter 2, I give some “simple” 2Dphysics that might help readers understand our work in this experiment. Also, I give a shortreview of previous experiments that studied the BKT transition in a continuous system.
In Chapter 3, I give some details of our experimental setup for creating a 2D system andcontrolling its temperature. I also discuss a couple of imaging artifacts that influence themeasurements. In the end of the chapter, I discuss two different imaging methods used
for measuring the coordinate-space and momentum-space density distributions. Chapter 4

4gives the model and the method we used to identify a very interesting regime that is beyond
mean-field description. I will detail the development of the model we call the semiclassical,
mean-field, local-density approximations. This model is used to compare with coordinatespace images. A couple interesting features of the model are also shown in this chapter. In
Chapter 5, I show the identification of the BKT temperature TBKT , below which the system
possesses a coherent fraction. The measurement of the coherent fraction is also presented
in this chapter. In Chapter 6, I summarize our experiment on exploring the interesting
physics around the BKT transition. Next, I present my personal outlook on the studiesof 2D ultracold-atom systems. Appendix A is a Physical Review Letter [2] I wrote on theinteractions between rotating optical lattices and vortex lattices. Appendix B is also a
Physical Review Letter [3] I co-authored on some precursor experiments to the work of thisthesis. There, we looked at the BKT transition in a superfluid on a 2D lattice. More details
about this experiment can be found in Ref. [16]. I am proud of both of these papers, butthey are so far removed from the topic of a 2D Bose gas in a continuous system, that I chose
not to include a discussion of these earlier works in the main body of my thesis.

Chapter 2
Interacting bosons in two dimensions

2.1

Introduction
It is known that the statistical behavior of systems can be very different in lower

dimensions from it is in three dimensions (3D). In 3D, many-body systems can possess along-range order (LRO) or undergo Bose-Einstein condensation. Take a Bose gas in 3D for
example. The system can undergo Bose-Einstein condensation at a critical temperature. For
this system at the critical temperature or lower, the atomic populations in excited states are
saturated, and atoms begin to accumulate in the ground state of the system, as shown in
Fig.2.1. For a 3D homogeneous system, the Bose-Einstein condensation happens when thephase-space density of the system reaches 2.612 [17]. For weakly interacting atoms in a 3D
homogenous system, the BEC transition still exists.
For a two-dimensional (2D) system, things change dramatically. In a 2D crystal, along-range order does not exist, as was pointed out by Peierls [18]. In a 2D crystal at a finite
temperature, the uncertainty in the relative position of two atoms is given byh(~u(~r) − ~u(0))2 i ∝ T ln

³r´
a

,

(2.1)

where ~u(~r) is the displacement vector of an atom from its equilibrium position at ~r, and a isthe lattice constant of the 2D crystal. Equation (2.1) diverges as the separation between two
atoms goes to infinity. This fact tells us that at any finite temperature, one cannot predict6

Figure 2.1: At a temperature T > Tc , atoms have a very small probability to populate any
state of the harmonic trap. However, at T < Tc , atoms begin to accumulate in the groundstate of the system. This phenomenon is known as Bose-Einstein condensation.

7

the position of an atom far away by knowing the position of another atom at the origin of thecoordinate. In contrast, in a 3D crystal, the uncertainty of the displacement of two atoms issmaller than the lattice constant a as long as the temperature of the system is smaller thanthe melting temperature. After Peierls’ finding, Mermin and Wagner [19] and Hohenberg[20] proved that the lack of long-range order (LRO) [21] at a nonzero temperature in a 2D
crystal is also true for other systems in lower dimensions.
For an interacting Bose gas in 2D, the above statement says there is no Bose-Einsteincondensation at any finite temperature. In a system with weakly interacting bosons, thedispersion of excitations of the system at small momentum is linear instead of quadratic[17, 22]. The linear part of the dispersion is usually referred to as phonons. The existenceof the low energy phonons causes phase fluctuations at finite temperatures. In a lower
dimensional system, the enhanced phase fluctuations can destroy the LRO or the condensate.
The absence of the BEC in a 2D system can be seen by examining the decay of its first-ordercorrelation functions.
2.2

First-order correlation function
The absence of Bose-Einstein condensation can be equivalently written, using first-

order correlation functions, as
lim

|~
r−~
r 0 |→∞

G(1) (~r, ~r 0 ) = 0,

(2.2)

where G(1) (~r, ~r 0 ) is the first-order correlation function defined as [23]G(1) (~r, ~r 0 ) ≡ hΨ̂† (~r)Ψ̂(~r 0 )i,(2.3)

where Ψ̂† (r) and Ψ̂(r) are, respectively, the creation and annihilation field operators thatcreate or annihilate a particle at position r. In the case of ~r 0 = ~r,G(1) (~r, ~r) = hΨ̂† (~r)Ψ̂(~r)i= n(~r).

(2.4)

8G(1) (~r, ~r) is also known as the one-body density matrix, which gives the density distributionn(~r) in the coordinate space.
In an infinite homogeneous system, an important relation exists between the momentum distribution of the system and the first-order correlation function. The momentum
distribution of the system can be expressed asn(~p) = hΨ̂† (~p)Ψ̂(~p)i,(2.5)

where Ψ̂† (~p) and Ψ̂(~p) are, respectively, the creation and annihilation operators in momentum
space. These operators areZ1
Ψ̂(~p) =d~re−i~p·~r/~ Ψ̂(~r),
2π~
Z
1
†
d~rei~p·~r/~ Ψ̂† (~r).
Ψ̂ (~p) =
2π~

(2.6)
(2.7)

Putting Eqs. (2.6) and (2.7) into Eq. (2.5), one finds thatAn(~p) =(2π~)2Zd~sG(1) (s)ei~p·~s/~ ,(2.8)

where ~s = ~r − ~r 0 , A is the area of the 2D system, and G(1) (s) = G(1) (|~s|). As one can see
from Eq. (2.8), the momentum distribution is the inverse Fourier transform of the firstorder correlation function. In other words, the first-order correlation function is the Fouriertransform of the momentum distribution, i.e.,1
G (~s) =A(1)Z
d2 p n(~p)e−i~p·~s/~ .

(2.9)

In the presence of a BEC, the momentum distribution would have a singularity at
p~ = 0, andn(p) = N0 δ(~p = 0) + nT (~p 6= 0),(2.10)

where nT is a smooth function of p~. One can see, after plugging Eq. (2.10) into Eq. (2.9), thefirst-order correlation function approaches a constant value of N0 /A, which is the condensate

9density, as s → ∞. In contrast, in the absence of the condensate, the first-order correlation
function goes to zero when s → ∞. As mentioned earlier, in 2D, there is no condensate at
any finite temperature. Therefore, one would expect the correlation function to approachzero at infinity for any 2D system with continuous symmetry.
2.2.1

The ideal Bose gas

How G(1) (s) decays depends on the behavior of nT . For a highly nondegenerate ideal
Bose gas in a 2D box, the distribution is well described by the Boltzmann distribution,~2 k 2

nk = e−( 2m −µ)/kb T ,(2.11)

where the chemical potential is determined by
Z ∞
2 k2
L2
−µ)/kb T−( ~2mN=dk2πke(2π)2 0L2
= 2 eµ/kb T ,λdb(2.12)

where λdb is the de-Broglie wavelength, L2 is the area of the 2D box, and N is the total
number of atoms. Using Eqs.(2.9), (2.11), and (2.12), one finds that in the high-temperature
regime, the first-order correlation function can be written as
G(1) (s) = ne−π( λs )2db.

(2.13)√In this regime, nλ2db ¿ 1, and G(1) (r) is a Gaussian with a width of λdb / 2π.
2 2

In the regime close to degeneracy, where nλ2db ∼ O(1), only the states with ~2mk ¿ kb T
are strongly occupied. In this case, the number of atoms in state k is given bykT
nk = ~2 k2.
−µ
2m(2.14)

Using Eqs. (2.9) and (2.14), one can get
G(1) (s) =2
K0 (r/lc ),λ2db(2.15)

10

where K0 is the zero order of the modified Bessel function of the second kind [24], and lc =p√2
λdb enλdb /2 / 4π. The asymptotic behavior (at large r) of K0 can be found as e−r/lc / r/lc .
As one can see from Eqs. (2.13) and (2.15), the correlation functions for 2D ideal Bose gases
eventually decay to zero at r → ∞. In a system close to degeneracy, the spatial extent ofcoherence lc can still be very long and increases exponentially according to its phase-space2

density i.e. lc ∝ enλdb /2 . Interestingly, in the 1D case, the spatial extent of coherence in thisregime increases only linearly according to nλ2db [25].
2.2.2

The weakly interacting Bose gas

Let’s consider the first-order correlation function for an interacting 2D Bose gas. At
high temperature, when nλ2db ¿ 1, the interactions have little effect. Therefore, one expects avery similar result with ideal-gas model. In the low temperature regime, Petrov et al. [26, 27]
found a similar result for the spatial extent of coherence for weakly interacting bosons in aharmonically trapped system. They used a method developed for 3D-trapped condensates
to calculate the excitation spectrum of the system. Then, using the excitation spectrum tocalculate the mean-square fluctuations of the phase, they found that at temperatures below
the critical temperature, the characteristic radius of the phase fluctuation was
Rφ ≈ λdb e

nλ2
db
2

,

(2.16)

where n is the density of the system. From the exponentially increasing spatial extent ofcoherence, one expects a slow decay in the first-order correlation function.
Calculating the first-order correlation function for a weakly interacting Bose gas requires some serious mathematics. However, at low temperature T ¿ TBKT , the calculation
can be simplified by a series of approximations. The approach I adopted here is a 2D versionof the 1D derivation in [25]. A similar 2D derivation can be also found in [28]. First, in the11framework of the second quantization, the Hamiltonian of the system can be written as
Z

2
† −~

U0
dr Ψ̂
∇ Ψ̂ +
2m
2

Ĥ =

Z

2

dr Ψ̂† Ψ̂† Ψ̂Ψ̂,(2.17)

where U0 is the interaction strength. Here, a contact interaction is adopted. If the fluctuationof the field is small, one can write the field operator asΨ̂(~r) =pn + δn(~r) eiθ(~r) ,(2.18)

where n is the the mean density, δn is the density fluctuation operator, and θ is the local
phase operator. The density and the phase operator obey the commutation relation[δn(~r), θ(~r 0 )] = −iδ(~r − ~r 0 ).

(2.19)

Plugging Eq. (2.18) into Eq. (2.17), one can find a new quadratic Hamiltonian Ĥ 0 in termsof the operator θ and δn,
~2 n
Ĥ =
2m
0

Z

¸
Z · 2
~ n ∇δn(~r) 2r) 2 22 δn(~(∇θ(~r)) d r +() + 2U0 n () d r.
2m2
2
2

2(2.20)

Next, by expanding the two operators in sinusoidal modes, one gets
´
X√ ³
~
~
θ=
2 θc,k cos(k · ~r) + θs,k sin(k · ~r) ,(2.21)

~k>0

δn =

´
X√ ³
2 δnc,k cos(~k · ~r) + δns,k sin(~k · ~r) ,(2.22)~k>0

where ~k = L2π2 (nx x̂ + ny ŷ), which are the eigenmodes of a 2D box with a size of L2 . Plugging
Eqs.(2.21) and (2.22) into Eq. (2.20), one can find for each wave vector k the corresponding
Hamiltonian is given byµĤk = L2¶
¢ n~2 k 2 ¡ 2¢
U0 ¡ 22
2
δnc,k + δns,k +θc,k + θs,k .
22m(2.23)

At a sufficiently low temperature, the contribution to Ĥk from the density fluctuation δn/n ismuch smaller than the contribution from the phase modulation θ. In this case, one can, to a12good approximation, keep only the phase modulation term in Eq. (2.23). If the temperature
is still much larger than U0 n, only the excitation and a smaller wave vector k are relevant,and classical statistics apply. Thus, for every quadratic degree of freedom the energy is
kb T /2, and one finds that22i=i = hθc,k
hθc,k

1
π
.
2
2
nλdb k L2(2.24)

Since the density fluctuation is neglected, the first-order correlation Eq.(2.3) can be
written as
G(1) (r) = nhe−i(θ(~r)−θ(0)) i=e−h(θ(~r)−θ(0))2 i/2(2.25).

The second equality holds because of the assumption of a Gaussian fluctuation in the phasefluctuation. The phase fluctuation is given byh(θ(~r) − θ(0))2 i =

X

³
´2 X
2
2
2hθc,ki sin(~k · ~r)2 .
2hθc,ki cos(~k · ~r) − 1 +(2.26)

k>0

k>0

Using Eq. (2.24), one obtainsh(θ(~r) − θ(0))2 i = 4(
Taking L → ∞, one can replaceP(1)

π X 1 − cos(~k · ~r)1).
nλ2db k 2 L2 k>0
k2

L 2
)
with ( 2π

R

−1
2πnλ2
db

(2.27)

d2 k. Then, one finds
R 1/ξ 1−cos(~k·~r) kd20k
G (r) = neµ ¶1/nλ2db
ξ
=n
,
r

(2.28)

where 1/ξ is the cutoff wave vector. Therefore, at a temperature of T ¿ TBKT , one expects
an algebraic decay in the first-order correlation function. At infinity, the spatial extent ofcoherence will eventually go to zero. However, for a finite system, the phase fluctuation
can still be very small. This algebraic decay of the coherence is the definition of a quasicondensate. This power-law decay of the first-order correlation function at large r is first13

obtained by Kane and Kadanoff [29]. Although there are other different definitions for aquasicondensate in recent literatures, I use the original definition throughout my thesis to
avoid confusion. Regarding the regime that is beyond mean-field behavior and lacks a longspatial extent of coherence, we call it “a very interesting regime.”
2.3

Berezinskii-Kosterlitz-Thouless (BKT) transition
In a 2D system, the long-wavelength excitations from phonons are strong enough todestroy the condensate or LRO, but not strong enough to completely suppress the superfluid.
As was first pointed out by Kane and Kadanoff [29] and then proved by Berezinskii [6],there exists a normal-to-superfluid phase transition at the critical temperature TBKT . Later,
Kosterlitz and Thouless [7] found that this phase transition is related to the unbinding of
the vortex-antivortex pairs that exist in a 2D system at a finite temperature, as shown in
Fig. 2.2. At a temperature lower than TBKT , the 2D system is a superfluid, but it doesn’t
possess a long-range order. This noncondensate superfluid is called a quasicondensate.
The critical temperature of the BKT transition can be estimated as follows. Consider
the energy of a free vortex Ev that appears in an infinite homogeneous 2D system. The
vortex considered here is a quantum object that appears in a superfluid. The atomic density
at the center of the vortex is zero, and the phase winding around the vortex is 2π. The
antivortex has the same properties, but has a phase winding of −2π. The velocity of thesuperfluid is
~vs =

~~
∇f (~r),m(2.29)

where f (~r) is the function representing the phase variation in a super system. Therefore,
the superfluid velocity around a vortex can be found to be
~vs =

~
θ̂,
mr

(2.30)

where r is the distance from the center of the vortex. The energy of the vortex can be written

14

Figure 2.2: A 2D system cannot have a long-range order or condensate, but it can still
undergo a BKT transition that is a normal-to-superfluid transition. At T > TBKT , thethermal fluctuation can proliferate free vortices or free antivortices in the atomic cloud.
Those vortices destroy the phase coherence of the system. However, at T < TBKT , vortex
and antivortex are bound. In this regime, the coherence of the system decays according toa power law. The phase is often referred to as quasicondensate.

15

as

Z R
Ev =
ξ

1
ns
2

µ

~
mr

¶2
d2 r,

(2.31)

where R is the size of the system, and ξ is the healing length of the vortex. Creating a freevortex increases the entropy of the system as
Sv = lnR2.
ξ2(2.32)

Consider the free energy F = Ev − T Sv , i.e.,
F = kb T (ns λ2db − 4) ln

R2
.
ξ2

(2.33)

One can see that the proliferation of free vortices would become energetically favorable
if the phase-space density was smaller than 4. One can also convert this critical phase-spacedensity to the critical temperature TBKT . Then, one gets
TBKT =

π~2 ns.
2mkb T(2.34)

When the temperature of the 2D system is lower than TBKT , the system is in a superfluid
phase, and the vortex-antivortex pairs are bound. However, when the temperature of thesystem is warmed up to TBKT , the vortex-antivortex pairs become unbound, and free vortices
proliferate. Then, the phase of the system becomes normal.
2.4

Experiments on 2D ultracold bosonic atoms
The first effort to create a 2D ultracold atom system was that of Ref. [30]. Next, aseries of papers from the Paris group [8, 9, 13, 31, 32] focused on the fluctuating regime ina 2D ultracold atom system. The Paris group observed the phase fluctuations of a 2D Bosegas. First, the group loaded a 3D cloud into a 1D optical lattice. The lattice sliced the 3D
cloud into many 2D layers. Second, the researchers used microwave pumping to get rid of
all but two central layers for interfering. From the interference pattern they obtained, they

16

extracted information about correlations in the 2D layers according to the method proposedin [33]. The exponent that described the coherence changed from 0.5 to 0.25, which was
expected as the system passed the BKT transition. The researchers saw dislocations (ordefects) show up in the interference pattern. They also found that the probability of seeing
dislocations had a temperature dependence. Both experimental results encouraged them torelate their observations to the BKT physics. Their observations also successfully generated
much interest in theoretical studying the fluctuating regime in a 2D ultracold atom system.
In a later paper [13], the same group identified the BKT transition with an observationof bimodality (a Gaussian + an inverse parabola) in a time-of-flight (TOF) image. Right
before the observation of the bimodal distribution, the researchers measured the criticalnumber by integrating the measured optical depth over the whole TOF image. Comparing
their measurements to the predictions of the ideal gas model, they saw a substantial deviationthat was, however, very different from what was expected in a 3D system.
In this paper [13], the researchers identified the BKT transition with the observationof a bimodal distribution. However, the appearance of bimodality strongly depended on thetwo models chosen for fittings. It is not clear to me that the two models, a Gaussian anda inverse parabola, are accurate enough to describe a 2D cloud at BKT temperature. Also,
since the ideal and mean-field models lack of the magic of phonons and vortices, using either
one to describe a system close to the BKT transition is, in principle, wrong. However, ina trapped system, the area of the fluctuating regime is small compared to the size of thesystem. Therefore, a global parameter like the critical number is not very different for amean-field model or for a beyond mean-field model.
In another experiment [14], the Gaithersburg group loaded sodium atoms into a lightsheet that created a single-layer 2D system. After releasing the atoms from the trap, athree-component (superfluid quasicondensate, non-superfluid quasicondensate, and normal
phase) density distribution was observed. The “non-superfluid quasicondensate” used in

17

the Gaithersburg group’s paper represented the regime where the density fluctuation of thesystem was suppressed, but without a long coherence length. Considering that an ideal
Bose-gas model would also predict a three-component like distribution, a description of athree-component distribution was just not very accurate for describing the lower-dimensiongas.

Chapter 3
Experimental details

3.1

Optical lattice
Optical potentials provide a variety of opportunities for studying ultracold atoms.

Laser beams can interfere with each other and create a rich variety of interference patterns.
Even before the creation of the Bose-Einstein condensate (BEC), scientists used this amazing
tool to study cold atoms in optical lattices. In the post-BEC era, superfluids in periodical
potentials have been the focus of studies involving ultracold atoms. Greiner et. al. [34]loaded a condensate into a 3D optical lattice and transferred the system from a superfluid
state into a Mott state. This work has triggered a new generation of experiments withthe possibility of realizing a quantum computer in an ultracold atom system. Also, theease of creating a high spatial frequency in an optical potential is conducive to creating
lower-dimensional systems.
An atom interacts with a light field via an induced-dipole moment. Considering theatom and light as a classical driven, damped oscillator, one can derive the interaction energy
and the scattering rate, respectively, as [35]3πc2 Γ
I,
2ω03 ∆µ ¶2
3πc2 Γ
Γsc =I,
2~ω03 ∆

Uint =

(3.1)
(3.2)

where c is the speed of light, ω0 is the natural frequency of the classical atomic oscillator,

19∆ is the detuning of the driven frequency to the natural frequency ω0 , Γ is the on-resonant
damping rate, and I is the field intensity. Equation (3.1) and (3.2) give a good approximation
to a real atom if Γ ¿ ∆, and there is no saturation effect.
A couple things can be easily seen from Eqs. (3.1) and .(3.2). First, if the light is red
detuned (∆ < 0), the shift of the interaction energy is negative. If there exists a spatial
dependence in the field, atoms would be pushed towards the field maximum (the nodes ofthe optical lattice). If the light is blue detuned (∆ > 0), atoms would be pushed towards thefield minimum (the antinodes of the optical lattice). Second, the interaction energy scales as(I/∆), whereas the photon scattering rate scales as I/∆2 . Therefore, using a large detuning
and a large intensity can keep the photon-scattering rate low while providing a large enoughpotential depth for trapping.
For a two-beam interference pattern, the field intensity can be written as
I(z) = I0 cos2

³π ´
z ,
d

(3.3)

where I0 is the maximum intensity, and d is the period of spatial modulation determined by
d=

λ
¡ ¢,
2sin 2θ

(3.4)

where θ is the angle between the two interfering beams, and λ is the wavelength of thelaser. If the interfering beams are elliptical Gaussian beams with waists of σx and σz and apower P , the maximum intensity I0 is given by 2P/πσx σz . When an atom placed in such afield, the interference pattern will create an array of trapping sites. This array of trappingsites is called a 1D optical lattice. For a single trapping site, the trapping frequency can beapproximated as

s
ωz =

α
∆

µ

2P
πσx σz

¶

π
,
d

(3.5)

where α = 3πc2 Γ/mRb ω03 .
In our experiments, the wavelength for the interfering laser is 532 nm. The angle
between the two interfering beams is 8◦ , which gives a lattice constant of 3.8 µm. The

20ωz(kH z)
2π

3

(a)

2
1
0

0

0.1

0.2

0.3

0.4

0.35

0.6

0.7

0.8

0.9

1

2.2

2.4

2.6

2.8

3

(b)

0.3

lz

0.5

P (W )

0.25
0.2

g̃

0.15

(c)

0.1

0.05

1

1.2

1.4

1.6

1.8

2

ω z /2π (kH z)

Figure 3.1: (a) Trapping frequency for a 87 Rb atom at the bottom of a single trapping site in
the 1D optical lattice as a function of the optical power per interfering beam. (b) size of the
2D Bose gas in the lattice direction as a function of trapping frequency (c) The dimensionless
2D interaction strength g̃ as a function of trapping frequency.

21

waist of the interfering beams, σx and σz , is 225 µm and 84 µm, respectively. Figure 3.1(a)shows the relation between the input power and the trapping frequency at the bottom ofa single trapping site. We measure the power of a single beam to be 300 mW. Using themodel above, one gets a trapping frequency of 1.55 kHz, which is in fine agreement with themeasured-trapping frequency of 1.4 kHz. Figure 3.1(b) gives the size, lz , of the 2D Bose gas
in the tight confinement (ẑ) direction. Figure 3.1(c) gives the dimensionless 2D interaction√
strength, g̃ ≡ 8π(as /lz ), where as is the 3D scattering length of the atom, and g̃ = 0.09
in our experiments, which compares to g̃ = 0.13 [9, 13] and 0.146 [32] in the Paris groupand g̃ = 0.02 in the Gaithersburg group [14]. As one can see from Fig. 3.1(c), g̃ has avery weak dependence on the power of the interfering beams. To reach g̃ = 1, which is thevalue for a strongly interacting system like a 4 He film, one has to increase the power by ahundredfold. Another way to increase g̃ is to reduce the lattice constant, d. As the lattice
constant becomes smaller and smaller, one has to worry about atoms tunneling between
different layers, however.
The experimental setup is shown in Fig. 3.2. The laser light is provided by a Verdi V10,
which has a wavelength of 532 nm. The light passes through an acousto-optic modulator(AOM), which is used to control and stabilize the intensity of the beam. The intensitystabilizing module is made of a microwave source, a microwave amplifier, a mixer, and a
feedback circuit that was designed by the JILA shop. After the AOM, the light is guided
through a single-mode polarization-maintaining (PM) fiber. When the lattice is not in use,
we keep the lattice continuously on at its close-to-maximum value such that the fiber andthe AOM can thermally equilibrate at their operating temperature. The light is blocked
by a mechanical shutter after the fiber. Before the shutter is open, the intensity of the
optical lattice is abruptly turned off by the intensity-stabilizing module. After the shutter
is opened, the optical lattice is gradually ramped up to the desired value in one second.
The two interfering beams are focused onto the atomic cloud by a cylindrical lens. The

Figure 3.2: Experimental setup for the 1D optical lattice. The total power we injected into the fiber was ∼1.8 W. The power
output from the fiber was ∼0.8 W. The laser light first passed through a polarizing beam-splitting (PBS) cube to select asingle polarization. A few percent of the laser light was split out by a pick-off window and sent to the intensity-stabilizing
module. After the laser beam was split into two interfering beams by the second PBS, they passed two mirrors (indicated by
the dashed line) that was mounted on the same post to increase rejection off common vibration modes.
22

23

resulting waists of the beams are 225 µm (long axis) and 84 µm (short axis). The laser lightis separated into two beams by a polarizing beam-splitting (PBS) cube. The power of thetwo beams is carefully balanced by rotating the λ/2 plate in front of the PBS. The λ/2 plateafter the PBS cube rotates the polarization back for maximal interference. The two mirrors
on the two arms of the interferometer are mounted to the same post to increase the rejection
off common vibration modes.
3.2

Time-averaged orbiting magnetic trap
In our experiments, the optical lattice is blue-detuned. Atoms are trapped in the nodesof the lattice, which are intensity minima. The optical lattice only provides confinement in
the z (or gravity) direction. The radial confinement comes from a time-averaged orbiting
magnetic trap, a so-called TOP trap [36]. The TOP trap was used for the creation of thefirst BEC in 1995 [37]. It is made of an anti-Helmholtz coil (quadrupole field) along thez-axis and a pair of Helmholtz coils (bias field) along the x- and y- axes. The bias field
generated by the Helmholtz coils can be made to rotate according to
Bb = B0x cos(ωT t)x̂ + B0y sin(ωT t + φ)ŷ.

(3.6)

The rotating bias field shifts the center of the quadrupole trap generated by the antiHelmholtz coil in a circular trajectory. To be able to confine atoms in the trap, the conditionω ¿ ωT ¿ ωL

(3.7)

must be met [38]. Here ω is the trap frequency, and ωL is the Larmor frequency splitting
between the spin states of an atom in a magnetic field. One of the advantages of the TOP
trap is its ability to become very round. By changing the value of B0x /B0y and φ separately,
we are able to minimize the difference of the trapping frequencies in two independent axesin the horizontal plane to within 3%. A 3% is far from the record roundness achieved in this24

lab [39]. However, it is good enough for the current experiment. Using a transverse roundtrap allows us to easily take an azimuthal average of our data to minimize statistical noise.
Because of gravity, the position of the cloud often sags from the trap center. The
trapping potential under the influence of gravity can be written asZ¯
1 T ¯¯¯~µ ¯BT OP (x, y, z, t)¯ dt + mgz,
U (x, y, z) =
T 0

(3.8)

where µ is the Zeeman-dependent magnetic moment of an atom, m is the mass of an atom,~ b , andg is the gravitational acceleration, T is the oscillating period of the bias field Br00~ T OP = ( Bz x + B0 cos(ωT t))2 + ( Bz x + B0 sin(ωT t))2 + (Bz0 z)2 .
B
(3.9)
2
2
The equilibrium position can be found by solving∂U (0, 0, z)= 0,∂z(3.10)−(B0 /B 0 )zmin = q 0 z .
z 2
( µB
) −1
mg

(3.11)

which gives

If the bias field is suddenly (1/ω À t À 1/ωL ) turned off, the magnetic field aroundthe equilibrium position, (0, 0, zmin ), can be expressed asr|B(x, y, z)| =(Bz0 zmin + Bz0 (z − zmin ))2 + (

Bz0 2
B0
x) + ( z y)2 .
2
2(3.12)

If zmin is large enough such that Bz0 zmin À (Bz0 /2)δR, where δR is the extent of the atomic
sample, the field can be approximated as~
B(x,y, z) ≈ (Bz0 zmin + Bz0 (z − zmin ))ẑ.

(3.13)In the experiment, an extra bias field of 24 G is added along the the z direction to further
reduce the magnetic-field curvature in the x and y directions. This magnetic gradient isused, in conjunction with a microwave pumping pulse, for selecting a single atomic layer for
imaging.

25

Figure 3.3: Instantaneous X-Z cross-section of the TOP trap. The contours represent thesame magnitude of the magnetic field at a particular time. The blue cross indicates thecenter of the quadrupole trap, and the red blobs indicate atoms loaded into a 1D opticallattice. The circle in the x-y plane with radius R is called the “circle of death.” Atoms canescape from here after a Majorana spin flip. Because of the gravity sag, the equilibrium
position is offset by zmin from the center of the TOP trap. With the current experimental
setup, zmin is larger than the size of the cloud. Therefore, the only atom-ejection mechanism
comes from the RF knife.

26

Since the radial confinement only comes from the TOP fields, it is possible to implement
rf evaporative cooling in the system. Evaporative cooling provides an easy and accurate
handle for manipulating the temperature of the system. The working principle of evaporative
cooling is selectively forcing hotter atoms in the atomic distribution away from the trap,
which lower the temperature of the system after thermal equilibration. This selective removal
comes from two mechanisms in a TOP trap. First, the applied radio-frequency field (RFknife) pumps atoms with higher energy from a trappable state to an untrappable state at theedge of the cloud. Second, atoms with higher energy can move to the center of the rotating
magnetic-quadrupole field where atoms can undergo a Majorana spin flip and escape fromthe trap. Both mechanisms help in the early stage of the evaporative process. In the final
stage of evaporation, only the selective removal from the RF knife is relevant, as shown in
Fig.3.3. The temperature of the central layer vs frequency of the RF knife is plotted in Fig.
3.4.
3.3

Imaging
Our detection of atoms is the absorption-imaging method. The advantage of absorption

imaging is that it is easy to set up compared to other imaging methods. One locks the probing
light to one of the atomic transitions and shines the light onto the atoms. The photons in the
probing beam will be scattered by the atoms and cast a shadow on the probing beam. Thisshadow is then imaged onto a CCD camera. Two separate images are often taken during
the imaging process, one with atoms and the other without atoms. A quantity called opticaldepth (OD) is often extracted after an imaging sequence. If the atomic transition is a cyclingtransition, a single atom can scatter many photons, which helps increase the signal-to-noise
ratio (SNR) of the images. Assume that the probing beam is fired for a long enough time,such that the photon scattering process reaches a steady state and can be faithfully described

27

240

T (nK)

200

160

120

80

40

0
5.7

5.8

5.9

6.0

6.1

6.2

6.3

6.4

f

RF

Figure 3.4: Temperature vs frequency of the RF knife. The temperature is obtained fromthe fit of the in-situ images of the atoms in the central layer of the 1D optical lattice. The
red line is the linear fit to the data.

28as
Is (x, y) = (1 − α)I0 e−OD(x,y) + αI0 ,(3.14)

where α represents the small fraction of the light in I0 that cannot interact with the atomic
sample because of polarization imperfection and tails in the laser-frequency distribution.
Here OD(x, y) is defined asZ
OD(x, y) =n(x, y, z)σ0 dz,(3.15)

where n(x, y, z) is the density of the atomic sample, and σ0 is the scattering cross-section.
As one can see, OD(x, y) can be found by taking
OD(x, y) = ln

(1 − α)I0.
Is (x, y) − αI0(3.16)

However, the CCD camera does not distinguish the light that can interact with the atomfrom light that cannot. The quantity we actually measure from the camera is ODmeas which
is defined as
ODmeas = ln

I0.
Is (x, y)(3.17)

Equation (3.16) can be rewrittten using ODmeas . One can show that
OD(x, y) = ODmeas (x, y) + ln (1 − α) − ln (1 − αe−ODmeas (x,y) ).

(3.18)

If α and αeODmeas are both small numbers (¿ 1), the difference between OD(x, y) andODmeas (x, y) can be approximated as
OD − ODmeas = α(eODsat − 1).

(3.19)

In our experiment, α ≈ 3%, which gives about a 10% correction if the maximum OD reaches
2. ODmeas reaches its saturation value ODsat , when (1 − αeODsat ) = 0. Thus, α can be
determined by saturating ODmeas according to
α = e−ODsat .

(3.20)

29

Although one can correct this saturation effect, it is not recommended because the ODmeas
is particularly sensitive near ODsat . Any imaging imperfection can cause large fluctuationsin the OD correction. Often times, α is also position dependent. If one can stay awayfrom the saturation value of ODmeas , one can still have a relatively linear mapping from
ODmeas to OD. However, a 3% non-interacting light gives a saturation value of 3.5. The
linear range (< 10% correction) for the ODmeas would just be 0 < ODmeas < 2. This lack ofa linear dynamic range is one of the drawbacks of the absorption-imaging method. In our
experiments, we limited ODmeas to a value smaller than 2 by putting only a small numberof atoms into the imaging state. We also allowed the atomic cloud to expand axially for 1
ms, which helps reduce the atomic density by a factor of 8, before the probing beam wasfired. This trick does not change the optical depth of atoms, but does help reduce imaging
artifacts. To convert the optical depth to a real atomic density, one has to know σ0 , whichis the scattering cross-section. However, the scattering cross-section depends strongly on
the polarization, the detuning of the probing beam, and the population of the atoms. This
makes accurate determination of the value of σ0 difficult. In our experiments, the scatteringcross-section, together with the transfer percentage, is the fitting parameter η. For a twolevel atom interacting with resonant light, the on-resonant scattering cross-section σ0 can beeasily calculated [40]; σ0 = 3λ2 /2π.
Our line-of-sight is along the axis of tight confinement, which means we do not needto do a deconvolution of our images to get the 2D density distribution. After the probing
sequence, we extract the OD image according to Eq. (3.17). As mentioned earlier, ourtrap is axially symmetric. Thus, the OD images also share the same symmetry. The axial
symmetry in the image of the atomic sample allows us to take azimuthal averages to minimize
the statistical noise in the absorption images, as shown in Fig. 4.3.
In our experiment, the intensity of the probing beam I is approximately 0.1 Is , where
Is is the saturation intensity for the atomic transition. The photon-scattering rate of the30atoms, Y , is given by
Y = Γρee ,(3.21)

where Γ is the decay rate of the excited state, which is (26 ns)−1 for the transition of |F =2i → |F 0 = 3i, and ρee is the population in the excited state, which is obtained by solving
the optical Bloch equation [40]. ρee is given by
ρee =

S0 /2
,
)2
1 + S0 + ( 2∆Γ

(3.22)

where S0 ≡ (probing beam intensity I / saturation intensity Is ), ∆ is the probing detuning,and Γ is the natural linewidth of the excited state. In the low-intensity regime, the scattering
rate is proportional to the intensity, i.e., Y ' Γ(S0 /2). In this case, the attenuation of thelight passing the atomic sample is well described by the Beer’s law, which leads to Eq. (3.14).
The intensity of the probe beam is 0.1 Is with a resonant frequency. The duration of
the probing time is set to 50 µs. During this time duration, an atom can scatter on average
100 photons. This scattering gives a 0.75 MHz Doppler shift and moves the atom 1.4 µm
away from its original position (the imaging plane). These are both negligible effects.
3.3.1

In-situ imaging method

In this method, the atoms are first prepared in the |F, mf i = |1, −1i state in the TOP
trap. Then, the 1D optical lattice is turned on. The 1D optical lattice slices the atomic
cloud into many layers. Tunneling between layers is negligible. Therefore, every single layer
can be viewed as an isolated 2D system. Atoms in individual layers are allowed to thermallyequilibrate before probing occurs.
We implement the in-situ imaging method as follows: First, we turn off the rotating
bias field such that the atomic cloud sees a transient magnetic-field gradient Bz0 . Second,
we turn on an extra bias field along the ẑ direction to reduce the magnetic-field curvatureacross the 2D cloud. Third, we apply a weak microwave pulse to pump a small fraction of

31

the atoms in one layer from |1, −1i to |2, −1i, as shown in Fig. 3.5(a). Fourth, we turn off
the the optical lattice to reduce the stray light on the camera. We also turn off most of themagnetic fields to increase the pumping efficiency for the probe light. Finally, we fire theprobing beam, which is resonant at the |F = 2i → |F 0 = 3i transition with σ + polarization.
Details of the experimental sequence are shown in Fig. 3.6.
By sweeping the frequency of the microwave pulse, we can see a layer structure on thenumber of atoms we transfer into the imaging state, as shown in Fig. 3.7. Given that themagnitude of Bz0 is 7.8 G/mm, a measured energy splitting of h ∗ 40 kHz between adjacentlayers gives a lattice constant of 3.8 µm. This measured lattice constant is in good agreementwith the one calculated from Eq. (3.4), which using the measured angle, θ = 8◦ , between
the two interfering beams.
3.3.2

Focusing-imaging method

Measuring momentum-space distribution is more difficult than it first sounds. Conceptually, it can be done by turning off atom-atom interactions, turning off the confinement,
allowing the cloud to expand ballistically for a very long time, and delaying the absorption image until the position of each atom is determined by its initial momentum (not by
its initial position in the cloud). Comparing this process to the well-known time-of-flightmethod, we anticipate two difficulties: turning off interactions and obtaining a sufficiently
long expansion time.
To remedy these problems, one needs to find a way to turn off the interactions and makethe required expansion time shorter. The method we use is a temporal-focusing techniquethat we refer to as the“focusing” imaging method. We first select a single layer via microwave
pulses. Then, we turn off the optical lattice and let the layer expand into a purely magnetic
trap (TOP). Because of the 140:1 aspect ratio of the cloud, the resulting expansion is initiallypurely axial, which very rapidly reduces the 3D density (and thus the strength of the atom-

Figure 3.5: Drawing of the atomic population during the imaging processes. (a) In the in-situ imaging method, a weak
microwave pulse pumps a small number of atoms in one layer from |1, −1i to |2, −1i for imaging. (b) In the focusing-imaging
method, a second strong microwave pulse, together with a time-varying magnetic field, is applied to swap the population
between |1, −1i and |2, −1i. The atoms in |1, −1i will still be trapped, while the atoms in the |2, −1i will be ejected fromthe trap. During focusing, after the atoms in |2, −1i are all gone, a pumping laser is fired to put the atoms in |1, −1i back to|2, −1i for imaging.

32

33

Figure 3.6: Qualitative schematic (not to scale) of the experimental timing of the in-situ
imaging method. (a) The magnitude of the rotating bias field B0 is 7.3 G. The rotating bias
field is turned off 1 ms before the microwave pulse fires. (b) The duration of the microwave
pulse is 50 µs. (c) The 1D optical lattice is turned off 1 ms before the probing occurs. (d)The magnitude of the extra bias field Bex is 24 G. This extra bias field is turned on at thesame time as the rotating bias field is turned off. (e) The magnitude of the magnetic-field
gradient is 7.8 G/mm. (f) The duration of the probing beam is 50 µs.

34

5

2.5x105

2.0x105

N (a.u.)1.5x105

1.0x104

5.0x10

0.0
50

100

150

200

250

300

f (6830.XXX MHz)

Figure 3.7: Atom number N vs microwave frequency f of the short pumping microwave pulse.
This plot shows a layer structure of the atomic sample. The data is taken continuously over
one hour. The lattice constant derived from the data is 3.8 µm.

35

atom interactions). This expansion does not affect the in-plane coordinate- or momentumspace distributions. After the near-instantaneous suppression of the (repulsive) atom-atominteractions, each atom undergoes essentially free motion in the radially symmetric harmonic
confinement in the x-y plane. At the focused time, tf = 41 ω2πr , the initial radial-momentum
distribution is mapped onto the final-coordinate space distribution, and vice versa. At
the point, we take an absorption image and scale the spatial coordinate by mωr to get
the momentum-space distribution that existed just as the focusing sequence began. Thisfocusing technique is an extension to 2D of a procedure originally developed for imagining
1D momentum distributions [11, 12].
The focusing-imaging method proceeds as follows: First, we turn off the rotating bias
field in the TOP trap. Then we use a short microwave pulse to pump a small fraction (5%)
of the atoms that are in one single layer from |1, −1i to |2, −1i. An extra bias field is added
to reduce the magnetic curvature across the 2D cloud. This step is identical to step 1 of
the in-situ imaging method. Second, we turn off the optical lattice, and fire the second
microwave pulse and ramp down the extra bias field. This step exchanges the populationin |1, −1i and |2, −1i through an adiabatic rapid passage (ARP), as shown in Fig. 3.5(b).
Third, we turn the rotating bias field back on and shift the extra bias field to a new value.
Fourth, we wait until the momentum-space distribution is mapped onto the coordinate space.
Finally, we turn off most of the magnetic fields before firing the probing beam. The detailsof the experimental sequence are shown in Fig. 3.8.
The mapping of the in-plane momentum-space distribution and the in-plane coordinatespace distribution at the focused time tf can be understood using a classical harmonic oscillator model. This model works well for describing the in-plane motion of a 2D Bose gas since
the interaction between atoms is greatly suppressed soon after the optical lattice is turned
off. The equations of motion of an atom with a initial velocity ṙ(0) at the initial position

36

Figure 3.8: Qualitative schematic (not to scale) of the experimental timing of the focusingimaging method. (a) The rotating bias field is switched back after the second microwave(ARP) pulse. The magnitude of the bias field is changed from its original value of 7.3
G to 35 G. Therefore, the trapping frequencies of the TOP trap change from (ωr , ωz ) =2π(10.2, 22.4)Hz to 2π(4.5, 10.8)Hz. We switch to a weaker trap to increase the focused time
tf , and thus, increase the sensitivity of the imaging method. (b) The second microwave(ARP) pulse is on for 400 µs. (c) The optical lattice is turned off after the first microwave
pulse. (d) The extra bias field is ramped down during the ARP pulse, then switched to 13 G
to bring the center of new TOP trap back to the imaging plane (its original position). The
magnetic fields (a), (d), and (e) are turned off before the probing occurs to avoid off-resonant
pumping. During the probing process, a pumping beam is added to pump atoms from the|1, −1i state to the |2, −1i state to interact with the probing beam. (f) Within 1 ms, after
most of the magnetic fields are turned off, the probing beam fires.

37

r(0) in a harmonic potential arer(t) = r(0)cos(ωr t) +ṙ(0)sin(ωr t),ωrṙ(t) = −r(0)ωr sin(ωr t) + ṙ(0)cos(ωr t),(3.23)(3.24)

where ωr is the trapping frequency of the confining harmonic potential. At the focused time,
tf = π/2ωr , one gets
mṙ(0) = mωr r(tf ),(3.25)

ωr r(0) = −ṙ(tf ).

(3.26)

From Eqs. (3.25) and (3.26), one can see the initial radial momentum-space distribution ismapped onto the final coordinate-space distribution at the focused time tf , and vice versa,
to within a scaling factor of mωr . At this point we can probe the cloud and reconstruct the
momentum distribution more accurately as compared to previous experiments [13, 14].
If the ratio of the trapping frequencies, ωz /ωr , is set to 2, at the focused time, theaxial motion of atoms makes a half oscillation and returns to its original size, while theradial motion of atoms bring them to the center of the trap. At the focused time, thedensity of the focused cloud might have an even higher density than it had before focusing.
The trapping frequency ratio in our experiment is set to a value slightly below 2. In thiscase, when the in-plane momentum-space distribution is mapped to coordinate space, thesize of the cloud in the axial direction is not yet at its minimum. Therefore, it lowers thedensity of the cloud at focused time such that the interactions do not influence the focussing
dynamics.

Chapter 4
Analysis of in-situ images

In our experiment, the thermal energy of the system, kb T , is on the same order with~ωz . Therefore, there will be populations on the axial excited state. The line of sight of ourimaging system aligns along the ẑ axis. Therefore, in the images we cannot distinguish those
populations from a ground-state population, which is the one that we are really interested
in. To extract a ground-state signal for BEC physics from our data, we calculate the axial
excited-state populations via a model, which is, in essence, the same model that the Parisgroup calls the Hartree-Fock, mean-field, local-density approximation [10].
4.1

The semiclassical mean-field modelFor bosons, the mean occupation of a single-particle state with wavevector k is given

1, where Ek is the energy of the state.
by the the Bose-Einstein distribution, Nk = e(Ek −µ)/kb T −1For our system, the thermal energy kB T À ~ωr , but kB T ∼ ~ωz . We treat the atomic motion semiclassically in the horizontal direction, while preserving discrete harmonic-oscillator
quantum levels in the z-direction. Thus, we get the two dimensional (2D) coordinate-spacedensity in the j th axial level1nj (~r) = 2hZZ1d2 p~ (ε(~p)+θj (~r)−µj )/k T
,
b
e
−1

(4.1)

39
2

p
where h is Planck’s constant, ε(~p) = 2m
is the free particle dispersion, and the local chemical

potential for the j th level is given by
X 4π~2
1µj = µglobal − mωr2 r2 − j~ωz −2(
as flj nl (~r)).
2
m
l6=j

(4.2)

The intralevel interaction energy is given byθj (~r) = 2(

4π~2 as)fjj nj ,m(4.3)

and the relevant mean-field interaction energies depend on fjl , which are the appropriately
normalized overlap integrals over the axial dimension between densities associated with theaxial quantum states j and l. With per-particle interaction energies comfortably less than
the axial energy spacing ~ωz , it is safe to treat interlevel axial interactions with what is, inessence, first-order perturbation theory. Evaluating the integral in Eq.(4.1), one getsnj (~r) = − ln(1 − e(−(θj (~r)−µj (~r))/kb T ) )/λ2db ,qwhere the deBroglie wavelength λdb =

2π~2
.
mkb T

(4.4)

For any given value of µ and T , nj aredetermined self-consistently. For kb T . ~ωz , the model converges in just a few iterations.
To get a better understanding of Eq. (4.4), one can break down the model into three
approximations: (1) the local-density approximation (2) the semiclassical approximation,and (3) the mean-field approximation. These approximations are discussed below:
4.1.1

Local-density approximation

The local-density approximation can be applied to a system where the density of theatomic cloud does not change significantly over the distance of the de Broglie wavelength of
an atom. This approximation works surprisingly well with magnetic traps, especially a weak
trap, like the one with ωr = 2π ∗ 10 Hz in our experiment . The local-density approximation
divides the system into many small boxes. Atoms inside each box are treated as if they werein a homogeneous potential with a potential offset determined by the external potential.

40

In this way, we can define a “local” chemical potential, a “local” density, and a “local”
phase-space density in each box.
4.1.2

Semiclassical approximation

The benefit of using the local density approximation is that it is very easy to find theenergy excitations for a free particle inside a box. For an ideal Bose gas confined in a box
potential, the total number of atoms in excited states as a function of the chemical potential
is given by
Nex =

X
i>0

1
e(Ei −µ)/kb T − 1

,

(4.5)

2

~
where Ei = 2m
( 2π~n )2 , and µ is the chemical potential. Replacing the summation in Eq. (4.5)
L iwith the integral, one getsLqNex =(2π~)qZdq p~1e(ε(p)−µ)/kb T − 1

,

(4.6)

where ε(p) = p2 /2m, and q is the dimensionality of the system. As one can see, the dispersion
relation is a classical single-particle dispersion. This is the “classical” part of the model. For
two dimensions, q = 2, and the integral in Eq. (4.6) can be carried out, leading to Eq. (??).
In a 2D harmonic trap, as a result of the local-density approximation, the chemical
potential in Eq. (??) is offset by the external potential; it should be replaced by the “local”
chemical potential, i.e.,1µ0 = µ − mωr2 r2 .
2

(4.7)

This substitution leads to the semiclassical model for an ideal Bose gas in a 2D harmonic
trap without mean fields, i.e.,
nex =

1
1
2 2
ln(1 − e( 2 mωr r −µ)/kb T ).
2
λdb

(4.8)

41

4.1.3

Mean-field approximation

Considering an interacting Bose gas, the easiest way to take into account the interactions might be via the introduction of mean fields. In the framework of mean-field theory,a many-body system with interactions is replaced with a free particle moving in an effective
potential. The mean field in our system consists of two parts: Intralevel and interlevel interactions. The “level” here means the energy levels of the axial-confining potential, which
we assume is quadratic.
Consider a quasi-2D system where the horizontal-confining potential is homogeneousand the axial-confining potential is quadratic. If Nj is the number of atoms in the j th axial
excited state, the density in the j th level isNj|ψj (z)|2 ,
L2

(4.9)

where L2 is the area of the quasi-2D system. An atom in the j th axial level feels a mean fieldcreated by the other atoms in the same level:Z
U 0 Njθj =2 2dz |ψj (z)|4L=2(U0 fjj )nj(4.10)=2ujj nj ,where nj is the 2D density and U0 = 4π~2 as /m. The factor of two arises from an implicit
assumption that the second-order correlation function at zero distance is two, as it would be
for fully fluctuating nondegenerate ideal bosons. The overlap integral represents the joint
probability of finding two atoms in the same level at the same position, z. As one can see
from Eq. (4.10), the 3D interaction strength U0 is renormalized by the overlap integral, fjj ,to get the effective 2D interaction strength ujj . The assumption that the scattering process
is still 3D is verified by the fact that the axial length of the quasi-2D system is still muchgreater than the 3D-scattering length, i.e., lz À as [41].

42Similarly, for an atom in the j th level, the mean field created by all the atoms in the
ith level is 2uij ni . uij is defined as
uij =U0 fij
Z=U0 dz |ψi (z)|2 |ψj (z)|2 .

(4.11)

Thus, the total interlevel mean field isX

2uij ni .

(4.12)i>j

For convenience, the interlevel mean-field energies are bundled with the chemical potential µ,radial-confining potential 12 mωr2 r2 , and the axial-excitation energy j~ωz such that the local
chemical potential isX
12uij ni (~r),µj = µ − mωr2 r2 − j~ωz −2i>j(4.13)

which is is identical to Eq. (4.2). With Eqs. (4.8), (4.10), and (4.13), one recovers Eq. (4.4).
4.1.4

The “no-condensate” model

Looking at Eq. (4.4), it’s clear why we call this the no-condensate model. For a given
value of T, we can choose a value of µ0 such that the model predicts an arbitrarily large n0 ,even one for which n0 λ2db À 1, without requiring a singularity at p = 0. More importantly,
however, this is a model from which all the many-body “magic” associated with degeneratebosons has been intentionally omitted. There’s no hint of the presence of anything like aphonon or other collective excitation. Still, the no-condensate model should do very wellwhere the phase-space density ρj ≡ nj λ2db < 1. As for the calculated value of n0 (~r), a
comparison with the naive no-condensate n0 allows us to quantify the tell-tale discrepancy.
Several interesting features can be found by solving Eq. (4.4). Figure 4.1 shows thesolutions of Eq. (4.4) graphically for five different chemical potentials. First, from Fig. 4.1,
one can see that as µj /kT reaches 0.24 and beyond, the density that satisfied Eq. (4.4) is43well approximated by (π/g̃)(µ0j /kb T )λ2db . With the fact that µ0j ∼ µ − 21 mωr2 r2 , the density
profile is an inverse parabola, as shown in Fig. 4.2. Since Eq. (4.4) is a no-condensate model,
the appearance of bimodality does not imply the existence of a condensate phase, as opposedto the 3D case. Second, one can see from Fig. 4.1 that there always exists a solution for Eq.
(4.4) for any chemical potential. This is another way of saying that this model predicts nocondensate.
4.2

Compressiblity
The isothermal compressibility is a quantity that measures a system’s ability to becompressed at fixed temperature. The definition of the isothermal compressibility is given
by
1
κ≡−
V

µ

∂V
∂P

¶
.

Using a thermodynamic relation, one can prove that
µ
¶
N̄ 2
∂ N̄
=
κ,
∂µ T,V
V
µ ¶
∂ n̄
or
= n̄2 κ,∂µ T,V(4.14)(4.15)

where N̄ is the average number of atoms in the system, and n̄ = N̄ /V . The quantity
µ

∂ n̄
∂µ

¶
,

(4.16)

T,V

which we call scaled compressibility, is extracted from our data and compared to the nocondensate model.
An analysis of a coordinate-space image proceeds as follows. The interesting physics
lie in n0 (~r), but we measure a density integrated through all the axial states, nmeas (~r) =Pnj (~r). Initially, we get an estimate of T , µ, and the OD scale (a multiplicative factor
that relates the observed OD to the nmeas ). In principle, the OD scale is calculable from thecrosssection of an atom absorbing light from a probing laser beam. However, drifts in the

44

30

25

y

20

15

10

5

0

0

5

10

15
20
phase space density (PSD)

25

30

Figure 4.1: The solution of Eq. (4.4) for five different chemical potentials. The blue curves
represent y = − ln(1 − e(−(θj (~r)−µj (~r))/kb T ) ) for five different chemical potentials µ0 /kb T , 0.14, 0 , 0.1, 0.24, and 0.5 (from left to right). The black curve represents y = nλ2db . Theintersection for each blue curve and the black curve gives the density that satisfies Eq. (4.4).
The three dashed lines indicate the values of the phase-space density where the last threeblue curves diverge.

45

25
20

PSD

15
10
5
0
−5
0

50

100

150

200

250

r (µm)

Figure 4.2: Coordinate-space distributions calculated from Eq. (4.4). The blue curve iscalculated with a chemical potential µ = 50 nK and a temperature T = 70 nK. The black
curve is calculated with µ = 6 nK, and T=180 nK. The dashed line is an inverse parabola
for comparison with the central part of the distribution from the blue curve.

46

efficiency of our microwave-selection pulse cause uncertainties in the OD scale. The estimate
comes from comparing the results of the no-condensate model to observed nmeas and fitting
the parameters T , µ and OD scale to the low phase-space regions of the cloud. With thechemical potential and temperature obtained from the fit, we can use Eq. (4.4) to evaluate
the excited-state distributions. Then, we find numerically a self-consistent solution to solveP
nj and n0 , where n0 = nmeas − nk , as shown in Fig. 4.3.
Once the corrected ground-state distribution, n0 , is extracted from nmeas , we calculate
the scaled compressibility, κs , at each imaging pixel, according to
¯
¯−1
dn0
dn0 ¯¯dµ0 ¯¯κs ==∗,dµlocaldr ¯r=i dr ¯r=i

(4.17)

where i is the number index for camera pixels. Although µ is not a quantity we can know
with great accuracy in a model-independent way, dµ0 /dr = −dVef f /dr ≈ −mωr r is readily
determined quite precisely, as the corrections to Vef f arising from the mean field of the axially
excited atoms are small and calculable. dn0 /dr is determined pixel by pixel from the values
extracted from our coordinate-space images. Combining the two spatial variations for eachpixel according to Eq. (4.17), we then get the scaled compressibility at every discrete radiusin our images. We plot the result vs the local ground-state phase-space density, ρ0 = n0 λ2db ,in Fig.4.4.
We compare our empirically scaled compressibility with the value given by the nocondensate model at the same density. For an observed value of n0 , we numerically solve
the no-condensate prediction, n0 = − ln(1 − exp(µ0 − θ0 (n0 ))/kb T )/λ2db , for µ0 . We thennumerically determine how n0 would change for a small change in µ0 . Thus, we extract ascaled no-condensate compressibility, which is plotted in Fig.4.4.
As a check of the local-density approximation that is central to our analysis, we show
images from two very different classes. The first, Fig. 4.4(a), for images with T = 171 nK
and a central ρ0 of about 9, and the second, Fig. 4.4(b), for images with T = 128 nK and acentral ρ0 of about 30. In both cases, the shape of the data is the same, and, in particular,phase space density ( ρ = n λ2)

47

12
10
8
6
4
2
0
0

50

100
r (µm)

150

200

Figure 4.3: Azimuthally averaged data from a coordinate-space image. The thin line shows
the measured density nmeas obtained from an in-situ image. The dashed lines are the calculated distributions from excited axial levels nj , j=1 to 5, from top to bottom. The thickline shows the pure axial ground-state distribution after correcting for nj>0 .

48

the value of ρ0 for which the empirical κ begins to deviate strongly from the no-condensate
model κ in both cases is about 4.
A note on the overall reliability of our data: when we allow our fitting routine toconsider only low phase-space regions of images, so that we can extract model-independent
fit values of T , µglobal , and the OD scale, we find that the fits are insufficiently orthogonal in
the latter two parameters. Large uncertainties in µglobal and OD scale result. We improve
the precision by fixing the OD scale using independent measurements of very low T clouds
in which the atoms are in a Thomas-Fermi inverted parabola with negligible noncondensed
wings. In this limit, we make the assumption that µglobal =u00 n0 (0). Fixing our OD scalewith this assumption means that our scaled values of κ, as shown in Fig. 4.4, are constrained
to saturate to 2.0 at very high ρ0 . In essence, we get high-precision measurements of the
behavior for exotic intermediate values of ρ0 , for the price of assuming a priori a goodunderstanding of low-degeneracy, mean-field behavior at high T , and long-coherence purecondensate behavior at low T .
4.2.1

Discussion

From the thermodynamic relations [42, 43], the number fluctuation for an ideal Bosegas can be expressed in terms of scaled compressibility (∂nk /∂µ), i.e.,hn2k i − hnk i2 = kb T

∂hnk i,∂µ(4.18)

where k represents the k th excited state for a 2D box. The mean occupation number for thestate k is given by Bose statistics as
nk =

1
2 2

k
−µ)
( ~2m

e

−1

.

(4.19)

Differentiating Eq. (4.19) with respect to µ, one can rewrite Eq. (4.18) as(∆nk )2 = hnk i + hnk i2 ,(4.20)

49

2.5

(a)

1

2.0

0.01

1.5
1.0

1

10

/ (1 / 2 u

00

)

1E-4
0.1

0.5
0.0
2.5
1

(b)
2.0

0.01

1.5

1.0

1E-4
0.1

1

10

0.5

0.0
0.1

1
0

n

0

10
2
db

Figure 4.4: Compressibility κ vs phase-space density ρ0 . (a) Measured κ extracted fromimages of samples at the same T as for the image in Fig. 5.2(c). Black circles are data
averaged over the values calculated from images of three separate clouds. The blue curves
are κnc derived from the no-condensate model. (b) Same but with κ extracted from imagesof lower-T samples. The paired vertical dotted lines indicate the location of the jump in
coherence discussed in the text.

50

where ∆nk =p

hn2k i − hnk i2 . One can see that the number fluctuation of state k comes fromtwo parts: The first term in Eq. (4.20) is atom shot noise, and the second term is excess
noise coming from the bunching effect [44, 45]. If there are N atoms occupying α states ina 2D box, one can show, using Eq. (4.20), that the number fluctuation of those N atoms inthe 2D box can be written ashN 2 i − hN i2 = hN i +α−1X
hnk i2 ,(4.21)

k=0

where N =

P

nk .

k

What we measured In Fig. 4.4 is the scaled compressibility of atoms in every camera
pixel projected in coordinate space. Therefore, the measured compressibility is
¯
∂N ¯¯
,
∂µ ¯i

(4.22)

where i is the number index of camera pixels. Using Eq. (4.18), (4.20), and (4.21), one canshow
X ∂hnk i
∂N
=
∂µ
∂µ
k
¢
1 X¡
=
hnk i + hnk i2
kb T k
=(4.23)

1(∆N )2 .
kb T

Therefore, one can that see our measured scaled compressibility is proportional to the number
fluctuation, and the proportionality constant is kb T . I would like to emphasize that the above
derivations come from considering an ideal Bose gas model. If interactions exist in the regime
close to degeneracy, the above relations will become inaccurate. The number fluctuation canbe directly measured from the experiment [46]. However, those measurements require carefulcalibration of the imaging system and a more stringent control over experimental conditions.
The plateau in Fig. 4.4(b) (ρ0 & 8) can be understood as follows. We know that thechemical potential is defined as the energy change due to adding one particle into the system

51

while keeping the entropy fixed, i. e.,
¯
∂E ¯¯
µ=
.
∂N ¯s

(4.24)The entropy of the system is kept constant if the atom is added into the ground state, andthere is only one ground state. In this case, the only energy associated with adding one
particle is the interaction energy with other atoms in the ground state, i.e.,µ = n0 u00 + 2

X

ui0 ni ,(4.25)

i>0

where uij is defined in Eq. (4.11). If the temperature of the system T is much smaller than
TBKT , one can ignore the contribution from the second term on the right-hand side of Eq.
(4.25). Then, the scaled compressibility is given by
∂n
1
∼
.
∂µ
u00

(4.26)

Therefore, for a pure condensate, one expects the value of the plateau in Fig. 4.4 to be 1/u00 .
Of course, we won’t reach that value at any finite temperature. However, the plateau indicates a macroscopic population in the ground state, thus, the existence of a quasicondensate
in our system.
The plateau can also be understood as the interactions suppressing the number fluctuations. Assume there are N particles confined in a 2D box, and those atoms equally occupy
the lowest α states. From Eq. (4.21), one can show that the number fluctuation of the box
is given byhN 2 i − hN i2 = hN i +

hN i2
.
α

(4.27)In the case of hN i À 1 and α ∼ O(1), one would expect the number fluctuation of thesystem to show a quadratic dependence on the mean number of atoms hN i. However,number (or density) fluctuations require energy when interactions become important. In
this case, these interactions will suppress the number (or density) fluctuation. Therefore,

52

the plateau indicates that the interactions in the system have become important. This also
signals the existence of a quasicondensate.

Chapter 5
Analysis of focusing images

For an infinite homogenous system, the first-order correlation function is the Fouriertransform of the momentum-space distribution. Thus, the inverse width of the momentumspace distribution gives the spatial extent of the coherent fraction in the system. At thetemperature T . TBKT , the first-order correlation function starts to decay algebraicallyinstead of exponentially. In a small system like ours, the system is still pretty coherent.
Therefore, one expects a sharp peak with a width of ∼ 1/R, where R is the size of coherent
part of the cloud, to appear in a momentum-space distribution, as shown in Fig. 5.1(a).
In Fig.5.2, we show the azimuthally averaged coordinate-space and momentum-space
distributions side by side. These images are taken back and forth between coordinate spaceand momentum space to ensure that the effects of temporal drifts in the experiment areminimized. As the temperature decreases, a sharp peak emerges from the center of the
momentum-space density distribution at T = 171 nK, as shown in Fig.5.2(g); however, down
to this temperature, there is no obvious change in the shape of the ground-state distribution
in coordinate space, as shown in Fig.5.2(a)-(c). A more peaked distribution in coordinate
space appears at a lower temperature, as shown in Fig.5.2(d). The width of the peak inmomentum distribution is inversely proportional to the spatial extent of coherence of thehigh density part of the cloud. One would expect a delta-function-like distribution for ahighly coherent system at T ¿ TBKT . However, because of finite imaging resolution, we are
not able to see that. Still, the width we observe from the images gives a spatial extent ofFigure 5.1: (a) The momentum-space density distribution at T . TBKT and (b) the in-trap density distribution under the same
experimental conditions. The azimuthally averaged data of both the momentum-space and the coordinate-space distributions
are shown in Fig. 5.2 (g) and 5.2 (c), respectively. At T = TBKT , a sharp spike emerges from the center of the momentum-space
distribution, signaling the appearance of a coherent fraction in the system. At T = TBKT , there is no obvious bimodality inthe coordinate-space density distribution.
54

55

coherence of about 9 µm in diameter, which is comparable to the spatial extent of the highdensity region in the cloud. From coordinate-space images taken under the same conditions
for which the coherence spike first appears in momentum space, we determine that coherence
happens when the central phase-space density ρ0 reaches a value of 8.0(0.7) [47] (this criticalvalue ρc is determined from looking at many more pairs of images than are presented in Fig.
5.2).
We emphasize that from the coordinate-space distribution alone, the identification ofa transition temperature would require model-dependent analysis of the smoothly varying
distribution. With access to both distributions at once, we readily see that a modest change
in the central phase-space density of < 15% causes the distribution at p = 0 to jump by afactor of three. Unambiguous, qualitative bimodality appears only in the coordinate-space
images at values of ρ0 that are 50% higher.
5.1

Mean-field model for momentum-space distribution
Can we model momentum-space density distributions? One reasonable starting point

for calculating the momentum-space density distribution is using the ideal-gas model. In
this model, there is no phonon, no vortex, no mean field, and no quasicondensate. Not onlyis this the simplest model, but also it gives us a good estimate of the temperature becauseinteractions have little effect on high-energy atoms. Using the semiclassical, local-density,and mean-field approximations, the number of atoms with momentum p at position r isgiven by
N (p, r) =1
e(p2 /2m+mωr2 r2 /2−µ)/kb T − 1

,

(5.1)

where ωr is the radial trapping frequency of the harmonic trap. Integrating Eq. (5.1) overthe coordinate space, one then gets the momentum-space density distribution, i.e.,
µ
¶
p2
1
−( 2m
−µ)/kb Tn(p) =
ln 1 − e.
(λmωr )2

(5.2)

56

15

15
T =

10

(e)

10
5

0
15

0
15
T =

176 nK

(b)

10

5
0
15
T =

10

171 nK(f)(c)

Density (a.u.)

db

2

(a)

5

10

= n

196 nK

5
0
15
(g)

10

resolution

5

5

0
15

0
15
T =

10

155 nK

(d)

5

0

0
50

100

r ( m)

(h)

10

5

0

imaging

150

0

1

2
p

3
-1

4

( m )

Figure 5.2: (a)-(d) Coordinate-space distributions and corresponding (e)-(h) momentumspace distributions. Two distributions in the same row are taken under near-identical conditions. The thin black curves give the azimuthal averages of nmeas from the raw images.
The thick blue curves in the coordinate-space distributions are the ground-state distribution
n0 after correcting for nk>0 . The spike in momentum that first appears in (g) has no corresponding dramatic change in coordinate space (c). The vertical dotted line in (c) represents
the inverse of the momentum resolution limit indicated in (g) and is thus a lower limit on
the coherence length of the population of low-p (high-coherence) atoms represented by thearea (about 1.4 % of total) under the spike in (g).

57

Figure 5.3 shows a fit to our data using Eq. (5.2). In the fitting routine, temperature T
is the only fitting parameter. The fitting points are weighted linearly according to theirdistance to the center of the system. Because Eq. (5.2) does not include contributions
from density-dependent mean fields, quasicondensates, or phonons, we therefore exclude thepcentral region, r < 2µ/mωr2 , from the fitting routine.
At large r, Eq. (5.2) shows the asymptotic behavior as
µ

1

2 2

ηe kb T e− 2 mωr r .

(5.3)

One can see in Eq. (5.3) that the OD scale η and the chemical potential µ are strongly
coupled. Therefore, it is very difficult for the fitting routine to get both the correct ηand the correct µ at the same time. However, one can at least get the temperature rightfrom fitting at large r. In this case, the signal-to-noise ratio in the tail of the distribution
becomes very important. We improve the precision by fixing the OD scale using independentmeasurements of very low T clouds in which the atoms are in a Thomas-Fermi inverted
parabola with negligible noncondensed wings.
To include the mean-field contribution, I first use Eq. (4.4) to calculate the in-trap
distributions for a trial set of T and µ. Once the in-trap distributions nj are found, the total
mean-field potential at the j th axial level, Mj , can be calculated, which gives
Mj (r) =X

2uij ni (r).

(5.4)

i

Then, according to Bose-Einstein statistics, the atoms with momentum p at the j th level are
given by1nj (p) =(2π~)2

Z ∞
0

d2 p

1
e

p2
( 2m
−µj )/kb T

−1

,

(5.5)

where µj = µ − 12 mωr2 r2 − Mj (r). Considering the fact that the image is taken in coordinate
space, one has to multiply the coordinates in real space by mωr to get the momentum
coordinates. Also, one has to consider that an image is taken pixel by pixel. With these

58

60

40

0
10

n

2D

(µm−2)

20

5

0
0

100

200

300

400

r (µm)
Figure 5.3: The ideal-gas model fit to the azimuthally averaged data that is shown in Fig.
5.2(g). The full-scale plot is shown in (a), while in (b), only 0 < n2D < 10 is shown toemphasize the departure of the measured data (Black curve) and the fitting curve (Bluecurve). The five red curves represent nj (p) with j = 0 to 5 from top to bottom. In this fit,
the chemical potential is fixed at zero, and a central area with a radius 80 µm is ignoredby the fitting routine. The boundary of this ignored area indicates the highest possible
momentum for a phonon in the system. The fitting curve shows a sharp peak at the centerof the distribution, which comes from the distribution of the lowest axial level, j = 0. As
one can see in the figure, the fitting is not satisfactory.

59

consideration, Eq. (5.5) can be written in the form ofnj (p) =1 X
1∗ 2πiδ 2 ,p2
12(2π~) i e( 2m −(µ− 2 mωr2 (iδ)2 −M (iδ)))/kb T − 1(5.6)

where i is the the number index of camera pixels, and δ is the pixel size in momentum
space. After summing up all the j th that are populated, one can compare the resulting total
momentum-space density distribution to the azimuthally averaged data pixel by pixel.
To summarize the fitting procedure, one has to (1) choose a trial set of µ, T , and OD
scale; (2) find a consistent solution for nj (r); (3) from nj (r), calculate M (r); (4) once M (r) is
found, calculate the j th -level momentum-space distribution; (5) summing contributions fromall j levels, compare the resulting value to the data and calculate the total deviation-squarevalue, and (6) stop the routine if a stable minimum of the total deviation-square value isfound; otherwise return to step 1.
Figure 5.4 shows the fit to a focusing image. As one can see from Fig. 5.4(b), thetwo curves from the data and the mean-field model start to resolve at the momentum whichcorresponds to the chemical potential. Since the mean-field model does not include any
phonon contribution, this deviation is understandable. Although this mean-field model inmomentum space can not account for contributions from phonons and vortices, this model
gives an accurate description on the nondegenerate part of the system. This part of thesystem is what we are least interested in. However, by subtracting the mean-field model
curve from a measured momentum-space density distribution, one can obtain the part inan image that is beyond the mean-field description. For an image taken at T . TBKT , theresidual distribution has a Gaussian-like shape. This residual density distribution includes
the contribution from a quasicondensate as well as atoms in the fluctuating region, as shownin Fig. 5.5. We then fit the residual distribution to a Gaussian, and determine the extent of
the coherence by the inverse of the fitted-Gaussian width. Considering the broadening fromimaging errors to the measured density distribution, the measured width of the sharp spikein the density distribution is overestimated. Therefore, the broadened width sets the lower

60

60

(a)

n2D (µm−2)

40
20
0
10

(b)
5

0
0

100

200
r (µm)

300

400

Figure 5.4: The mean-field model fit to the azimuthally averaged data shown in Fig. 5.2(g).
In the fitting routine, chemical potential µ and temperature T are the fit parameters. The
full-scale plot is shown in (a), while in (b), only 0 < n2D < 10 is shown to emphasize thedeparture of the measured data (black curve) and the fitting curve (blue curve). The five
red curves represent nj (p) with j = 0 to 5 from top to bottom. One can see the sharp central
peak in the fitting curve of Fig. 5.3 is suppressed by mean fields. After including mean-field
contributions, the fit is now satisfactory.

61bound of the real spatial extent of coherence. Considering our imaging resolution, the minimum detectable width of a Gaussian distribution in momentum space implies that coherence
extends over a central disk of radius 4.5 µm. Figure 5.6 shows the fitting temperatures ofseveral focusing images. The temperatures from fitting momentum-space images shows agood agreement with the temperatures obtained from fitting of in-situ images, of course, as
one expects.

62

50

Residual density ( µm−2)

40

30

20

10

0
0

20

40

60

80

100

r (µm)
Figure 5.5: Residual density distribution. The density distribution is obtained by subtracting
the mean-field model fit (blue line) from the measured density distribution (black line) in
Fig. 5.4. The distribution includes the atoms that are beyond mean-field description. The
width of the distribution is limited by the imaging resolution, not by the coherence of thesystem.

63

240

200

T (nK)

160

120

80

40

0
5.7

5.8

5.9

6.0

f

RF

6.1

6.2

6.3

6.4

(MHz)

Figure 5.6: Temperatures obtained from fitting focusing images vs final rf-knife frequency.
The temperatures fit by the mean-field model in momentum space show a good agreementwith the temperature obtained from fitting in-situ images in coordinate space. The red line
comes from the fitting of the temperatures shown in Fig. 3.4

Chapter 6
Conclusion and outlook

6.1

Conclusion
In Fig. 4.4, we see that the deviation of the measured scaled compressibility κ fromthe mean-field predicted value κnc starts at a phase-space density of ρ0 ∼ 3. From Fig. 5.2,
we see that the long coherence does not appear in an atomic cloud until κ/κnc = 1.7, or
the phase-space density ρ0 = 8. When 3 < ρ0 < 8, a very interesting regime is identified.
A 2D Bose gas in this regime is beyond the mean-field description, but it does not yet
possess a coherent fraction, or quasicondensate. This regime, which does not exist in three
dimensions, is a product of the enhanced interactions associated with reduced dimensionality.
The appearance of a long coherence in a 2D Bose gas starts at ρ = 8(0.7), which is in a very
good agreement with the prediction of the BKT transition from the Monte-Carlo simulation[15].
We tried very hard to present our data in a model-free fashion. Of course, there areassumptions in our analyses. We assumed that the semiclassical, mean-field, and localdensity approximations are good for describing the atomic populations in the axial-excitedstates. We don’t assume any model for the ground-state population. Nor as we assume thephase-space density at which the BKT transition occurs. All our measurements are done
locally, pixel by pixel. We do not measure global quantities, e.g., the critical number Nc of thesystem for the BKT transition. I would like to emphasize that the goal of our measurements

65

was not to prove the existence of the BKT transition in a 2D system, which had alreadybeen proven first in 2D systems [48, 49, 50] and then in a 2D ultracold-atom system [9].
Our goal for the experiment described in this thesis was to provide better measurements tocompare with theory.
6.2

Outlook
How interactions participate in the statistical behavior of a many-body system is always

an interesting topic. In lower dimensional systems, enhanced interactions come from reduced
dimensionality. In an extreme 2D case, the microscopic collision behavior of atoms willbe also influenced by reduced dimensionality. These properties make a 2D ultracold-atomsystem an extremely interesting playground for physicists. I think there are still many
interesting experiments to be done and many interesting questions to be answered. From
my perspective, I list some of these experiments and questions here.
(1) The first possible extension of this work would be to go to a more strongerly interacting regime, where g̃ & 1. In this regime the interaction between atoms is so largethat any quasicondensate would be seriously depleted. How the quasicondensate
becomes depleted as a function of g̃ is an interesting question.
(2) The second possible extension would be to measure directly the first-order correlation
at a temperature T < TBKT . The theory tells us the first-order correlation function
has an algebraic decay on the separation between two points in the system. Can wemeasure it? How does the decay change as a function of temperature? How does
the finite-size effect influence the first-order correlation function?
(3) How does a 2D quasicondensate establish its phase coherence? By controlling how
quickly the system enters the BKT regime, one can study couple interesting questions: How much time does the system need to establish its coherence? How many

66

excitations will this process induce into the system as a function of the ramping
speed?
(4) What are the properties of a vortex state in a 2D system? A vortex lattice shouldbe melted by the increasing effect of thermal fluctuations. Can one see the melting
process? How does a vortex lattice melt? There are also predictions that at anextremely high rotation rate and at extremely low temperature, the system willshow a quantum Hall-like behaviors. Can we realize this regime?

Bibliography[1] Both appendix A and B have already appeared in print in Phys. Rev. Lett .
[2] S. Tung, V. Schweikhard, and E. A. Cornell, Observation of vortex pinning in BoseEinstein condensates, Phys. Rev. Lett. 97, 240402 (2006).
[3] V. Schweikhard, S. Tung, and E. A. Cornell, Vortex proliferation in the BerezinskiiKosterlitz-Thouless regime on a two-dimensional lattice of Bose-Einstein condensates,
Phys. Rev. Lett. 99, 030401 (2007).
[4] V. Schweikhard, I. Coddington, P. Engels, S. Tung, and E. A. Cornell, Vortex-latticedynamics in rotating spinor Bose-Einstein condensates, Phys. Rev. Lett. 93, 210403(2004).
[5] I. Coddington, P. C. Haljan, P. Engels, V. Schweikhard, S. Tung, and E. A. Cornell,
Experimental studies of the equilibrium vortex properties in Bose-condensed gas, Phys.
Rev. A 70, 063607 (2004).
[6] V. L. Berezinskii, Destruction of long-range order in one-dimensional and twodimensional system possessing a continuous symmetry group-II Quantum systems, Soviet Physics JETP 34, 610 (1971).
[7] J. Kosterlitz, and D. J. Thouless, Ordering, metastability and phase transitions in two
dimensional systems, J. Phys. C: Solid State Physics 6, 1181 (1973).
[8] Z. Hadzibabic, S. Stock, B. Battelier, V. Bretin, and J. Dalibard, Interference of anarray of independent bose-einstein condensates, Phys. Rev. Lett. 93, 180403 (2004).
[9] Z. Hadzibabic, P. Krüger, M. Cheneau, B. Battelier, and J. Dalibard, Berezinskiikosterlitz-thouless crossover in a trapped atomic gas, Nature 441, 1118 (2006).
[10] Z. Hadzibabic, P. Krüger, M. Cheneau, S. P. Rath, and J. Dalibard, The trapped
two-dimensional bose gas: from Bose-Einstein condensation to Berezinskii-KosterlitzThouless physics, New J of Phys 10, 045006 (2008).

68[11] I. Shvarchuck, C. Buggle, D. S. Petrov, K. Dieckmann, M. Zielonkoski, M. Kemmann,
T. G. Tiecke, W. von Klitzing, G. V. Shlyapnikov, and J. T. M. Walraven, Bose-Einsteincondensation into nonequilibrium states studied by condensate focusing, Phys. Rev. Lett.
89, 270404 (2002).
[12] A. H. van Amerongen, J. J. van Es, P. Wicke, K. V. Kheruntsyan, and N. J. van Druten,
Yang-Yang thermodynamics on an atom chip, Phys. Rev. Lett. 100, 090402 (2008).
[13] P. Krüger, Z. Hadzibabic, and J. Dalibard, Critical point of an interacting twodimensional atomic bose gas, Phys. Rev. Lett. 99, 040402 (2007).
[14] P. Cladé, C. Ryu, A. Ramanathan, K. Helmerson, and W. D. Phillips, Observation of
a 2D bose gas: From thermal to quasicondensate to superfluid, Phys. Rev. Lett. 102,
170401 (2009).
[15] N. V. Prokofév, O. Ruebenacker, and B. V. Svistunov, Critical point of a weakly interacting two-dimensional Bose gas, Phys. Rev. Lett. 87, 270402 (2001).
[16] V. Schweikhard, Ultracold Bose gases under rotation, in lattice potentials, and both,
Ph.D. thesis, University of Colorado at Boulder (2008).
[17] C. J. Pethick, and H. Smith, Bose-Einstein condensation in dilute gases (Cambridge
university press, Cambridge, 2002).
[18] R. E. Peierls, Quelques propriétés typiques des corps solides, Ann. Inst. Henri Poincaré5, 177 (1935).
[19] N. D. Mermin, and H. Wagner, Absence of ferromagnetism or antiferromagnetism inone- or two-dimensional isotropic Heisenberg models, Phys. Rev. Lett. 17, 1307 (1966).
[20] P. C. Hohenberg, Existence of long-range order in one and two dimensions, Phys. Rev.
158, 383 (1967).
[21] O. Penrose, and L. Onsager, Bose-Einstein condensation and liquid Helium, Phys. Rev.
104, 576 (1956).
[22] L. Pitaveskii, and S. Stringari, Bose-Einstein condensation (Oxford university press,
Oxford, 2003).
[23] M. Naraschewski, and R. J. Glauber, Spatial coherence and density correlations oftrapped Bose gases, Phys. Rev. A 59, 4595 (1999).
[24] G. B. Arfken, and H. J. Webber, Mathematical methods for physicists (Academic press,
MA, USA, 2005).
[25] I. Bouchoule, N. J. van Druten, and C. I. Westbrook, Atoms chips and one-dimensioanl
bose gases, ArXiv:0901.3303 .

69[26] D. S. Petrov, M. Holzmann, and G. V. Shlyapnikov, Bose-Einstein condensation inquasi-2D trapped gases, Phys. Rev. Lett. 84, 2551 (2000).
[27] D. S. Petrov, D. M. Gangardt, and G. V. Shlyapnikov, Low-dimensional trapped gases,
J. Phys. IV 116, 5 (2004).
[28] Z. Hadzibabic, and J. Dalibard, Two dimensional bose fluids: An atomic physics perspective, ArXiv:0912.1490 .
[29] J. W. Kane, and L. P. Kadanoff, Long-range order in superfluid Helium, Phys. Rev.
155, 80 (1967).
[30] A. Gölitz, J. M. Vogels, A. E. Leanhardt, C. Raman, T. L. Gustavson, J. R. Abo-Shaeer,
A. P. Chikkatur, S. Gupta, S. Inouye, T. Rosenband, and W. Ketterle, Realization of
Bose-Einstein condensates in lower dimensions., Phys. Rev. Lett. 87, 130402 (2001).
[31] S. Stock, Z. Hadzibabic, B. Battelier, M. Cheneau, and J. Dalibard, Observation ofphase defects in quasi-two-dimensional Bose-Einstein condensates, Phys. Rev. Lett. 95,
190403 (2005).
[32] S. P. Rath, T. Yefsah, K. J. Günter, M. Cheneau, R. Desbuquois, M. Holzmann,
W. Krauth, and J. Dalibard, Equilibrium state of a trapped two-dimensional Bosegas, Phys. Rev. A 82, 013609 (2010).
[33] E. A. A. Polkovnikov, and E. Demler, Interference between independent fluctuatingcondensates, Proc. Natl. Acad. Sci. USA 103, 6125 (2006).
[34] M. Greiner, O. Mandel, T. Esslinger, T. W. Hänsch, and I. Bloch, Quantum phasetransition from a superfluid to a Mott insulator in a gas of ultracold atoms, Nature
415, 39 (2002).
[35] R. Grimm, and M. Weidemüller, Optical dipole traps for neutral atoms, Adv. At. Mol.
Opt. Phys. 42, 95 (2000).
[36] W. Petrich, M. H. Anderson, J. R. Ensher, and E. A. Cornell, Stable, tightly confining
magnetic trap for evaporative cooling of neutral atoms, Phys. Rev. Lett. 74, 3352 (1995).
[37] M. H. Anderson, J. R. Ensher, M. R. Matthews, C. E. Wieman, and E. A. Cornell,
Observation of Bose-Einstin condensation in a dilute atomic vapor, Science 269, 198(1995).
[38] J. R. Ensher, The first experiments with Bose-Einstein condensation of 87 Rb, Ph.D.
thesis, University of Colorado at Boulder (1998).
[39] I. R. Coddington, Vortices in a highly rotating bose-condensed gas, Ph.D. thesis, University of Colorado at Boulder (2004).

70[40] H. J. Metcalf, and P. van der Straten, Laser cooling and trapping (Springer-Verlag, New
York, 1999).
[41] In our experiment, lz /as ∼ 60 .
[42] L. D. Landau, and E. M. Lifshitz, Statistical physics, part 1 (Butterworth-Heinemann,
Oxford, 1976), third edn.
[43] K. Huang, Statistical mechanics (Jonh Wiley & Sons, Inc, New York, 1987).
[44] M. Schellekens, R. Hoppeler, A. Perrin, J. V. Gomes, D. Boiron, A. Aspect, and C. I.
Westbrook, Hanbury Brown Twiss effect for ultracold quantum gases, Science 310, 648(2005).
[45] T. Jeltes, J. M. McNamara, W. Hogervorst, W. Vassen, V. Krachmalnicoff,
M. Schellekens, A. Perrin, H. Chang, D. Boiron, A. Aspect, and C. I. Westbrook,
Comparison of the Hanbury BrownTwiss effect for bosons and fermions, Nature 445,
402 (2007).
[46] In the preparation of this thesis, I learned that Cheng Chin’s group at University of
Chicago were studying universal behavior of a 2D Bose gas. In their experiment, theymeasured the number fluctuation directly from the noise of their images. .
[47] Our observed value ρc = 8.0(0.7) is consistent with predictions from [15], which giveρc = 8.3. our observation that κ/κnc = 1.7 at ρ0 = ρc could be compared to theory via
e.g. the slopes of lines in fig. 3 in ref. [28] .
[48] G. Agnolet, D. F. McQueeney, and J. D. Reppy, Kosterlitz-Thouless transition in helium
films, Phys. Rev. B 39, 8934 (1989).
[49] A. I. Safonov, S. A. Vasilyev, I. S. Yasnikov, I. I. Lukashevich, and S. Jaakkola, Observation of quasicondensate in two-dimensional atomic hydrogen, Phys. Rev. Lett. 81,4545 (1998).
[50] D. Resnick, J. C. Garland, J. T. Boyd, S. Shoemaker, and R. S. Newrock, KosterlitzThouless transition in proximity-coupled superconducting arrays, Phys. Rev. Lett. 47,
1542 (1981).
[51] D. S. Petrov, G. V. Shlyapnikov, and J. T. M. Walraven, Phase-fluctuating 3D BoseEinstein condensates in elongated traps, Phys. Rev. Lett. 87, 050404 (2001).
[52] T. Kinoshita, T. Wenger, and D. S. Weiss, Observation of a one-dimensional TonksGirardeau gas, Science 305, 1125 (2004).
[53] T. Stöferle, H. Moritz, C. Schori, M. Köhl, and T. Esslinger, Transition from a stronglyinteracting 1d superfluid to a Mott insulator, Phys. Rev. Lett. 92, 130403 (2004).

71[54] K. W. Madison, F. Chevy, W. Wohlleben, and J. Dalibard, Vortex formation in a stirred
bose-einstein condensate, Phys. Rev. Lett. 84, 806 (2000).
[55] V. Schweikhard, I. Coddington, P. Engels, V. P. Mogendorff, and E. A. Cornell, Rapidly
rotating Bose-Einstein condensates in and near the lowest landau level, Phys. Rev. Lett.
92, 040404 (2004).
[56] J. R. Abo-Shaeer, C. Raman, J. M. Vogels, and W. Ketterle, Observation of vortex
lattices in Bose-Einstein condensates, Science 292, 476 (2001).
[57] N. R. Cooper, N. K. Wilkin, and J. M. F. Gunn, Phys. Rev. Lett. 87, 120405 (2001).
[58] D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and P. Zoller, Phys. Rev. Lett. 81,
3108 (1998).
[59] P. W. Stephens, P. A. Heiney, R. J. Birgeneau, P. M. Horn, D. E. Moncton, and G. S.
Brown, High-resolution x-ray-scattering study of the commensurate-incommensuratetransition of monolayer kr on graphite, Phys. Rev. B 29, 3512 (1984).
[60] P. Bak, Commensurate phases, incommensurate phases and the devil’s staircase, Rep.
Prog. Phys. 45, 587 (1982).
[61] T. Goto, T. Kimura, G. Lawes, A. P. Ramirez, and Y. Tokura, Ferroelectricity and giant
magnetocapacitance in perovskite rare-earth manganites, Phys. Rev. Lett. 92, 257201(2004).
[62] M. Baert, V. V. Metlushko, R. Jonckheere, V. V. Moshchalkov, and Y. Bruynseraede,
Composite flux-line lattices stabilized in superconducting films by a regular array ofartificial defects, Phys. Rev. Lett. 74, 3269 (1995).
[63] A. N. Grigorenko, S. J. Bending, M. J. V. Bael, M. Lange, V. V. Moshchalkov, H. Fangohr, and P. A. J. de Groot, Symmetry locking and commensurate vortex domain
formation in periodic pinning arrays, Phys. Rev. Lett. 90, 237001 (2003).
[64] J. E. Villegas, S. Savel’ev, F. Nori, E. M. Gonzalez, J. V. Anguita, R. Garcia, and J. L.
Vicent, A superconducting reversible rectifier that controls the motion of magnetic flux
quanta, science 3025, 1188 (2003).
[65] C. Reichhardt, C. J. Olson, R. T. Scalettar, and G. T. Zimányi, Commensurate andincommensurate vortex lattice melting in periodic pinning arrays, Phys. Rev. B 64,
144509 (2001).
[66] J. W. Reijnders, and R. A. Duine, Pinning of vortices in a bose-einstein condensate byan optical lattice, Phys. Rev. Lett. 93, 060401 (2004).
[67] H. Pu, L. O. Baksmaty, S. Yi, and N. P. Bigelow, Structural phase transitions of vortex
matter in an optical lattice, Phys. Rev. Lett. 94, 190401 (2005).

72[68] C. Wu, H. D. Chen, J. P. Hu, and S. C. Zhang, Vortex configurations of bosons in an
optical lattice, Phys. Rev. A 69, 043609 (2004).
[69] R. Bhat, L. D. Carr, and M. J. Holland, Bose-Einstein condensates in rotating lattices,
Phys. Rev. Lett. 96, 060405 (2006).
[70] Mask (b): φ1 = 2.5 mm and d1 = 11.5 mm. The lattice constant of the triangular opticallattice is 7.8 µm, giving ωc = 0.806 ωr .; mask (c): φ2 = 2.5 mm and d2 = 13.5 mm. The
lattice constant of the square lattice is 7.0 µm, giving ωc = 0.866 ωr. The focal lens of
two lenses are 300 mm and 250 mm, respectively. .
[71] I. Coddington, P. C. Haljan, P. Engels, V. Schweikhard, S. Tung, and E. A. Cornell,
Experimental studies of equilibrium vortex properties in a bose-condensed gas, Phys.
Rev. A 70, 063607 (2004).
[72] D. E. Sheehy, and L. Radzihovsky, Vortex lattice inhomogeneity in spatially inhomogeneous superfluids, Phys. Rev. A 70, 051602 (2004).
[73] The chemical potential is determined from the axial Thomas-Fermi radius of the condensate obtained from in-trap images, and is proportional to the peak density of thecondensate .
[74] There is a constant measurement offset to the values of θOL − θV L in this data set.
we have used in-trap images of vortex lattices to confirm that the absolute value ofθOL − θV L is zero for locked lattices .
[75] |s(ksq )|=a2 +A1 −A2
U −U ; A1 , A2 , U0 , and δU are fit to the data. .
1+Exp[ δU 0 ][76] L. J. Campbell, M. M. Doria, and J. B. Kadtke, Energy of infinite vortex lattices, Phys.
Rev. A 39, 5436 (1989).
[77] A. Leggett, Bose-Einstein condensation in the alkali gases: Some fundamental concepts,
Rev. Mod. Phys. 73, 307 (2001).
[78] D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and P. Zoller, Cold bosonic atoms
in optical lattices, Phys. Rev. Lett. 81, 003108 (1998).
[79] S. Richardr, F. Gerbier, J. H. Thywissen, M. Hugbart, P. Bouyer, and A. Aspect,
Momentum spectroscopy of 1d phase fluctuations in Bose-Einstein condensates, Phys.
Rev. Lett. 91, 010405 (2003).
[80] D. Hellweg, L. Cacciapuoti, M. Kottke, T. Schulte, K. Sengstock, W. Ertmer, and J. J.
Arlt, Measurement of the spatial correlation function of phase fluctuating bose-einsteincondensates, Phys. Rev. Lett. 91, 010406 (2003).
[81] R. Gati, B. Hemmerling, J. Fölling, M. Albiez, and M. K. Oberthaler, Noise thermometry with two weakly coupled Bose-Einstein condensates, Phys. Rev. Lett. 96, 130404(2006).

73[82] M. Tinkham, Introduction to Superconductivity (McGraw-Hill, Inc., New York, 1996).
[83] Within a dataset, the ramp-down rate is kept fixed, tr = τ × vOL /1.3 kHz, τ = 18 ms,
if not otherwise indicated .
[84] T. Simula, and P. Blakie, Thermal activation of vortex-antivortex pairs in quasi-twodimensional Bose-Einstein condensates, Phys. Rev. Lett. 96, 020404 (2006).
[85] A. Trombettoni, A. Smerzi, and P. Sodano, Observable signature of the BerezinskiiKosterlitz-Thouless transition in a planar lattice of Bose-Einstein condensates, New J.
Phys. 7, 57 (2005).
[86] Three circularly polarized laser beams (λ = 810 nm) intersect in a tripodlike configuration, with θ = 6.6◦ angles to the z-axis. Calculation of the optical dipole potential [35]includes counterrotating terms and interaction with both the d1 and d2 lines, as wellas the fictitious magnetic field due to the circular polarization. The tilted bias field of
the TOP trap makes P ≈ 0.5[35] .
[87] J. Javanainen, Phonon approach to an array of traps containing Bose-Einstein condensates, Phys. Rev. A 60, 4902 (1999).
[88] K. Burnett, M. Edwards, C. W. Clark, and M. Shotter, The Bogoliubov approach tonumber squeezing of atoms in an optical lattice, J. Phys. B 35, 1671 (2002).
[89] To avoid radial flows during vOL ramp-up, rT F is kept constant by balancing the latticeenhanced mean-field pressure with radial confinement due to the optical lattice envelope,by the choice of a 67µm 1/e2 intensity waist. Axial (z) confinement is due to the magnetic
trap alone .
[90] In the axial condensate region between z = −rz /3, +rz /3, where according to our 3D
GPE simulations 85% of the tunnel current is localized and hence the relative phase ismeasured, axial phase fluctuations [51] vary between ≈ 600 mrad (“cold” data in Fig.
B) and 800 mrad (“hot”) in the regime j/t ≈ 1 .
[91] j is averaged over junctions within the central 11% of the array area, from which allquantitative experimental results are extracted. .
[92] D. Ananikian, and T. Bergeman, Gross-Pitaevskii equation for Bose particles in a
double-well potential: Two-mode models and beyond, Phys. Rev. A 73, 013604 (2006).
[93] M. Naraschewski, and D. Stamper-Kurn, Analytical description of a trapped semi-ideal
Bose gas at finite temperature, Phys. Rev. A 58, 2423 (1998).
[94] Following [37], in simulations we count a vortex if all three phase differences in anelemental triangle of junctions are ∈ (0, π), or if all are ∈ (−π, 0) .

74[95] The very short annihilation time
√ of configuration i pairs [ τ < 5 ms in fig. B(c) ] is not
unexpected: their spacing d/ 3 = 2.8 µm is comparable to the diameter of a vortexcore
in the bulk condensate (after vOL rampdown) .
[96] L. Pitaevskii, and S. Stringari, Thermal vs quantum decoherence in double well trapped
Bose-Einstein condensates, Phys. Rev. Lett. 87, 180402 (2001).

Appendix A
Observation of vortex pinning in Bose-Einstein condensates [1, 2]

We report the observation of vortex pinning in rotating gaseous Bose-Einstein condensates(BEC). Vortices are pinned to columnar pinning sites created by a co-rotating optical latticesuperimposed on the rotating BEC. We study the effects of two types of optical lattice,triangular and square. In both geometries we see an orientation locking between the vortex
and the optical lattices. At sufficient intensity the square optical lattice induces a structural
cross-over in the vortex lattice.
Some of the most appealing results from recent work in superfluid gases have had to
do with lattices, either optical lattices [34, 52, 53] or vortex lattices [54, 55, 56]. These twokinds of lattices could hardly be more different! The former is an externally imposed periodic
potential arising from the interference of laser beams, while the latter is the self-organized
natural response of a superfluid to rotation. As distinct as these two periodic structures
may be, there are reasons for trying to marry them in the same experiment. For one thing,
the extreme limits of rapid rotation (in the case of vortex lattices) [57] and deep potentials(in the case of optical lattices) [58] both lead to the same thing: correlated many-bodystates. For another, there is considerable precedent, from various subdisciplines of physics,
for interesting effects arising from the interplay between competing lattices [59, 60, 61].
Moreover, the pinning of superconducting flux vortices to an array of pinning sites in solids
is an area of very active research as well [62, 63, 64, 65]. With these considerations in mind,
we undertook a preliminary experimental study of the effects of a rotating optical lattice

76
on a vortex lattice in a Bose-condensed sample of 87 Rb. The density of the superfluid issuppressed at the antinodes of the two-dimensional standing wave pattern of the opticallattice. These antinodes then become pinning sites, regions of low potential energy, for thesuperfluid vortices. Vortices can lower their interaction energy by arranging themselves tobe as far apart as possible from one another. The competition between these effects has
been examined in several theoretical works [66, 67]. Also [68, 69] discuss similar systems inthe strong interacting region
The setup for creating a rotating optical lattice is shown in Fig. A.1(a). A mask witha set of holes is mounted onto a motor-driven rotary stage, and a laser beam (532 nm) is
expanded, collimated, and passed through the mask. After the mask the resulting three
beams are focused onto the BEC. The interference pattern at the focus constructs a quasi2D optical lattice. The geometry and spatial extent of the triangular or the square opticallattice is determined by the size and layout of the holes and the focal length of the second
lens. For the pinning sites to appear static in the frame of a rotating BEC, the rotationof the two lattices must be concentric, and mechanical instabilities and optical aberrations(which lead to epicyclic motion of the pinning sites) must be particularly minimized. Evenso, residual undesired motion is such that the strength of the optical lattice must be kept
at less than 30% of the condensate’s chemical potential or unacceptable heating results overthe experiment duration of tens of seconds. We work perforce in the weak pinning regime.
The experiments begin with condensates containing ∼ 3 × 106 87 Rb atoms, held in the
Zeeman state |F = 1, mf = −1i by an axial symmetric magnetic trap with trapping frequencies {ωr , ωz } = 2π{8.5, 5.5}Hz. Before the optical lattice, rotating at angular frequencyΩOL , is ramped on, the BEC is spun up [55] close to ΩOL . This leads, before application ofan optical lattice, to the formation of a near perfect triangular vortex lattice with a random
initial angular orientation in inertial space. Through dissipation a vortex lattice can come toequilibrium with an optical lattice, with their rotation rates and angular orientations locked.

77

Figure A.1: (a) Schematic diagram of our setup for the rotating quasi-2D optical lattice.
Layouts of the masks for a triangular (b) and square (c) optical lattices. (d) and (e) are
pictures of triangular and square optical lattices, respectively. For details of the layouts see
[70].

78

In the absence of pinning sites, a vortex lattice with areal density of vortices nv will rotate
at (approximately [71, 72]) Ω = ( ~π)nv . This suggests that for an optical lattice with anmareal density of pinning sites nOL , locking between the two lattices will be facilitated if the
optical lattice rotates at the commensurate frequency Ωc ≡ ( ~π)nOL .
m
We measure the angular difference θOL − θV L between the orientation of the opticaland vortex lattice in reciprocal space (see Figs. 2(a)–2(b)). Fig. A.2(c) shows θOL − θV L asa function of the pinning strength with an optical lattice rotation rate ΩOL = 1.133Ωc =0.913ωr . The strength of pinning is characterized by the ratio Upin /µ (µ is the chemical
potential of the condensate [73]), which gives the relative suppression of the superfluid density
at pinning sites. We can see the initially random angular difference between the two lattices
becomes smaller as the pinning strength Upin /µ increases. For Upin /µ & 0.08 , the angulardifferences become very close to the locked value. Figure A.2(d) shows the phase diagram.
The data points and error bars mark the minimum pinning strength (Upin /µ)min above which
the lattices lock. We observe two distinct regimes. First, for small rotation-rate mismatch,(Upin /µ)min is rather independent of the rotation-rate mismatch. Second, for rotation-ratemismatch beyond the range indicated by the dashed line in Fig. A.2(d), angular orientation
locking becomes very difficult for any Upin /µ in our experiment. Instead, an ordered vortex
lattice with random overall angular orientation observed at low Upin /µ transforms into adisordered vortex arrangement at high Upin /µ.
This box–like shape of the locked region in Upin − ΩOL space is worth considering. In
a simple model, vortex motion in our system is governed by a balance of the pinning force→−
and the Magnus force. The pinning force is F pin (x) ∝ Upin /d, where Upin and d are thestrength of the pinning potential and its period, respectively. The Magnus force, acting on→−−→
a vortex moving with velocity →v vortex in a superfluid with velocity −v f luid is F mag (x) ∝
→
−
→
→
h
n(x) (−
v vortex − →v f luid ) × −κ where −κ = (m)ẑ, and n(x) is the superfluid density. A locked→−→→
vortex lattice will co-rotate with the pinning potential, giving −v vortex (r) = −v OL (r) = Ω OL ×

79

Figure A.2: (a) Triangular optical lattice and (b) vortex lattice in reciprocal space. Inset are
the corresponding original real-space CCD-camera images. (c) The difference in orientation
θOL −θV L versus the strength of pinning Upin /µ ( the peak of the optical potential normalized
by the condensate’s chemical potential) for the rotation rates ΩOL = 1.133Ωc = 0.913ωr .
With increasing pinning strength, θOL −θV L tends towards its locked value [74]. (d) Minimum
pinning strength needed for orientation locking between two lattices as a function of therotation rate of the optical lattice. The dashed and dotted lines are discussed in the text.

80
−
→
→
r , whereas the superfluid velocity in a solid-body approximation is −v f luid (r) = ~πn rθ̂ =m v−→→
Ω f luid × −r . Comparing the magnitudes of both forces at r = R(Ω)/2, where R(Ω) is thecentrifugal-force modified Thomas-Fermi radius, we obtain a minimum strength for pinning(Upin /µ)min ≈ ( 2√1 3 R(Ω)/d) × (ΩOL − Ωf luid )/Ωc .
What will be the fluid rotation rate Ωf luid in the presence of the pinning potential? On
the one hand, if vortices are tightly locked to the optical lattice sites, we have Ωf luid = Ωc .
The minimum strength (Upin /µ)min resulting from this assumption is plotted as solid line in
Fig. A.2(d). The lack of predicted decrease of (Upin /µ)min to zero around Ωc may be due
to long equilibration times in a very shallow pinning potential, as well as slight mismatches
in alignment and initial rotation rate of the BEC and the pinning potential. The ease oforientation locking with increasing rotation rate mismatch is less easy to explain in thismodel. On the other hand, in the weak-pinning regime, the vortex lattice can accommodate
a rotation rate mismatch by stretching/compressing away from the pinning sites. This allows
the fluid to co-rotate with the optical lattice (Ωf luid ≈ ΩOL ) and reduce the Magnus force.
This leads to a very low minimum pinning strength, as suggested by our data. However, thevortex lattice’s gain in pinning energy decreases rapidly in the locked orientation when themismatch between vortex spacing and optical lattice constant increases to the point where
the outermost vortices fall radially in between two pinning sites. Then the preference for the
locked angular orientation vanishes. This predicted limit is indicated by the vertical dottedlines in Fig. A.2(d).
In the absence of a pinning potential, the interaction energy of a square vortex lattice
is calculated to exceed that of a triangular lattice by less than 1% [76], thus it is predicted[66, 67] that the influence of even a relatively weak square optical lattice will be sufficient
to induce a structural transition in the vortex lattice. This structural cross-over of a vortex
lattice is observed in our experiment. Figure A.3 shows how the vortex lattice evolves
from triangular to square as the pinning strength increases. Over a wide range of pinning

81

Figure A.3: Images of rotating condensates pinned to an optical lattice at ΩOL = Ωc =0.866ωr with pinning strength Upin /µ=(a) 0.049 (b) 0.084 (c) 0.143, showing the structural
cross-over of the vortex lattice. (a)–(c) are the absorption images of the vortex lattices
after expansion. (d)–(f) are the Fourier transforms of the images in (a)–(c). ko is taken by
convention to be the strongest peak; ktr1 , ksq , and ktr2 are at 60◦ , 90◦ , and 120◦ , respectively,from ko .

82

Figure A.4: Structure factors (a) |S(ksq )| (¥), (b) |S(ktr )| (N) (average of |S(ktr1 )| and|S(ktr2 )|), and (c) |S(ko )| (F) versus the strength of the pinning lattice at the commensurate
rotation rate Ωc . |S(ksq )| is fitted by [75]. The fitting leads to a maximum value 0.707 of|S(ksq )|. An ideal square vortex lattice would have |S(ksq )|=1.

83

strengths, we observe that there is always at least one lattice peak in reciprocal space thatremains very strong. We define this peak to be ko . Lattice peaks at 60◦ and 120◦ from ko arereferred to as ktr , and, together with ko , their strength is a measure of the continued presenceof a triangular lattice. A peak at 90◦ , referred to as ksq , is instead a signal for the squarenessof the vortex lattice. With increasing pinning strength (Fig. A.3(a–c), or (d–f)), we see thetriangle to square crossover evolve. At intermediate strengths (Upin /µ = 0.084), a family of
zigzag vortex rows emerges, indicated by the dotted lines in Fig. A.3(b); in reciprocal space
we see the presence of structure at ktr and ksq .
We quantify the crossover by means of an image-processing routine that locates eachvortex core, replaces it with a point with unit strength, Fourier transforms the resultant
pattern, and calculates structure factors |S| [67] based on the strength of the images atlattice vectors ksq , ktr , and ko . In Fig. A.4, we see with increasing optical potential theturn-on of |S(ksq )| balanced by the turn-off of |S(ktr )|. We use a fitting function to smooththe noisy data of |S(ksq )|. The structure crossover takes place around Upin /µ ≈ 8%, in roughagreement with predictions of Upin /µ ≈ 5% from numerical simulations [67] and Upin /µ ≈ 1%from analytic theory for infinite lattices [66]. The fact that one lattice peak remains strong
for all pinning strengths (the stars (F) in Fig. 4) suggests that as the pinning strength isincreased, one family of vortex rows represented by ko in Fourier space locks to the squarepinning lattice and remains locked as the shape cross-over distorts the other two families ofvortex rows into a square geometry. The effects of various rotation rates and optical potential
strengths on the squareness of the vortex lattice is summarized in Fig. A.5. We surmise that
there are a number of effects at play. When ΩOL differs from Ωc , pinning strength is required
not only to deform the shape of the vortex lattice from triangular to square, but also tocompress or expand it to match the density of the optical lattice sites. At higher opticalintensities, we know from separate observations that imperfections in the rotation of the
optical lattice lead to heating of the condensate, which may limit the obtainable strength of

84

the square lattice.
A dumbbell-shape lattice defect (Fig. A.6) is sometimes observed in the early stages ofthe square vortex lattice formation when ΩOL > Ωc . In the weak-pinning regime, the defectwill relax towards the equilibrium configuration by pushing extra vortices at the edge of thecondensate outside the system. Defects of this nature, involving extra (or missing) vortices,
are the exception and not the rule in our observations, even for ΩOL 6= Ωc . In an infinite
system, the physics of the lattice-lattice interaction would likely be dominated by thesepoint defects. In our finite system, would-be incommensurate lattices can accommodate bystretching or compressing.

85

Figure A.5: Effect of square pinning lattice. (a) Contours of |S(ksq )| are plotted, showing theeffect of the rotation rate and pinning strength on the squareness of the vortex lattice. (b)The maximum observed squareness. In (a) and (b), for each rotation rate, the data points
are extracted from fits such as that shown in Fig. A.4 for ΩOL = Ωc . The vertical dottedline plus arrow shows the possible range of Ωc consistent with the uncertainty in nOL .

86

Figure A.6: Image of a dumbbell-shape defect consisting of two vortices locked to one pinningsite during the formation of the square vortex lattice. Dotted lines are to guide the eye.

Appendix B
Vortex proliferation in the Berezinskii-Kosterlitz-Thouless regime on a
two-dimensional lattice of Bose-Einstein condensates [1, 3]

Appendices We observe the proliferation of vortices in the Berezinskii-KosterlitzThouless regime on a two-dimensional array of Josephson-coupled Bose-Einstein condensates.
As long as the Josephson (tunneling) energy J exceeds the thermal energy T , the array isvortex-free. With decreasing J/T , vortices appear in the system in ever greater numbers.
We confirm thermal activation as the vortex formation mechanism and obtain informationon the size of bound vortex pairs as J/T is varied.
One of the defining characteristics of superfluids is long-range phase coherence [77],which may be destroyed by quantum fluctuations, as in the Mott-insulator transition [34, 78],or thermal fluctuations, e.g. in one-dimensional Bose gases [79, 80] and in a double-wellsystem [81]. In two dimensions (2D), Berezinskii [6], Kosterlitz and Thouless [7] (BKT)
developed an elegant description of thermal phase fluctuations based on the unbinding ofvortex-antivortex pairs, i.e. pairs of vortices of opposite circulation. The BKT picture
applies to a wide variety of 2D systems, among them Josephson junction arrays (JJA), i.e.
arrays of superfluids in which phase coherence is mediated via a tunnel coupling J between
adjacent sites. Placing an isolated (free) vortex into a JJA is thermodynamically favored if
its free energy F = E − T S ≤ 0. In an array of period d the vortex energy diverges witharray size R as E ≈ J log(R/d) [82], but may be offset by an entropy gain S ≈ log(R/d) due
to the available ≈ R2 /d2 sites. This leads to a critical condition (J/T )crit ≈ 1 independent

88of system size, below which free vortices will proliferate. In contrast, tightly bound vortexantivortex pairs are less energetically costly and show up even above (J/T )crit . The overall
vortex density is thus expected to grow smoothly with decreasing J/T in the BKT crossover
regime.
Transport measurements, both in continuous superfluids [48, 49] and superconducting
JJA [50] have confirmed the predictions of BKT, without however directly observing its
microscopic mechanism, vortex-antivortex unbinding. A recent beautiful experiment [9] in acontinuous 2D Bose gas measured the phase-phase decay function through the BKT crossover, and saw evidence for thermal vortex formation. For related theoretical studies see
e.g. [84]. In this work we present more detailed vortex-formation data, collected in a 2D
array of BECs with experimentally controllable Joephson couplings. The system was studied
theoretically in [85].
Our experiment starts with production of a partially Bose-condensed sample of 87 Rb
atoms in a harmonic, axially symmetric magnetic trap with oscillation frequencies {ωx,y , ωz } =2π{6.95, 15.0} Hz. The number of condensed atoms is kept fixed around 6 × 105 as the temperature is varied. We then transform this system into a Josephson junction array, as
illustrated in Fig. B. In a 10 s linear ramp, we raise the intensity of a 2D hexagonal opticallattice [86] of period d = 4.7µm in the x-y plane. The resulting potential barriers of height
VOL between adjacent sites [Fig.B(b)] rise above the condensate’s chemical potential aroundVOL ≈ 250 − 300 Hz, splitting it into an array of condensates which now communicate only
through tunneling. This procedure is adiabatic even with respect to the longest-wavelength
phonon modes of the array [87, 88] over the full range of VOL in our experiments. Each ofthe ≈ 190 occupied sites (15 sites across the BEC diameter 2 × RT F ≈ 68 µm [89]) now
contains a macroscopic BEC, with Nwell ≈ 7000 condensed atoms in each of the central wells
at a temperature T that can be adjusted between 30 − 70 nK. By varying VOL in a range
between 500 Hz and 2 kHz we tune J between 1.5 µK and 5 nK, whereas the “charging”

89

Figure B.1: Experimental system. (a) 2D optical lattice intensity profile. A lattice ofJosephson-coupled BECs is created in the white-shaded area. The central box marks thebasic building block of our system, the double-well potential shown in (b). The barrier height
VOL and the number of condensed atoms per well, Nwell , control the Josephson coupling J,
which acts to lock the relative phase ∆φ. A cloud of uncondensed atoms at temperature T
induces thermal fluctuations and phase defects in the array when J < T . (c) Experimental
sequence: A BEC (i) is loaded into the optical lattice over 10 s, suppressing J to values
around T. We allow 2 s for thermalization. To probe the system, we ramp off the lattice ona faster timescale tr [83] and take images of the recombined condensate. When J is reduced
below T (ii)-(vii), vortices (dark spots) appear as remnants of the thermal fluctuations in
the array.

90

energy Ec , defined in [77], is on the order of a few pK, much smaller than both J and T . Inpthis regime, thermal fluctuations of the relative phases ∆φT h ≈ T /J are expected, whilequantum fluctuations ∆φQ ≈ (Ec /4J)1/4 are negligible [77].
The suppression of the Josephson coupling greatly suppresses the energy cost of phasefluctuations in the x-y plane, between condensates, J[1 − cos(∆φ)], compared to the cost of
axial (z) phase fluctuations inside the condensates [51]. As a result, axial phase fluctuations
remain relatively small, and each condensate can be approximated as a single-phase object[90].
After allowing 2 s for thermalization, we initiate our probe sequence. We first takea nondestructive thermometry image in the x-z plane, from which the temperature T and,from the axial condensate size Rz , the number of condensed particles per well, Nwell , isobtained (see below). To observe the phase fluctuations we then turn down the opticallattice on a time-scale tr [83], which is fast enough to trap phase winding defects, but slow
enough to allow neighboring condensates to merge, provided their phase difference is small.
Phase fluctuations are thus converted to vortices in the reconnected condensate, as has been
observed in the experiments of Scherer et al. [37]. We then expand the condensate by afactor of 6 and take a destructive image in the x-y plane.
Figure B(c) illustrates our observations: (ii)-(vii) is a sequence of images at successively
smaller J/T (measured in the center of the array [91]). Vortices, with their cores visible asdark “spots” in (iii)-(vii), occur in the BEC center around J/T = 1. Vortices at the BEC edge
appear earlier, as here the magnetic trap potential adds to the tunnel barrier, suppressing
the local J/T below the quoted value. That the observed “spots” are indeed circulationcarrying vortices and antivortices is inferred from their slow ≈ 100 ms decay after the opticallattice ramp-down, presumably dominated by vortex-antivortex annihilation. From extensive
experiments on vortices in our system we know that circulation-free “holes” fill so quickly
due to positive mean field pressure, that they do not survive the pre-imaging expansion.

91

Vortices with identical circulation would decay by dissipative motion to the BEC edge, in
our trap over & 10 s.
To investigate the thermal nature of phase fluctuations, we study vortex activation
while varying J at different temperatures. For a quantitative study, accurate parameter
estimates are required. The Josephson-coupling energy J is obtained from 3D numerical simulations of the Gross-Pitaevskii equation (GPE) for the central double-well system [92, 81] [Fig.B(b)], self-consistently including mean-field interactions of both condensed
and uncondensed atoms [93]. A useful approximation for J in our experiments is [91]:J(VOL , Nwell , T ) ≈ Nwell × 0.315 nK exp[Nwell /3950 − VOL /244Hz](1 + 0.59 T /100nK). The
finite-T correction to J arises from both the lifting-up of the BEC’s chemical potential andthe axial compression by the thermal cloud’s repulsive mean field, but does not take into
account the effects of phase fluctuations on J (in condensed-matter language, we calculate
the bare J). Nwell is determined by comparison of the experimentally measured Rz , toRz (VOL , Nwell , T ) obtained from GPE simulations. Both experimental and simulated Rz are
obtained from a fit to the distribution of condensed and uncondensed atoms, to a ThomasFermi profile plus mean-field-modified Bose function [93]. In determination of all J values,
there is an overall systematic multiplicative uncertainty ∆J/J = ×
÷ 1.6, dominated by uncertainties in the optical lattice modulation contrast, the absolute intensity calibration, andmagnification in the image used to determine Nwell . In comparing J for “hot” and “cold”
clouds (see Fig. B) there is a relative systematic error of 15% associated with image fitting
and theory uncertainties in the thermal-cloud mean-field correction to J.
The qualitative results of our work are consistent whether we use an automated vortexcounting routine or count vortices by hand, but the former shows signs of saturation error
at high vortex density, and the latter is vulnerable to subjective bias. As a robust vortexdensity surrogate we therefore use the “roughness” D of the condensate image caused by
the vortex cores. Precisely, we define D as the normalized variance of the measured col-

92

Figure B.2: Quantitative study of vortices. The areal density of vortices is quantified by
the plotted D defined in the text. D is extracted only from the central 11% of the condensate region [circle in inset (a)] to minimize effects of spatial inhomogeneity. (b) D vs J
for two datasets with distinct “cold” and “hot” temperatures. Each point represents one
experimental cycle. The increase in D with decreasing J . 100 nK signals the spontaneous
appearance of vortices, while the “background” D . 0.01 for J & 200 nK is not associatedwith vortices. Vortices clearly proliferate at larger J for the “hot” data, indicating thermal
activation as the underlying mechanism. The large scatter in D at low J is due to shot noiseon the small average number of vortices in the central condensate region. (c) same data
as in (b), but averaged within bins of size ∆[log(J)] = 0.15. Error bars of D are standard
errors. (d) same data as (b), but plotted vs J/T . “Cold” and “hot” datasets almost overlap
on what appears to be a universal vortex activation curve, as confirmed by averaging [inset(e)], clearly revealing the underlying competition of J and T .

93

umn density profile from a fit to a smooth finite-T Bose profile [93], with a small constantoffset subtracted to account e.g. for imaging noise. To limit spatial inhomogeneity in J,
caused by spatially varying condensate density and optical lattice intensity, to < 10%, D isextracted only from the central 11% of the condensate area which contains 20 lattice sites[Fig. B(a)]. Comparison to automated vortex-counts shows that D is roughly linear in theobserved number of vortices, irrespective of the sign of their circulation, with a sensitivityof ≈ 0.01/vortex.
Figure B shows results of our quantitative study. In Fig. B(b), we plot D vs J for
two datasets with distinct temperatures. At large J & 200 nK a background D . 0.01 is
observed, that is not associated with vortices, but due to residual density ripples remaining
after the optical lattice ramp-down. Vortex proliferation, signaled by a rise of D above≈ 0.01, occurs around J ≈ 100 nK for “hot” BECs and at a distinctly lower J ≈ 50 nK
for “cold” BECs [confirmed by the averaged data shown in Fig. B(c)], indicating thermal
activation as the vortex formation mechanism. Plotting the same data vs J/T in Fig. B(d)shows collapse onto a universal vortex activation curve, providing strong evidence for thermal
activation. A slight residual difference becomes visible in the averaged “cold” vs “hot” data
[Fig.B(e)], perhaps because of systematic differences in our determination of J at different
temperatures.
The vortex density D by itself provides no distinction between bound vortex-antivortex
pairs and free vortices. In the following we exploit the flexibility of optical potentials to
distinguish free or loosely bound vortices from tightly bound vortex-antivortex pairs. A“slow” optical lattice ramp-down allows time for tightly bound pairs to annihilate before
they can be imaged. By slowing down the ramp-down duration τ [inset of Fig. B (a)], we
therefore selectively probe vortices at increasing spatial scales. This represents an attempt
to approach the “true” BKT vortex unbinding crossover that is complementary to transportmeasurements employed so successfully in superconductive and liquid Helium systems.

94

Figure B.3: (a) Vortex density D probed at different optical lattice ramp-down timescales τ .
A slow ramp provides time for tightly bound vortex-antivortex pairs to annihilate, allowing
selective counting of loosely bound or free vortices only, whereas a fast ramp probes bothfree and tightly bound vortices. A fit to the vortex activation curve determines its midpoint(J/T )50% , its 27% − 73% width ∆(J/T )27−73 , and the limiting values D< (D> ) well below(above) (J/T )50% . (b) A downshift in (J/T )50% is seen for slow ramps, consistent with theoccurrence of loosely bound or free vortices at lower J/T only. (c) Mapping between rampdown timescale τ and estimated size of the smallest pairs surviving the ramp (upper axis).
The difference D< − D> measures the number of vortices surviving the ramp (right axis).
Comparison to simulated vortex distributions yields a size estimate of the smallest √surviving
3, II: d,pairs (upperaxis).
Inset:
smallestpossible
pair
sizes
in
a
hexagonal
array,
I:d/√
III: 2d/ 3.

95

Figure B(a) shows vortex activation curves, probed with two different ramp-downtimes. Two points are worth noticing: First, a slow ramp compared to a fast one showsa reduction of the vortex density D< in arrays with fully randomized phases at low J/T .
The difference directly shows the fraction of tightly bound pairs that annihilate on the longramp. Second, a slower ramp shows vortex activation at lower (J/T )50% , confirming that freeor very loosely bound vortices occur only at higher T (lower J). Specifically, the data clearlyshow a range around J/T ≈ 1.4 where only tightly bound pairs exist. Figure B(b) quantitatively shows the shift of (J/T )50% from 1.4 to 1.0 with slower ramp time. We can makea crude mapping of the experimental ramp-down time-scale to theoretically more accessible
vortex-antivortex pair sizes as follows: In Fig. B(c), we see the decrease of the saturated(low-J/T ) vortex density D< with increasing ramp timescale τ . The right axis shows theinferred number of vortices that survived the ramp. We compare this number of surviving vortices to simulations [94] of a 20-site hexagonal array with random phases. In thesesimulations we find, on average, a total of 10 vortices, 6 of which occur in nearest-neighbor
vortex-antivortex pairs [configuration I in Fig. B(c)], 1.7 (0.4) occur in configuration II (III)respectively, and 1.9 occur in larger pairs or as free vortices. Experimentally ≈ 11 vortices
are observed for the fastest ramps, in good agreement with the expected total number ofvortices. For just somewhat slower ramps of τ ≈ 5 ms, only 3 vortices survive, consistent
with only vortices in configuration II & III or larger remaining (indicated in Fig. B, topaxis) [95]. For τ & 30 ms ramps less than 2 vortices remain, according to our simulations√
spaced by more than 2d/ 3. Thus we infer that ramps of τ ≈ 30 ms or longer allow time√
for bound pairs of spacing . 2d/ 3 to decay before we observe them. The downward shift√of (J/T )50% in Fig. B(b) thus tells us that loosely bound pairs of size larger than 2d/ 3, orindeed free vortices, do not appear in quantity until J/T ≤ 1.0, whereas more tightly boundvortex pairs appear in large number already for J/T ≤ 1.4.
A further interesting observation concerns the width of the vortex activation curve.

96

The relative width, determined from fits to data such as the ones shown in Fig. B(a), is∆(J/T )27−73 /(J/T )50% ≈ 0.3, independent of ramp-down duration. This width is neither
as broad as in a double-well system [96, 81], where the coherence factor rises over a range∆(J/T )27−73 /(J/T )50% ≈ 1.4, nor as broad as expected from our simulations [94] of an arraypof uncoupled phases, each fluctuating independently with ∆φRM S = T /J, for which wefind ∆(J/T )27−73 /(J/T )50% ≈ 0.85. Presumably collective effects in the highly multiply
connected lattice narrow the curve. On the other hand, the width is 3 times larger thanthe limit due to spatial inhomogeneity in J, suggesting contributions to the width due tofinite-size effects or perhaps revealing the intrinsically smooth behavior of vortex activationin the BKT regime.
In conclusion, we have probed vortex proliferation in the BKT regime on a 2D latticeof Josephson-coupled BECs. Allowing variable time for vortex-antivortex pair annihilation
before probing the system provides a time-to-length mapping, which reveals information on
the size of pairs with varying J/T .