Precise Calibrations of Few-Body Physics in Potassium-39:
Experiment and Theory
by
Xin Xie
B.S., Physics, Fudan University, Shanghai 2013
M.S., Physics, University of Colorado, Boulder 2016

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
2020

DocuSign Envelope ID: D1745571-062F-48B4-B709-E9A00E883057

This thesis entitled:
Precise Calibrations of Few-Body Physics in Potassium-39:
Experiment and Theory
written by
Xin Xie
has been approved for the

masters

doctoral degree program in:

Physics
5/4/2020

Eric Cornell
Committee Chair Name
5/4/2020

Jun YeCommittee Member Name

The final copy of this thesis has been examined by the signatories, and wefind that both the content and the form meet acceptable presentation standardsof scholarly work in the above mentioned discipline.
Student ID number
IRB protocol #IACUC protocol #102616589

iii
Xie, Xin (Ph.D., Physics)
Precise Calibrations of Few-Body Physics in Potassium-39: Experiment and Theory
Thesis directed by Prof. Eric A. Cornell

Ultracold atoms provide a versatile toolbox for the study of diverse quantum phenomena.
The minimal platform to present quantum mechanics is at single particle level. Two-particle interference effects such as Hanbury Brown-Twiss are a consequence of the particle statistics. Exotic
quantum effects in which statistics and interactions are combined appear in three-particle systems
as elucidated by e.g. Efimov in the early 1970s. Efimov’s work can be further extended to four- or
five-particle systems and applied to general dilute quantum gases in the weakly interacting regime.
In the many-body context, it may also assist the analysis of quench dynamics in unitary gasses
where perturbation theory breaks down.
A systematic study of few-body physics is made possible by the use of magnetic Feshbachresonances. One can greatly benefit from the improved experimental precision when comparing theobservables to generic few-body models. In this thesis, I begin by describing the efforts we spenton achieving overall robustness and accuracy on our potassium-39 quantum gas machine, including
the hardware upgrades we made for future scientific projects. Then following a brief review on
the realization of Efimov physics in cold atoms, I expand on a series of scattering experiments weperformed to understand few-body physics with effectively repulsive inter-particle interactions, asa complement to our previous investigations of the attractive-interaction regime. At the end, I
combine our knowledge from both regimes together to discuss the global picture of Efimov physics
in homonuclear systems at a beyond proof-of-principle level.

Dedication“ Blind to all fault, destiny can be ruthless at one’s slightest distraction. ”— Jorge Luis Borges, Ficciones

To my family, for your unconditional support.
To my friends, for your company and tolerance.
To healthcare workers, for protecting my loved ones.

v

Acknowledgements

There are many people I would like to express my gratitude to in this thesis. My late advisor
Deborah Jin taught me various expertise in the lab. As time passed by those knowledge inevitably
started to fade out. What lingers in my memory about her is the short walk that we took together
along the cliff when we were at a conference at Rhode Island. I could see the Atlantic Ocean, andsmell the seaweed. Everything is still vivid when I close my eyes.
My other advisor Eric Cornell has been the most important person throughout my graduate
school study. His infinite patience and playful personality make him the greatest teacher and friend.
Eric is the smartest, most considerate and financially supportive advisor. The statistical probability
is extremely low to find all these amazing qualities in a single person who needs also be good with
one-handed softball and Mandarin Chinese.
Since 2016, Eric and Jun Ye have been co-advising our project on potassium-39. Jun introduces to us his hard-core, single-minded, second-to-none instrumentation techniques. He influenced
me with his strong attitude toward numerous things such as error bars and theoretical models.
When facing technical difficulties, Jun is usually the only one who keeps the faith even everyone
else loses it. And this is because he has programmed himself to do so.
I cannot finish my thesis without the help from JILA personnel. I began my first machining
project under the training of Kim Hagen. I am grateful for all the crafty and friendly instrument
shop staff especially Todd Asnicar, Hans Green, and James Urich. Research can be frustrating,but machining is always full of fun. I would like to say thank you to the electronic staff Terry
Brown, James Fung-A-Fat, and Carl Sauer. Thank all of you for chatting with me and sharing

vi
your magic tricks. Time really flies when the solder iron is in hand (and Jan Hall takes the next
seat at the counter). I appreciate all the computer rescue provided by J. R. Raith, Corey Keaslingand Jim McKown, who have been frequently dragged down to the lab after hours. I appreciate
the administrative work done by Krista Beck, Amy Allison, Cynthia Torres and Agnieszka Lynch.
Your warmhearted assistance is essential to me accomplishing my degree. I would also like to thankBrian Lynch for tracking my purchases, and Jason Ketcherside for building maintenance. Your help
is super important and allows me to focus on my research.
My labmates are those who I spent the longest time with in the past seven years. I havegained endless power from their encouragement. Cathy Klauss is the senior graduate student who
introduced me to many fundamental skill sets in AMO experiment. We spent tremendous efforts on
repairing our old rubidium-85 machine. And eventually we were able to make the biggest rubidium85 BEC worldwide and studied strongly interacting dynamics with it. Cathy is gifted with stunningpresentation and writing skills. I later work with graduate students Roman Chapurin and Michael
Van de Graaff (a.k.a. Vandy) on the potassium-39 experiment. We converted this new generation
fermi gas machine into a bose gas machine in 2016-2017. Roman is highly motivated and productive.
We went through the optimization of every single step in the whole cooling process of potassium.
This became the most valuable experience for my independent research later on. Vandy and I
took lots of data on the scattering experiment for repulsive inter-particle interactions. Vandy is amaster of optics. He is meticulous about details and full of creativity. I believe he will discover a
lot of interesting physics in the days to come. I miss the time I spent with graduate students Philip
Makotyn, Rabin Paudel, Noah Schlossberger, undergraduate students David Goldberger, Carlos
Lopez-Abadia, Joshua Giles, Jared Popowski, Bjorn Sumner. It was the most peaceful memorythat one day in the afternoon everyone in the lab was so immersed in the work on their hands.
Hours of time elapsed in this quietness before anyone realized it. In addition, I will always remember
the relaxing and informative conversations with the folks from the EDM, KRb, strontium, cavity,
tweezer, cold molecule, and membrane labs.
I am often impressed by how startlingly well theory and experiment agree with each other invii
the field of cold atom scattering and interactions. I feel really fortunate to be acquainted with anumber of brilliant theorists in the past years. Paul Julienne (a.k.a. “the wizard”) has repeatedly
amazed me with his encyclopedic knowledge on the scattering properties of any atomic species.
He explains physics in the most precise language and always points to the right literature for mein the blink of an eye. José D’Incao is our expert on three-body physics. His numerical models
are the most powerful theoretical toolbox that we have. Besides, I appreciate all the enlightening
discussions I had with John Bohn, Doerte Blume, Bo Gao, and many others who I exchanged
e-mails with but have not had a chance to meet in person.
Lastly, I would like to thank my parents for always encouraging me to pursue my dreams.
With your guidance, I am able to enjoy my hobbies and choose the career that I love. I would also
like to thank Sudi Chen and Ruqi Yi, for helping me develop confidence as well as an independent
personality, and for being there for me on the rainy days.

viii

Contents

Chapter
1 Introduction1

1.1

What is Universality and Why is it Important? . . . . . . . . . . . . . . . . . . . . .

1

1.2

Thesis Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

52 Apparatus
2.1

2.2

2.3

8

Apparatus in Few-Body Physics Project . . . . . . . . . . . . . . . . . . . . . . . . .

8

2.1.1

Laser Cooling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

9

2.1.2

Evaporative Cooling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

2.1.3

Crossed-Dipole Trap Calibrations . . . . . . . . . . . . . . . . . . . . . . . . . 16

2.1.4

Laser Linewidth Measurements . . . . . . . . . . . . . . . . . . . . . . . . . . 19

2.1.5

Absorption Imaging Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

2.1.6

Magnetic Field Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

High-Resolution Top Imaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
2.2.1

Design of a Compact Imaging System of NA 0.4 . . . . . . . . . . . . . . . . 36

2.2.2

Bench Test Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38

2.2.3

Analysis of Optical Aberrations . . . . . . . . . . . . . . . . . . . . . . . . . . 45

2.2.4

Top Imaging Calibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

2.2.5

Science Camera Calibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53

Fast Magnetic Field Coils Setup. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53

ix
2.3.1

Fast-B Servo and Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55

2.3.2

Compensation for Induced and Eddy Current . . . . . . . . . . . . . . . . . . 57

3 Efimov Physics: A Brief Review
3.1

3.2

3.3

From Two to Three. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62

3.1.1

Two-Particle Scattering in Zero-Energy Limit . . . . . . . . . . . . . . . . . . 62

3.1.2

Finite-Energy Effects and Effective Range . . . . . . . . . . . . . . . . . . . . 67

3.1.3

A Pedagogical Illustration of Three-Particle Scattering . . . . . . . . . . . . . 69

Universality in Efimov Physics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
3.2.1

Observation through Loss . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72

3.2.2

Van der Waals Universality . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75

3.2.3

Zero-Range Universality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76

An Outline of Our Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82

4 Experiments on Feshbach Molecules
4.1

4.2

4.3

62

83

Feshbach Molecule One-Body Lifetimes . . . . . . . . . . . . . . . . . . . . . . . . . 83
4.1.1

Measurements of Feshbach Molecule Lifetimes . . . . . . . . . . . . . . . . . . 83

4.1.2

A Coupled-Channel Model for the Lifetime of the Quasibound State . . . . . 88

Finite Temperature Effects in Atom-Dimer Reactions . . . . . . . . . . . . . . . . . . 91
4.2.1

Discovery of an Atom-Dimer Resonance . . . . . . . . . . . . . . . . . . . . . 91

4.2.2

Zero-Temperature Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96

4.2.3

Finite-Temperature Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98

4.2.4

Temperature Dependency of Atom-Dimer Resonance . . . . . . . . . . . . . . 101

Analysis of Systematic Effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102
4.3.1

Atom-Dimer Thermal Equilibrium . . . . . . . . . . . . . . . . . . . . . . . . 102

4.3.2

Atom-Dimer Chemical Equilibrium . . . . . . . . . . . . . . . . . . . . . . . . 105

4.3.3

Dimer-Dimer Inelastic Collisions . . . . . . . . . . . . . . . . . . . . . . . . . 105

4.3.4

Residual Breathing Motion of the Atom Cloud . . . . . . . . . . . . . . . . . 111

x
4.3.5
4.4

In Situ Detection of Molecules . . . . . . . . . . . . . . . . . . . . . . . . . . 112

Bound-to-Bound Transitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114

5 Experiments on Free Atoms

118

5.1

Experimental Sequence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118

5.2

Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120

5.3

Analysis of Raw Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124

5.4

Finite Temperature Model of Three-body Recombination . . . . . . . . . . . . . . . 127

5.5

Analysis of Systematic Effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133

6 Summary and Outlook

140

6.1

Fifty Years of Efimov Physics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140

6.2

Future Directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143

Bibliography

145

Appendix
A Algorithms: Absorption Imaging Simulation

154

A.1 Absorption Imaging Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 154
A.2 OD Correction as Function of I/Isat . . . . . . . . . . . . . . . . . . . . . . . . . . . 157
A.3 Fit Function for Isat . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 158
B Algorithms: Analysis of Optical Aberrations

160

B.1 PSF and MTF for an Ideal Pupil . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 160
B.2 Determination of Zernike Polynomial Coefficients . . . . . . . . . . . . . . . . . . . . 161
C Algorithms: Finite-Temperature Effects in Atom-Dimer Relaxation

164

C.1 Parametrization of Atom-Dimer Scattering Amplitude . . . . . . . . . . . . . . . . . 164

xi
C.2 Atom-Dimer Chemical Equilibrium: a Toy Model . . . . . . . . . . . . . . . . . . . . 166
C.3 Fit Function for Dimer-Dimer Collision Rate Constant . . . . . . . . . . . . . . . . . 167
D Algorithms: Finite-Temperature Effects in Three-Body Recombination

169

D.1 L3 Fitting Protocol I . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 169
D.2 L3 Fitting Protocol II (Simultaneous Fit) . . . . . . . . . . . . . . . . . . . . . . . . 172
D.3 S-Matrix Formalism . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173
E Data and Fits: Three-Body Recombination at Finite Temperatures

179

(0)

E.1 Extract a+ from 460 nK Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 179(0)E.2 Extract ap from 410 nK and 230 nK Data . . . . . . . . . . . . . . . . . . . . . . . 180(1)E.3 Exploration of a+. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 181

xii

TablesTable
4.1

Atom-dimer resonance measurement conditions and fit results . . . . . . . . . . . . . 102

5.1

Three-body recombination conditions and fit results . . . . . . . . . . . . . . . . . . 122

6.1

Summary of three-body observables characterizing the Efimov spectrum in 39 K acrossthe Feshbach resonance at 33.58 G . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142

xiii

FiguresFigure
2.1

Photographs of the laboratory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8

2.2

The vacuum system and optical transitions of 39 K . . . . . . . . . . . . . . . . . . . 10

2.3

Timeline of laser cooling process . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12

2.4

Λ-type configuration for gray molasses and phase lock signal . . . . . . . . . . . . . . 14

2.5

Timeline of evaporative cooling process

2.6

Background heating and number loss in optical dipole trap . . . . . . . . . . . . . . 18

2.7

Calibrations of trapping frequencies for optical dipole trap . . . . . . . . . . . . . . . 20

2.8

Self-heterodyne interferometer optical layout and beat note spectrum . . . . . . . . . 21

2.9

Comparison of different fitting functions for the self-heterodyne beat note spectrum. . . . . . . . . . . . . . . . . . . . . . . . . 17

24

2.10 Coarse-grained method for absorption imaging simulation . . . . . . . . . . . . . . . 26
2.11 Absorption imaging simulation results . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.12 Comparison of two fitting models for Isat . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.13 Application of the recoil model to the measured optical depth . . . . . . . . . . . . . 31
2.14 Photographs of rf/microwave antennae surrounding the science cell and antennacomponents. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.15 Return loss measurement on a network analyzer. . . . . . . . . . . . . . . . . . . . . 34
2.16 Plug-in low-pass filter for bias field servo inputs . . . . . . . . . . . . . . . . . . . . . 35
2.17 Rf spectra taken with various rf radiation power being present . . . . . . . . . . . . 37
2.18 Design of our compact top imaging system . . . . . . . . . . . . . . . . . . . . . . . . 39

xiv
2.19 Preliminary assessment of the aspheric lens with shearing interferometer . . . . . . . 39
2.20 Optics layout for the top imaging bench test setup . . . . . . . . . . . . . . . . . . . 40
2.21 Photograph of the point source setup and SEM image of the pinhole . . . . . . . . . 41
2.22 Point spread function (PSF) of the top imaging system generated from a 0.5 µm
pinhole . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
2.23 Study the pixelation effect by numerically binning the PSF . . . . . . . . . . . . . . 44
2.24 Explore the finite depth of focus of the PSF . . . . . . . . . . . . . . . . . . . . . . . 46
2.25 The PSF and MTF associated with an ideal pupil for fluorescence and absorption
imaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
2.26 PSF fit with the wavefront distortion expanded in Zernike polynomials . . . . . . . . 50
2.27 Optics layout for the top imaging setup on the science table . . . . . . . . . . . . . . 51
2.28 Calibration of top imaging camera . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
2.29 Schematic of the fast-B coils feedback servo setup . . . . . . . . . . . . . . . . . . . . 56
2.30 Tuning up the fast-B coils feedback servo . . . . . . . . . . . . . . . . . . . . . . . . 58
2.31 Induced current in bias coils due to fast-B field jump . . . . . . . . . . . . . . . . . . 59
2.32 Compensating the induced current and eddy current caused by the fast-B field jump

61

3.1

Two-channel model of Feshbach resonances . . . . . . . . . . . . . . . . . . . . . . . 63

3.2

Feshbach resonances in 39 K for the F = 1 spin-polarized gases . . . . . . . . . . . . . 65

3.3

Two shallowest bound states in 39 K for the cc channel . . . . . . . . . . . . . . . . . 65

3.4

Zeeman energies of the five spin configurations in the multi-channel model for 39 K,
total mF = −2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66

3.5

Fraction of each channel in Feshbach molecule as a function of magnetic field . . . . 66

3.6

A coupled-channel calculation of bound states near threshold for 39 K . . . . . . . . . 67

3.7

Effective range re as a function of sres . . . . . . . . . . . . . . . . . . . . . . . . . . 68

3.8

A coupled-channel calculation of re for 39 K . . . . . . . . . . . . . . . . . . . . . . . 69

3.9

Jacobi coordinates describing the relative positions of three identical particles . . . . 70

xv
3.10 Efimov attraction and scaling invariance . . . . . . . . . . . . . . . . . . . . . . . . . 71
3.11 Schematic plot of Efimov spectrum . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
(0)3.12 Summary of experimentally determined values of a− in homonuclear systems asfunction of sres . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
3.13 Van der Waals and zero-range universality in Efimov physics . . . . . . . . . . . . . 79
3.14 Summary of experimental and theoretical results for Efimov’s ratios between threebody observables for a < 0 and a > 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
4.1

Zeeman energies of the d-wave decay channels for the incoming cc channel . . . . . . 84

4.2

Magneto association of Feshbach molecules and atom purification . . . . . . . . . . . 86

4.3

Timeline of Feshbach molecule lifetime experiments . . . . . . . . . . . . . . . . . . . 87

4.4

Refinement of the d-wave sector of the two-body coupled-channel model . . . . . . . 89

4.5

Bound state calculations for 39 K in sd-wave basis . . . . . . . . . . . . . . . . . . . . 92

4.6

Real and imaginary parts of a in vicinity of the narrow d-wave resonance

4.7

Time evolution of dimers measured at 1200 a0 with atoms being the majority in the. . . . . . 93

sample . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95
4.8

Temperature dependence of atom-dimer relaxation coefficient βAD as a function of
scattering length a . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97

4.9

Finite-temperature model predictions for atom-dimer reaction rate coefficient βAD . 100

4.10 Summary of the three-body observables extracted from finite-temperature and zerotemperature models respectively . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
4.11 A close look at the atom-dimer thermal equilibrium . . . . . . . . . . . . . . . . . . 104
4.12 Atom-dimer chemical equilibrium condition as function of scattering length . . . . . 106
4.13 A two-step ramp to large a while staying adiabatic with respect to Eb . . . . . . . . 108
4.14 Atom formation as function of a in dimer-dimer experiment . . . . . . . . . . . . . . 109
4.15 Time evolution of dimers and atoms with dimers initially being the majority in thesample . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110

xvi
4.16 Dimer-dimer relaxation coefficient into deeply bound molecules . . . . . . . . . . . . 111
4.17 Residual breathing motion of atoms in atom-dimer relaxation experiment . . . . . . 112
4.18 Expansion of post-dissociation molecule cloud . . . . . . . . . . . . . . . . . . . . . . 113
4.19 Projections of shallow dimers onto asymptotic atomic basis set . . . . . . . . . . . . 115
4.20 Bound-to-bound transition spectra . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
5.1

Timeline of three-body recombination experiments . . . . . . . . . . . . . . . . . . . 119

5.2

Dynamic range of individual condition in three-body recombination experiments . . 120

5.3

Compilation of L3 measurements for a < 0 and a > 0 . . . . . . . . . . . . . . . . . . 122

5.4

Three-body recombination coefficient L3 scaled by a4 . . . . . . . . . . . . . . . . . . 123

5.5

Two- and three-body loss coefficients extracted from atom samples of various densities at 3000 a0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125

5.6

Illustration of L3 simultaneous fit that combines samples of two different densities
but same temperature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126

5.7

Two-body loss rate coefficient L2 as function of B and a . . . . . . . . . . . . . . . . 128

5.8

Using the S-matrix method to calculate L3 at finite temperatures . . . . . . . . . . . 130

5.9

Finite-temperature model predictions for three-body recombination rate coefficient
(a>0)

L3

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131

5.10 Contributions to the three-body recombination rate into the shallow dimer fromangular momenta J = 1 − 6 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132
5.11 A coupled-channel calculation of −Im(a) in sdg-wave basis near 34 G

. . . . . . . . 135

5.12 Rf lineshape for the |1, −1i → |2, 0i transition at around 2000 a0 . . . . . . . . . . . 135
5.13 Vary the trap depth in three-body recombination experiment at large a . . . . . . . 136
5.14 Frequency and intensity modulation on the ODT beams . . . . . . . . . . . . . . . . 137
5.15 Comparison of three-body recombination data with the predictions by a toy modelthat includes avalanche effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138
E.1 Finite-temperature fit for recombination minimum, 460 nK condition, J = 0 . . . . . 179

xvii
E.2 Finite-temperature fit for recombination maximum, 410 nK and 230 nK conditions,
J = 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 180
E.3 Finite-temperature fit for recombination minimum, 230 nK and 210 nK conditions,
J = 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 181
E.4 Finite-temperature fit for recombination minimum, 150 nK and 60 nK conditions,
J = 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 182

Chapter 1

Introduction

1.1

What is Universality and Why is it Important?
The concept of universality has been broadly discussed in the context of critical phenomena

in condensed matter. Physics systems with completely different properties at the atomic level
obey universal power-law thermodynamics when they undergo phase transitions such as liquid/gas
transition, paramagnetic/ferromagnetic transition and so on so forth. In the critical regions which
can be defined by the correlation length l much longer than the intrinsic length scale r0 , different
systems can be divided into a small number of classes based on the same critical exponent of
the power laws [1]. Deeply rooted in these universal phenomena is the so-called renormalizationgroup (RG) theory by Kenneth Wilson [2]. By integrating out the ultraviolet degrees of freedom(or performing a renormalization group transformation) in the space of coupling constants, people
obtain the quantity called RG flow which leads in the infrared limit to interesting scaling-invariant
physics described only by long-distance degrees of freedom such as temperature and pressure. The
above mentioned continuous scaling invariance in phase transitions represents the simplest topologyof the RG flow that is termed a fixed point. The next complicated structure is a limit cycle whichheralds discrete scaling invariant behaviors at the macroscopic level.
Universal physics in neutral atoms with large scattering length can be derived in the same
framework yet exhibit richer meanings of multiple layers. This is explained from a general dimensional argument in this Introduction section and illustrated with more microscopic details in
Ch. 3.

2
Neutral atoms interact with each other through van der Waals (vdW) force at long distance
whose range is usually parametrized by rvdW . When atoms stay in their electronic ground state,
rvdW is only s few nanometers. In the low-collision-energy limit, the characteristic length scaleof scattering processes is the s-wave scattering length a. Intuitively this is because the scattering
wavefunction is so delocalized at low temperatures that the majority of it lies outside the interaction
potential thus simply follows a free-particle-like sinusoid. The tiny fraction of the wavefunctionthat overlaps with the interaction potential effectively sets the short-range phase of the sinusoid
which turns into an offset from the origin in the center-of-mass frame. This offset is the so-called
scattering length a. When looking from long distances, the two particles are repulsive when there
is a positive offset and attractive when there is a negative offset regardless of the sign of theunderlying interaction potential [3]. Sometimes the scattering length a can acquire a much larger
absolute value compared to rvdW either due to some fine tuning of interaction potential or theuse of a Feshbach resonance. It then becomes the only relevant length scale that people need tounderstand the scattering- and bound-state observables in two-body problems. For instance, thebinding energy Eb of Feshbach molecules depends only on a in the large a limit: Eb = ~2 /ma2 .
This is the first layer of universality that may be called the “scattering-length universality” [4].
Low-energy two-body observables display continuous scaling invariance under the transformation
on scattering length a, the coordinate r, and the time t by a positive number ζ
a −→ ζa,r −→ ζr,t −→ ζ −2 t.

(1.1)

Violations of this universality are suppressed by powers of rvdW /|a|. In this regime, a also determines some macroscopic properties in mean-field problems for a collection of particles.
If we crank up the scattering length with the help of a Feshbach resonance, eventually we
can reach the limit of |a0 | → ∞ which is called the unitary limit; for two particles colliding with
momenta ±~k in the center-of-mass frame, the s-wave contribution to the scattering cross section
saturates to σ(a→±∞) = 4π/k 2 (the phase shift δl=0 (k) modulo π approaches π/2), which is thelargest it can be consistent with the constraint that the scattering S-matrix is unitary as a result

3
of conservation of probability. The relevant length scale in this case is played by 1/k which can be
related to the de Broglie wavelength λ for bulk gases. This can be understood by the fact that theinteraction between neutral atoms is short-range in nature no matter how large a is — it can only
take place when the interacting entities are in contact with each other. |a0 | → ∞ is also called theresonant limit and it has been a breeding ground for new physics. A well-known example is thatdifferent strongly interacting spin- 12 fermions share a universal expression for the chemical potentialµ that is purely set by the number density nµ=ξ

~2(6π 2 n)2/3 ,
2m

(1.2)

where ξ is the Bertsch parameter that has been determined both theoretically and experimentally[5]. Similarly, in the two-body sector of strongly interacting degenerate bosons or Bose-Einsteincondensates (BECs), the divergent λ gets replaced by the inter-particle spacing d which scales asn−1/3 and sets the energy and time scales for the development of quantum dynamics [6, 7, 8].
The third layer of unversality is unique to the three-body sector and only applies to bosons
and heteronuclear fermions in the absence of Pauli blocking effect. In addition to the aforementioned
two-body length scales, low-energy observables may also have dependence on another scale which
is sometimes named the three-body parameter (3BP). Vitaly Efimov first discovered that in theresonant limit there are an infinite number of arbitrarily-shallow three-body bound states (named
trimers) subject to a binding energy spectrum of the form(n)ET−→ (e−2π/s0 )n−n0~2 κ20m(n, n0 ∈ Z),(1.3)

where s0 = 1.00624 for three identical bosons [9]. In this expression, a wave number κ0 acts as the
3BP and fixes the absolute locations of the trimer states. This spectrum displays discrete scaling
invariance under the transformation
κ0 −→ κ0 ,r −→ ζ0n r,

t −→ ζ0−2n t,(1.4)

where ζ0 = eπ/s0 . Efimov trimers are peculiar in that their spatial extents can be infinitely large as
n → +∞. Stepping away from the resonant limit, the above symmetry also holds at finite a under4
the discrete transformation
κ0 −→ κ0 ,a −→ ζ0n a,r −→ ζ0n r,

t −→ ζ0−2n t.

(1.5)

This result gives rise to the famous Efimov’s three-body bound state spectrum in two-dimension
space as functions of E 1/2 and 1/a, which is regarded as a remarkable example of self-similarity. The
original derivation of Efimov’s spectrum is based on a zero-range (δ function) two-body interaction
potential, which is approximately satisfied in real atomic systems at finite a when |a|  r0 . Thisapproach brings up another powerful concept of the zero-range limit which has been frequently
used in analytic methods by introducing 3BP as a boundary condition at r → 0.
Last but not least, the recent discovery of a modern unversality has really brought the
study of few-body physics with ultacold atoms into focus. This is the well-known van der Waals
universality. As discussed above, people use zero-range theory to predict the geometry and scalingof trimer spectrum, then fix the absolute values of three-body observables with the help of a threebody parameter (3BP). The necessity of introducing a 3BP has to do with a shortcoming of thetwo-body zero-range potential, namely that it leads to an effective three-body potential that can
support an infinite number of trimer states with infinitely small size and infinitely high energy(as n → −∞). This unbounded-from-below spectrum is referred to as Thomas collapse or fall of
the particles to the center [10]. One would expect this behavior to break down below a distancecomparable with the range of real potentials such as rvdW . The function of 3BP is to encapsulate
the effect of two-body and possibly three-body interactions at short distances hence pins down thetrue ground trimer state. It can take the form of a wave number, a length, an energy or a phasedepending on the formalism used. Strictly speaking, only in the large n limit does it represent thetrue 3BP since the universal scaling factor only approximately holds for the low-lying trimer states
in real systems. Given that to access the excited states with n ≥ 1 in practice requires extremely
low temperature, trimer branches with n = 0, 1 have been serving as the 3BP as compromise.
A convenient quantity frequently measured is where the ground trimer state intersects the free(n=0)

atom threshold in the a < 0 region, denoted as a−. In the early days this empirical 3BP5
was not anticipated to be known a priori due to its sensitivity to the detailed form of the shortrange forces [11]. As people collected extensive experimental evidences on various atomic species,(n=0)

they discovered a universal position of a−which is roughly ten times rvdW [12]. This rather

unexpected result has later been traced back to an origin related to the 1/r6 form of the van der
Waals potential between two neutral atoms [13, 14, 15, 16]. This original version of van der Waals
universality assumes a virtually pure spin configuration during the collision processes (sometimes
termed single-channel model ) and applies to the broadest Feshbach resonances in species like 133 Cs.
It has later on been generalized to two-, three- and even more scattering channels to emulate more
complicated spin configurations [17, 18]. Notably the remarkable predicting power of this multichannel van der Waals theory suggests that additive pair-wise interaction (specifically van der
Waals) plays the dominant role in the three-body problem while short-range forces that explicitly
involve three particles at the same time (e.g. Axilrod-Teller) are somehow always negligible [19].
Much of contents in this thesis is along the line of this last kind of universality to explore the“imperfections” of the 3BPs associated with the low-lying trimer branches.

1.2

Thesis Contents
There are two major technical challenges in the cooling processes of 39 K. One is the small

hyperfine splitting energy of the electronic excited state 2 P3/2 on D2 transition that is normally
exploited for laser cooling alkali metals [20, 21]. Such small energy level splittings are comparable
to the natural linewidth of the excited states therefore spoil closed optical transitions. The second
challenge is the small elastic scattering cross section between atoms at hundreds of microKelvinand zero magnetic field due to the presence of a Ramsauer-Townsend minimum [22, 23]. This
suppresses the efficiency of evaporative cooling thus necessitates the use of an all-optical trapping
scheme combined with a magnetic-field-tunable scattering length.
In Chapter 2, I first introduce the apparatus that we use to perform all few-body physics
experiments. In order to repeatedly produce 39 K gas samples in a well-controlled quantum stateincluding temperature, number density, and spin state, we start off with a laser cooling stage that

6consists of Magneto-Optical Trap (MOT) collection, compressed-MOT, and sub-Doppler coolingwith Λ-enhanced gray molasses. We then move on to evaporative cooling in a crossed-dipole trapwith well-calibrated trapping frequencies. I describe an interferometric method we use to track thespectral linewidths e.g. of the laser that probes the atoms. Another thing to focus on is how weunderstand our absorption imaging process to obtain trustworthy number and spatial distribution
measurements of the atoms. I also describe how we achieve stable magnetic-field control as wellas its robustness against strong rf radiation. Next I detail the efforts we made to upgrade thehardware of our 39 K machine for future envisioned experiments. This includes the implementationof a high resolution imaging setup of numerical aperture about 0.4 and a fast magnetic-field control
system that can be used to initiate quench dynamics on a BEC .
In Chapter 3, I begin with briefly introducing major previous experimental and theoreticalresults on Efimov physics to illustrate the evolution of understanding in this field. As to the global
picture of Efimov’s spectrum as a function of the coupling strength, I refer to a two-channel field
theory [15, 24] which qualitatively predicts some interesting trends. When it comes to the absolutevalues of three-body observables, I emphasize the quantitative contributions from a numericalapproach done in hyperspherical coordinate [25]. I then discuss the missing pieces of the puzzle andwhat makes 39 K a promising candidate in filling in the unexplored territory.
In Chapter 4 and 5, I expand on the main body of our few-body physics project. We
measured a series of two- and three-body observables including the lifetime of Feshbach molecules, dwave resonance at 63 G, atom-dimer reaction rate, dimer-dimer reaction rate, molecule-to-molecule
transition by microwave coupling, and three-atom recombination rate. The study of each observable
is divided into the following aspects: design of experimental sequence, data analysis, systematic
analysis of experimental uncertainties, theoretical model or numerical simulation. In this thesis, I
focus on the repulsive interaction region with the s-wave scattering length a > 0. This is a regime
that people frequently encounter in all kinds of cold-atom experiments that involve collisions.
Nevertheless the seemingly straightforward dynamics can sometimes lead to difficult-to-interpret
results if certain details of the experiment are not properly handled.

7
In Chapter 6, I conclude the work presented in this thesis by discussing the connections
between different three-body observables. In the end I briefly cover ongoing projects and future
directions.

Chapter 2

Apparatus

2.1

Apparatus in Few-Body Physics Project

Figure 2.1: Photographs of the laboratory as of early 2020. Top panel: the laser table; bottom
panel: the science table. Our lab is located in X1B30 in the basement of JILA X-wing. The labtemperature is stabilized at 20 ± 0.5 ◦ C. The lab humidity is controlled at 37 ± 3 %.

Although I was not involved in the basic infrastructure of the 40 K experiment dating back toaround 2012, I participated in the entire second half of the construction process ca. 2016, with thegoal of modifying it into a bose gas machine with large atom number and relatively fast repetition
rate. Properties of 39 K and overviews of our quantum gas machine can also be found in Rabin

9
Paudel’s and Roman Chapurin’s these [26, 27] as complementary materials to this thesis.

2.1.1

Laser Cooling

Our experiment consists of three vacuum chambers (see the left panel in Fig. 2.2 for a top
view of the whole setup), two for 3D-MOTs and one smaller glass cell for doing science. The first
MOT acts as a source MOT where we collect and pre-cool the atoms from ambient gas. It is a Pyrex
cell with 9 windows and six dispensers with potassium and rubidium sources in them (see Fig. 4.1
in [26] for a picture of the bare cell). It typically takes less than 15 minutes for the dispensers tosaturate 39 K vapor on a daily-operation basis. The temperature of MOT1 is constantly maintainedat around 65°C to prevent the vapor from depositing on the walls. To enhance the MOT1 loading
rate (which is found to affect the MOT2 loading rate), we achieved our best performance at the
maximal possible gradient of MOT1 at around 10 G/cm. The MOT1 cell has been equipped withthree UVA LEDs for the purpose of light-induced atom absorption (LIAD), whose efficiency turns
out to degrade over several months’ time. So in practice we are not using this approach any more.
The 1/e loading time constant of MOT1 is 1-2 seconds. The lifetime of the atoms in MOT1 is on
the same scale.
The atoms are continuously pushed from the first MOT to the second MOT by a laser beamduring the MOT load stage. The push efficiency appears to be the limiting factor in MOT2 loading
rate whose 1/e time constant is about 10 seconds. The linearly-polarized push beam is pickedup from the zeroth-order of an AOM in the repumper beam path on D2 transition, and carefully
focused down to its 1 mm waist offset to the side of MOT2 cloud to avoid blasting the atoms in
MOT2 away. As the atoms travel from MOT1 to MOT2, they are further guided and supported
against the gravity by a permanent hexapole magnetic field transverse to the transfer tube. The
5.5 mm inner radius and 47 cm length of the transfer tube require a sensitive alignment of thepush beam. We pulse on the push beam repeatedly when walking the mirrors and optimize thefluorescence signal from the MOT2 cloud on the fly.
MOT2 is a stainless steel chamber with 7 viewports. Its quadrupole field is provided by the

10

2P
3/2

F’=3 (14.4MHz)
F’=2 (-6.7)
F’=1 (-16.1)
F’=0 (-19.4)

2P
1/2

F’=2 (20.8)

39K

push beam
F’=1 (-34.7)

2

767nm

770nm
F=2 (173.1)

2S

3/2

cart transfer
quadrupole oruniform bias field

F=1 (-288.6)

compress 100 G/cm

Γ/2π = 6 MHz

Figure 2.2: The vacuum system and optical transitions of 39 K. Left panel: top view of the vacuum
system and the functions of individual components. Right panel: optical transitions of the D1 -line(red) and D2 -line (cyan) of 39 K.

11water-cooled coils mounted on a moving cart. These coils confine the atoms as they are transferredfrom MOT2 to the science cell. We purposely use low quadrupole gradient 6 G/cm (in the strongconfinement direction) to increase the MOT2 capacity at the price of slow filling rate. We employ
both D2 - and D1 transitions to perform Doppler cooling and sub-Doppler cooling stages in MOT2.
The timeline for the entire laser cooling sequence is summarized in Fig. 2.3. Information on thelayouts of MOT lasers, part numbers of the optics and frequency stabilization systems can be foundin [27]. MOT2 collection (or Doppler cooling process) is optimized iteratively since the push beam
detuning, the magnetic field gradient, and MOT D2 laser detunings are all more or less coupled.
We have the conventional three orthogonal retro-reflecting beam configuration. We suspect thatthe retro-reflecting beams may limit the MOT number due to the dark shadows created by thedense cloud on the retro-reflected beams. We operate the MOT D2 lasers at high I/Isat to havesufficient optical pumping while being able to stay far-red-detuned at the same time. The MOT2
cloud is relatively hot and puffy due to the strong re-absorption of photons. We implement a short
end-MOT stage (equivalent to D2 compressed-MOT) by further red-detuning the laser frequencies
in order to reduce the cloud size for the purpose of fluorescence imaging. It is not completely trivial
to calculate atom number directly from the measured fluorescence intensity due to the complicatedenergy levels (a detailed six-level analysis on 40 K MOT can be found in [26]). Therefore, the MOT2
numbers quoted in this thesis are estimated based on the transfer efficiencies between successive
stages plus the final atom number measured in science quadrupole trap using the more reliable
absorption imaging method.
To further lower the temperature of the cloud, we borrowed the idea of an unconventional
D2 and D1 combined compressed-MOT from [28, 29]. D2 repump beam and D1 trap beam work
together to cool down the atoms to about 200 µK with a number transfer efficiency of 80 % from
the MOT. We then load the atoms to the Λ-enhanced gray molasses stage [30, 28, 31] with anefficiency of over 90 %. This cooling scheme employs the D1 transitions |F = 2i → |F 0 = 2i (“trap”
laser) and phase-coherent |F = 1i → |F 0 = 2i (“repump” laser) to form a Λ-type configuration (seeFig. 2.4 (a)). At the two-photon resonance condition δ = 0, a velocity-selective dark state can be12

D2 trap (F=2 to F’=3)
Power mW
Detuning Γ

D2 repump (F=1 to F’=2)
Power mW
Detuning Γ

Push (F=1 to F’=2)
Power mW
Detuning Γ

200
-7.6

0

-7.8

200

-1.8

-3.4

-3.6

10

0

-5
10

0

D1 trap (F=2 to F’=2)
Power mW
Detuning Γ

110

0.5

0

2.3

10
0

D1 repump (F=1 to F’=2)
Power mW
Detuning Γ

40
2.3

0

4

0

B field gradient
Gauss/cm

100
10

0

MOT1&2,

end-MOT,

N:10E10,
10-20sec, 6mK

20ms

D2+D1
c-MOT,

30

D1QT load and compress
0.5sec, N: 10E9
molasses

10ms, 220μK 8ms, 10μK

start: 23μK, end: 86μK

Figure 2.3: Typical timeline of laser cooling for 39 K. MOT laser powers, detunings and magnetic
field gradient are shown for each stage. Note that the time intervals for different stages are nottrue to scale, and the numbers can vary a bit every time we re-optimize the sequence. The 1/e2
beam waist of MOT1 lasers is about 1.4 cm, the 1/e2 beam waist of MOT2 lasers is about 1.4 cm.
The push beam waist is about 1 mm.

13created (sometimes called velocity-selective coherent population trapping (VSCPT) in literature).
Benefiting from the narrow spectral width of this dark state, a temperature much lower than withthe regular polarization-gradient cooling (PGC) can be achieved. We ramp down the powers ofboth lasers to 1/10 of the initial value over 8 ms to gradually reduce the capture velocity. We areable to constantly reach below 10 µK at the end of the gray molasses. The duration of molasses
is kept as short as possible to minimize the diffusion of the cloud. Our trap laser is locked to39

K saturated-absorption spectroscopy via a modified Pound-Drever-Hall (PDH) lock [27]. We use

a double-pass AOM to bridge the frequency gap and also to adjust the global detuning ∆. The
repump laser is locked to the trap laser via a heterodyne beat note lock (see Fig. 2.4 (b)(c) for
its spectrum). The two lasers are combined on a 70:30 non-polarizing beam splitter to inject a
1W Tapered-Amplifier (TA). A home-built optical cavity is used to monitor the power ratio of
the two tones. Note that as the Distributed Bragg Reflector (DBR) lasers age over a few years’
time, the contrast of phase lock peak will slowly decay. However, this does not impact the overall
performance of our laser cooling process because even without the Λ-enhanced cooling power westill have low enough temperature in molasses to load the atoms into the quadrupole trap.
Following the laser cooling sequence, we suddenly turn on (in less than 0.5 ms) the quadrupole
trap (QT) at about 30 G/cm in the strong-confinement direction. About 1/3 of the atoms, thosein the low-field-seeking spin state |F = 1, mF = −1i, are loaded into the QT. In principle, doing
optical pumping (OP) prior to the QT load could greatly improve this efficiency. However due tothe high density of our cloud, we could not perform an efficient OP without significantly heating
up the cloud. Following the QT load, we adiabatically ramp up the cart magnetic field gradient
to 100 G/cm in 100 ms. This step is to reduce the size of the cloud as initially loaded from themolasses so that it can travel through the half-meter-long transfer tube that is only 1.5 cm in innerdiameter. About 90 % of the atoms can successfully arrive at the science cell in 3.5 seconds with aposition accuracy of a few microns. We then ramp down (up) the cart (science) QT simultaneously
to prepare for optical trap load. The science QT is further ramped up to 193 G/cm as the cart
moves away from the science cell back to its MOT2 position. We consistently have 3 × 108 atoms

14

(a)

(b)

3

1

δ

Δ

Ω1

(c)

Ω2
Γ
2

Figure 2.4: Λ-type configuration for gray molasses and phase lock signal. (a) The three-level Λtype configuration for gray molasses cooling. A strong standing wave with Rabi frequency Ω2and a weaker standing wave with Rabi frequency Ω1 illuminate the atoms at the same time. The
global detuning ∆ is blue-detuned from the D1 transition. The single-photon transition linewidthΓ = 6 MHz. The narrower two-photon transition linewidth leads to enhanced cooling efficiency thatlies within 0.5 Γ near the Raman detuning δ ≈ 0 (see Fig. 2.26 in [27]). (b) The beat note signal
measured with a photo detector is monitored on a spectrum analyzer. The 3dB point represents
approximately twice the laser linewidth which reads 300 kHz in this case. (c) Zoom in on the phaselock peak whose contrast is on the order of 20dB. The span of the phase lock peak between themarkers is about 17 kHz which is limited by the jittering of the VCO box.

15
at 180 µK in science QT.
Our science cell is a fused-silica octagon with seven 0.625” windows on the side and two 1.7”
windows on top and bottom (see Fig. 4.3 in [26] for a picture of the bare cell). Characterizationsof the transmissions can be found in [26] and [32]. We stock additional sample windows for the
purpose of off-line test on our imaging systems.

2.1.2

Evaporative Cooling

Unlike species such as 87 Rb, it is difficult to directly evaporate 39 K in a quadrupole trap due toits small elastic collision rate at low magnetic field. One solution is to perform sympathetic coolingwith species like 41 K or 87 Rb [33, 34]. In order to reduce the number of lasers used in the experimentthus improve its overall stability and reliability, we decided to take the advantage of high-power
infrared fiber amplifier technique and do evaporation in a optical dipole trap (ODT). An ODT helps
us tightly confine the atoms and allow for a uniform DC magnetic field that can be used to adjust theelastic scattering cross-section of the atoms freely. We use a 500mW CW 1064 nm single-frequency
Nd:YAG laser (Coherent, Mephisto 500NE) to seed a 50W fiber amplifier (Coherent, NUA-1064PC-0050-D0) with a customized output beam size of 1 mm. A 25W 28 × 28 µm waist ODT beam(named “H1 beam”) is shone onto the atoms for 8 seconds as they are being held in science QT. The
trap depth is calculated to be 2.7 mK. An H-bridge consisting of four mechanical relays switches
the polarity of one of the QT coils to switch from the original anti-Helmholtz configuration into
Helmholtz configuration during the ODT load. Within the 200 ms switching time, a shim fieldthat is perpendicular to the bias field is temporarily turned on to maintain the quantization axis
for the spin-polarized gas. Omitting this crucial step will lead to spin impurities. The timeline of
QT/ODT load sequence is summarized in Fig. 2.5. Working with a high power laser comes withchallenges involving long-term mechanical and thermal stability from month to month (see [27] forpart numbers of the optics, detailed layouts and how we handled numerous technical issues). Dueto the relatively low atomic density of 1010 cm−3 and small collision rate about 1.3 Hz in the QT,
we end up with 6×106 (up to 1×107 ) atoms at 300−400 µK being transferred to ODT. To improve

16
this loading efficiency, one idea might be to recycle the H1 beam to form a crossed-ODT [34] inorder to gain in the trap depth; another idea could be to borrow the idea of in-trap (sub-)Doppler
cooling from the optical-lattice- or optical-tweezer-based experiments [35, 36, 37]).
We ramp down the power of H1 beam over about 3 seconds time for the stage-1 evaporation.
Meanwhile the atoms are gradually caught by our science crossed-ODT (x-ODT). The science xODT is up to 10 µK deep and formed by a sheet beam (“H2”) along the horizontal direction (roughly45◦ to the H1 beam) plus a cylindrical beam (“V”) along the vertical direction. The trap centerof x-ODT lies spatially below the H1 trap center due to the different gravity sags. The stage-2
evaporation is done in this weakly-confined x-ODT for the purpose of generating low density clouds
in the end. Five orders of magnitude increase in the phase space density (PSD) can be achieved
with this two-stage evaporation scheme. The temperature of the cloud is mainly set by the trap
depth of H2 beam, – fast atoms are leaked out primarily from the bottom of the sheet trap undergravity; the density of the cloud can be adjusted by changing the power on the V beam which
provides confinement along two horizontal directions. The evaporation steps inside the magenta
box of Fig. 2.5 are empirically fine-tuned to prepare the samples at desired initial temperature andatom number. In equilibrium condition, the temperature, kB T of our cloud is found to be about
1/8 of the calculated values of the trap depth. Upon finishing the evaporation, we adiabatiallyramp up the trap depth to be about 10 times the final temperature of the cloud in order to prepare
redundancies for any heating effects when measuring the lifetimes of the atoms. Omitting this step
will lead to extra number loss due to ongoing evaporation.

2.1.3

Crossed-Dipole Trap Calibrations

Huge amount of effort is spent on suppressing the mechanical, thermal, electronic noise on
the ODT system. Fig. 2.6 presents a systematic characterization of the lifetime and heating rate
for different initial conditions. We are able to achieve very low heating rate on the order of nK/s.
The background lifetime (50 - 100 seconds) of the atoms is limited by the vacuum quality rather
than the off-resonance scattering from the ODT photons.

17

H1 beam(28x28μm)

25

τ = 0.4s

Power W

V beam(140x140μm)
Power Wτ = 1s2

H2 beam(80x800μm)

10Power Wτ = 1sScattering length
as a 0

100

-20

200

B field DC
Bias Gauss
Shim Gauss

B field gradient
Cart G/cm
Science QT G/cm

5
100

40

200

MOT2 to Cart toScience sciencetransfer QT transfer
3.8sec43

OT load,mainly into H1
8sec

QT to bias, Evap1 & pancake
Evap2 &
390ms,
trap load,
adiabatic ramp up
PSD: 10E-5

3.5sec, PSD: 10E-3

PSD: 0.1 to >1

400ms

Figure 2.5: Typical timeline of transferring the atoms from cart quadruple trap (QT) into science
QT, optical dipole trap (ODT) load and evaporative cooling. Shim field is orthogonal to bias fieldand provides a transient quantization field when bias coils switch from anti-Helmholtz to Helmholtz
configuration. Note that the steps inside magenta box are to be adjusted in order to obtain a target
initial condition (e.g. temperature, number density). The time intervals for different stages are notshown to scale.

18

(b)

0.026

background loss rate (s-1)

measured at a=63a0 on |1,-1> state

20nK

0.024

background heating rate (nK/s)

(a)

0.022

40nK
0.020
0.018

100-150nK

0.016

24nK

300nK

0.014

370nK
H1+V

0.012

10

300nK
40nK

370nK
H1+V

100-150nK

1

24nK
0.1

20nK

0.010

measured at a=63a0 on |1,-1> state

0.01

0.008
1

2

3

4

5

6

7

8

total ODT intensity at atoms (107W/m2)

9

10

1

2

3

4

5

6

7

8

9

10

total ODT intensity at atoms (107W/m2)

Figure 2.6: Background loss (a) and heating (b) rates measured at small positive positive scattering
length for H2+V dipole trap except for the red circles which are for H1+V trap. The initial
temperatures of individual data points are labeled on the plots.

19
Another important calibration of ODT is the trapping frequencies. We measure the trapping
frequencies by using the slosh motion of the atoms which is initialized either with a short gravity
pull or a side kick from a small magnetic field gradient. The parameters for these steps are carefully
adjusted to give a small slosh amplitude compared to the cloud sizes. The initial condition of theclouds are chosen to give relatively long slosh time without significant dephasing. In the presenceof the decay of oscillation contrast, we fit the slosh motion y(t) into a damped sine curve withfrequency f , damping time constant τ , time offset t0 as below
y(t) = y0 + Ae−t/τ sin[2πf (t − t0 )],and derive the natural frequency f0 = f 2 + (1/τ )2

1/2

(2.1)

of the oscillation according to [38]. The

trapping frequencies are compared to the calculated values based on measured beam waists andbeam powers. Our trap model and calibrations are summarized in Fig. 2.7.

2.1.4

Laser Linewidth Measurements

We employ a self-heterodyne interferometer to detect the spectral linewidths of our Distributed Bragg Reflector (DBR) lasers which we use to inject Tapered-Amplifiers (TAs) for coolingand probing the atoms. The setup is shown in the schematic in Fig. 2.8. Our DBR laser diodes
are fabricated by Photodigm (PH767DBR080T8, PH770DBR080T8). By virtue of the very shortcavities, DBR lasers are highly robust against vibrations and can stay locked for weeks or evenmonths. The drawback is that their linewidths are relatively broad, usually 150 − 600 kHz. Alongwith other imperfections in the imaging system (e.g. non-perfect polarization), this can broaden
the probe lineshape and lead to underestimated number density. These systematic effects usuallyare not an issue since the D2 -lines of alkali metals are as wide as 5 − 10 MHz. But from time to timethey do compromise the precision of bulk gas experiments. In particular, for the three-body recombination rate measurement (ref, chap4), 10 % uncertainty on atomic density calibration will cause
the recombination rate coefficient to be off by 20 %. Therefore it is necessary to run diagnostics ona regularly basis.

20

(a)

(b)

(c)

(d)

(e)
160

80

140
60

(Hz)

(Hz)

120
40
Model
Equation
Plot
y0
xc
A
P
Reduced Chi-Sqr
R-Square (COD)
Adj. R-Square

20

100

Power0y = y0+A*(abs(x-xc)^P);\q(f_x)0±00±037.99831 ± 0.28608
0.5 ± 015.56837
0.998070.99807

Model
Equation
Plot
a
b
xc
Reduced Chi-Sqr
R-Square (COD)
Adj. R-Square

800

freqwithSag2 (User)a*lambertW(b*(x-xc),0,1)\q(f_z)93.7417 ± 12.39145
5.14709 ± 1.79925
0.66789 ± 0.06929
11.12735
0.99763
0.99668

60
0.0

0.5

1.0

1.5

2.0

V power set voltage (V)

2.5

3.0

1.0

1.2

1.4

1.6

1.8

2.0

2.2

2.4

H2 power set voltage (V)Figure 2.7: Calibrations of trapping frequencies for optical dipole trap. (a)(b)(c) Top panels:
representative ODT potentials (solid curves) and parabolic approximation (dashed curves) at thetrap center along three spatial dimensions; bottom panels: corresponding first derivatives, e.g.
restoring forces along three dimensions. (d)(e) Calibrations of horizontal and vertical trapping
frequencies with the slosh motion of the cloud. Fitting uncertainties are within 0.5 %. Red curves
represent fitting functions as shown in the inset tables. Gravity sag is accounted for in the fitting
function for vertical (ẑ) direction.

21

flipper mirror

(a)

(b)
100

λ/2 PBS

2km delay length
150m delay length

λ/2

AOM

amp

amp
spectrum
analyzer

long
bare
fiber

VCO
80MHz

voltage (mV)ND filter

10

PD
1

76

78

80

82

84

frequency (MHz)BS

Figure 2.8: Self-heterodyne interferometer optical layout and beat note spectrum. (a) Schematic
optical arrangement for the self-heterodyne interferometer. The fiber component inside the reddashed-line box that acts as the delay length can be swapped out. (b) Typical spectra downloaded
from the spectrum analyzer. Two different delay lengths (blue and red lines) are plotted.

22
When there is no access to a high-finesse optical cavity or a reference laser as the local
oscillator, self-heterodyne provides a convenient method to characterize a diode laser’s linewidth.
The basic idea of this technique is to convert the phase or frequency noise of the laser into thefluctuations of the light intensity via a Mach-Zehnder type interferometer [39, 40]. Specifically, theoptical field to be calibrated is mixed with a delayed replica of itself (usually by an optical fiber)
and then the AC component of the resulting beat note signal is detected on a fast photodiode(PD). One question people frequently ask is that how long the delay length needs to be to meet therequirement of being much longer than the coherence length of the laser. We tested the viability ofthis interferometric method by doing a systematic study on the delay fiber length. Low-cost fiberrings of various lengths (Fiber XP, singlemode FC-FC 9 µm) are used in the setup. In the ideal case,the phase dynamics of the interfering optical fields are completely uncorrelated on the PD. This
leads to an easy-to-interpret power spectrum which consists of a broadband Lorentzian pedestal(coming from the white noise) and a small Gaussian portion (associated with the 1/f “pink” noise)near the peak (as depicted by the blue line in Fig. 2.8 (b)). For short delay lengths on the order ofthe coherence length of the laser, the phase dynamics of the interfering fields are partially correlated.
This gives rise to a narrow δ-function peak sitting on top of a “scalloped-Lorentzian” (as depictedby the red line in Fig. 2.8 (b)) and complicates the analysis of the measured spectrum.
Generically the beat note signal on PD is associated with the Fourier transform of the intensity
correlation function. By assuming a white noise spectrum on the laser frequency, we can explicitly
write down the beat note signal as the following [39],V (f ) = V1 V2 ,

(2.2)

where
A4 ∆f
,
∆f 2 + (f − f0 )2

∆f sin (2π (f − f0 ) τd )−2π∆f τd
V2 (f ) = 1 − ecos (2π (f − f0 ) τd ) ++ Voffset .
f − f0V1 (f ) =Here the delay time τd = L/c, L is the difference in two beam path lengths which is approximately

23
equal to the length of the delay fiber, and c is the speed of light. ∆f represents the laser linewidth.
Note that the δ-function part is omitted for it introduces instability in the fitting process. When
we use the above Eqn. 2.2 to fit our beat note signal, the delay length is a free parameter that is
extracted from the period of the scallops (a good sanity check). In the limit of ∆f τd  1, Eqn. 2.2
reduces to a Lorentzian shape. As a comparison, we also fit the data with a Lorentzian function.
Note that the laser linewidth equals one half of the full-width-half-maximum (FWHM) of the power
spectrum in the case of a Lorentzian lineshape. The two fitting functions are presented in Fig. 2.9(a) and (b) respectively.
As shown in Fig. 2.9 panel (c), the two functional forms give essentially the same result fordelay lengths over 1km. For short delay lengths, Lorentzian function is no longer valid; however thefull equation of the power spectrum can still predict a more or less correct linewidth. This simple
study justifies the empirical criterion that the laser phases of the two arms can be considered
uncorrelated so long as the delay length is at least 6-7 times the coherence length. We notice thatthe contrast of the scallops predicted by this white noise model is slightly off. By accounting forboth white and pink noise [41, 40] in the frequency noise spectrum, one will be able to get a better
fit over a larger dynamic range.

2.1.5

Absorption Imaging Analysis

We rely on absorption imaging method to obtain information of the temperature and numberof the atoms. In the few-body physics studies, we typically work in the low-probe-intensity regime
to adapt to the clouds’ low density. Detailed explanation on the side-imaging setup, calibrations of
the imaging resolution and the camera can be found in [27]. In this thesis I focus on the calibrationsof the saturation intensity especially how we account for the Doppler shift induced by the atomic
recoil momentum.
Absorption imaging detects the shadow cast by an atom cloud in a probe laser beam that
illuminates the atoms. The atoms interact with the laser field and scatter photons into the vacuum
modes. We use the cycling transition |2, −2i → |3, −3i to probe the atoms, which is a good two-level

24

(a)

voltage (mV)

15

(c)

10

800

5

0
76

78

80

82

84

frequency (MHz)(b)
15

fit linewidth (kHz)700

Coherence Length

full equation of the noise spectrum
Lorentzian fit

600

500

voltage (mV)

400

10
300
100

1000

delay length (meter)

5

0
76

78

80

82

84

frequency (MHz)Figure 2.9: Comparison of the self-heterodyne beat note spectrum (black dots in panel (a)(b)) fit
with Lorenztian (red line in (a), red open diamond in (c)) and Eqn. 2.2 (blue line in (b), blue
solid circles in (c)). The vertical dotted line in panel (c) indicates the coherence length of the laser
extracted from the fitting result in the long-delay-length limit, which is consistent with the value
given by Photodigm.

116 m

25
approximation. The scattering rate of photons is [21]
Γsc =

Γ
I˜
,
2 1 + 4δ̃ 2 + I˜

(2.3)

where I˜ = I/Isat is probe intensity in unit of the saturation intensity, δ̃ = δ/Γ is the probe frequency
detuning in unit of the natural linewidth Γ. According to the modified Beer-Lambert law [42], theRcolumn density of the atoms ρ = n(z)dz along the imaging line of sight ẑ isσ0 ρ = − 1 + 4δ̃ 2 ln (If /Ii ) + I˜i − I˜f ,(2.4)

where σ0 = 3λ2 /(2π) is the resonant scattering cross section for optical wavelength λ. In the limit
of low saturation intensity I˜i  1 and at zero detuning, the right-hand-side of the above equation
reduces to the optical depth expression that people normally use to obtain the atomic column
density σ0 ρ ≈ OD ≡ −ln (If /Ii ), with Ii and If being the light intensities in “light” (initial) and“shadow” (final) frames. In this case, OD is independent of probe intensity and probe duration. In
practice, the signal detected on a camera is the total number of photons (or fluence, usually called“light counts”) that arrive at the CCD within a certain probe duration which is gated by an AOM.
A third fame taken with the camera shutter closed (called “dark” frame) is always subtracted out
from the light and shadow frames so that the CCD offset is removed.
In the cases of light-weight atom species or a long imaging time, the atoms can be easily
accelerated during the imaging process. An atom absorbs momentum kick of ~kr = h/λ in ẑ
direction from scattering a single photon. The associated recoil velocity is vr = ~kr /m for atomic
mass m. In addition to gaining a net momentum pz in the longitudinal direction, the atoms also
spread diffusively in the transverse directions due to a non-zero variance of pz . Note that the probe
intensity’s effect on pz implicitly depends on the displacement ẑ as the light field decays whentraveling through the cloud. We implement a numerical model to simulate the longitudinal velocity
and displacement of the atoms in coarse-grained space-time. The basic idea of the method is shownin Fig. 2.10. For practical imaging times of tens of µs, the atoms have gained velocity but not
moved appreciably compared to the size of the cloud. Therefore the atoms can be presumed to bestationary during the imaging process.

26①

orientations ofpropagation

②

I1
δ1
z1

I2
δ2
z2

I3
δ3
z3

…

I1e3
δ1e3
z1e3

t2

I1
δ1
z1

I2
δ2
z2

I3
δ3
z3

…

I1e3
δ1e3
z1e3

I1
δ1
z1

I2
δ2
z2

I3
δ3
z3

…

I1e3
δ1e3
z1e3

…

t1

t1e3

Figure 2.10: Coarse-grained method for absorption imaging simulation. The probe duration is 60µs, I/Isat = 0.1. The imaging time of an atomic cloud with density n(z) is coarse-grained into a
103 × 103 grid. The absorption process of the probe light is propagated first along z then alongt. The cross-section of the atomic density distribution along z is plotted before (top panel) andafter (bottom panel) the imaging process. The initial distribution of the atoms is assumed to be
Gaussian. See Appx. A for the algorithm of this model.

27
To visualize the effects of recoil momentum, we compare the predictions of OD vs. I/Isat by
Eqn. 2.4 and by the numerical model as presented in Fig. 2.11 panel (a). For short detection time of
10 µs, the modified Beer’s law is a good approximation up to I/Isat = 1. For the detection time of
40 µs that we normally use on our side imaging, the two models start to deviate at I/Isat as low as
0.1. Eqn. 2.4 overestimates the OD by assuming that the atoms scatter the light constantly at the
maximal rate for detuning δ = 0. As shown in Fig. 2.11 panel (b), the fictitious “sheets” of atoms
at the front, middle and back of the cloud transport by different distances due to the different
effective intensities seen by the atoms. For I/Isat = 0.1 and probe duration longer than 120 µs, thedifferent atom sheets start to reorder. This will further complicate the dynamics and distort theoriginal spatial distribution of the atomic density (i.e. Gaussian). Due to a limited depth of focus(DoF) of the imaging optics, this translation motion will also potentially compromise the imaging
resolution.
A key parameter in absorption imaging technique is the saturation intensity Isat . The lack of
a good knowledge of Isat (in unit of light counts on a CCD camera) will introduce large systematic
errors on the absolute values of atom number and density. This problem will become more serious
when operating in the high-probe-intensity regime. One approach to measure Isat is by directly
calibrating the efficiency of the CCD camera (i.e. the product of quantum efficiency and analogto-digital conversion gain), and translate the theoretical Isat in unit of mW/cm2 to camera counts.
The accuracy of this method sometimes can be limited by not knowing the exact total transmissionof the numerous optics in the imaging beam path (in particular in between the atoms and the CCD,
including the science cell window, imaging optics, band-pass filters etc.). For instance, we foundthat some camera models in the lab have been installed an uncoated window for the CCD vacuum.
For this reason we often choose to use atoms to calibrate Isat . We fix the probe duration in eachscan and image the atoms with a wide range of probe intensities I. Ii and If are extracted byaveraging the light counts within a small box at the center of the cloud where the column density
is relatively homogeneous. We then fit the data of OD vs. Ii into the modified Beer’s law and ournumerical model respectively (see Fig. 2.12 panel (a)). Note that the fitting parameter Isat is in28

Probe Duration 40μs

−𝑙𝑛 𝐼𝑓 Τ𝐼𝑖

Probe Duration 10μs

−𝑙𝑛 𝐼𝑓 Τ𝐼𝑖(a)

𝐼/𝐼𝑠𝑎𝑡(b)

𝐼/𝐼𝑠𝑎𝑡
𝐼
= 0.1
𝐼𝑠𝑎𝑡

displacement in 𝑧Ƹ (μm)

Atom Displacement at

probe duration (μs)

Figure 2.11: Absorption imaging simulation results. (a) Optical depth OD = −ln (If /Ii ) as afunction of the probe intensity I scaled by Isat . Only in the limit of I/Isat  1 does OD/σ0
represent the column density of the cloud ρ. Blue line: model using Eqn. 2.4; red open circles: ournumerical simulation. (b) Displacements of atom sheets at the front, middle, and back of the cloud
as a function of the probe duration. The simulation apparently breaks down when the different
sheets of atoms start to reorder at around 120 µs.

29fact the saturation intensity integrated over a certain probe duration.
As shown in Fig. 2.12 panel (b), fitting scheme using Eqn.2.4 yields an Isat that depends on
the probe duration in a non-linear fashion. On the contrary, our numerical model gives the expected
linear behavior. This difference in Isat qualitatively agrees with the observed difference in OD as
discussed above (remember Fig. 2.11 (a)), – ignoring the recoil momentum leads to overestimatedOD therefore a smaller Isat value is needed to correct for the OD. As a sanity check, our Isat for40 µs probe duration is 3300(300) counts according to the numerical fitting model that includes
the recoil effect. From an independent calibration of the CCD camera efficiency, Isat turns out tobe roughly 5008 counts for the same probe duration. The Isat is only 1700(200) from the fitting
function that uses the recoil-free model.
In the end, we apply our numerical model with calibrated Isat to the raw images of the atomic
cloud in a pixel-by-pixel manner. Fig. 2.13 panel (a) provides this conversion factor for certain probe
durations in the standard low intensity regime of I/Isat = 0.1. Based on this conversion factor,
we can estimate how large systematic error there will be if using Beer’s law OD = −ln (If /Ii ) toapproximate the column density of the cloud. For the lowest density clouds we have, this error isalready on the order of 10 %. For imaging extremely low density clouds with relatively long probe
pulses such as in [43], this correction factor can be even larger and becomes critical in obtaining
the correct shape of the cloud.

2.1.6

Magnetic Field Control

Rf and microwave spectroscopy provides an essential component of the toolkit for quantum
gas experiments especially in quantum state preparation and precision metrology. To achieve highaccuracy with this technique, a stable control of the DC magnetic field plays is critical. The
shot-to-shot magnetic field fluctuation in our experiment is determined to be 2 − 5 mG (standard
deviation). And benefiting from the long-term stability of its absolute value, the uncertainty on
the magnetic field can be averaged down to within 1 mG by taking only 20 shots. We calibrate ourmagnetic field noise by driving different hyperfine transitions on the atoms (i.e. |1, −1i → |2, 0i,

30

(a)

−𝑙𝑛 𝐼𝑓 Τ𝐼𝑖

Recoil Model Fit (40μs Probe)−𝑙𝑛 𝐼𝑓 Τ𝐼𝑖

Recoil-Free Model Fit

Comparison of Two Fitting Methods

(b)

4000
18000
16000

3000

14000
12000
2000

10000
8000

1000

Model

Line

Equationy = A + B*x

ReducedChi-Sqr0.094546000
4000

0.99734

Adj. R-Square

Value
Isat with recoil

Isat from recoil model (cts)Isat from recoil-free model (cts)

20000Standard Error

A

0

0

B

81.99087

1.59754

0

2000
0

0

100

200

probe duration (us)Figure 2.12: Comparison of two fitting models for Isat . (a) Two fitting protocols using Eqn.2.4 (leftpanel) and our numerical model that takes into the recoil-momentum (right panel) respectively.
The data is taken with 40 µs probe duration. (b) Comparison of the results from two fitting models(black circles: left y-axis; blue diamonds: right y-axis) as a function of probe duration. The solid
blue line is a linear fit to the blue diamonds with fitting parameters listed in the inset table.

31

(a)

OD/𝜎0 𝜌Conversion Factor for Different Probe Durations

𝜎0 𝜌

(b)

𝝈𝟎 𝝆 for Column Density(conversion
factorx5exaggerated)exaggerated)
OD corrected(5 times

OpticalDepth
OD uncorrected

1.2

1.2

20

20
40

1

40

1

60

60
0.8

80

0.8

80
100

100

0.6

0.6

120

120
140

0.4

140

0.4

160

160
0.2

180
200

0.2

180
200

20

40

60

80

100

120

140

160

180

200

20

40

60

80

100

120

140

160

180

200

Figure 2.13: Application of the recoil model to the measured optical depth. (a) OD to column
density conversion factor for different probe durations at I/Isat = 0.1. (b) Application of theconversion factor to a raw image (left penal) pixel by pixel. The conversion factor is enhanced by
5 times to visualize its effect on the density distribution of the cloud (right panel).

32|1, −1i → |2, −2i, whose ∆µ differ by a factor of 10) and measuring the variance of the transition
fraction P = N|2,xi /N|1,−1i (x = 0, 1) at the microwave detuning where P ≈ 0.5. The shot-to-shot
atom number fluctuation is independently checked at zero detuning where P = 1 and turned out
to be constantly below 2%.
We put a lot of efforts into improving the stability of the bias field control. The schematic
diagram for our science coils control system can be found in Fig. 2.43 in [27]. This setup is built
based on a JILA coil driver (JD013A2) designed by Brian Demarco. This PI feedback servo includes
only a minimal number of op-amps hence suffers less from the internal noise. Two copies of this coildriver with different servo gains are combined at the servo output in order to broaden the dynamic
range. High current up to over 300 A in science quadrupole trap and low current less than 15 A inbias field can be switched back and forth. An Insulated-Gate Bipolar Transistor (IGBT, Semikron)serves as the current actuator. Fluxgate transducers (Danisense for in-loop, LEM for out-of-loop)
along with high stability sensing resistors (Vishay VPR221ZT) are implemented to monitor thecurrent.
As to the rf and microwave sources, we employ a programmable frequency generator (versatilefrequency generator, VFG-150) to produce fast pulses with the frequency, amplitude and phaseedited in a step-by-step manner. It can output frequencies up to 150 MHz and features by a step
switching time of 5 ns with approximately 1000 step number. We control the VFG via USB 2.0
which typically takes a couple of seconds to upload each sequence 1 . This data transport happens
during MOT filling stage so as not to take extra experiment time. The signal from VFG is mixed up
to half a GHz with two frequency synthesizers (from Keysight) to generate the hyperfine transition
frequencies. To do so we picked a Mini-Circuits mixer (part number ZX05-10+) that has the weakest
harmonics among various models. Analog frequency modulation has also been implemented on thesynthesizers with external DAC for modulating slower pulses. The synthesizers are programmed
via SCPI commands. Air-cooled amplifiers from Mini-Circuits (20W, 20 − 1000 MHz, part number
1

The hardware driver and Matlab scripts for VFG-150 are provided by Photonics Technologies Ltd and Department Physik, University of Siegen.

33
ZHL-20W-13+) are installed to broadcast the signal. Note that we use dedicated Topward power
supplies and AC outlets for these big rf amplifiers because each one of them can draw huge amountof power from the wall transiently during the pulses and also that the rf signal can easily propagate
backwards into our electronic network through any poorly-shielded connectors or cables.
We use home-made loop antennae to deliver the rf and microwave to the atoms. The atoms arein the near-field regime of the dipole radiation hence mainly sense the magnetic field that oscillates
in the axial direction. Our science cell is currently equipped with 5 antennae as shown in Fig. 2.14
panel (a). A photograph of a one-loop antenna can be found in the inset of Fig. 2.14 panel (b). The
head of our antennae are made of 16 gauge magnetic wire and semi-rigid co-ax. The co-ax supports
the antenna in a robust way so that the antenna is not in contact with nearby dielectric materials
and also its coupling to the environment is maintained at a constant level. All these help to keepa steady radiation efficiency and thus reproducible Rabi frequency from day to day. The downsideof a co-ax is its relatively high capacitance of about 30pF/ft. This effective capacitance together
with the number of loops determine the maximal frequency we can operate. The antenna headis mounted on a impedance matching box through non-magnetic SMA connectors. High-voltage
trimmer capacitors up to 100 pF are used to tune the resonance frequency as well as to match the
impedance. By carefully tuning the antenna after it is mounted in place, return loss is limited to
be on the order of a few percent (see Fig. 2.15 for some reflection measurements).
After we installed the 105 MHz longitudinal antennae which are co-linear with the bias coils,
we found that their radiated powers are strongly coupled into the bias coils and influence the bias
field servo system. In the far field, this rf radiation mainly propagates in the horizontal plane. It
can unlock all lasers in the lab if pulsed on for a few seconds at full power. RFI rectification happens
when rf signal gets picked up by conductors at the inputs, outputs and power supply pins on opamps or in-amps through capacitive, inductive or radiation coupling. We made a series of trials toscreen out this spurious signal on our bias-field control system. For instance we implemented lowpass filters on the inputs of the bias servo, used shielded twisted-pair cable on the servo output,
enclose all feedback components in a shielding box, and relocated the setpoint signal from an

34
(a)

(b)
bias coils

stub
tunersfast-B coils

Figure 2.14: Photographs of rf/microwave antennae mounted surrounding the science cell andantenna components. (a) Five antennae surrounding the science cell can be seen on this photograph.
They are mounted as close to the atoms as possible. We have three transverse antennae (labeled
1,2,3) that address spin-flip transitions and two longitudinal antennae (labeled 4,5) that address
non-spin-flip transitions. The bias field along vertical direction is provided mainly by bias coils andsometimes by fast-B coils for fast switching time. (b) Top left: photograph of a single-loop antenna
head. The diameter of the loop of about 2.5 cm is made to be slightly larger than the side windows
on the science cell. Top right: small trimmer capacitor that we use for tuning the resonancefrequency of the antenna. Bottom: big trimmer capacitor for impedance matching enclosed insidea rf-shielded enclosure.

(a)

(b)

Figure 2.15: (a) Reflection (on linear scale) of a microwave antenna measured on a network analyzer
as the in-place tuning of the antenna is finished. The Q-factor is about 100 in our case. (b)Reflection (top trace, channel 1) and transmission (bottom trace, channel 2) received by a pick-up
coil. We always use a small-diameter (a few mm) pick-up antenna to detect the radiation as we dopreliminary tunings. This is to ensure that the observed drop in reflection actually corresponds toa peak in power radiation. Sometimes there are electronic resonances that occur mainly within thetuning box (e.g. the signal to the right on the figure) usually lying at higher frequencies. These
spurious signals are misleading since they can also be received by the pick-up antenna.

35
FPGA board in the control computer to an independent DAC board inside the shielding box thatis controlled only through data and digital lines. We found that the rf radiation enters into ourbias servo system mainly through the inputs (i.e. signal from the current transducer) and also gets
picked up on the servo output to the gate of IGBT. We made some plug-in low-pass filters thatcan be easily installed (as shown in Fig. 2.16). Most rf noise can be rejected by these filter boxes.
Doubling the number of these filter boxes on each port does not further suppress the noise. We also
replaced the servo output cable with a short well-shielded twisted-pair cable (silver plated copper).
These two remedies have already greatly improved the robustness of the bias field control system.
We found that adding clamp-on ferrite chokes or Mini-Circuits in-line low-pass filters have very
minor effects in our case.
(a)

(b)

Figure 2.16: Photograph (a) and schematic (b) of the plug-in low-pass filter. We use bolt-in π-filters(5500 pF, cutoff frequency a few MHz) to shunt rf noise to the chassis.

We always use atoms as the ultimate caliper for rf rectification effects considering that magnetic field or current sensors potentially also suffer from the rf radiation and can shift their measurement in an unpredictable way. A 2 ms rf signal at about 105 MHz is pulsed on at its full power 20
W to perturb the bias servo system; meanwhile a second Gaussian π-pulse is sent to the atoms forfield detection by driving a hyperfine transition. The Gaussian pulse peak power is only 0.1% thatof the rf square pulse. As presented in Fig. 2.17 panel (a), with the original bias servo setup, therf radiation can cause a 0.25 G/W shift on the DC level of the bias field. With the modified setup,

36
we are able to mitigate the shift down to 0.15 mG/W (150 ppm/W). This level of stability on themagnetic field plus the high power rf signal allow us to drive some molecule-to-molecule transitions
at certain frequencies. However unfortunately due to the extremely small dipole moment of thesetransitions we ended up still being limited by the Rabi frequency. And adding more power to theantennae or shifting the resonance frequency to a few MHz lower value resulted in higher noise levelon the bias field, causing us to to lose our precision on the spectroscopy.

2.2

High-Resolution Top Imaging
We take the advantage of the large optical access on the top and bottom windows on thescience cell to do high-resolution in situ imaging of the density distributions of the atoms. This top
imaging scheme is compatible with our pancake-shaped dipole trap and ideal for 2D or lattice-based
cloud geometry.

2.2.1

Design of a Compact Imaging System of NA 0.4

The main idea of our high-resolution imaging system design is to achieve diffraction-limited
performance with minimal pieces of optics in the objective lens so that we have enough clearance forour moving cart coils. We end up using a Thorlabs CNC-polished asphere (part number AL4532-B)
as the main component that sets the numerical aperture (NA), together with a thin plano-convex(LA1978-B) to compensate for the optical aberrations induced by the cell window. The two pieces
are mounted in an Ultem holder which is eddy-current-free. The structure of our objective lensis shown in Fig. 2.18. The specifications of the objective according to the design are listed as thefollowing: effective focal length 31.03 mm; NA (defined as the acceptance cone) 0.384; exit pupil
diameter (defined as the image of the stop) 26.15 mm; limiting aperture (limited by the spacer
ring in between the two pieces) 22 mm; Airy disk radius at 767 nm 1.13 µm. The eyepiece we use
is a 1500 mm plano-convex from CVI optics (PLCX-50.8-772.6-C). This combination of objective
and eyepiece gives a paraxial magnification of 48.70. As a comparison, our side imaging that we
mainly use for time-of-flight measurements has NA about 0.17 and Airy disk radius of 3.49 µm.

37(b) 2.0x105

2.0x105

no rf
no rf
0.1W
0.2W
0.4W
0.8W

1.5x105transfered atom number

transfered atom number(a)

1.0x105

5.0x104

0.0
-14

-12

-10

-8

-6

-4

-2

0

2

|1,-1> to |2,0> transition (+446565 kHz)

4

6

GaussianPulse (Equation

A*exp(-(x-xc)^2*(

ReducedChi-Sqr9.68768E7
0.96821

Adj. R-Square

Valu

no rf
no rf
6W
12W
25W
50W

1.5x105A
tau
N

1.0x105

5.0x104

-10

-8

-6

-4

-2

0

offset
xc

0.0
-16

Model

2

|1,-1> to |2,0> transition (+446565 kHz)

Figure 2.17: Rf spectra taken with various rf radiation power being present. The magnetic momentum of the hyperfine transition used here is about 0.15 kHz/mG. (a) Without good shieldingin the servo system, the magnetic field is shifted significantly by the rf signal. (b) With propershielding and filtering, the rf signal neither affects the center nor the width of the rf lineshapes. We
took microwave spectra without rf radiation being present before and after all the scans in order
to trace the background magnetic field fluctuation that takes place during the experiment. Solidlines are fitting curves with Fourier-limited Gaussian function.

1

38
Its magnification is only 6.7. The top imaging objective requires more delicate alignment than theside imaging lenses. In order to have all the degrees of freedom on the lens mount, we installed alinear translation stage from Newport for XYZ axes, a goniometer from Edmund (66-534) for pitchaxis, and a rotary stage from Edmund (66-516) for roll axis.
The big asphere in our objective is trimmed down to fit in the small empty space encircled
by our magnetic coils. As a quick assessment of the quality of this asphere, we generate an 1
inch diameter collimated beam with it and then use a lateral-shearing interferometer to check the
wavefront errors qualitatively. We did this test on the two big aspheres that are cut by JILA shop
and Thorlabs respectively. This method turned out to be limited by the imperfect point sourcethat we make by focusing down a laser beam with a small molded asphere (8 mm focal length, NA
0.5, part number C240TME-B). As shown in Fig. 2.19, the distortions on the interference fringes
for the two lenses are very alike. The fact that this distortion pattern is independent of the rotation
angle of the big aspheres indicates that the aberrations are mostly from the molded asphere. In
the end, we picked the shop-cut one that has relatively less distortions.

2.2.2

Bench Test Setup

To calibrate the resolution of the top imaging system, we built a testing setup with a small
pinhole as the point source, and a sample window to emulate the cell window. The optics layout
is presented in Fig. 2.20. The pinhole is mounted on an indexed repeatability mount (part number
NX1NF) so that it can be swapped out easily without too serious misalignment. We also preparedlarger-diameter pinholes (1 µm, 10 µm) for the purpose of preliminary alignment, – larger pinholes
allow higher transmitted power therefore assist the concentricity and back-reflection alignments ofthe downstream optics. A photograph of the point source setup is shown in Fig. 2.21 panel (a).
We purchased the pinholes from National Aperture. The smallest pinhole we have (part number
1-0.5+HS) is 0.5 µm in diameter and the tunnel is about 3 µm thick. An SEM image of the pinhole
is shown in Fig. 2.21 panel (b). It could take many hours to locate the real pinhole among all thedust particles under the SEM microscope. We found it easier to look for the perforation structure

39
(a)

(b)

16.7

0.57

2.2 6.59 6.35

10.46

Figure 2.18: Design of our compact top imaging system in (a) Solidworks and (b) Zemax. The
numbers in (b) are central thichnesses in unit of mm.

(a)

(b)

Figure 2.19: The JILA shop-cut (a) and Thorlabs customized (b) aspheres display similar fringes
under a shearing interferometer. This pattern is actually associated with the 8 mm molded aspherethat we use to generate the point source.

40
on the back of the substrate. A Vortex laser at 767 nm with about 5 mW output power serves as
the light source. By focusing down the beam on the pinhole with a 30 mm achromatic doublet, weare able to get a sufficient amount of 100 nW light traveling through the pinhole.

0.5um
pinholesamplewindow

767nm source

repeatability
platform

2x expander

objective
f1 = 32mm
eyepiece
f2 = 1500mm

DataRay CCD

Figure 2.20: Optics layout for the top imaging bench test setup.

The image of a small pinhole is called the point spread function (PSF) of an imaging system.
As been plotted in Fig. 2.22, the intensity of a pinhole image can be approximated by the Airy diskfunctionI(r) =

2J1 (r)
r

2
,

(2.5)

where J1 is the Bessel function of the first kind of order one, r is the radial coordinate in the Fourier
plane of the objective. We can relate r to the Cartesian coordinates in the image plane throughr = 2π(x2 + y 2 )1/2 /(M λN.A.), where x and y are in unit of pixels, λ is the wavelength of the light,
M is the magnification of the imaging system, and N.A. stands for numerical aperture. By fitting
the PSF with the Airy disk function, we obtain an NA of about 0.4, which meets our expectation.
The Strehl ratio is estimated to be over 90%, well within the diffraction limit.
After finalizing the alignment, we study the relative sensitivity of each component in thissetup. We found that the tilt degree of freedom on the objective is critical in attaining the smallest
spot size. It is coupled to the translation degree of freedom so both of them need to be adjusted
at the same time. We found that it is actually easier to see the tilt by purposely defocusing the

41

(a)

(b)
sample window
30mm focus lensxyz 3 axis stage

10/1/0.5um
pinholeUltem objective
holderFigure 2.21: (a) A photograph of the point source setup with part numbers labeled on the picture.
(b) An SEM image of the 0.5 µm pinhole (top) and dimensions of its substrate (bottom).

42
objective and using the higher order annular structure as a reference. We use micrometers on
all five knobs on the objective lens mount including translations and rotations. The longitudinaltranslation knob consists of fine and coarse micrometers for easy focusing. The concentricity of
the asphere and plane-convex is another crucial factor. Assembling the lenses is a meticulous job.
Machining errors on the lens mount and the spacer ring can lead to extra aberrations especiallyanisotropic portions. We found that the tilt degrees of freedom on the sample window and eyepiecelens do not matter as much. Regular knobs on these two optics elements would suffice. The
translation of the sample window is also insensitive for this kind of intermediate NA.
On the bench test setup, we use DataRay beam profiling cameras to image the PSF. The
two models we use have pixel sizes of 3.2 µm and 4.4 µm. For such small pixels, discreteness ofthe sampling is not a concern. For our science camera, however, the pixel size is as large as 16 µm.
We chose the magnification to be 48 so that our Airy disk radius is a couple of times larger thana single pixel. In general, the Nyquist criterion can be satisfied when the Abbe radius is at least
twice as big as the pixel size 2 . To study how pixelation can affect the imaging resolution, we
numerically bin the adjacent pixels into groups to mimic the effect of sampling with larger pixels.
The borders of the groups are selected in a random fashion. As shown in Fig 2.23, larger binning
numbers lead to an effectively larger spot size as well as a slightly wider spread in the spot sizes.
Images taken with our science camera are equivalent to a binning number of 5 on this DataRay
camera whose pixel size is 3.2 µm. We verified that for a 5 × 5 binning the resolution stays as goodas the no binning case.
Finally we investigate the finite depth of focus (DoF) effect which is of practical importance tous due to the thickness of the cloud. The DoF is characterized by the Rayleigh length of a Gaussian
beam and is inversely proportional to NA squared. In the test setup, we spatially translate theobjective and observe the evolution of PSF as a function of the longitudinal displacement. As shownin Fig. 2.24 panel (a), the central peak of the PSF diminishes very quickly and also asymmetricallyon different sides of the focal plane. This will lead to a complicated resolution function on the clouds
2

The Abbe radius is defined as rAbbe = λ/(2N.A.). See Section 2.2.3 for detailed analysis of the optical aberrations.

43

Figure 2.22: Point spread function (PSF) of the top imaging system generated from a 0.5 µm
pinhole. The region inside the yellow square is fit into Airy disk function and the fitting parameters
are shown in the top right corner. Slices of data cut through the center of the cloud along thevertical (red) and horizontal (green) lines are shown in the right and top panels respectively. The
corresponding fitting functions are plotted as black lines in each panel. Both panels zoom in on
the wiggles around the main peak.

44

(a)
1x1

(b)

10x10

5x5

2x2

Pixelation Simulation
0.60

sx

0.55

sy

mag = 24,rAbbe = 1.4*px

s from fit (mm)

0.50
0.45
0.40
0.35
0.30

mag = 48,
rAbbe = 2.8*px

0.25
0

2

4

6

8

10

12

14

16

binning number

Figure 2.23: Study the pixelation effect by numerically binning the PSF. (a) Original image (leftmost) and the original image binned by different numbers. The third image from the left (5 × 5binning) is what the PSF would look like on our science camera. (b) By doing various binning
numbers on the original image and fitting the binned images into a 2D Gaussian function, one cansee how the artificial large pixels compromise the imaging resolution through undersampling. The
Gaussian r.m.s. width σ is related to Airy disk radius through σ ≈ 0.4rAiry . The difference in σxand σy indicates anisotropic optical aberration terms.

45that are thicker than 10 µm along the vertical direction. We machined two threaded apertures outof Ultem so as to iris down the objective to as small as half inch in its diameter. This gives an
elongated DoF up to 26 µm as measured in Fig. 2.24 panel (b). The spot size at focus in this caseis correspondingly twice as big as that of the full aperture.

2.2.3

Analysis of Optical Aberrations

As can be seen from the PSF fitting in Fig. 2.22, although the Airy function yields a reasonablevalue of NA, it does not describe the amplitudes of wiggles perfectly. This deviation is due tothe presence of optical aberrations, i.e. the objective lens imprints a spatially dependent phase
distortion on the wavefront of the light field in addition to its primary effect which is to cut off highspatial frequencies in the momentum space. We played with the measured PSF and determined
the type of aberrations in terms of Zernike polynomials.
An optical imaging system transforms radiation sources into blurred pictures. The pictures
can be mathematically expressed as a convolution of the original distributions of the radiation withthe point response function of the imaging system, – this is where the PSF comes into play. The
PSF in momentum space is named modulation transfer function (MTF) (or optical transfer function(OTF)) which is basically a low pass function that only allows information with spatial frequencies
below the so-called Abbe frequency to be transmitted. Both PSF and MTF are related to the exit
pupil function, which describes the physical process that actually happens when a radiation field
passes through the objective lens. For an ideal lens, it is a box function with unity transmission.
We now derive the PSF and MTF associated with an ideal pupil for fluorescence and absorption
imaging respectively.
In an imaging process, the incident field E0 is shone on an atom (a realistic point source) andscattered by the atom in the object plane. The scattered field in the image plane can be written as∆E ∝ E0 eiδ p̃(k), where δ is phase shift due to probe laser detuning, p̃(k) is Fourier transform of
the exit pupil function p(r). Here, k = 2π(N.A.)r/(λra ) relates the position r in the object planeto the position R = M r in the image plane, where ra is the radius of the pupil aperture and M is

46

(a)

z=30 um

40 um

48 um

54 um

62 um

σ=0.37um(b)
aperture and adapter

z=30 um

40 um

50 um

70 um

90 um

σ=0.64um

Figure 2.24: Explore the finite depth of focus of the PSF. (a) By translating the objective lensalong longitudinal direction ẑ, the central peak on the PSF diminishes very quickly. This shows
the sensitivity of the PSF to defocus for this full aperture of about 23 mm n diameter. The
corresponding Rayleigh length is only 5.5 µm. (b) We made apertures in different diameters (threequarter inch, half inch) to iris down the objective. As shown in the bottom row, an aperture of 12mm in its diameter leads to a longer Rayleigh length of about 26 µm at the price of a twice bigger
spot size at focus.

47
the magnification [44]. For fluorescence imaging, the intensity of the scattered light is collected ona camera. The PSF of this imaging process is therefore
P (r) = |∆E|2 ∝ |eiδ p̃(k)|2 .

(2.6)

For absorption imaging, the image of an atom is instead formed by interfering the scattered fieldwith the incident field. In the limit of a small scattering amplitude, the transmission of the incident
beam is t2 = |E0 +∆E|2 /|E0 |2 ≈ 1+2Re[∆E/E0 ], where Re[.] refers to the real part [45]. Therefore
the PSF in absorption imaging is
P (r) ∝ Re[∆E/E0 ] = Re[eiδ p̃(k)].

(2.7)

In both imaging schemes, the MTF is the Fourier transform of the PSF, i.e. M (k) ≡ |P̃ (k)|. In
Fig. 2.25, the PSF and MTF for an ideal exit pupil are plotted in radial coordinate since both
functions demonstrate azimuthal symmetry.
In practice, an optical lens can introduce spatially inhomogeneous transmission as well asphase distortions. This leads to a complex exit pupil functionp(r, θ) = eiΦ(r/ra ,θ) .

(2.8)

We assume a homogeneous unity transmission across the pupil here. By plugging Eqn. 2.8 into
Eqn. 2.6 or Eqn. 2.7, we can obtain the resolution functions for fluorescence and absorption imagingsin the presence of optical aberrations. In general, Φ(r/ra , θ) can be expanded in Zernike polynomials, among which the leading terms are known as tilts, astigmatism, defocus, trefoil, coma,
spherical, etc. [46].
We fit our measured PSF from the pinhole test with the wavefront distortion that consists ofthe first eight terms of Zernike polynomials 3 . We expand the phase of the pupil function as the3

There are different conventions for the prefactors in these terms. We simply followed the basis functions from
Wikipedia.

48

(a)

r

r

r

r

r

(b)

r

Figure 2.25: The PSF and MTF associated with an ideal pupil for fluorescence (a) and absorption
imaging (b) as a function of radial coordinate. The NA and λ used here are the same as in our top
imaging system.

49
following series
!
 4
 2
√
Φ(r/ra , θ)rr
= S0 5 6
−6
+1
2π
ra
ra
 2
 2
√
√
r
r
+ α1 6cos 2θ + α2 6
sin 2θ
ra
ra
!
 2
√
r
+β 3 2
−1
ra
 !
 !
 3
 3
√
√
r
r
r
r
−2
sin θ + δ2 8 3−2cos θ+ δ1 8 3
ra
ra
ra
ra
 3
 3
√
√
r
r
+ 1 8sin 3θ + 2 8
cos 3θ.
rara
The values of the polynomial coefficients in the above equation are found to be S0 = −0.0978,
α1 = −0.0065, α2 = 0.0008, β = 0.0027, δ1 = 0.0053, δ2 = 0.0136, 1 = −0.0021, 2 = −0.0008 for
an NA of 0.384. The fitting results are summarized in Fig. 2.26. The dominant contribution comesfrom the spherical aberration. The fitting program can be found in Appx. B.

2.2.4

Top Imaging Calibration

We implemented this top imaging setup to our apparatus. The optics components andarrangement are shown in Fig. 2.27. Before the objective lens has been placed in the beam path,
the probe laser beam is guaranteed to propagate rigorously along the holes on the optical table sothat we can quickly estimate the total beam path length. The major challenge in the top imaging
alignment comes from the large magnification. We could easily miss the image of our cloud on
the camera if any one of the reflecting mirrors is not in the right angle. Plus we are not able todo the alignments on-the-fly due to the long duty cycle of the experiment. We did a step-by-steppreliminary alignment with the probe beam first so as to stay close to the correct beam path asbest as we can. These steps turned out to work great for us. Once everything has been set up on
the table, it only took us a day or two to see some signal from the atoms, first in a 1D dipole trap,
then in a crossed dipole trap.
 Optimize vertical dipole trap (VDT) beam alignment. The VDT beam co-propagates with

50

(a)

(b)

Figure 2.26: PSF fit with the wavefront distortion expanded in Zernike polynomials. (a) Top-leftpanel: measured PSF on log scale; top-right panel: fitting function with the first eight terms in
Zernike polynomials; bottom panels: horizontal and vertical cuts through the center of the PSF onlog scale to visualize the wiggles. (b) The wavefront distortion (i.e. the phase component of theexit pupil) determined from the PSF fit. In principle, a compensating phase plate with the same
phase pattern but opposite amplitude can be used to completely suppress the aberrations.

767nm fiber
launcher(w/8mm lens)

1st polarizer(LPVIS100)

2x telescope λ/21064nm VDT and lattice
latticefocusing lens

caged system

51

band pass filter(FF01-775/140-25)

low pass
dichroic
EMCCD
170x170um field of view

2nd polarizer(LPVIS100)

objective

high passdichroic

side view top view

caged systemeyepiece

Figure 2.27: Optics layout for the top imaging setup on the science table. The distance between
objective and eyepiece is set to equal f1 + f2 .

52
the top probe beam in between the two dichroics therefore serves as a rough guidance.
Note that our VDT beam is not perfectly perpendicular to the cell top/bottom windows
for prevention of back-reflection.
 Collimate the top probe beam and overlap it on its back-reflection from the cell top window.

Use a printed target to center the probe beam onto the window.
 Iris down the optical pumping (OP) beam to be about less than 1 mm. This OP beam

travels through the same fiber as the top probe beam. We translate the OP beam while
imaging the atoms from the side to make sure that the OP beam is hitting on the atoms andtransferring them into the imaging spin state. The OP beam serves as a better guidance,
although its beam size is still larger than the field of view of our objective lens, which isestimated to be about 0.1 mm in radius.
 Use the top probe/OP beam to align the optics downstream (e.g. dichroic, mirrors, eye-piece, camera). We also use an auxiliary diode laser (Thorlabs compact laser module) to
propagate backwards from the camera as another guidance. This helps handle the extremely
long and snaky beam path.
 The objective is inserted in the end. Three back-reflected beam spots can be observed,

two small nearly-collimated ones and a large diverging one. They are from different optics
surfaces in the objective. Centering them on top of each other fixes the tilting angles ofthe objective. Now the objective should be parallel with the cell top/bottom windows.
 We use a home-made lens tube adapter (see the inset of Fig. 2.24) to connect an iris tothe back of the objective at a distance of around its focal length. We align the translationof the objective by closing down this aperture while observing the focused beam being
isotropically cut in from the edge. This was the same method that we used when aligning
the side imaging objective.
We also paid great attention to the polarization of the probe beam. On the output of the

53
caged fiber launcher, the extinction ratio of the polarization is better than 105 : 1 (set by thefilm polarizer LPVIS100 from Thorlabs, damage threshold 1W/cm2 continuous block, 5W/cm2
continuous pass). Birefringence on cell windows, dichroics and mirrors reduces the extinction ratio
to be > 103 : 1. This number should be sufficient for us to do non-demolition detections such as
Faraday rotation imaging [47, 48, 49].

2.2.5

Science Camera Calibration

Our science camera for the top imaging is a Princeton Instruments ProEM 512B. The quantum efficiency of the CCD is around 77 % at 767 nm. We use the shot noise of dark frames tocalibrate its noise level and ADC gain. As shown in Fig. 2.28 panel (a), we took two dark frames
at different exposure times and extracted the mean dark counts from them. Assuming that all the
pixels are equivalent, the standard deviation of the dark count per pixel is equal to the standard√
deviation of the two frames’ difference divided by 2. Each frame contains a statistically large
number of 512 × 512 pixels. In the short exposure time limit, the noise is dominated by the readnoise. For long exposure times, we can extract the gain of the CCD from the shot noise on theq
2 + δ 2 , anddark current. In Fig. 2.28 panel (b), we fit the total noise per pixel into δtotal = δcts
read√
δcts = GδNe = G Ne , where G is the ADC gain. From this fit, the read noise turns out to be
7.3(8) cts. (counts, or least significant bits), and the gain is 0.98(4) cts. per electron. From themean dark count in either frame, the bias reads 603 cts., and the dark current is estimated to be
33 cts. per pixel per second in the linear regime. The above measurements are taken at −40 ◦ C,
read out rate is 1 MHz and controller gain set to be 1.

2.3

Fast Magnetic Field Coils Setup
Our fast magnetic coils (named fast-B coils) are a pair of 11 turn coils made of 20 gauge

magnetic wires. Each coil is 0.9” in inner radius and 1.0035” in outer radius. The inductance percoil is about 20-30 µH. The magnetic field from the fast-B coils is calibrated to be 2.25 G/A. So weneed about 4.44 A current in order to generate a magnetic field change of 10 Gauss. This is enough

54(a) 10,000(b) 100
Model
Equation
Plot
A
B
Reduced Chi-Sqr
R-Square(COD)
Adj. R-Squaremean (count)standard deviation (count)frame 1frame 2

1,000

1E-6

1E-5

1E-4 0.001

0.01

0.1

1

exposure time (sec)

10

100

1,000

sqrt2_user (User)((A*sqrt(x))^2+B^2)^0.5
std_diff
0.99272 Â± 0.02255
7.33274 Â± 0.81754
6.63836
0.985120.98429

10

1

10

100

1,000

10,000mean - offset (counts)

Figure 2.28: Calibration of top imaging camera. (a) Mean total dark current in count per pixel vs.
camera exposure time. The total dark current consists of readout bias plus dark current. For eachexposure time we took two frames with camera shutter closed. (b) Standard deviation of the darkcount per pixel vs. mean dark count per pixel. The solid line shows a fitting function for the total
dark noise which is the quadrature sum of read noise and shot noise.

55
for us to jump the scattering length of the atoms from +200 a0 (at 38 G) to −200 a0 (at 28 G).
Note that the polarity of the fast-B coils is set to be opposite to that of the bias coils. The fast-B
coils are air-cooled and in close proximity to the the science cell, hence any heating of the coils
should be avoided. It is worth noting that due to the strong inductive couplings between the bias
coils and fast-B coils (they are very similar in diameters as can be seen from Fig. 2.14 panel (a)),
closing the loop of the fast-B coils can affect the bias feedback system if the proportional gain forthe latter is aggressive. One should make sure the noise level of the bias current stays the samewith and without the fast-B coils being present. For the bias field, we use Keithley multimeter
to measure the voltage drop across the sensing resistor on the out-of-loop current transducer. By
using fast mode on the multimeter, ppm level noise can be detected at a sampling rate of 10/PLC(power line cycle). For the fast-B field, we use a fast oscilloscope to calibrate the noise level on theout-of-loop Hall probe. The Hall probe is clamped onto the one of the current leads and placed 1
meter away from the coils. It is a quite sensitive device and prone to picking up noise from random
sources.

2.3.1

Fast-B Servo and Setup

The fast-B feedback servo setup is shown in Fig. 2.29 with model numbers of individualcomponents labeled on it. We employ the JILA coil drive designed by T. Brown, and a small
resistor 0.25 Ω to sense the current. A clamp-on fast current probe is used for monitoring thecurrent. For safety reason, the servo is first tuned up with a pair of test coils which are of the same
dimensions as the fast-B coils. A fast-blow fuse is installed to break the loop after 0.5 seconds when
the current exceeds 3 A. We also use the analog control of the fast-B coil power supply to gate its
output. For a control voltage of 3.33 V, the power supply can output 100 V. The response time ofthe power supply to the analog signal is about 80 ms.
For tuning up the coil driver, we start off with all PID trimpots fully off. Then the gain andPI frequency are increased proportionally until oscillations in the current show up. The lead spanand bandwidth are engaged in the end with the gain backed off and PI frequency increased when

56

Fast-B coils, ~27𝜇𝐻 per coil+ ~90V (300V max)
XHR 300-3.5-

MOSFET
IRF830servo out Coil Driver
QQ116A3
error in

0.25ΩVCS331

Hall probe monitor
TCP202A, 30MHz bw

control V(high Z)

Function Generator
33522A(use DAC for slow jump)

Figure 2.29: Schematic of the fast-B coils feedback servo setup. The inset photograph shows a test
coil that we wound on a lathe.

57
necessary. The above adjustments are done iteratively in small steps. We are able to obtain a 1090% rise time of about 3 µs on the fast-B current for a 10 Gauss field jump. The setpoint voltage isprovided by an arbitrary function generator in this test. The results are shown by the scope tracesin Fig. 2.30 panel (a)(b). We applied about 94 V to the coils so that the current changing rate isless limited by the coil inductance according to V = −L (dI/dt). Due to the non-linearity of the
MOSFET, a smaller jump distance requires a slightly different tuning of the servo. We also test
the fast-B performance with a DAC line as the setpoint input as shown in Fig. 2.30 panel (c)(d).
The DAC has a finite slew rate of about 1 V/µs. However, this does not seem to be limiting the
fast-B current rise time in our case.

2.3.2

Compensation for Induced and Eddy Current

When the fast-B field is being changed very rapidly, it induces current in our bias coils andnearby metals. The induced current in bias coils oscillates at the resonance frequency which isdetermined by the coil inductance and also the bias coil power supply voltage. As presented inFig. 2.31, the peak-peak amplitude of this oscillation amounts up to 1.5 G. And the damping time
constant is about 120 µs. Since we do not trap atoms by magnetic field, the bias coils are the only
ones that respond to the fast-B field change. The eddy current in our setup has a more minor effecton the total field. It curves the edge on the fast-B pulse and slows down its rise time from a fewµs to about 200 µs.

We tried some RC snubbers on the bias coils but could not suppress the induced current
oscillations very well. In the end as the fast-B field is pulsing on, we jump the bias current to ahigher value (corresponding to a larger magnetic field and a smaller scattering length) at about 100µs prior to the fast-B leading edge. As shown in Fig. 2.32 panel (a), this jump on bias coils generates

a current oscillation with an opposite phase and exactly the same amplitude as the induced current
from the fast-B pulse, therefore the two effects cancel each other. With the induced current on bias
coils being removed, we can visualize the eddy current effect as presented by the rf spectroscopy inFig. 2.32 panel (b) (black circles). This eddy current is fairly smooth therefore can be easily fixed

58

(a)

(b)

2% overhang, but only 4𝜇𝑠
𝑡𝑟𝑖𝑠𝑒 (10-90%) = 3𝜇𝑠
4.4𝐴 jump = 10𝐺

(c)

(d)

Figure 2.30: Tuning up the fast-B coils feedback servo. Blue traces: setpoint signal; yellow traces:clamp-on hall probe monitor. (a) The 10-90% rise time of the current is measured to be about 3µs for a 10 Gauss jump. We use 94 V on the power supply for the coils. We use extremely low
duty cycles (< 0.5%) on the function generator that gives the set voltage when tuning this servo
in order to minimize the heating. (b) Zooming in on the small overshoot of the current. Notethat the delay between the setpoint and the actual current can be removed by pre-turning on the
MOSFET gate with a small set voltage like -10 mV. (c) Use a DAC line as the set voltage insteadof the arbitrary function generator. The rise time of the fast-B current is minimally affected by
the finite slew rate of the FPGA DAC line. (d) The full fast-B pulse of about 1 ms duration. A
small pre- and post-pulse positive setpoint completely shuts off the fast-B current to prevent itfrom oscillating at such high power supply voltage.

59

(b)

Fast-B monitor

33.0

32.5

magnetic field (G)

(a)

32.0

31.5

AC current transducer for bias coils
31.0
0

10

20

30

40

50

60

70

80

90 100 110 120 130

thold at pulse center (ms)Figure 2.31: Induced current in bias coils due to a fast-B field jump of 8 G. (a) Induced current
measured on an out-of-loop AC wide-band current transducer (Pearson Electronics, Inc). (b) The
actual magnetic field calibrated with rf spectroscopy. The horizontal error bars represent the pulse
width of the microwave for a hyperfine transition. The solid line is a fit with damped sine function.
Note that the two different detection methods are consistent in the oscillation period and amplitude,
which indicates that the bias coils are the only source of this induced current since the mechanicalrelays on all other coils are open at this moment.

60
by feeding-forward to the setpoint on the fast-B field. The data in Fig. 2.32 panel (b) show one
such example in which we overshoot the fast-B setpoint by about 1 G and then ramp it back to
the target field by two 100 µs steps. Some iterations of this optimization process might be needed
if one desires a even more fine control of the scattering length especially at short hold times.

61

(a)

(b)
28.6

DC monitor for bias coils

magnetic field (G)AC current transducer
for bias coils

fast-B square pulse
eddy compensated

28.4

28.2

28.0

27.8

DC monitor for bias coils (induced current by fast-B)
27.6
10

100

1000

thold at pulse center (ms)Figure 2.32: Compensating the induced current and eddy current caused by the fast-B field jump.
(a) The current wiggles on bias coils induced by a fast-B rapid change (yellow trace) can be canceledout by jumping the bias coil current (gray trace) with an opposite phase and a properly-adjusted
amplitude. Note that the noise on the bias current is mainly from the picked-up signal along thelong BNC cable on this out-of-loop monitor. (b) Black circles: by canceling out the induced current
in bias coils, we can see the eddy current effect that manifests as the corner of the square pulse
being smoothed over. The solid line indicates a fit with exponential function. Green triangles:by “overshooting” the setpoint on the fast-B field by about 1 G and slowly ramp it back we canremove most of this eddy current.

Chapter 3

Efimov Physics: A Brief Review

Comprehensive discussions on Efimov physics can be found in review papers such as [50, 51,
52, 53]. In this thesis, I focus on some fundamental aspects of the research highlights in this field,
including both experimental and theoretical works.

3.1

From Two to Three

3.1.1

Two-Particle Scattering in Zero-Energy Limit

In ultra-cold scattering theories, a channel can be specified by the spin quantum numbers of
the two colliding atoms, e.g. α = {f1 mf1 f2 mf2 lml }, where l and ml specify the angular momentumof the collision 1 . A channel is referred to as “open” when the asymptotic channel energy Eα =Ef1 mf1 + Ef2 mf2 is lower than the system total energy E; and “closed” if otherwise. The scattering
potentials for different channels are shifted with respect to each other due to the Zeeman effect.
In the two-channel picture of a magnetic Fano-Feshbach resonance (FR), a pair of atoms in the
open channel (or entrance channel) are coupled to a shallow bound state in the closed channel
via virtual spin-exchange process. As shown in Fig. 3.1 panel (a), the energy difference between
the two channels can be conveniently tuned by an external magnetic field. As a result, the neardegeneracy between the scattering state and shallow bound state leads to a scattering resonance inthe s-wave scattering length as as shown in Fig. 3.1 panel (b). In the zero collision energy limit,
1

We closely follow the notations introduced in [54] in this section.

63
as can be parameterized asas (B) = abg 1 −∆B − B0,(3.1)

where ∆ represents the width of the FR in units of Gauss. The a > 0 and a < 0 regions correspondto effectively repulsive and attractive interactions between atoms. Scattering processes in general
are contributed by the background scattering (the first term in Eqn. 3.1, given by abg which inmost cases is on the order of rvdW ) plus the resonant part (the second term in Eqn. 3.1).

(a)

(b)
𝛿𝜇 ∙ 𝐵asclosed channel
open channel0𝑎>0
repulsive

EnergyΔ𝑎<0attractive

A

0

Atomic separation r
B

B-B0
Eb

𝑎𝑠 ≫ 𝑟𝑣𝑑𝑊

Figure 3.1: (a) Two-channel model of Feshbach resonances, and (b) the s-wave scattering length
as as a function of magnetic field B. The cartoon picture shows that the delocalized wave packetsof two atoms interact with each other at a distance of as which can be much longer than rvdW inthe vicinity of a Feshbach resonance.

The shallow bound state in dressed state picture becomes the Feshbach molecule and its
binding energy is Eb . The wavefunction of a Feshbach molecule can be expressed as a superpositionof the bare bound state and the background scattering state as|ψb i =√

Zφc |ci + χbg |bgi,(3.2)

where Z ∈ [0, 1] indicates the weight of the closed-channel component and depends on the magnetic
field detuning (B − B0 ). Both Z and Eb vanish in the limit of B → B0 .
One important parameter that characterizes the strength of the coupling between the open

64
and closed channels is sres . In a field theory formalism, sres can be defined assres =

ā
,
r∗

(3.3)

where ā ≈ 0.956 rvdW and r∗ is sometimes called the range of the FR [55]. For 39 K, rvdW =1
1
2 4 = 64.61 a . We can relate r ∗ to the coupling amplitude g between the atomic and0
2 mC6 /~
molecular field through g 2 = 32π/r∗ and determine its value by matching the experimentallymeasured profile as (B). Another equivalent approach to calculate sres is to use
sres =abg ∆ δµ
,
ā Ē

(3.4)

where δµ stands for the differential magnetic moment between the open and closed channel, andĒ = ~2 /(mā2 ) is the energy scale associated with the van der Waals force. Feshbach resonanceswith sres  1 ( 1) are termed entrance- (closed-) channel dominated resonances in that thenear-threshold scattering and bound states have the spin character of the open channel over a large(small) fraction of the width ∆. Open-channel dominated resonances have been widely studied
with 133 Cs and 85 Rb, etc. A single-channel approximation is sometimes sufficient in the case of133

Cs. While closed-channel dominated resonances are relatively trickier to access due to theirnarrow widths ∆. A multi-channel description is usually needed to model such resonances.
In 39 K, a number of FRs are of intermediate coupling strengths. The FRs on all three F = 1
spin states at low magnetic field are plotted in Fig. 3.2. We focus on the |F = 1, mF = −1ispin-polarized gas (called cc channel) in this thesis. The shallow bound states associated with thisscattering channel are shown in Fig. 3.3. For the FR at 33.5820 G, five different spin configurations
need to be involved in a full coupled-channel calculation. As shown in Fig. 3.4, the relative hyperfine
energies of these channels are plotted with respect to the mean energy value. The fraction of eachchannel in the Feshbach molecule state is also plotted in order to visualize their weights as a functionof magnetic field. All above calculations are done by P. S. Julienne, et al.
Although the strength parameter sres is originally introduced for a two-channel model, we
can define a multi-channel analogy of it for our FR in 39 K. As depicted in Fig. 3.6 and discussed

65
39

K

200

33.58Gs basis only

aa
bb
cc

a (a0)

100

0

-100

-200
0

50

100
B (Gauss)

150

200

Figure 3.2: Feshbach resonances in 39 K for the F = 1 spin-polarized gases, 0 − 200 G, s-wave basis;a → |F = 1, mF = +1i, b → |F = 1, mF = 0i, c → |F = 1, mF = −1i. The FR we focus on is theone at 33.5820 G in the cc channel (red lines). (Fig. 3.2 to 3.5 are all generated by P. S. Julienne.)0E/h (MHz)

-50
-100
-150

cc

-200
-250
-300
0

50

100
B (Gauss)

150

200

Figure 3.3: Two shallowest bound states in 39 K for the cc channel, 0 − 200 G, which are associatedwith the two FRs at 33.5820 G and 162.34 G respectively.

66
39

K mf=-2 channel energies

E=0 represents nonrelativistic mean (no hyperfine)3

E/h (GHz)

2
1
ee
df

0

ce
bd
cc

-1
-2
-3
0

200

400

600

800

B (G)

1000

Figure 3.4: Hyperfine energies of the five spin configurations in the multi-channel model for 39 K,
total mF = −2, 0 − 1000 G; b → |F = 1, mF = 0i, c → |F = 1, mF = −1i, d → |F = 2, mF = −2i,e → |F = 2, mF = −1i, f → |F = 2, mF = 0i.

1
1
2
3
4
5
Sum 2+3+4+5
ratio 2:3

Fraction

0.8
0.6
0.4
0.2
0
30

40

50

60

70

80

B (G)

Figure 3.5: Fraction of each channel in Feshbach molecule as a function of magnetic field, 30 − 80 G;
1 → cc, 2 → bd, 3 → ce, 4 → df, 5 → ee. The cc entrance channel plays a predominant role near
the pole of the resonance at 33.5820 G.

67
in [18], the dashed lines indicate the magnetic moments of two free atoms and the bare boundstate in their linear regimes respectively. The differential magnetic moment is calculated to beδµ/h ≈ 3.43 MHz/G. The value of sres is determined to be about 2.45 in this approach.
-0.5

-0.6

-0.7

-0.8

-0.9

-1

0

50

100

150

200

Figure 3.6: Atomic energy (solid black line) and coupled-channel calculation of bound state energies(red lines) near threshold for 39 K, 0 − 220 G. The figure is borrowed from [18]. The dashed lines
represent the atomic and molecular magnetic moments near and away from B0 respectively. Theirdifference is δµ/h ≈ 3.43 MHz/G.

3.1.2

Finite-Energy Effects and Effective Range

For incoming atoms of finite momenta ±k in the center of mass frame, the low-energy expansion of the scattering amplitude can be expressed in effective range approximationf (k) =

1.
−1/a + re k 2 /2 − ik

(3.5)

The effective range re is another useful parameter which accounts for the finite energy effects up
to order k 2 . For a two-channel model, it is a function of magnetic field


abg 2
re (B) = −2r∗ 1 −+ re(1) ,
a(B)

(3.6)

68
where


Γ( 14 )4ā  ā 2(1)
re (a) =ā 1 − 2 +
2
6π

a

(3.7)

a

is the single-channel approximation for large sres values (see supplemental materials of [18]). The fullfunction of re Eqn. 3.6 is plotted in Fig. 3.7 for various a. As a → ∞, re (B0 ) = −2ā/sres +2.91791ā.
For sres  1, re → 3ā; for sres  1, re → −2r∗ ; at sres ≈ 1, re crosses zero as a result of continuity.
Although the exact crossing point is model-dependent [56], this behavior implies a suppressed role

𝑟𝑒 /𝑎ത

of finite range in scattering and bound state problems.

𝑎→∞
10a
𝑎=
10𝑎ത
a𝑎 =
5a5𝑎
ത
a
a

𝑠𝑟𝑒𝑠

Figure 3.7: Effective range as a function of sres evaluated at three different values of a. abg is set
to 0 in this plot to represent species such as 39 K and 7 Li.

A physical interpretation of the change of sign in re can be found in E. P. Wigner’s old
concept of time delay in two-body scattering [57]). As depicted in Fig. 1 of [58], the s-wave phaseshift at finite wave vector leads to a spatial shift ∆r which is equal to re according to
δ0 (k) =π re k
−
+ ...
2
2

⇒

∆r = −2dδ= re + ...
dk

(3.8)

at a → ∞. The scattered wavepacket from far field appears as if it originates at ∆r rather thanr = 0. For open-channel dominated resonances, the situation is the same as in as single channel
potential with a van der Waals tail, where ∆r ' 3ā. In this case the particles can not get closer tothe origin (in the center-of-mass frame) than about 3 ā. The situation is completely different for

69
closed-channel dominated resonances, where the two colliding particles form an intermediate boundstate. The associated long time delay translates into a negative value of ∆r, i.e. the scattered wave
appears to come from negative relative distances.
A coupled-channel result on re as a function of B for 39 K cc channel is shown in Fig. 3.8 andcompared to a single-channel calculation. By matching Eqn. 3.6 to the coupled-channel model, we
arrive at sres = 2.68 which is consistent with the alternative approach given in 3.1.1.
200
180
160
140
120
100
80
31

32

33

34

35

36

Figure 3.8: A full coupled-channel calculation of re for 39 K (red line) near the pole B = B0 . The
figure is borrowed from [18]. The cyan line is single-channel formula Eqn. 3.7.

3.1.3

A Pedagogical Illustration of Three-Particle Scattering

For three particles located at ~x1 , ~x2 and ~x3 , eliminating the center of mass, the system canbe described in Jacobi coordinates:~rij = ~xj − ~xi ,~xi + ~xj
2
√
~xk −
ρ
~ij,k =
2
3

(3.9)

where (i, j, k) are chosen among (1, 2, 3). Three-body problems can be greatly simplified by doinga frame transformation into the so-called hyper-spherical coordinates:r = R sin α,

ρ = R cosα(3.10)

70
where R =p
r2 + ρ2 is the hyper-radius that sets the overall spatial extent of the system, andα = arctan(r/ρ) ∈ [−π/2, π/2] is the Delves hyper-angle. In this representation, the motions of
three or more particles can be quasi-separated into hyper-radial and hyper-angular wavefunctions.

𝑟Ԧ12

1

2

3
𝜌Ԧ
2 12,33𝜌Ԧ
2 13,2

1

1
𝑟Ԧ13

2

3

3

3
𝜌Ԧ
2 23,12
𝑟Ԧ23

3

Figure 3.9: Jacobi coordinates describing the relative positions of three identical particles.

In 1970, V. Efimov discovered an effective long-range attractive interaction in three-particle
systems with short-range inter-particle interaction near a scattering resonance [59, 60]. In thenear-resonance condition, the details of the interactions can be ignored as the scattering length
as sets the only length scale in two-particle scattering processes. Efimov proved that with twobody resonant interactions, the ground potential of three-body hyper-radial motion has an effective
attraction of the form
V0 (R) = −|s0 |2 + 1/4,
R2(3.11)

where s0 is a constant that depends on the mass ratio. For three identical bosons, s0 ≈ 1.00624.
Heuristically, the long-range interaction arises from two particles interacting with each other at along distance by effectively exchanging the third one [52]. The idea is shown by the cartoon in
Fig. 3.10 panel (a). This process becomes clearer in the Born–Oppenheimer picture (e.g. covalentbond between atoms) especially for mass-imbalanced systems with two heavy particles and one
light particle, in which the Efimov effect is in fact strengthened (i.e. featured by larger |s0 | values).
The effective long-range attraction is sometimes said to be of kinetic origin since it happens withthe aid of low-energy diffused wave packets as its media [61].

71

(a)
A

A

C
B

B

~𝑎𝑠

C
~𝑎𝑠

(b)

~−

0

1
𝑚𝑅2

Hyper radius R

Energy0

Energy

0

~−

0

1
𝑚𝑅2

Hyper radius R

Figure 3.10: Efimov attraction and scaling invariance. (a) Effective long-range interaction arises
in a three-particle system in which each pair of atoms are resonantly interacting. Particle A and
B can interact at a long distance by exchanging a “fast” particle C. Note that real particles areindistinguishable bosons in homonuclear case. (b) A sketch of this long-range attraction as afunction of hyper-radius R. Dragging the potential from the dashed black line to solid black line
accomplishes a scaling transformation on R. The bound states which are indicated by the coloredsolid lines remain in the same place after this transformation.

72
A consequence of this effective long-range attraction is scale invariance. A ∝ 1/R2 formof potential energy scales as the kinetic energy ∝ d2 /dR2 under scaling transformation R → λR.
As a result, the quantized bound state spectrum in a potential such as the one give by Eqn. 3.11
displays discrete scale invariance with a scaling factor of e2π/|s0 | . Namely, an infinite numberof three-body bound states form a geometric series of energies accumulating at the threshold at
unitarity. This scaling invariant spectrum associated with the effective long range potential isqualitatively explained in Fig. 3.10 panel (b). Only three of the infinite number of bound states
are drawn in the picture.
Remarkably, Efimov has shown that this discrete scale invariance can be extended to finite
scattering lengths. The bound state spectrum can be plotted in a two dimensional space of 1/a andκ, where κ is the wave number of the bound state. A sketch of this well-known Efimov’s scenario
can be found in Fig. 3.11 panel (a). On the a > 0 side, the trimer states co-exist with the Feshbachdimer state. On the a < 0 side, the trimers can stay bound even when the inter-particle potential
is not strong enough to support a two-body bound state. This is the so-called Borromean state.
For three identical bosons, the scaling factor between adjacent trimer states in this (1/a, κ) basis
is eπ/|s0 | ≈ 22.7. Efimov spectrum is famous for its self-similarity appearance. Familiar examplesof self-similarity (or fractality) include the Matryoshka dolls and Romanesco broccoli as shown in
Fig. 3.11 panel (b).

3.2

Universality in Efimov Physics

3.2.1

Observation through Loss

In his original work, Efimov considered nuclear systems such as triton (the nucleus of tritium
3H), and the Hoyle state of 12 C (three α particles). Yet these nuclear systems do not fully satisfy

the short-range interaction condition due to the presence of competing Coulomb repulsion in thesesystems. Neutral atoms, more than thirty years after Efimov proposed his ideas, have become anovel and ideal platform for the study of Efimov physics by virtue of an adjustable interaction

73

𝜅

(a)

a-0

22.7

22.7

Et1
Et0

1
𝑎

Et2
a*1

(b)

Figure 3.11: (a) Schematic plot of Efimov spectrum as function of 1/a in unit of wave number κ.
The area above the x-axis indicates three-free-atom continuum. The blue shaded area indicates
two-atom bound state plus a free atom. The gray shaded area indicates three-atom bound states.
Three trimer states are drawn in colored solid lines with an artificially reduced scaling factor.
They emerge at the three-atom threshold for a < 0 and connect to the atom-dimer threshold for
a > 0. An infinite number of trimer states that accumulated at the resonance pole are omitted.
A two-body potential in realistic systems also supports deeply bound dimer states whose binding
2
energies are on the order of ~2 /(mrvdW) and beyond. Those dimers are indicated conceptually as
the dashed line near the bottom of the plot. (b) Nesting dolls and broccoli, two other examples
that also present self-similar feature. Pictures are clipped from Wikipedia.

74
strength in the vicinity of a Feshbach resonance.
In contrast to 4 He for which the interaction potential only hosts one dimer state, potentials
for alkali metals are deep enough to support a large number of ro-vibrational states whose binding
2
energies are on the order of van der Waals energy ~2 /(mrvdW) and beyond. An Efimov trimer inalkali metals is more of a quasi-bound state (i.e. resonances embedded in the continua formed by
one atom and one diatomic molecule) and can dissociate into a deeply bound dimer plus a free atom.
This process leads to a very limited lifetime of the Efimov trimers and disfavors direct observationsof the trimers by spectroscopic methods. Due to this reason, Efimov trimers are usually detected
through enhanced losses of the atoms (e.g. three-body recombination). As encircled by the dottedlines in Fig. 3.11 (a), a trimer state intersects the free-atom threshold (a < 0) and atom-dimer
threshold (a > 0) at certain a’s. The trimers at these locations can greatly enhance the losses byacting as intermediate decay products [62] and manifest as resonance-like peaks in the loss rate.
(0)

One such feature, the tri-atomic resonance at a = a− , occurs where the ground trimer statecrosses the free atom threshold. It is first measured in 133 Cs and later has been repeated withother atomic species [63]. Although it is quite difficult to access higher order trimer branches inorder to verify the scaling invariance aspect of the Efimov spectrum, it is still convincing that thisloss enhancement is linked to Efimov physics because of the Borromean nature of the trimer state.
(0)

Physically a− sets the smallest size of the Efimov series for a < 0. This parameter in fact isa missing piece of information in the Efimov physics based on zero-range theories. As has beenshown by L. H. Thomas, a zero-range potential can lead to three-body bound states with arbitrarily
large binding energies therefore a lower-bound for the interaction range has to be introduced. The(0)experimentally determined a− provides such a lower-bound. It had been long anticipated to be
unknown a priori due to its close relationship with the details of the inter-particle forces at shortdistances [64].

75
3.2.2

Van der Waals Universality

Three-body observables such as the above mentioned loss enhancements are associated withthe short-range boundary condition that regularizes the zero-range theory of three particles. It issometimes termed three-body parameter or “3BP” as a shorthand. This boundary condition in
general depends on the detailed form of the short-range potential, in analogy to the two-body case.
For instance we know that for an attractive potential with a deep well, the scattering wavefunction accumulates a large phase inside the well and a minute change in the potential can modify
the scattering length drastically. The difficulty in predicting the 3BP for three-body attractive
potentials can be understood in an analogous way. Generic three-body potentials for neutral atoms
can be written in terms of a pairwise sum of van der Waals interaction, plus genuine non-additive
three-body forces such as Axilrod-Teller interaction [53].
The discovery of van der Waals universality shed the first light on the 3BP problem. Ex(0)perimental groups measured a constant value of a− in 133 Cs after sampling a series of FRs thatare located at from a few Gauss to over 800 Gauss [12]. This experimental result suggests that thethree-body short-range physics is not substantially altered by the absolute value of the magnetic(0)

field. Extensive experimental evidence later on confirms that a− is often about 9 times rvdW regardless of the atomic species being used (see Fig. 21 in [52] for a summary plot). The microscopic
mechanism of this universal location is revealed by a numerical calculation done in hyper-sphericalcoordinates. Wang et al found that the ground hyper-radial potential exhibits a steep repulsive
barrier at a hyper-radius R ≈ 2 rvdW instead of a strong attraction as one would naively expect fromvan der Waals potentials [13]. This repulsive barrier prevents the hyper-radial wavefunction from
penetrating into the short-range part of the hyper-spherical potential where the situation becomes
more species specific. They also proved that this effective three-body repulsion arises from theacceleration of the two-body relative motion due to the attractive van der Waals potential. This
acceleration makes the two atoms spend little time at short separations and equivalently leads to a
suppression of the two-body probability. As a consequence, a deformation of the three-body system

76
takes place as it changes from its configuration at large separations to that at small hyper-radii.
At larger hyper-radius, the most probable configurations of hyper-angular wavefunction have an
elongated-triangle geometry. While at smaller hyper-radii on the order of rvdW , the suppressionof two-body probability forces the system to adopt an equilateral configuration. This deformationof the three-body system gives rise to an increase of the kinetic energy cost and manifests as arepulsive barrier in the ground hyper-radial potential. An illustration of this change in three-bodygeometry can be found in Fig. 22 of [52]. Briefly, the van der Waals universality in 3BP has its
origin in the two-body van der Waals universality in the r > rvdW region, because the two-bodyradial probability density is suppressed in the r < rvdW region as if there were a repulsive core in
the potential.
The above conceptual explanation is based on a single-channel van der Waals potential.
This is a good approximation near Feshbach resonances characterized by sres  1. One would
expect deviations from this theory when multiple channels are coupled by hyperfine interaction
and participate in the scattering processes together. 39 K is such a species with an intermediatesres value. Indeed, in our lab we made the first definitive observation of a deviation from vdW
universality among homonuclear systems [18]. A compilation of all experimental values are presented in Fig. 3.12. The gray band indicates the prediction of single-channel vdW universal value(0)

a− /rvdW = −9.73±15%. The uncertainty is determined by the variance of details of the potentials
used in individual theories. For the FR at 33.5820 G in 39 K, we experimentally nailed down the(0)value of a− to be −14.05(17)rvdW with unprecedented precision. Our result has been reproduced
by a coupled-channel three-body model done by J. P. D’Incao with relatively good agreement. In
this model, we incorporate realistic spinor physics that we learned from our Feshbach resonancecalibration project.

3.2.3

Zero-Range Universality

Single-channel van der Waals universality is unique to the interaction potential of a 1/r6
form between two neutral atoms. In the original, more general version of Efimov physics, discrete

77

(a)

(b)
18
16

15

our work

14

10

12
10

5

8

1

2

3

4

5

6

7

6
0.1

1

10

100

1000

(0)

Figure 3.12: (a) A summary of experimentally determined values of a− in homonuclear systems
as function of sres . 7 Li: Rice University, Bar-Ilan University, ENS; 39 K: LENS, Aarhus University,
JILA; 85 Rb: JILA; 133 Cs: University of Innsbruck. The red star indicates our reported value of(0)(0)
a− for 39 K. The open blue circles indicate previous results on the same a− . The dashed line andgray band represent the range of van der Waals universality prediction. (b) Numerical calculations(0)of a− based on a single van der Waals potential (green squares) and a multi-channel model (black
triangles). N represents the number of s-wave diatomic bound states in the model. The red band
indicates our experimental value. Both figures are cited from [18].

78
scaling invariance and a universal shape of the trimer levels in the (1/a, E) plane are the twocentral results of a zero-range theory. In this zero-range universality framework, single-channel
vdW universality plays the role of finite-range corrections to the low-lying trimers on the Efimov
spectrum. In other words, the benefit of knowing the lower bound of the spectrum comes at the(n)cost of severe distortions to the ground and first excited trimer energies Et(n = 0, 1). This

also shifts the corresponding locations of the scattering resonances where these trimers meet the
(0)

(1)

(1)

(0)

scattering continua (e.g. a− and a∗ ), as well as the ratios between these locations (e.g. a∗ /a− ).
It is worth noting that the Feshbach molecule state is also modified by the finite range effect in the(n)small a region. The atom-dimer resonances a∗

are then related to the distortions of both dimer

and trimer energy levels on the a > 0 side. In particular, the ground trimer state actually never
meets the atom-dimer scattering threshold [65].
The above descriptions on vdW and zero-range universality are illustrated in Fig. 3.13. Twofrequently measured observables, three-body recombination rate coefficient L3 and atom-dimer
relaxation rate coefficient βAD are plotted as a function of rvdW /a in this figure. The scattering
(n)

(n)

resonances a− , a∗

(n)

(n)

as well as the interference minima a+ (and maxima ap(n)≡ a+ eπ/(2s0 ) , witheπ/(2s0 ) being one half of the Efimov period [66]) are associated with corresponding crossing points
on the trimer spectrum through the green dotted lines. The portions of the trimer energies andthree-body observables plotted in solid black lines indicate predictions based on zero-range theory.
The relative locations of these trimer branches are fixed by zero-range universality. The dashedblack lines indicate the portions which are subject to modifications due to finite-range effects inreal systems. The locations of these low-lying resonances relative to the more excited ones dependon the details of the interaction potentials.
The basic idea of Fig. 3.13 is the following. In real atomic systems, the experimental observation of a complete Efimov spectrum is thwarted by two factors. One is the finite-range effects
which distort the Efimov energy levels at small a’s. These regions are indicated by the faint yellow
squares whose upper bounds are approximately set by ±re multiplied by a prefactor. The otherfactor comes from the systematic errors when a is large compared to the de Broglie wavelength λ

79

(0)
𝑎−
/𝑟vdW = absolute scale
(𝑛+1)
(𝑛)
𝑎∗
/𝑎−
= −1.065(𝑛)

𝑎𝑝 /𝑎−𝑛 = −1(𝑛+1)𝑎−,∗,+,𝑝(𝑛)𝑎−,∗,+,𝑝= 22.7

re

T,n
re

(𝑛)

(2)

(1)

𝑎∗

𝑎∗

(1)
𝑎−

−1

𝑎<0
(0)
𝑎−

𝑎 = ±∞

(1)

𝑎𝑝

𝑎>0

(1)
𝑎−

𝑎+
(2)

𝐸𝑡

(0)

𝑎+

(0)

(0)

𝑎+

(2)
𝑎∗

/𝑎4

𝑎𝑝

(𝑎>0)

+ +

(0)

(1)

𝑎+

𝐿3

𝐿3

(𝑎<0)

/𝑎4

𝛽AD /𝑎

(0)
𝑎−

(𝑛)

𝑎𝑝 /𝑎+ = 22.7

T,n

+ +

1

𝑟vdW
𝑎

(1)

(1)

𝑎∗

𝐸𝑡

+

(0)

𝐸𝑡

𝐸𝑑

(n)

Figure 3.13: Efimov universality. The bottom panel shows trimer (Et ) and dimer (Ed ) energies
as a function of a−1 , with energy levels distorted to make clearer the resonance locations. The nth(n)(n)
trimer is resonant with three free atoms at a = a− , with the dimer at a = a∗ . The middle panel
depicts the rate constants for inelastic collisions. The three-body recombination constant L3 peaks(n)(n)
at each value of a− ; while the atom-dimer relaxation constant βAD peaks at each value of a∗ .
The entire structure is log-periodic with period 22.7. For a > 0, L3 /a4 shows quantum interference,
√
(n)
(n)
with each local maximum ap spaced from the corresponding local minimum a+ by 22.7. Each(n)(n)
value of ap is related to the corresponding a− by a factor of exactly −1, but is offset from the(n)nearest value of a∗ by 6.5%. Going from relative to absolute values of a requires a single absolute
scale indicated by the oval circle. All these ratios are for the zero-range limit; the dashed curves
suggest the possibility of perturbations as |a| enters the region, indicated by yellow shading, where
it is not large compared to the effective range re . The range of the yellow zone may be “adjusted”
by choosing atomic species with different re . The purple shading represents large |a| regions thatare prone to systematic effects such as those caused by finite temperature and density.

80or the inter-particle distance n−1/3 of the sample. This region is indicated by the purple square
whose lower bounds are set by the temperature and number density of the cloud 2 . The above
two factors leave us very limited windows for an observation of Efimov’s original predictions basedon zero-range interaction. These windows are indicated by the white gaps on both sides of the FRpole.
One signature of the zero-range universality is the universal ratios between three-body observables for a < 0 (i.e. a− ) and for a > 0 (i.e. a∗ , a+ and ap ). We collect all experimental resultson these ratios in homonulear systems as summarized in Fig. 3.14 and compare them to a zero-range
theory [50] and a two-channel theory that incorporates the FR coupling strength sres [15]. One of
the conclusions from the latter theory is the dependence of the ratios on sres , which can be traced
to the close relationship between sres and the effective range re (re is in general a function of a andsres as discussed in Section 3.1.2). As pointed out by the two-channel theory, a suppressed re isassociated with an enlarged universal window, meaning that the ratios obey the zero-range results
even on the low-lying trimer state with n = 0; whereas for open-channel dominated resonanceswith relatively large re , these ratios slowly approach its zero-range values as n increases from n = 0
to n > 1 (see Table 1 in [15]). This predicted trend is evidenced by the experimental values of
(1)

(0)

a∗ /a− in panel (a) of Fig. 3.14. In the zero-range limit, the three-body recombination maxima
(n)

ap

(n)

turns out to be about equal to a∗within a few percent [50, 53]. Therefore we also add thisratio to this summary plot.
(0)

(0)

In panel (b) of Fig. 3.14, we compile results on another ratio a+ /a− . The scattering(0)interference a+ is not within the scope of the two-channel predictions in [15], so we only compare(0)

the experimental values to the zero-range theory. The observable a+ lies deeply inside the small aregion for all available species. Therefore it is not surprising to see no particular dependence on re
(0)

(0)

for the ratio a+ /a− . This ratio hence does not serve as a promising candidate for the realizationof zero-range universality.
2

These technical issues can be conquered by generating extremely cold samples with low densities for instancein a microgravity environment. We collaborate with the Cold Atom Laboratory (CAL) at NASA’s Jet Propulsion
Laboratory (JPL) and a group at Washinton State University (WSU) on this project

81

(a)
1.2

1.0
39

K, JILA| ratio |0.8

7Li, Bar-Ilan

0.6

133

Cs, Innsbruck(experiment)

0.4(experiment)(zero-range)(zero-range)(Zwerger&Schmidt)(Zwerger&Schmidt)

0.2

(b)

| ratio |0.4

39K, LENS

7

Li, Rice

39

K, JILA

0.2

133

Cs, Innsbruck(zero-range)

0.0

-1

0

1

2

3

Figure 3.14: Summary of experimental and theoretical results for Efimov’s ratios between threebody observables for a < 0 and a > 0 as function of re . re is evaluated at unitarity a → ∞ using
the model given in the supplemental materials of [18]. Definitions of the three-body observables areillustrated in Fig. 3.13. The zero-range results are predicted by [50, 53]. Zwerger&Schmidt refers
to a two-channel model that incorporates the FR coupling strength and provided by [67]. Notethat the atom-dimer threshold takes the universal functional form in this theory. The accuracy ofthis model could be readily improved by using more realistic values of the dimer binding energy
[67].

82

3.3

An Outline of Our WorkIn our previous work [18], we studied Efimov physics in the attractive interaction regime, orthe so-called Borromean state as defined in Section 3.1.3. The main results are the following: We calibrate a FR in 39 K with unprecedented precision by doing Feshbach molecule disso-

ciation spectroscopy; this helps refine the two-body coupled-channel model in the s-wave
sector.
 We precisely locate an Efimov tri-atomic resonance for a < 0, which is a first definitive

deviation from the single-channel vdW universal value in homonuclear systems.
 We develop a three-body coupled-channel model based on realistic two-body vdW poten-

tials, and reproduce the measured location and width of the tri-atomic resonance withrelatively good agreement.
In this thesis, we extend our precise calibrations of few-body physics to the a > 0 region.
The main results are the following: By mapping out the lifetimes of Feshbach molecules as function of magnetic field, wecalibrate the d-wave sector of the two-body coupled-channel model.
 We discover an atom-dimer resonance in 39 K and associated finite temperature effects.
 We locate the loss rate interference minimum and maximum for a > 0.
 We observe signatures of zero-range universality in 39 K by relating the three-body observ-

ables for a < 0 and a > 0.

Chapter 4

Experiments on Feshbach Molecules

In contrast to the attractive interaction regime, Efimov physics in the repulsive interaction
regime is relatively less explored and not as well understood. To achieve a global picture of Efimov
physics in 39 K, we perform a series of two- and three-body collision experiments in the a > 0 region.
In this chapter, we focus on the experiments on Feshbach molecules (also known as the shallow
dimers, or sometimes dimers for short).

4.1

Feshbach Molecule One-Body Lifetimes
Our science state |F = 1, mF = −1i is not the ground spin state in the Zeeman manifold,meaning that atoms or molecules in the |1, −1i + |1, −1i configuration (or cc channel) can spontaneously relax into lower Zeeman states accompanied by an energy release whose magnitude isa function of magnetic field. A Feshbach molecule is a quasibound state in nature in our case.
As depicted in Fig. 4.1, for the incoming cc channel, the decay channels are bc, ac and bb spin
configurations. And the energy differences between the exit and entrance channels are on the orderof tens of MHz at a magnetic field around 34 G.

4.1.1

Measurements of Feshbach Molecule Lifetimes

The experimental timeline for molecule lifetime measurements is shown in Fig. 4.2 andFig. 4.3. We use magneto-association to produce the shallow dimers in a standard way. We
scan the magnetic field to a negative scattering length and then hold it there for 10 ms to damp

84

39

K+

E/h (MHz)

50

0

39

K spin channels

cc,-2

ac,0

bc,-1

ab,1
bb,0

-50
aa,2

-100
0

50

100

150

B (G)

Figure 4.1: Zeeman energies of the d-wave decay channels (bc, ac, and bb) for the incoming cc
channel, 0 − 150 G. a → |F = 1, mF = +1i, b → |F = 1, mF = 0i, c → |F = 1, mF = −1i.

85
out the cloud’s breathing motion which is excited by the rapid magnetic field change. We thenramp the field back to the positive side at about 500 a0 as our home field. At this field, we adjust(or eliminate) the number density of the residual atoms. We apply a Gaussian-enveloped adiabatic
rapid passage (ARP) pulse followed by an optical blasting pulse to eject the atoms either thoroughly or partially. For our measurement of molecular one-body decay, we remove over 99.9% of
the atoms. We optimize this purification process such that it only takes about 0.5 ms and is highly
reproducible from day to day. The binding energy of the molecules at this home field is about 500
kHz, so that the molecules are unaffected by the ∆ = ±150 kHz span of the ARP pulse which flips
the free atoms into |2, −2i state with the molecules intact.
To produce molecule samples with different initial conditions (e.g. temperature and number
density), we adjust the final stages in the evaporation process prior to the association step. Afterall these preparations, we jump (within 0.5 ms) the magnetic field to a target science field andlet the dimer clouds evolve. At the end of the interrogation, we jump the field back to 500 a0 for
detection. An ideal detection scheme for the molecules would be in situ absorption imaging inwhich the spatial distribution of molecules can be preserved. However since we deal with molecule
samples of extremely low densities in this spontaneous decay measurement, we choose an alternative
method which is to convert them into free atoms by sweeping the magnetic field back to the negative
scattering length side so that both parent atoms can be imaged. This method enhances the signal
and applies to experiment in which the spatial distribution information of the molecules are notcritical.
We map out the molecule population as a function of interrogation time at various magnetic
fields. The population decay data is fit into an exponential function to extract the 1/e lifetime of
the molecules. In survey of the large a regions (at a couple of thousand a0 ), we only fit the tailsof the data so that the initial decay which is potentially contaminated by the collision-induced loss
processes can be excluded. We can comfortably work with low atom numbers from 500 up to 3000
(peak OD from 0.05 up to 0.3) by virtue of the low background noise in our imaging system. The
molecule lifetime τbg as function of field is presented in Fig. 4.4 panel (a) together with various

86

(a)

(b)
imaging light|2, −2ۧmolecules
blasting light|1, −1ۧ①. magneto-association of③A. Dissociate molecules w/
microwave to imaging state|2, −2ۧ

Feshbach moleculesor|1, −1ۧ−𝜇 ∙ 𝐵|1, −1ۧ

②. Adjust/eliminate
atomic density

OD

atoms
100 µmimaging light|1, −1ۧ

optical pumping③B. Remove atoms, thendissociate molecules w/ B fieldsweep to free atoms for imaging

Figure 4.2: Magneto association of Feshbach molecules and atom purification processes (a) and anabsorption image of the molecule cloud which is spatially separated from the free atom cloud after
the atoms absorbing the blasting photons (b). The step inside the dashed square is skipped for
preparation of a full atom number in the atom-dimer collision experiment. The ARP pulse used inthis step is made fast and diabatic in order to partially remove the atoms.

87−𝜇 ∙ 𝐵|1, −1ۧ

|1,0ۧ

|1,1ۧ

atom-dimer relaxation

one-body decay

Vary hold at
science field

Bias before science(preparation)

at 500𝑎0 for |1, −1ۧ
population control
at 500𝑎0 for |1,0ۧ
cleaning rf pulses
10ms

Hold 10ms at− 100𝑎0 to damp
out the breathing

Optimizemolecule
production

A. in situ
detection ofmolecules

B. Convert
molecules into
atoms for detection

Figure 4.3: Timeline of Feshbach molecule lifetime experiments, B vs. time. The double-headed
arrows indicate stages that are to be adjusted. The ramps to and from the science field are donein a diabatic way in order to reduce the loss of molecules during the magnetic field sweeps. The
blue arrows indicate that at science field two major collision processes can happen. One is molecule
one-body decay which is introduced in this section. The other is atom-dimer relaxation which we
discuss in 4.2.1.

88
theoretical models.

4.1.2

A Coupled-Channel Model for the Lifetime of the Quasibound StateFor incoming scattering channel with orbital angular momentum L = 0 for the relative motionof the two atoms, the bound state responsible for the resonance is also of predominantly L = 0
character. However, there is weak coupling from L = 0 to the channels with L = 2. This coupling is
entirely responsible for the inelastic processes involving two free atoms or two atoms in the shallow
bound state. The coupling operator isV̂ d (R) = λ(R)[ŝ1 · ŝ2 − 3(ŝ1 · ~eR )(ŝ2 · ~eR )],(4.1)

where ŝ1 and ŝ2 are the electron spin operators for the two atoms, and ~eR is a unit vector along
the internuclear vector of length R, and the amplitude has the functional formλ(R) = Eh α

2


ASO e

−γR/a0
gS2 a30
+
,
4R3

(4.2)

where α is the fine structure constant, gS is the g-factor of electron. The second term in theabove equation is the magnetic dipolar interaction between the two electronic spins. The firstterm represents the second-order spin–orbit interaction [68], which has been previously ignored in
calculations on light species such as K, Na and Li 1 .
In Fig. 4.4 the measured molecule lifetime is compared to a coupled-channel calculation
converged with an automated procedure by M. D. Frye and J. M. Hutson. The calculation (dashedline) which neglects second-order spin–orbit coupling shows qualitatively the correct behavior, witha strong peak in the lifetime close to the pole of the resonance and a weaker peak around 46 G.
The peak close to the resonance occurs because in this region the quasibound state is dominatedby the component in the incoming channel, whose spatial extent scales as a. The subsidiary peakoccurs because of interference between the decay amplitudes from the closed-channel and incoming
components.
1

The calculations in this section are to be published in [69] and contributed by M. D. Frye and J. M. Hutson.

89

5 0 0

2 0 0

1 0 0

6 0

4 0

2 5

1 .5

(a )

1 0 2

2 4

2 3

2 2

6 3

6 4

6 5

(b )

1 .0

0 .5

1 0 1

3 5

4 0

4 5

5 0

5 5

0 .0

6 2

Figure 4.4: Refinement of the d-wave sector of the two-body coupled-channel (cc) model. (a)Lifetime of Feshbach molecules (shallow dimers) measured after the free atoms being removed
from the trap. Two different dimer conditions hnD i = 2 × 109 cm−3 , T = 70 nK (red circles) and
hnD i = 6 × 1010 cm−3 , T = 280 nK (green diamonds) are presented. Solid line: cc-model that
includes spin-spin dipolar interaction plus second-order spin–orbit coupling of optimal magnitude;
dashed line: cc-model that includes only spin-spin dipolar interaction; dotted line: model from[70]. (b) Comparison of the experimental free-atom spin-relaxation rate constant (teal diamonds)
with cc-model calculations. The gray region indicates the experimental systematic uncertainty on
atomic density. Solid and dashed curves are produced with the same sets of parameters as being
used for the corresponding curves in (a).

90
The position of the subsidiary peak depends on the ratio of the decay amplitudes from theclosed-channel and incoming/open components of the quasibound state. A short-range term in the
coupling ASO is introduced to vary this ratio. The results are not very sensitive to the exponentγ, so we fix it at the value 0.7196 obtained for Rb [71]. The solid black line in Fig. 4.4 shows the
lifetime with the value of ASO fitted to the experimental lifetimes. The resulting value, ASO = −0.8,
gives a second-order spin–orbit interaction of −0.015 cm−1 at the zero-energy inner turning pointof the triplet state of K2 .
Köhler et al. [70] developed a model for the lifetime of the quasibound state near a resonanceand obtained the expression
τ=

4πa3
,
L2

(4.3)

where L2 is the zero-energy rate constant for two-body inelastic collision from the incoming channel.
Physically speaking, the relaxation of Feshbach molecule in the B → B0 limit can be viewed as
half-collision process in a pre-assigned atomic configuration (in particular encoded in the volume
4πa3 ). As shown by the dotted line in Fig. 4.4, this model is a good approximation close to thepole of the resonance. The fact that this mode fails to reproduce the height and location of thesubsidiary peak is not surprising since the incoming closed-channel components play a crucial rolein τ but not in L2 .
The second-order spin–orbit interaction partially cancels the magnetic dipolar coupling. For
the optimized value of ASO , λ(R) crosses zero inside the classically forbidden region for the triplet
state. In 39 K the overall coupling is thus dominated by the magnetic dipolar term, and the effect of
the second-order spin–orbit term is to reduce slightly the total decay rates due to spin relaxation.
Fig. 4.4 panel (b) shows the rate coefficient for two-body inelastic loss across the narrow d-wave
resonance near 63.609(4) G with full-width-at-half-maximum of 0.567(12) G. There is a small
discrepancy in the calculated position of the resonance, due to remaining deficiencies in the singlet
and/or triplet potential curves. Without ASO the prediction overestimates the height of the peak(dashed line), whereas the prediction including ASO agrees with the height and width within the91experimental uncertainties. This provides an independent sanity check on the accuracy of theoverall coupling function λ(R).
As illustrated in Fig. 4.5, this narrow d-wave resonance at 63.6 G is associated with a ramping
d-wave bound state that lies at higher magnetic field. This d-wave resonance is extremely weak, asevidenced by the barely noticeable change in the scattering length going across the pole. As shownin Fig. 4.6, the singularity is almost totally extinguished by the strong resonance decay.

4.2

Finite Temperature Effects in Atom-Dimer Reactions
A through understanding of two-body physics in 39 K provides a foundation of the experimentson three-body physics. In this section I discuss in particular how Feshbach molecules interact withfree atoms, which is essentially a three-body process. Two atom-dimer resonances have beenpreviously located in 133 Cs [72, 73]. The Cs2 molecules are relatively long-lived in their ground
Zeeman sublevel. The locations of the atom-dimer resonances are found to agree with a numericalmodel based on van der Waals interaction [73]. A lack of discernible temperature dependency on
the resonance peaks in these works contradicts the predictions by [74, 75].

4.2.1

Discovery of an Atom-Dimer ResonanceIn the atom-dimer reaction experiment, we measure the number and size of the dimers in situ
as a function of interrogation time as shown in Fig. 4.7. We apply a Gaussian-enveloped rf pulse todissociate the dimers and spin-flip one of the parent atoms into the imaging state. The frequencyof this rf pulse is chosen to be at the empirical peak response of the molecule dissociation spectrum.
This rf pulse is of full power so that it can deplete the molecule population in a short amount oftime about 120 µs. We then take absorption images on |2, −2i spin state at 500 a0 field with thequantization axis provided by the bias coils. Only one parent atom from each molecule is imaged inthis method. We repeat the measurements with various translational temperatures of the molecules.
For the lowest temperature 60 nK condition, we only have enough imaging signal by detecting bothparent atoms. This is realized in the alternative detection scheme in which the atoms are removed

92

39

K+

39

K

cc channel

E/h (MHz) (relative to Ecc=0)

0

-100

-200
s basis63.58(4)G d-wave resoance
Atomic open channel E
d -2
d -1
d0

-300

-400
0

20

40

60

80 100 120 140 160 180 200
B (G)Figure 4.5: Bound state calculations for 39 K in sd-wave basis, 0 − 200 G, plotted relative to the
cc channel energy as zero. The three bold black lines are the open channel thresholds (bc, ac, bb
from top to bottom). Above such thresholds there are the quasi-continua, and below are the truebound states. Bound states are shown for (L = 0, Ms = −2) and for L = 2 with Ms = −2, -1,
and 0. A bound state in the (L = 2, Ms = −1) channel that lies below the bc threshold and above
the bb and ac thresholds (represented by the green diamonds) projects to near where the 63.6 G
resonance is observed (indicated by the green arrow). Figure is generated by P. S. Julienne.

93

39

K2 cc channelRe(a) (a0)2423Im(a) (a0)

0
-0.1
-0.2
63

63.2

63.4

63.6
63.8
B (Gauss)

64

64.2

Figure 4.6: Real and imaginary parts of a in vicinity of the narrow d-wave resonance.

94
at end of the hold time and molecules are swept into free atoms by a magneto dissociation process(Fig. 4.2 step 3B). The temperature information therefore may not be accurately maintained in thiscase. The width of the molecules along the tightest confinement direction is about 5 µm, which islarge compared to our imaging resolution of 1.5 µm. Both molecule width and imaging resolution
are given in root-mean-square radius.
In the previous section, I describe how we measure the intrinsic lifetime of the dimers τbgwith atoms completely ejected from the trap. By skipping the atom population control step, we
let the molecules interact with atoms in addition to the molecules decaying by themselves. The
total lifetime of molecules is denoted as τtot . Since the atom cloud causing the inelastic decay
is 10 to 30 times denser than the dimers, the reaction processes mainly happen between atoms
and dimers rather than within dimer themselves. We again fit the dimer population ND with an
exponential form and extract its 1/e lifetime as a function of scattering length a. We then extract
the atom-dimer reaction coefficient from the rate equation:
dND
ND
ND
≡−
= −βAD hnA iND −
,
dt
τtot
τbg

1  −1−1,βAD =τtot − τbg
hnA i

(4.4)


3/2
where hnA i = NA mω̄ 2 / (4πkB TA )is the mean density of atoms in a harmonic trap, where
NA is the atom number, m is the atomic mass, ω̄ is the geometric mean trap frequency and TA isthe atom temperature. It is worth noting that in our optical trap, the polarizability of dimers is
twice that of free atoms. The trap frequencies of the atoms and the dimers are hence the same.
And because of the unequal mass, the density distribution of atoms and dimers are different. The
spatial overlap between atoms and dimers as a result leads to an overestimation of the relaxation
harmonic in a harmonic trap in the above equation by a factor of √8 ≈ 1.5 comparedcoefficient βAD27homogeneous
to βADfor homogeneous gases,8
homogeneous
harmonicβAD= √ βAD
.
27

(4.5)

As can be told from Fig. 4.8 panel (a), the molecule lifetime gets significantly reduced in the

95

N u m b e r

2 x 1 0 3
1 x 1 0 3(a )0

T e m p e ra tu re (n K )

3 0 0

(b )

2 0 0
1 0 0
0
1

2

H o ld tim e ( m s )3

4

Figure 4.7: Time evolution of dimers measured at 1200 a0 with atoms being the majority in thesample. (a) Dimer number decay (blue squares) fitted with an exponential function (solid line)with a 1/e lifetime τtot . (b) Temperature of dimers extracted from the horizontal (x̂, red circles)
and vertical (ŷ, gray diamonds) widths of the cloud from in situ images of the molecules. The
heating rate on the order of a few nK/s is minimal during the short hold time.

96
presence of free atoms. By following Eqn. 4.4 we subtract out the one-body decay rate from thedimer total decay rate and obtain the atom-dimer relaxation coefficient βAD at different temperatures as shown in Fig. 4.8 panel (b). A resonance-like peak at around 900 a0 is quite pronounced
for all temperatures. With the highest temperature condition, we vary the atom density by ejecting
two thirds of the atoms from the trap to confirm that the atom density in this experiment acts asa scaling factor that is divided out to get βAD .
Throughout this thesis, the magnetic field is calibrated with rf spectroscopy. We use the most
up-to-date Feshbach resonance calibration given by [18] for the conversion from magnetic field toa. The uncertainty on scattering length amounts to about 5% at around 20, 000 a0 and is negligible
below 10, 000 a0 .

4.2.2

Zero-Temperature ModelIn the zero collision energy limit, the interaction between atoms and dimers can be describedin an effective field theory (EFT) formalism [50]. In the presence of deep dimers, the scattering
length between an atom and a dimer isaAD = (1.46 + 2.15 cot[s0 ln(a/a∗ ) + iη∗ ]) a.

(4.6)

At the atom-dimer scattering threshold E = −Eb , the s-wave differential cross section is |aAD |2 .
The effects of deep dimers encoded in η∗ is to fill in the zeroes of the cross section and prevent
it from saturating the unitarity bound. As to inelastic collision, the dimer relaxation constant isdeduced from the imaginary part of aAD and given byβAD =20.3 sinh(2η∗ )~a,
2sin [s0 ln(a/a∗ )] + sinh η∗ m
2

(4.7)

which is a log-periodic function of a∗ .
The scaling behavior βAD ∝ ~a/m can be understood as the following argument. For areaction event between a free atom and a dimer of size a in a unit volume, three things need tohappen. 1) The atom and dimer approach each other at a distance on the order of a. 2) The three

97

(a)

Dimer lifetime (ms)100dimer one-body lifetime
dimer react with atoms:

10

1

100

1,000

100

1,000

(b)
10-8

10-9

10-10

Figure 4.8: (a) Lifetimes of dimers with and without atoms being present at various temperatures.
Error bars on dimer lifetimes are extracted from the fitting routine and include only the statisticalnoise on dimer number. Gray dots represent the intrinsic lifetimes of the dimers. Gray dashed line
is a guide to the eye. (b) Temperature dependence of atom-dimer relaxation coefficient βAD as afunction of scattering length a. Only atomic density is varied between the two highest temperature
data sets. Error bars on βAD stand for 1σ propagated uncertainty involving the statistical errorof atom density as well as the uncertainty of dimer lifetimes. Solid lines are fitting curves witha zero-range finite-temperature model [74] at corresponding temperatures. The inset of (b) is a
cartoon picture that illustrates the role of a trimer in enhancing the atom-dimer relaxation process
as an intermediate decay product.

98
atoms together come from distance a to R0 which represents the recombination region. 3) Oncethey are at R0 they recombine at a certain rate. The probability for the first two steps are a3 and(R0 /a)2 (as a property of the scale-invariant hyperradial wavefunction in the R0 < R < a region([52], etc.), and the recombination rate of the final step is ~/(mR02 ). In the end, the total rate
coefficient turns out to be proportional to ~a/m. The reason that the short-range length scale R0
drops out of the picture can also be viewed as a result of the “mode matching” between three free
atoms within a distance of a and the Efimov hyperradial wavefunction [76]. In other words, theEfimov channel in the R0 < R < a region can transport atoms to the short distances without any
barriers.

4.2.3

Finite-Temperature Model

We follow the model introduced in [74] to incorporate finite temperature effects in the fitting
function for the atom-dimer relaxation data. The decay rate of dimer number is expressed as athermal average of the rate constant,
dND
1
=− 3
dt
2π

Z ∞

2

Z ∞

r dr
0

2

Z kD

Z 1

k dk ×

P dP
0

2

0

dxnA (pA , r)nD (pD , r)g(k),

(4.8)

−1

where P = pA + pD is total wave number of atom and dimer, k = 32 pA − 31 pD is relative wave
number of atom and dimer, and x is the cosine of the angle between P and k. We only evaluate the√integral of k up to the dimer break-up wave number kD = 2/( 3a) for the model to be valid. We
assume Boltzmann distribution on nA (pA , r) and nD (pD , r) since the temperature here is always
higher than 2 × Tc for atoms. The different density distributions of atoms and dimers in a harmonic
trap come into play in this step.
The function g(k) is related to the atom-dimer inelastic scattering cross section,g(k) ≡

3~k  total
3~k inelastic
elastic
σAD(k) =σAD (k) − σAD
(k) ,
2m
2m

(4.9)

total (k) and σ elastic (k) from the atom-dimer scattering
where m is the atomic mass. We derive σADADamplitude fAD (k) that is parametrized by series of the product ka [74]. The short-range properties

99
of atom-dimer scattering are encoded in fAD (k) by a∗ and η∗ which serve as the fitting parametersin our case. We then evaluate the integrals in Eqn. 4.8 to fit βAD from 300 a0 up to 2500 a0 . The
upper bound guarantees that the temperature of the data always stays below one tenth of the dimer
break-up energy. The lower bound is due to the assumption that this model breaks down outside
the universal regime.
We use the Boltzmann distribution function for a harmonic trap whose geometric mean trap
frequency is denoted as ω̄,ni (pi , r) =

~ω̄
kB T

3

with the normalization condition being

R

 2 2


~ pi
mω̄ 2 r2
1+,
Ni exp −kB T 2mi
2(4.10)

R
d3 r (2π)−3 d3 pi ni (pi , r) = Ni (i = A, D). By combining

Eqn.4.8 with Eqn.4.4, we obtain the following expression for βAD (T, a)4~3
βAD (T, a) = √π (mkB T )3/2

Z kD
0


3~2 k 2k dk exp −g(k).
4mkB T2

(4.11)

The algorithms for the numerical evaluation of this formula can be found in Appx. C. Some predictions from this model are plotted in Fig. 4.9 at various temperatures. In the T → 0 limit thisfinite-temperature model perfectly reproduces the zero-temperature EFT result multiplied by afactor of 1.5 which accounts for the spatial overlap between atom and molecule clouds. It is worth
noting that in the large a limit, this model predicts a suppressed βAD . This trend is consistent withan independent model explained by [77]. To understand this behavior, consider that for ka  1
atoms do not interact beyond the distance of 1/k. The three-step argument in this regime is now
modified as the following. 1) The two atoms in the dimer approach each other at a distance 1/kwith a probability of 1/ka. This can be calculated from the universal dimer wavefunction. 2) The
third atom also gets to the distance 1/k with a probability of 1/k 3 assuming a unit volume. 3)
The three atoms come to R0 from 1/k with a probability of (kR0 )2 . The recombination rate at
the short range is again ~/(mR02 ). In the end, the total rate becomes βAD ∝ ~/(ak 2 ). Briefly, thesuppression comes from the mode-matching that happens only at distances 1/k instead of a [76].
The above description presumes a well-defined atom-dimer channel for ka  1.

100

𝛽𝐴𝐷

(a)

𝑎 𝑎0

𝛽𝐴𝐷

(b)

𝑎 𝑎0

Figure 4.9: Finite-temperature model predictions for atom-dimer reaction rate coefficient βAD . (a)Compare the finite-temperature predictions at 0.01 nK (red line) and 0.1 nK (blue line) to the zerotemperature EFT prediction (black line). Inset: zoom in on the atom-dimer resonance peak on
linear scale. The zero-temperature theory is multiplied by a factor of about 1.5 to account for the
spatial overlap between the atom and dimer clouds in thermal equilibrium in a harmonic trap. This
comparison justifies the validity of the finite-temperature theory in the T → 0 limit. (b) βAD vs.
a at various temperatures. The three highest temperatures represent the conditions used in our
experiment. Their peak height is highly suppressed due to the finite temperature effects.

101
In the fitting process, the temperature of the sample is an input parameter to the fitting
function and measured with absorption images on atoms after a long time-of-flight. As depictedby the solid lines in Fig. 4.8 panel (b), this finite temperature model captures the shape of theatom-dimer resonance across the whole temperature range accessed in our experiment.

4.2.4

Temperature Dependency of Atom-Dimer ResonanceIn the above mentioned two fitting protocols, short-range physics is lumped into three fitting
parameters, a∗ , η∗ and magnitude A∗ . Similar to the definitions in the a < 0 region, a∗ denotes theatom-dimer resonance position at zero temperature and describes how close the three particles can
get together in the presence of any repulsive inner wall on the potential; η∗ denotes the width of theresonance and quantifies how absorptive the wall is; a∗ and η∗ bundled together serve to fix the realand imaginary part of the three-body parameter (3BP). A∗ indicates the global numerical factor
in front of βAD . It is associated with the total effect of deeply bound molecules on the inelastic
exit channel, which can be a model-dependent number. It also handles experimental errors on theabsolute density of atoms.
The fitting parameters are summarized in Fig. 4.10 as well as in Table 4.2.4. In contrast to
the phenomenological fit results from a zero temperature model [78], the finite-temperature model(1)

reveals an energy independent three-body observable a∗that we associate with the first excited(1)

Efimov trimer state. We determine a∗ to be 884(14) a0 = 13.7(2) rvdW from the weighted mean of(1)

the four experimental points, indicated as the solid line in Fig. 4.10 panel (a). Combining a∗ with
(0)

(1)

(0)

the three-body recombination resonance a− = −908(11) a0 [18], we obtain a∗ / a−

= 0.97(2).

This number is only 8% lower than the universal ratio of 1.06522 obtained from the zero-range(0)

theory [50, 53]. Considering the fact that our a− = −14.05(17) rvdW dramatically deviates fromthe value according to the single-channel van der Waals universality, the ratio reported in thiswork implies that short-range physics shifts the absolute positions of the two resonances in thesame direction while leaving their relative positions unchanged. Interestingly enough, the widthparameter η∗ overlaps with the previously measured η− for a < 0 within experimental uncertainties

102
T̄atomnK
218(13)
213(12)105(5)
62(3)hnA i
1010 cm−3
8.6(7)
25.9(1.4)
7.5(6)
3.7(2)hnD i
109 cm−3
7.8
7.8
6.4
1.4

ω̄/2π
Hz
42.7
42.7
24.4
17.1

(1)
a∗

η∗

A∗

a0
901(28)853(43)862(28)896(24)0.31(3)0.37(6)0.26(4)0.26(3)0.85(9)0.67(9)1.14(19)
1.55(15)

Table 4.1: Atom-dimer resonance measurement conditions and fit results. The magnitude A∗
includes corrections for the density overlap between atom and dimer clouds according to the temperature measurements on those two components respectively. The error bars on temperature anddensity represent the standard deviation of statistical noise.

(see Fig. 4.10 panel (b)). With the above observations, we report a 3BP that continuously connects
a < 0 and a > 0 at short range, although it is approached along disparate scattering channels inthese two regions in the long range.

4.3

Analysis of Systematic Effects

4.3.1

Atom-Dimer Thermal Equilibrium

The temperature of the samples used as an input parameter in the fitting function is measured
with the atom component. Since the molecule creation efficiency is closely correlated with thelocal properties of the atom cloud in a harmonic trap, the molecule component may be createdsomewhat out of equilibrium with the atom component. This causes uncertainty on the spatial
overlap between the atom and molecule clouds that enters the fitted global magnitude of βAD .
This issue is inevitable in our case since the decay of molecules is relatively fast compared to thethermalization rate. Extracted from their in-trap sizes, the effective temperature of molecules canbe up to 30% higher than the temperature of atoms (Fig. 4.11). The magnitude of βAD whichscales as T −3/2 thus can be at most 70% off from its correct value.

103

P e a k lo c a tio n ()

1 ,0 0 0
8 0 0
6 0 0
0 .6

W id th

(a )

(b )

0 .4
0 .2
1 0 0
2 0 0
T e m p e ra tu re (n K )3 0 0

Figure 4.10: Summary of the three-body observables extracted from finite-temperature and zerotemperature models. (a) Whilst the zero-temperature model (black diamonds) gives a phenomenological peak location of the atom-dimer resonance that shifts up as temperature decreases [78], thefinite-temperature model (magenta circles) reveals the true short-range atom-dimer 3BP that is
independent of temperature. (b) Temperature effect plays a remarkable role in the width of theresonance. It broadens the resonance feature hence leads to an incorrectly large η∗ in the zerotemperature fit. The η∗ from the finite-temperature model is free of this broadening effect andreflects the intrinsic width of the Efimov trimer state at zero energy. It is consistent with η− (openteal circles) reported in our previous work [18].

104

M o le c u le te m p e r a tu r e ( n K )

5 0 0
4 0 0
3 0 0
2 0 0
1 0 0
0
0

1 0 0

A to m

2 0 0

te m p e ra tu re (n K )3 0 0

Figure 4.11: Comparison of atom temperature and molecule temperatures in the atom-dimer reaction experiment. The out-of-equilibrium of the two components stems from two processes. One ismagneto association of molecules which has a spatially selective efficiency. The other is the rapid
inelastic collisions between atoms and molecules subsequent to the magneto association. The latterprocess preferably consume the molecules near the center where the atoms are the densest. Notethat the sagging between atom and molecule clouds due to gravity is negligible in our case.

105
4.3.2

Atom-Dimer Chemical Equilibrium

Near the pole of the resonance at positive a, the temperature of the cloud is no longer low
compared to the binding energy of the shallow dimers. In the absence of any dissipative processes,
atom and dimer relative populations are fixed by the chemical equilibrium condition. The dimer
density is related to the atomic density by√
nD = n2A λ3 2 2 eEd /(kB T ) .

(4.12)

This equilibrium condition is described by the Saha-Boltzmann equation [79]. Two mutually
time-reverse processes – three-atom recombination into shallow dimers and atom-dimer collisional
breakup – balance each other to form a dynamic equilibrium [80, 77]. We numerically estimate thedensity of created atoms as function of scattering length at various dimer sample initial conditions(see Appex. C for details of the program). A box trap is used here for simplicity. As depicted
in Fig. 4.12, atoms start to show up gradually around 1500 − 3000 a0 for the temperature rangeof 60 − 300 nK. In experiments we indeed observe free atoms emerging from initially pure dimer
samples at large a’s (see e.g. Fig. 4.14).

4.3.3

Dimer-Dimer Inelastic Collisions

One qualitative way to tell if dimer-dimer relaxation is competing with atom-dimer relaxation
is to see if the dimer decay data deviates from the exponential functional form. According to thiscriterion, inelastic collisions between dimers is not a concern in our experiment on atom-dimer
relaxation. A small deviation from exponential function can only be seen by going beyond 2000 a0
and also cranking up the dimer density at the same time. To quantify the dimer-dimer relaxation
rate in this large a region, we prepare pure dimer samples of relatively high density and extract
the dimer-dimer inelastic collision coefficient γDD . The rate equations turn out not trivial in thislarge a region for two reasons. The first reason is that the dynamical field ramp can projectmolecule and atom states onto each other when ramping beyond roughly 2500 a0 . The second
reason is that chemical equilibrium in the large a region requires a finite population of atoms at

106

𝐸𝑏𝑖𝑛𝑑 (𝑘𝐻𝑧)𝑛𝑎𝑓 (1010 𝑐𝑚−3 )𝑇𝑓𝑖𝑛𝑎𝑙 (𝑛𝐾)

𝑛𝑑𝑓 (1010 𝑐𝑚−3 )(a)𝑎 (𝑎0 )

atom

dimer

𝑎 (𝑎0 )𝑛𝑎𝑓 (1010 𝑐𝑚−3 )(b)𝑎 (𝑎0 )

Figure 4.12: Atom-dimer chemical equilibrium condition as function of scattering length for (a)dimer initial temperature 100 nK and density 2.1 × 1010 cm−3 . Eb is assumed to take the universal
form. It can be seen that atoms start to show up at around 2000 a0 . (b) The same calculation isrepeated with various initial dimer temperatures.

107
the experimentally available temperatures as discussed previously. In our case, these created atoms
can consume molecules and lead to additional losses.
One solution to the first issue is to pick a new home field at 3000 a0 to purify the molecules aswell as for detection. Molecules are not energetically separated from atoms very well due to its smallbinding energy of 12 kHz at this field. We take the advantage of the magnetic field insensitivity
on the |F = 1, mF = −1i to |F = 2, mF = 0i transition and apply Gaussian π-pulses with Fourier
width of ±3 kHz to spin-flip the atoms. We then use two cascaded ARPs to further bring atoms
into the |F = 2, mF = 2i state to receive the blasting light. Detection of atoms or molecules is doneat the same field with the final spin state being |F = 2, mF = −2i instead. And another solution
is to design a two-step magnetic field ramp to large a region so as to stay adiabatic with respect
to Eb throughout the ramp (Fig. 4.13).
To address the second issue, we image atoms and dimers individually to track both populations as a function of interrogation time. We take in situ images on both components to get thebest signal-to-noise and to directly obtain their density distributions. We use the following coupled
rate equations to describe the evolution of dimer density nD and atom density nA :
dnD
nD
0
=−
− βAD nA nD − γDD n2D − δAD n2D + δAD
n2A nD + higher order terms in nA ,
dt
τbg
dnA
0
= − βAD nA nD − 2δADn2A nD + 2δAD n2D + higher order terms in nA ,dt

(4.13)

0
are thewhere γDD denotes dimer-dimer relaxation coefficient into deeply bound states, δAD and δAD

dimer dissociation/formation coefficients for atom-dimer conversion. For atom-dimer conversion, weonly consider the leading order binary process and its inverse ternary process K2 +K2 ↔ K2 +K+K
given that atom density is small compared to dimer density [80]. Molecules containing more than
two atoms are ignored here for simplicity. Keeping only the linear terms in nA and using theempirical βAD values, we numerically fit the coupled rate equations and extract γDD . This fitting
process is presented in Fig. 4.15.
We map out γDD as a function of a in Fig. 4.16. The range in a is too limited to display any
scattering length dependence. Below 2000 a0 , we hardly observe any deviation from exponential

108

(a)

(b)

Figure 4.13: A two-step ramp to large a while staying adiabatic with respect to Eb . (a) Due tothe non-linear forms of the Eb and a as function of B, a single linear magnetic field ramp can
become highly diabatic near the end. (b) By doing two linear ramps with the latter step being
slow, we can reduce the dimer-to-atom conversion associated with the dynamics induced by theramp. Adiabaticity is qualitatively characterized by the dimensionless quantities (dEb /dt)/Eb2 and(da/dt)/a/Eb here.

109

(a)

(b)
0.1ms hold
10.1ms hold

14000

12000

12000

10000

10000

Nmolecule

Natom

14000

8000
6000

8000
6000

4000

4000

2000

2000

0

0.1ms hold
10.1ms hold
10.1ms hold, w/ 1b decay subtracted

01000

science field (a0)1000

science field (a0)Figure 4.14: Atom formation as function of a in dimer-dimer experiment. A home field at 3000 a0
is used to kill off the atoms. The data at short hold time 0.1 ms reflects the effect from the field
ramp. The data at long hold time 10.1 ms reflects the dynamics at a certain science field. By
subtracting out the dimer one-body decay rate, we obtain the dimer-dimer loss as presented by theopen blue diamonds in (b). Part of the lost dimers are converted into free atoms as shown in (a).
The rest of the lost dimers decay into deeply bound states. In the former process, the chemical
equilibrium is in fact never reached due to the dissipative nature of the dimers. The latter process
is what we are interested in and can be associated with the loss constant γDD .

110

( µm

-3

)

2 .0 x 1 0 -2
1 .0 x 1 0 -2

(a )

0 .0

(b )

( µm

-3

)

4 .0 x 1 0 -3

2 .0 x 1 0 -3

0 .0

0 .1

1

H o ld tim e ( m s )

1 0

1 0 0

Figure 4.15: Time evolution of dimers and atoms at 150 nK, 3500 a0 with dimers initially being the
majority in the sample. (a) Dimer number (orange circle) decay can be associated with multiple
inelastic processes. (b) Formation of a small number of atoms (purple diamonds) is driven by theatom-dimer conversion towards chemical equilibrium. This is similar to the process reported in [81]except that our dimer state is highly dissipative in nature. Dimer and atom data are simultaneously
fitted with Eqn. 4.13 as indicated by the orange and purple lines respectively.

111
decay form even with our highest dimer density 2.5 × 1010 cm−3 , neither can we observe any atoms
being created. With the dimer densities used in the atom-dimer experiment (1.4 − 7.8) × 109 cm−3 ,the lower bound on the dimer-dimer relaxation time is on the order of hundreds of milliseconds
to seconds. This is significantly longer than the time scales of the atom-dimer relaxation in the
scattering length range of interest. Therefore we can safely ignore the collisions between dimers in
the data analysis in Section 4.2.1.
1 .0 x 1 0 -8
8 .0 x 1 0 -9
6 .0 x 1 0 -9
4 .0 x 1 0 -9
2 .0 x 1 0 -9
0 .0
3 ,0 0 0

4 ,0 0 0

1 5 0 n K

3 0 0 n K

5 ,0 0 0

6 ,0 0 0a (a 0)Figure 4.16: Dimer-dimer relaxation coefficient into deeply bound molecules measured above
3000 a0 . γDD is extracted from the fit with the coupled rate equations 4.13. A slight difference inmagnitude can be resolved between the two temperatures. This measured γDD has the same orderof magnitude as a theoretical simulation done for 133 Cs in [82].

4.3.4Residual Breathing Motion of the Atom Cloud

The atom clouds present residual breathing motion during the atom-dimer relaxation process
in spite of the 10 ms hold time added to damp it out. The breathing motion is at twice the trapping
frequencies and is relatively slow compared to the dimer lifetimes. In order to estimate the atom
density uncertainty introduced by this breathing mode, the sizes of the atom cloud are being traced
for at least one period along the slow direction as presented in Fig. 4.17. The amplitude of the112
fractional size change is up to about 5%. The actual correction to the atomic density is less than
3% during the short dimer lifetimes. They are calculated based on a sinusoidal fit whose phaseis extracted from in situ sizes of the atoms. The corrections to atom densities are included inindividual conditions listed in Table 4.2.4.
(a)

(b)
30

Model

sloshWithDecay (User)

Equationy0+A*sin(2*pi*f*(x-x0)/1000)*exp(-x/tau)

0.71112

Adj. R-Square

Value

28

Standard Error

A

1.04483

0.1719y025.77855

0.06025x0-1.13831

0.18565f

156.32974

2.04981

tau

40.86276

21.1332844sx (px)sz (px)sigz

46

0.12342

Reduced Chi-Sqr

26

42
sloshWithDecay (User)

Model

y0+A*sin(2*pi*f*(x-x0)/1000)*exp(-x/tau)

Equation

0.56204

Reduced Chi-Sqr

240.69464

Adj. R-Square

40

Value

sigx

1.533960.17491y042.76140.15243x0

14.13261

0.93996f

25.7836

2.28667tau

22Standard Error

A

1E7

0

38
0

5

10

15

20

25

30

35

hold time (ms)

0

5

10

15

20

25

30

35

hold time (ms)Figure 4.17: Residual breathing motion of atoms in atom-dimer relaxation experiment along thevertical (a) and horizontal (b) directions measured at a long time-of-flight. The breathing frequencies are twice the trapping frequencies along corresponding directions. The inset tables show thefitting parameters in a damped sine function.

4.3.5

In Situ Detection of Molecules

As a Feshbach molecule is dissociated by microwave field, the parent atoms acquire excess
kinetic energies from the blue-detuned microwave photons. This leads to a fast expansion of thecloud as can be seen from Fig. 4.18. We make sure that the time gap in between the dissociation
pulse and imaging pulse is of minimal duration so that the parent atom cloud can largely preserve
the in situ density distribution of the molecule cloud.

113

(b)
20

36

18

34

16

32

14

30

sx (px)sz (px)

(a)

12
10

28
26

8

24

6

22

4

20

2

18
100

1000

actual expansion time (us)

100

1000

actual expansion time (us)Figure 4.18: Expansion of post-dissociation molecule cloud along the vertical (a) and horizontal(b) directions. The x-axis expansion time is defined as the interval between the microwave dissociation pulse center and imaging light pulse center. The fast expansion of the dissociating atoms is
consistent with the excess energy of the blue-detuned microwave photons.

114

4.4

Bound-to-Bound Transitions
For the remainder of this chapter, we explore the viability of using bound-to-bound transitions

to detect the real component of the atom-dimer and/or dimer-dimer scattering length. The idea is
to take the advantage of a narrow transition lineshape and a stable final state to probe small meanfield energy shifts in the loosely-bound upper level. One potential challenge is the tiny transition
dipole moment that is limited by the Franck-Condon overlap between the initial and final states.
The three bound states of interest can be found in Fig. 3.3 and Fig. 3.6.
We calculate the off-diagonal matrix element between two bound states to estimate how small
the transition dipole is. The coupled-channel representation of the wavefunction of a certain boundstate can be expressed asΨα (r) =5Xfν(α) (r) |νi ,(4.14)

ν=1

where α = −1, −2, −3 indexing the three shallowest bound states, ν = 1, 2, 3, 4, 5 indicating the(α)

five channels on the basis of free atom asymptotes, fν (r) represents the projections onto thesechannels. Fig. 4.19 gives the amplitude of projections as function of inter-nuclear distance r. The
transition dipole can then be calculated by evaluating the off-diagonal matrix elements
Z ∞hΨα | µB gS Sz + µN gI Iz |Ψα0 i =
0

(α0 )

hν| µB gS Sz + µN gI Iz |ν 0 i fν(α) (r) · fν 0 (r)r2 dr,(4.15)

where µB , µN are the Bohr magneton and nuclear magneton, gS and gI are the g-factors of theelectrons and nuclei, respectively, and Sz , Iz are the corresponding spin operator components inthe direction of the magnetic field. A longitudinal magnetic rf drive is used since the total mF isconserved in this kind of transition. The dipole moments turn out to be really small, only 0.05 µB
for the (−1) → (−2) transition and 0.01 µB for the (−1) → (−3) transition at 35 G. And it decreases
even more as we go to larger a.
The small dipole moment severely limits the Rabi frequencies on these transitions. We
improve our magnetic field control system by suppressing the rf rectification as described in Section
2.1.6 so that we can send rf or microwaves at their maximal powers to the atoms. In the end, we are115

Projections onto Asymptotic Atomic Basis Set
(a)

0.20
0.15

1 -1
2 -1
3 -1
4 -1
5 -1
1 -2
2 -2
3 -2
4 -2
5 -2

projections

0.10
0.05
0.00
-0.05
-0.10
-0.15
-0.20
10

100

1000

r (a0)

(b)

0.20
0.15

1 -1
2 -1
3 -1
4 -1
5 -1
1 -3
2 -3
3 -3
4 -3
5 -3

projections

0.10
0.05
0.00
-0.05
-0.10
-0.15
-0.20
10

100

1000

r (a0)

Figure 4.19: Projections of shallow dimers onto asymptotic atomic basis set as function of internuclear distance r, evaluated at 35 G. 1 → cc channel, 2 → db channel, 3 → ec channel, 4 → fd
channel, 5 → ee channel. (-1) indicates the shallowest bound state (Feshbach molecule), (-2) and(-3) indicate the second and third shallowest bound states. The binding energies of the second and
(−2)
(−3)
third bound states are calculated to be about Eb= 105 MHz and Eb= 363 MHz respectively
at 35 G field. It is worth noting that the cc (open) channel component of Feshbach molecule state
has a large weight in the classically forbidden region, which is as expected.

116
able to shovel a decent fraction of Feshbach molecules into the (−2) bound state by square pulses.
The molecule population transferred to (−2) state manifests as losses in the (−1) state as shownin Fig. 4.20. Our pulse duration is limited by the magnetic field sensitivity of this transition about
0.5 µB . The uncertainty on the center frequency is 100-200 Hz which is unfortunately not goodenough to resolve the small energy shift due to atom-dimer or dimer-dimer interactions. For the(−1) → (−3) transition which is less sensitive to the magnetic field but requires even higher pulse
power, we are only able to see some number losses from the (−1) state by using a frequency-scanning
ARP.

117

(b)

20000

20000

15000

15000

Nmol

Nmol

(a)

10000
Gaussian_user (User)

Model

Gaussian_user (User)

Equationy0+A*exp(-(x-xc)^2/(2*w^2))

Equationy0+A*exp(-(x-xc)^2/(2*w^2))

ReducedChi-Sqr

5189377.7
141

ReducedChi-Sqr

3542532.
08948

Value

N

0
0.000

0.80278

Adj. R-Squar

5000

0.005

10000

Model

0.010

rf freq (105.5+MHz)

Standard Erro

0.71981Adj. R-Square

5000

Value

Standard ErroA

-12436.46021

909.0882

A-8495.66482

579.30587

w

0.00159

1.53019E-4w

0.00261

2.56135E-4

y017589.78589

552.07207

y0

18597.95897

348.58046

xc

0.009641.26019E-4

xc

0.00908

2.10994E-4

0.015

0.020

N

0
0.000

0.005

0.010

0.015

0.020

rf freq (105.5+MHz)

Figure 4.20: Bound-to-bound transition spectra for the (−1) → (−2) transition at about 500 a0
field. Nmol stands for the remaining number of Feshbach molecules. A square pulse of 200 µs
duration is applied in (a). Its Fourier width (r.m.s.) is about 1.7 kHz. A square pulse of 150 µs
duration is applied in (b). Its Fourier width (r.m.s.) is about 2.3 kHz. The fitting functions andparameters are listed in the inset tables. Both spectra have Fourier-limited widths.

Chapter 5

Experiments on Free Atoms

In Ch. 4, we elaborated on three-body Efimov physics involving Feshbach molecules as an
essential component. The other set of three-body observables in the a > 0 region can be locatedby scanning along the free atom scattering threshold instead of along the atom-dimer scattering(a>0)

threshold. Unlike on the a < 0 side, the three-body recombination rate constant L3

is free of

scattering resonances because the trimer states never meet the free atom scattering continuum on(a>0)

the positive side. Instead, L3

is a sinusoidal function of a, featuring a series of interference
(n)

(n)

minima and maxima, which are denoted by a+ and ap

respectively in this thesis. The parentage

of associated Efimov trimer state is indexed by the integer n, whose physical meaning is the number(a>0)

of nodes in the three-body hyperradial wavefunction. Knowledge of L3

are complementary tothat on βAD since the atom-dimer scattering potential which acts as the entrance channel in atom(a>0)

dimer reaction now becomes the exit channel in three-atom recombination. Although L3

hasbeen previously measured on a variety of atomic species, we are the first group to report all values
(n)

(n)

of a+ , ap(n)and a∗

that can be used to characterize the Efimov spectrum for a homonuclear

system.

5.1

Experimental Sequence
The experimental sequence for the three-body recombination measurements is relatively

straightforward. At the end of the evaporation sequence, we ramp up the trap depth adiabatically to suppress residual evaporation. A typical timeline of the experiment can be found in

119
Fig. 5.1. The final trap depth is set to be 10 times the temperature of the cloud whereas the trap
depth is found to be 8 times the temperature of the cloud in equilibrium as presented in 2.1.3. The
peak value of the phase-space density (PSD) of the clouds is always restricted to below 1 to ensure(a>0)

Boltzmann statistics. Due to the small magnitude of L3 on the positive side (L3(a<0)

orders of magnitude smaller than L3

is in fact threeat the same value of |a|!), this peak PSD is made to be(a<0)

a few times higher than that used in the L3

measurements in order to avoid too long inter-rogation times. We avoid using BECs because they introduce density-related many-body effects
which can be difficult to interpret. In order to ensure that the atoms are in thermal equilibrium, wemust keep the three-body recombination rate small compared to the two-body elastic collision rate,i.e. Γ3,inel /Γ2,el  1. For our case, Γ3,inel /Γ2,el = (0.05 − 0.21)λna2 for finite scattering lengths
and Γ3,inel /Γ2,el = 0.04nλ3 for unitary gases [83], where λ is the de Broglie wavelength and n is
the number density. Therefore, a small PSD is always desired to satisfy the thermal equilibrium
condition.

Bias b4 science(preparation)

Hold 100ms at 60𝑎0
for final thermalizationb4 science

Hold 10ms at 500𝑎0 for
1,0 cleaning rf pulses

ramp to60𝑎0

Vary hold at
science field

Hold 6ms at
500𝑎0 for uwave
Trap/bias off
for imaging

Figure 5.1: Timeline of three-body recombination experiments, B vs. time. The arrows indicate
stages that are to be adjusted.

In order to cover a wide range dynamic in a, we stitch various initial conditions together.
As illustrated in Fig. 5.2 panel (a), the dynamic range of a certain condition is limited on the low
a side due to the competing timescales from the two-body loss and background loss, and on thehigh a side by the lack of thermal equilibrium. The timescales of the background loss (i.e. vacuum

120
lifetime) are given in 2.1.3. All conditions are generated in the H2+V crossed dipole trap setup
except for the highest temperature one which is generated in the H1+V setup for the purpose of
maximizing the number density. As shown in Fig. 5.2 panel (b), at small scattering lengths, theatom loss becomes dominated by the background loss. This limits our lowest sampled a to be ataround 100 a0 . The competing timescales in this small a region were not scrutinized in previous
works.
(a)

(b)

10-20

10-22

•
•
•

loss is too fast (tens of ms)lack of equilibrium
capped by unitary limit

10-26
10-27
10-28

10-24
230(26) nK

10-29

at 460(30) nK

10-26

•
•
10-28
100

10-30

loss is too slow (tens of seconds)
bg loss and 2b loss are overwhelming
1,000

zero-T theory

10,000

100

200

300

400

Figure 5.2: (a) Dynamic range of individual condition in three-body recombination experiments.
Here the 230 nK data is used as an example. (b) Comparison of the three-body recombination
rate coefficient L3 (squares) with theoretical two-body inelastic collision rate L2 (circles) as well asbackground loss rate α measured on the |1, 1i ground spin state (diamonds). L2 and α are rescaledby empirical hn2 i/hni and hn2 i respectively to give the proper dimension.

5.2

Experimental Results(a>0)

The measured L3

as function of a at various temperatures is presented in Fig. 5.3 together(a<0)

with our previously reported L3

in [18]. Six different conditions are combined together to coveran eight order-of-magnitude range on L3 . The good overlaps between adjacent conditions indicate minimal density-dependent systematic effects (see 2.1.5 for discussions on absorption imaging
method). Since the overall scaling of L3 obeys an a4 law as anticipated for pair-wise interactions

121
(a>0)

[84], this scaling factor is divided out from L3

in Fig. 5.4 in order to emphasize the logarithmic

periodic modulation on L3 due to Efimov physics.
(n)

(n)

The destructive interference minima a+ (or constructive interference maxima ap ) can beelucidated in the picture of a Fabry-Pérot cavity. As illustrated by the inset cartoon in Fig. 5.4 rightpanel, the couplings between different hyperspherical potentials behave effectively as “mirrors” inthe three-body recombination pathways. The location of the inner mirror R1 is set by the range
of interaction r0 , and the location of the outer mirror R2 is set by the scattering length a in low
energy limit (R2 = min(|a|, 1/k) in general) [85]. The two decay pathways are associated withweakly and deeply bound molecular channels respectively. They are coherent and can interferewith each other. By scanning the scattering length a (or equivalently the “length” of the cavity),
interference minima and maxima show up with a period that is determined by the Efimov scaling
factor eπ/|s0 | defined in Ch. 3. A more rigorous expalanation of generic recombination pathways
can be found in Fig. 9 of [53].
(0)In order to extract a+ , we numerically fit the data set of 460 nK with an S-matrix formalismthat is introduced in [86]. The fitting curves are represented by the red solid lines in Fig. 5.4 and(0)Fig, 5.3. In this model, short-range physics is lumped into three fitting parameters, a+ and η+
accounting for the location and contrast of the interference minimum, and a scaling factor A+ foradjusting the global magnitude. Only collisional process of zero total angular momentum has been
(0)

(0)

included in the model. We obtain a+ = 246(6) a0 , η+ = 0.20(2) and scaling factor of unity.
This result agrees with the fit with a zero-energy model within uncertainties [50], indicating the(0)

negligible effects from finite temperature in this small a region. To extract ap , we then fix thecontrast parameter to 0.20 in the same finite temperature model and fit the data sets of 410 nK and
230 nK with their empirical temperatures as inputs to the fitting function. These fitting results are
summarized in Table 5.2. All data and fits are plotted individually in Appx. E. As we move on to(1)

larger a’s, severe temperature effects lead to noise and obscure the first excited minimum a+ . A
better understanding of the finite temperature effects and also systematic effects is needed to truly(1)

nail down the location of a+ .

122

10-20

10-20

𝐸=0

𝑟0

𝑎

10-22

10-22

10-24

10-24

60(6) nK150(6) nK210(23) nK230(26) nK
410(40) nK
460(30) nK

10-26

10-28

1,000

10,000

10-26

10,000

10-28
100

1,000

Figure 5.3: Compilation of L3 measurements for a < 0 (left panel) and a > 0 (right panel)
respectively. Solid lines represent finite-temperature model fits. The dotted line on the left panel
represents the a4 scaling in the absence of Efimov physics. The inset cartoon picture illustrates
the role of a trimer in enhancing the three-body recombination process as an intermediate decay
product. The dashed line on the right panel represents recombination into deeply bound dimers
with η = 0.20. The inset cartoon picture illustrates the coherent interference between different
decay paths in three-body recombination into shallow dimers.

Ti
nK
370
300
150

T̄
nK
460(30)410(40)230(26)

hnii
1011 cm−3
51.2(4.2)
26.4(2.8)
9.9(9)

ω̄/2π
Hz
260.1
79.4
54.6

(0)

a+
a0
246(6)
–
–

(0)

ap
a0
–
893(28)860(28)

(0)

η+

A+

0.20(2)0.20(0)0.20(0)0.99(8)0.74(2)0.77(2)(0)

Table 5.1: Three-body recombination conditions and fit results. a+ and ap are fit in a local rangeof scattering lengths to avoid presumptions about the absolute value of Efimov period. The η+with fitting errors of zero are fixed parameters in the fit. Ti and T̄ denote initial temperature andtemperature averaged over the decay time respectively. The error bar on T̄ represents the varianceduring the interrogation time. hnii is the initial density-weighted density. The error bar on hnii
represents the standard deviation of statistical noise. ω̄ is the geometric mean of trap frequency.

123

2 .5 x 1 0 -3 6

6 0 (6 ) n K
1 5 0 (6 ) n K
2 1 0 (2 3 ) n K
2 3 0 (2 6 ) n K
4 1 0 (4 0 ) n K
4 6 0 (3 0 ) n K

2 .0 x 1 0 -3 6
1 .5 x 1 0 -3 6
1 .0 x 1 0 -3 6
5 .0 x 1 0 -3 7
0 .0

1 0 0

1 ,0 0 0

1 0 ,0 0 0

Figure 5.4: Three-body recombination coefficient L3 scaled by a4 , measured at 460(30) nK (redsquares), 410(40) nK (yellow left-pointing triangles), 230(26) nK (green circles), 210(23) nK (purple down-pointing triangles), 150(6) nK (blue right-pointing triangles), 60(6) nK (cyan diamonds).
Colored solid lines represent finite-temperature model fits either to the neighborhood of the maximum or the minimum. The gray dotted line represents the recombination into deeply bound dimer
states with η = 0.20.

124

5.3

Analysis of Raw Data
In the atom number loss measurements, absorption images are taken after a long time-of-flight to track both population decay and heating as a function of interrogation time. To describe
the three-atom recombination decay process, we use the following rate equation for atom numberN:
1 dN= −αbg − L2 hni − L3 hn2 i,
(5.1)
N dt

3/2
where hni = N mω̄ 2 / (4πkB T )is the spatially averaged density-weighted density, hn2 i =√3is the averaged density-weighted square density, and αbg is the backgroundN 2 mω̄ 2 / 2 3πkB T
loss rate which is calibrated for each individual condition on the lowest spin state |F = 1, mF = 1i
at about 60 a0 . The temperature of the sample T is being tracked at the same time and fitted
into polynomial functions which are then plugged into the above rate equation in the form of clouddensity. The d wave resonance at 63.6 G introduced at the end of 4.1.2 is measured in a similarmanner with the L3 term being neglected in the small a region.
Ideally, we use relatively high density cloud while keeping the phase space density low sothat the leading contribution to the number loss is from the L3 term. This is the case for the data
below 2000 a0 . We fix L2 to its theoretical value in our fitting function and float the value of L3 inthis regime. We call this method fitting protocol I.
The density dependence of the losses in large a region turns out to be not so straightforward.
To figure out the actual decay mechanisms that are involved, we vary the density of the cloud byten orders of magnitude and float both L2 and L3 in the fitting process. The obtained L2 andL3 values as shown in Fig. 5.5 are relatively constant. This implies that these two terms alone
can deplete the density dependence and no four-body or higher order terms are needed in the rateequation.
Although we can always extract L2 and L3 by floating them in the rate equation, the fitting
uncertainty is fairly large due to the inseparable time scales of two- and three-body losses fora single sample condition. To suppress the fitting error, we combine decay data sets that are

125

1 0 -1 1

1 0 -1 2
1 0 101 0 -2 2A v e r a g e d e n s ity ( c m-3)

1 0 11

1 0 -2 3

Figure 5.5: Two- and three-body loss coefficients extracted from atom samples of various densities
at 210 nK, 3000 a0 . The dotted line shows the theoretical prediction of L2 at this field.

126
taken on two samples of drastically different densities but identical temperatures and fit themsimultaneously. This fitting method is presented in Fig. 5.6. In this fitting protocol II, L2 and L3
are shared parameters between data sets while atom numbers and temperatures from each data
set are independent. All the data points above 1000 a0 in Fig. 5.3 and Fig. 5.4 come out of fitting
protocol II. The fitting algorithms are documented in Appx. D.
(a )

3 0

2 0

(d )

6 x 1 0 3

3 0

2 0

3 x 1 0 3

N u m b e r

5 x 1 0 3
4 x 1 0 3
3 x 1 0 3

4 0

T e m p e ra tu re (n K )T e m p e ra tu re (n K )N u m b e r

(b )

(c )

4 0

2 x 1 0 3

0 .0 1

H o ld tim e ( s )0 .1

0 .0 10 .1
H o ld tim e ( s )1Figure 5.6: Illustration of fitting protocol II that combines atom samples of two different initial
densities but same time-averaged temperatures. (a)-(d) Time evolution of number (blue squares)
and temperature (red circles) for the high density condition (solid symbols) and low density condition (open symbols) respectively. The data traces are taken at 3000 a0. Red solid lines are fitting
curves for temperatures with polynomial functions. Blue solid lines are fitting curves for atom
numbers with Eqn. 5.1.

Above 2000 a0 , we discover some enhanced two-body losses that can not be explained by thetwo-body coupled-channel theory. The abnormally high L2 is plotted in Fig. 5.7. This inconsistency
is rather unexpected considering the perfect overlap between experiment and theory on the dimer
spontaneous dissociation lifetime in the same region (see Section 4.1.2). Our detailed two-bodymodel of 39 K does not predict a resonance in the vicinity of a = 3000a0 . Our analysis protocol II
yields not only the unexplained excess two-body decay at a = 3000a0 , but also scattered temperature dependence of the three-body rate for a > 1500a0 . Taken together these results suggest thatthe overall explanation may come from some combination of uncharacterized systematic errors. We
note that as a becomes larger in magnitude, the binding energy of the dimer is no longer much

127
greater than kB T , opening up other collisional channels. Also the energy released by three-bodydecay into the shallow dimer becomes no longer large enough to eject the decay products withabsolute certainty [87, 77]. Finally, the cloud begins to become collisionally opaque, opening thepossibility of decay products undergoing some collisional cascade. Our simple attempts to develop
a model for our observations could not establish any of these mechanisms as leading to a definitive
explanation. A more sophisticated model that incorporates effects such as secondary collisions will
be needed to track down the kinetics in this large positive a regime.
Below in Section 5.5, hypothesized sources of systematic effects in atom loss rate experimentare discussed. However, none of them are well supported by the auxiliary measurements that we
made.

5.4

Finite Temperature Model of Three-body Recombination
We follow closely the S-matrix method introduced in [86] to fit L3 as a function of a inorder to account for the finite temperature effects. The three-body recombination coefficient isformulated as:
Z ∞L3 = Nlost K3 (T ) = Nlost

0



(J=0)
(J=0)
dEE 2 e−E/(kB T ) Kshallow (E) + Kdeep (E)Z ∞,3!
dEE 2 e−E/(kB T )

(5.2)

0

where Nlost = 3 denotes the number of atoms being lost per collision event, the factor of 3! in the
denominator accounts for the indistinguishable nature of three identical particles, and E denotes the
(J=0)

(J=0)

total energy of the three particles. Kshallow (E) and Kdeep (E) represent the decay rate constants
into the shallow dimer and deep dimers respectively at finite energies . They can be expressed interms of universal scaling functions s12 (ka), s22 (ka) and s11 (ka):!
√s212 e2iθ−2η 2(1 − e−4η )|s12 |2 ~a4
144 3π 2(J=0)−
Kshallow =1 − s22 +,(ka)4m
1 − s11 e2iθ−2η
|1 − s11 e2iθ−2η |2
√144 3π 2 (1 − e−4η ) 1 − |s11 |2 − |s12 |2 ~a4(J=0),
Kdeep =(ka)4m
|1 − s11 e2iθ−2η |2(5.3)

where the momentum k is related to E through E = ~2 k 2 /m with the atomic mass m, and theangle θ is defined as s0 ln(a/a+ ) with s0 ≈ 1.00624 for three identical bosons. We only include the128

8 ,0 0 0

3 ,0 0 0

2 ,0 0 0

1 ,5 0 0

1 ,2 0 0

1 ,0 0 0

1 0 3

L 2 th y

1 0 2

2 1 0 (2 3 ) n K
1 5 0 (6 ) n K
6 0 (6 ) n K
3 0 (3 ) n K

1 0 1
1 0 0
1 0 -13 3 .6

3 3 .83 4 .03 4 .23 4 .43 4 .6

3 4 .8Figure 5.7: Two-body loss rate coefficient L2 as function of B and a. The open symbols aremeasurements taken with pure atom samples and extracted by using fitting protocol II. Theseexperimental points present a strong enhancement around 33.9 G (or 3000 a0 ) and show little temperature dependence. This enhancement gradually diminishes for a below 2000 a0 . Solid red circles
are deduced from the measured shallow dimer lifetimes τbg in the large a limit by using a modelfrom [70]. The dashed line indicates the prediction from our refined two-body cc model on sd wave
basis.

129
contribution of zero total orbital angular momentum J = 0 in the above equations. The numericvalues of sij (i, j = 1, 2) matrices are taken from [86] in the region of 0 < ka < 1. Prior to doing
(J=0)

(J=0)

the thermal averages, Kshallow and Kdeep

are smoothly extrapolated to their threshold values in

the limit of ka → 0 where this numerical method loses its accuracy [88]:√3) sin2 θ + sinh2 η ~a4,msinh (πs0 + η) + cos2 θ√(3/2) × 128π 2 (4π − 3 3) sinh(2πs0 )sinh(2η) ~a4(J=0)
Kdeep (T = 0) =.
sinh2 (πs0 )
sinh2 (πs0 + η) + cos2 θ m

768π 2 (4π − 3(J=0)
Kshallow (T = 0) =
2

(5.4)

We use the measured temperatures of each condition as a fixed input parameter in the fit.
(0)

Two free fitting parameters a+ and η+ which account for the location and width of the three-body(0)

recombination minimum can then be extracted. When fitting for the maximum ap , θ is substituted
with π/2 + θ. The details of the algorithms can be found in Appx. D.
We also play with this model to further explore the role of angular momentum terms with
J > 1 in three-body recombination into the shallow dimer. The contribution from arbitrary angularmomentum J takes the form
K

(J)

√
144 3π 2 (2J + 1)fJ (ka) ~a4
(E) =
,
(ka)4
m

(5.5)

where fJ (ka) is related to the imaginary part of the atom-dimer scattering phase shift through(J)fJ (ka) = 1 − exp −4 Im[δAD (E)] . The numerical values of fJ (ka) are provided by D. Kang. For
three identical bosons, the leading term is J = 2. We extrapolate fJ (ka) to small values of ka byusing their asymptotic values according to the threshold behaviors,fJ (ka) = aJ (ka)2λJ +4 + bJ (ka)2λJ +6 ,(5.6)

where λ1 = 3 and λJ = J for J = 2 to 6. As presented in Fig. 5.10 panel (b), the contribution
from J > 1 terms can even overtake the J = 0 term at large enough a according to this model.
This observation could qualitatively explain why the second recombination minimum is completely
smeared out by the finite temperature effects. This observation is consistent with an independent
model given in [89].

130

(a)

(b)

Figure 5.8: Using the S-matrix method to calculate L3 at finite temperatures. (a) Universal scaling
functions s12 (left) s22 (middle), and s11 (right) as functions of the product ka. The modulus |sij |and the phase arg[sij ] of the functions are shown in the top and the bottom panels, respectively.
(b) The rescaled integrand Kshallow and Kdeep in the J = 0 sector as defined in Eqn. 5.4. The
different curves are evaluated for a = 1 − 1000 a0 as function of ka. The variance on individual
Kshallow curves is related to the modulation due to Efimov physics. The variance on individual
Kdeep curves is much smaller and cannot be resolved on this scale.

131

(b)

𝐿3 (𝑐𝑚3 /𝑠𝑒𝑐)

𝐿3 /𝑎4 (𝑐𝑚6 /𝑠/𝑎0 )(a)𝑎 (𝑎0 )𝑎 (𝑎0 )

Figure 5.9: Finite-temperature model predictions for three-body recombination rate coefficient (a)(a>0)(a>0) 4L3(b) and L3/a in the total angular momentum J = 0 sector. It is worth noting that
here the finite temperature effects shift the recombination minimum to higher a; contrarily the(a<0)
resonance peaks on L3
are shifted to lower a.

132

(c)

𝐿3 /𝑎4 (𝑐𝑚3 /𝑠𝑒𝑐)(a)𝐸/𝐸𝑏 or 𝑘𝑎 2𝐿3 (𝑐𝑚3 /𝑠𝑒𝑐)

𝐿3 (𝑐𝑚3 /𝑠𝑒𝑐)(b)𝑎 (𝑎0 )𝑎 (𝑎0 )

Figure 5.10: Contributions to the three-body recombination rate into the shallow dimer fromangular momenta J = 1 − 6 (a) The rescaled integrands for J = 1 to 6 as function of (ka)2 . (b)Comparison of contributions from J = 0 (blue line) and J > 1 (red line) at 20 nK (bottom left) and
100 nK (bottom right). The yellow lines represent total loss rate with J = 0 to 6. (c) The effects
of J > 1 terms tend to shift the recombination minimum to smaller a. This trend is intriguingly
opposite to the finite temperature effects from the J = 0 term as shown in Fig. 5.9.

133

5.5

Analysis of Systematic Effects
We spent a tremendous effort to understand the anomalously high values of L2 around 3000 a0 .

Does it actually come from technical artifacts or real physics? According to a coupled-channel
calculation on sdg-wave basis, there is no resonance-like features in this region at all as depicted
in Fig. 5.11. And we already show in Section 4.1.2 that this model successfully reproduces the
Feshbach molecule lifetime data without any scaling fudge factor. Both two-body loss rate L2 andmolecule one-body lifetime τbg have their roots in the d-wave sector in the two-body Hamiltonian.
The excellent agreement between experiment and theory on τbg suggests that the enhancement on
L2 could be caused by systematic effects rather than real resonances. To summarize, the systematiceffects that have been ruled out are listed as the following: The atom loss rate is uncorrelated with the magnetic field ramping time to/from the sciencefield in the range of 0.5 to 50 ms.
 There is no detectable Feshbach molecule population during the loss processes.
 There is no abnormal magnetic field noise at this scattering length, as evidenced by the rf

spectroscopy in Fig. 5.12.
 Increasing the trap depth by a factor of 2 does not affect the loss rate (see Fig. 5.13 (a)).

The initial temperature and density of the cloud generated in these two traps are delicately
kept the same (Fig. 5.13 (b)).
 The frequency difference between two ODT beams (vertical and horizontal) does not playa role. The frequency difference is to prevent interference pattern between the two beams.
 We do frequency and amplitude modulations on the two ODT beams to purposely en-large the technical noises on the ODT frequency and intensity. This experiment shownin Fig. 5.14 rules out any Raman-type two-photon transitions in which two free atoms
are associated into Feshbach molecules which are then lost from the trap via subsequentcollisions with other atoms.

134Interestingly, none of the above technical concerns appears to be correlated to the two-body loss
peak that we have observed around 3000 a0 .
Another systematic effect that has been frequently suspected in previous experiments [90,91, 92] but never resolved in a quantitative way [93, 94, 87, 95] is the so-called avalanche effect.
An avalanche effect is hypothesized to happen at large a’s where the Eb of the shallow dimers issmall compared to the trap depth. The released energy on one of the recombination products,
shallow dimers, is not high enough for them to escape the trap. The remaining shallow dimers
will then collide with free atoms with a relatively large (in)elastic cross section at a small relative
momentum. These secondary inelastic collisions can lead to extra loss of atoms hence modify the
population decay vs. hold time in a sophisticated way. These shallow dimers could in fact only
exist in the trap for a transient amount of time thus evade regular detection methods. Or their
population could be just too small to be observed.
We introduce a toy model to emulate the effects of secondary collisions between shallow
dimers and atoms. We use empirical values on the atom-dimer relaxation rate constant βAD andtheoretical values on K2 in the following coupled rate equations
dnA= −3 K3deep + K3shallow n3A − 2K2theory n2A − βAD nD nA − αnA ,
dt
nD
dnD
= 3K3shallow n3A − βAD nD nA −
.
dt
τbg

(5.7)

The predicted decay curves by these rate equations are plotted in Fig. 5.15 and compared to ourdata. This model can indeed grasp the formation of shallow dimers at the early stage of the threebody recombination process. And these shallow dimers decay away very quickly by reacting withthe atoms. However, the population of these shallow dimers is not large enough to lead to persistentatom number losses at long hold times.
Although the toy model seems to suggest a negligible role of the existence of shallow dimers
in the trap, it is too primitive to completely rule out their effects in the kinetics of the atoms.
A Monte Carlo simulation such as in [94, 87] will be needed to take into full account the energy
dependent atom-dimer (in)elastic collision cross section as well as the trajectories of the atoms and

135

39

K coupled channels -Im(a) for cc channel-Im(a) (a0)sd basis
sdg basis

1

33.7

33.8

33.9

34

B (G)

Figure 5.11: A coupled-channel calculation of −Im(a) in sdg-wave basis near 34 G. The numerics
are provided by P. S. Julienne.

4000

Model

GaussianPulse (User)

Equation

A*exp(-(x-xc)^2*(2*pi*tau)^2)+offset

ReducedChi-Sqr

100832.7709
10.84648

Adj. R-SquareValue

B = 34.10405+/-0.28 mG

A
N

0.067

offsetxc

3000

atom number

Standard Error

3213.13174

tau

62.60809
0000.02350.04687

2000

1000

0
-5

-4

-3

-2

-1

0

1

2

3

4

5

mwave frequency (446.810 MHz + kHz)

Figure 5.12: Rf lineshape for the |1, −1i → |2, 0i transition at around 2000 a0 . The fitting function
is shown in the inset table. The r.m.s. width of the Gaussian pulse is fixed to the set value of 67µs in temporal domain.

136

(a)

(b)

1,050 1,300

460nK trap depth
760nK trap depth

2.5x1010

1E-10
2.0x1010

1.5x1010

1E-11

80
70

1E-12

1E-13
1,000

460nK trap depth
760nK trap depth
theory10,000

60
50
1,000

10,000Figure 5.13: Vary the trap depth in three-body recombination experiment at large a. (a) Two
different trap depths give essentially the same result on L2 . The vertical lines indicate the scattering
lengths at which Eb is equal to the trap depth. (b) The time-averaged conditions of the samples
generated in two trap depths are controlled to be about the same.

137

(a)

no FM
FM at 12.2kHz, sideband amp ~ carrier amp
FM at 12.2kHz, FM signal is modulated w/ white noise

3.0x104

0.25

0.20

Tz (mK)

atom number

2.5x104

2.0x104

0.15
1.5x104

1.0x104

0.10
0.0

0.1

0.2

0.3

0.4

0.5

0.0

0.1

0.2

thold (s)

0.3

0.4

0.5

0.3

0.4

0.5

thold (s)

(b)
0.35

no amplitude modulation
10% AM on V beam (0.3W setpoint)
7% AM on H2 beam (4W setpoint)

4

0.30

2.0x104

Tz (mK)

atom number

2.5x101.5x104

0.25

0.20

0.15
1.0x104

0.10
0.0

0.1

0.2

0.3

thold (s)

0.4

0.5

0.0

0.1

0.2

thold (s)Figure 5.14: Frequency (a) and intensity (b) modulation on the ODT beams, initial temperature
120 nK, science field 3000 a0 . In (a), we modulate the frequency of the ODT at a single frequency of
Eb = 12.2 kHz, as well as in a broad-band fashion respectively. In (b), we modulation the intensityof two ODT beams at frequency of Eb = 12.2 kHz up to a modulation index of 10%. Black circles
represent loss data without any modulations. None of these modulations enhances the number
losses or heating rates.

138

(a)

(b)

Figure 5.15: Comparison of three-body recombination data with the predictions by a toy modelthat includes avalanche effects. The time-averaged temperature is 200 nK and the science field
is 3000 a0 . The black line represents three-body recombination assuming that the shallow dimers
all leave the trap immediately. The blues lines represent atom population (a) and shallow dimer
population (b) assuming that the shallow dimers all remain in the trap after being created.

139
dimers in the trap with a finite depth.

Chapter 6

Summary and Outlook

6.1

Fifty Years of Efimov Physics
Recent decades have witnessed tremendous progress in the field of few-body physics. This

burgeoning field has been meanwhile motivating both technical and computational developments
in ultracold atoms. We dedicate our work to improving the experimental precision on two- andthree-body observables as well as to characterizing the systematic effects that have been previously
neglected. Our goal is to set a technical standard for the studies of few-body physics with ultracold
atoms. Comparison of our work with generic few-body models can be made along two independent
axes of magnetic field and collision energy respectively. As regards physics, we emphasize theglobal picture of Efimov’s scenario by connecting the three-body observables in the a < 0 and
a > 0 regions together. Within a very limited universal window as illustrated in Fig. 3.13, we areable to observe signatures of Efimovs’s original zero-range results manifested in 39 K in the vicinityof an intermediate-strength Feshbach resonance.
All the measured three-body observables in our work are summarized in Table 6.1 and com(1)

(0)

(1)

(0)

pared to theoretical models of various layers of complexity. The two ratios a∗ /|a− | and a∗ /(apwhich lie inside the universal window agree with the zero-range results within the experimental
(0)

(1)

(0)

uncertainties. The width parameters associated with each location – η− , η∗ and η+ – are foundto be consistent with each other. This discovery for the first time verifies another hypothesis inzero-range theory which is the continuity of a three-body short-range parameter across the pole of
the Feshbach resonance [96]. As discussed in 3.2.3, these zero-range universal ratios come from a

141
(0)

combined effect of 1) an enlarged universal window and 2) an especially large value of |a− |/rrdW
that is unique to 39 K. The former effect is a consequence of a small effective range parameter re ina two-channel effective field theory [15]. As to the absolute locations of the three-body observables,
we demonstrate in 39 K that two-body van der Waals potentials alone can describe the three-body
physics in the intermediate coupling strength regime if hyperfine interactions are properly incorporated.

−626(9)−846(19)−908(11)–0.21(1)0.25(1)–0.28(2)884(14)–246(6)

91
0.20(2)–

Observables for a > 0
(1)
(0)
(0)
η∗
a+ /a0η+

215(1)a∗ /a0876(28)–(0)ap

0.97(2)0.271(7)0.96(3)

Efimov’s ratios
(0)
(0)
(0)
(0)
a+ /|a− | ap /|a− |1.065220.209914
1
0.34
0.14
–
(1)
(0)
a∗ /|a− |Table 6.1: Summary of three-body observables characterizing the Efimov spectrum in 39 K across the Feshbach resonance at 33.58 G. The
locations of the three-body observables and their ratios are compared to theoretical models of various layers of complexity.

Model
zero-range EFT [53]sc vdW [97]mc vdW [18] [69]Experimental

Observables for a < 0(0)(0)
a− /a0
η−

142

143

6.2

Future Directions
If one contrasts the summary figures Fig. 3.12 and Fig. 3.14, the three-body observables are

noticeably not as well sampled for a > 0 as in the a < 0 region. This relatively less-explored
region opens up a lot of opportunities for future few-body physics studies in both homonuclear andheteronuclear systems. Interplay between atoms and Feshbach molecules as well as between two
Feshbach molecules in this region present technical challenges in quantum state preparation andcontrol of noise.
Feshbach resonances with intermediate and small coupling strengths provide another interesting direction. In [18], we demonstrate a tri-atomic resonance in 39 K that deviates from van der
Waals universality. This observed trend seems to agree with a two-channel field theory in which
three-body observables are predicted to be linked with sres (or equivalently the two-body scattering effective range re ) [15, 24]. Nevertheless, two independent groups have reported that in 7 Li,(0)which also features an intermediate sres , the tri-atomic resonances |a− | ' 9 rvdW are in perfectagreement with van der Waals universality [98, 99]. This brings up the question of what otherelements are contributing to three-body physics in addition to the effective range. A three-bodymulti-channel calculation may shed light on this interesting case. For sufficiently narrow Feshbachresonances with sres  1, it has been proposed that a novel kind of universality exists in which
the three-body observables are related to the Feshbach resonance range parameter r∗ [100]. Suchnarrow resonances can be found in species such as 23 Na and 87 Rb. It will also be enlightening tosee if the multi-channel calculations based on van der Waals potentials can still do a good job inthis regime.
Last but not least, ongoing cold atom projects in microgravity environment provide a promising and novel tool to unveil the few-body dynamics near the pole of the Feshbach resonance. Explorations on this regime have been inhibited by the finite density and temperature effects for theexperiments on earth. And a large magnetic field gradient is usually needed to levitate the cloud ina shallow trap [101]. As I am concluding my thesis, tremendous progress has been made in gener-

144
ating atom samples at record-breaking low temperatures in space station. There is vision that one
could make use of these low-temperature and low-density samples to make accurate observations(1)of the first excited Efimov tri-atomic resonance and locate a− .

Bibliography[1] Michael E. Fisher. Renormalization group theory: Its basis and formulation in statisticalphysics. Rev. Mod. Phys., 70:653–681, Apr 1998.
[2] Kenneth G. Wilson. The renormalization group and critical phenomena. Rev. Mod. Phys.,
55:583–600, Jul 1983.
[3] Jun John Sakurai. Modern quantum mechanics; rev. ed. Addison-Wesley, Reading, MA,
1994.
[4] Private communications with Paul S. Julienne.
[5] George A. Baker Jr. The MBX challenge competition: a neutron matter model. Int. J. Mod.
Phys. B, 15:1314–1320, 1999.
[6] Tin-Lun Ho. Universal thermodynamics of degenerate quantum gases in the unitarity limit.
Phys. Rev. Lett., 92:090402, Mar 2004.
[7] P. Makotyn, C. E. Klauss, D. L. Goldberger, E. A. Cornell, and D. S. Jin. Universal dynamics
of a degenerate unitary bose gas. Nature Physics, 10(2):116–119, 2014.
[8] John Corson. Dynamics of Quenched Ultracold Quantum Gases. PhD thesis, Boulder, CO,
2016-07 2016.
[9] V. Efimov. Low-energy properties of three resonantly interacting particles. Sov. J. Nucl.
Phys, 29(4):1058–1069, 1979.
[10] L. H. Thomas. The interaction between a neutron and a proton and the structure of H3 .
Phys. Rev., 47:903–909, Jun 1935.
[11] José P D’Incao, Chris H Greene, and B D Esry. The short-range three-body phase and otherissues impacting the observation of efimov physics in ultracold quantum gases. Journal of
Physics B: Atomic, Molecular and Optical Physics, 42(4):044016, feb 2009.
[12] M. Berninger, A. Zenesini, B. Huang, W. Harm, H.-C. Nägerl, F. Ferlaino, R. Grimm, P. S.
Julienne, and J. M. Hutson. Universality of the three-body parameter for Efimov states inultracold cesium. Phys. Rev. Lett., 107:120401, Sep 2011.
[13] Jia Wang, J. P. D’Incao, B. D. Esry, and Chris H. Greene. Origin of the three-body parameter
universality in Efimov physics. Phys. Rev. Lett., 108:263001, Jun 2012.

146[14] C. Chin.
Universal scaling of Efimov resonance positions in cold atom systems.
arXiv:1111.1484, 2011.
[15] R. Schmidt, S.P. Rath, and W. Zwerger. Efimov physics beyond universality. The European
Physical Journal B, 85(11):386, Nov 2012.
[16] José P. D’Incao, J. Wang, Brett Esry, and Chris Greene. The universality of the Efimov
three-body parameter. Few-Body Systems, 54, Nov 2012.
[17] Yujun Wang and Paul S. Julienne. Universal van der Waals physics for three cold atoms near
Feshbach resonances. Nature Physics, 10:768, Aug 2014.
[18] Roman Chapurin, Xin Xie, Michael J. Van de Graaff, Jared S. Popowski, José P. D’Incao,
Paul S. Julienne, Jun Ye, and Eric A. Cornell. Precision test of the limits to universality infew-body physics. Phys. Rev. Lett., 123:233402, Dec 2019.
[19] B. M. Axilrod and E. Teller. Interaction of the van der Waals type between three atoms. J.
Chem. Phys., 11:299, Jun 1943.
[20] William D. Phillips. Nobel lecture: Laser cooling and trapping of neutral atoms. Rev. Mod.
Phys., 70:721–741, Jul 1998.
[21] Harold J. Metcalf and Peter van der Straten. Laser Cooling and Trapping. Springer-Verlag,
New York, 1999.
[22] N. F. Mott and Harrie Stewart Wilson Massey. The theory of atomic collisions. 1949.
[23] John L. Bohn, James P. Burke, Chris H. Greene, H. Wang, P. L. Gould, and W. C. Stwalley.
Collisional properties of ultracold potassium: Consequences for degenerate bose and fermi
gases. Phys. Rev. A, 59:3660–3664, May 1999.
[24] Christian Langmack, Richard Schmidt, and Wilhelm Zwerger. Efimov states near a Feshbachresonance and the limits of van der Waals universality at finite background scattering length.
Phys. Rev. A, 97:033623, Mar 2018.
[25] E. Nielsen, D.V. Fedorov, A.S. Jensen, and E. Garrido. The three-body problem with shortrange interactions. Physics Reports, 347(5):373 – 459, 2001.
[26] Rabin Paudel. Probing local quantities in a strongly interacting Fermi gas and theconstruction of an ultracold Fermi gas apparatus. PhD thesis, Boulder, CO, 2017-05 2017.
[27] Roman Chapurin. Precise Measurements of Few-Body Physics in Ultracold 39 K Bose Gas.
PhD thesis, Boulder, 2019-06 2019.
[28] G. Salomon, L. Fouché, P. Wang, A. Aspect, P. Bouyer, and T. Bourdel. Gray-molasses
cooling of 39 K to a high phase-space density. EPL (Europhysics Letters), 104(6):63002, Dec
2013.
[29] G. Salomon, L. Fouché, S. Lepoutre, A. Aspect, and T. Bourdel. All-optical cooling of 39 K
to Bose-Einstein condensation. Phys. Rev. A, 90:033405, Sep 2014.

147[30] Dipankar Nath, R Kollengode Easwaran, G. Rajalakshmi, and C. S. Unnikrishnan. Quantuminterference-enhanced deep sub-Doppler cooling of 39 K atoms in gray molasses. Phys. Rev.
A, 88:053407, Nov 2013.
[31] M. Gröbner, P. Weinmann, F. Meinert, K. Lauber, E. Kirilov, and H.-C. Nägerl. A newquantum gas apparatus for ultracold mixtures of K and Cs and KCs ground-state molecules.
Journal of Modern Optics, 63(18):1829–1839, 2016.
[32] Jacob Covey. Enhanced optical and electric manipulation of a quantum gas of KRb molecules.
PhD thesis, Boulder, 2017-07 2017.
[33] L. De Sarlo, P. Maioli, G. Barontini, J. Catani, F. Minardi, and M. Inguscio. Collisional
properties of sympathetically cooled 39 K. Phys. Rev. A, 75:022715, Feb 2007.
[34] Robert L. D. Campbell, Robert P. Smith, Naaman Tammuz, Scott Beattie, Stuart Moulder,
and Zoran Hadzibabic. Efficient production of large 39 K Bose-Einstein condensates. Phys.
Rev. A, 82:063611, Dec 2010.
[35] S. J. M. Kuppens, K. L. Corwin, K. W. Miller, T. E. Chupp, and C. E. Wieman. Loading
an optical dipole trap. Phys. Rev. A, 62:013406, Jun 2000.
[36] Nicolas Schlosser, Georges Reymond, Igor Protsenko, and Philippe Grangier. Sub-poissonianloading of single atoms in a microscopic dipole trap. Nature, 411(6841):1024–1027, 2001.
[37] T. Grünzweig, A. Hilliard, M. McGovern, and M. F. Andersen. Near-deterministic preparationof a single atom in an optical microtrap. Nature Physics, 6(12):951–954, 2010.
[38] Bruce H. Karnopp. Introduction to Dynamics. Addison-Wesley, Reading, MA, 1974.
[39] L. Richter, H. Mandelberg, M. Kruger, and P. McGrath. Linewidth determination fromself-heterodyne measurements with subcoherence delay times. IEEE Journal of Quantum
Electronics, 22(11):2070–2074, 1986.
[40] Hanne Ludvigsen, Mika Tossavainen, and Matti Kaivola. Laser linewidth measurements using
self-homodyne detection with short delay. Optics Communications, 155(1):180 – 186, 1998.
[41] L. B. Mercer. 1/f frequency noise effects on self-heterodyne linewidth measurements. Journalof Lightwave Technology, 9(4):485–493, 1991.
[42] G. Reinaudi, T. Lahaye, Z. Wang, and D. Guéry-Odelin. Strong saturation absorption imaging
of dense clouds of ultracold atoms. Opt. Lett., 32(21):3143–3145, Nov 2007.
[43] D Genkina, L M Aycock, B K Stuhl, H-I Lu, R A Williams, and I B Spielman. Feshbachenhanceds-wave scattering of fermions: direct observation with optimized absorption imaging.
New Journal of Physics, 18(1):013001, Dec 2015.
[44] Joseph W Goodman. Introduction to Fourier optics. Roberts & Co. Publishers, Englewood,CO, 2005.
[45] Chen-Lung Hung and Cheng Chin. In Situ Imaging of Atomic Quantum Gases, chapter
Chapter 6, pages 101–120.

148[46] Max Born and Emil Wolf. Principles of Optics: Electromagnetic Theory of Propagation,
Interference and Diffraction of Light (7th Edition). Cambridge University Press, 7th edition,
1999.
[47] Ryuta Yamamoto, Jun Kobayashi, Kohei Kato, Takuma Kuno, Yuto Sakura, and Yoshiro
Takahashi. Site-resolved imaging of single atoms with a Faraday quantum gas microscope.
Phys. Rev. A, 96:033610, Sep 2017.
[48] S. Kadlecek, J. Sebby, R. Newell, and T. G. Walker. Nondestructive spatial heterodyne
imaging of cold atoms. Opt. Lett., 26(3):137–139, Feb 2001.
[49] Erin Marshall Seroka, Ana Valdés Curiel, Dimitrios Trypogeorgos, Nathan Lundblad, and
Ian B. Spielman. Repeated measurements with minimally destructive partial-transfer absorption imaging. Opt. Express, 27(25):36611–36624, Dec 2019.
[50] Eric Braaten and H.-W. Hammer. Universality in few-body systems with large scattering
length. Physics Reports, 428(5):259 – 390, 2006.
[51] Francesca Ferlaino and Rudolf Grimm. Forty years of efimov physics: How a bizarre prediction
turned into a hot topic. Physics, 3(9), 2010.
[52] Pascal Naidon and Shimpei Endo. Efimov physics: a review. Reports on Progress in Physics,
80(5):056001, Mar 2017.
[53] José P D’Incao. Few-body physics in resonantly interacting ultracold quantum gases. Journalof Physics B: Atomic, Molecular and Optical Physics, 51(4):043001, Jan 2018.
[54] Cheng Chin, Rudolf Grimm, Paul Julienne, and Eite Tiesinga. Feshbach resonances in ultracold gases. Rev. Mod. Phys., 82(2):1225, 2010.
[55] Richard Schmidt. From few- to many-body physics with ultracold atoms. PhD thesis,
München, 2013.
[56] Private communications with Wilhelm Zwerger.
[57] Eugene P. Wigner. Lower limit for the energy derivative of the scattering phase shift. Phys.
Rev., 98:145–147, Apr 1955.
[58] H.-W. Hammer and Dean Lee. Causality and the effective range expansion. Annals of Physics,
325(10):2212 – 2233, 2010.
[59] V. N. Efimov. Weakly bound states of three resonantly interacting particles. Yadern. Fiz.,
12:1080–91, Nov 1970.
[60] V. Efimov. Energy levels arising from resonant two-body forces in a three-body system.
Physics Letters B, 33(8):563 – 564, 1970.
[61] Private communications with Pascal Naidon.
[62] B. D. Esry, Chris H. Greene, and James P. Burke. Recombination of three atoms in theultracold limit. Phys. Rev. Lett., 83:1751–1754, Aug 1999.

149[63] T. Kraemer, M. Mark, P. Waldburger, J. G. Danzl, C. Chin, B. Engeser, A. D. Lange,
K. Pilch, A. Jaakkola, H.-C. Nägerl, and R. Grimm. Evidence for Efimov quantum states inan ultracold gas of caesium atoms. Nature, 440:315, Mar 2006.
[64] P. F. Bedaque, H.-W. Hammer, and U. van Kolck. The three-boson system with short-range
interactions. Nuclear Physics A, 646(4):444 – 466, 1999.
[65] Ludwig W. Bruch and Katuro Sawada. Inequality relating the ground-state energies of two
and three bosons. Phys. Rev. Lett., 30:25–27, Jan 1973.
[66] K. Helfrich, H.-W. Hammer, and D. S. Petrov. Three-body problem in heteronuclear mixtureswith resonant interspecies interaction. Phys. Rev. A, 81:042715, Apr 2010.
[67] Private communications with Richard Schmidt.
[68] Frederick H. Mies, Carl J. Williams, Paul S. Julienne, and Morris Krauss. Estimating bounds
on collisional relaxation rates of spin-polarized 87 Rb atoms at ultracold temperatures. Journalof research of the National Institute of Standards and Technology, 101(4):521–535, 1996.
[69] Xin Xie et al, in preparation.
[70] Thorsten Köhler, Eite Tiesinga, and Paul S. Julienne. Spontaneous dissociation of long-range
Feshbach molecules. Phys. Rev. Lett., 94:020402, Jan 2005.
[71] C. Strauss, T. Takekoshi, F. Lang, K. Winkler, R. Grimm, J. Hecker Denschlag, and E. Tie87 Rb . Phys.
mann. Hyperfine, rotational, and vibrational structure of the a3 Σ+2
u state ofRev. A, 82:052514, Nov 2010.
[72] S. Knoop, F. Ferlaino, M. Mark, M. Berninger, H. Schöbel, H.-C. Nägerl, and R. Grimm.
Observation of an Efimov-like trimer resonance in ultracold atom-dimer scattering. Nature
Physics, 5:227, Feb 2009.
[73] A. Zenesini, B. Huang, M. Berninger, H.-C. Nägerl, F. Ferlaino, and R. Grimm. Resonant
atom-dimer collisions in cesium: Testing universality at positive scattering lengths. Phys.
Rev. A, 90:022704, Aug 2014.
[74] K. Helfrich and H.-W. Hammer. Resonant atom-dimer relaxation in ultracold atoms. EPL(Europhysics Letters), 86(5):53003, Jun 2009.
[75] K. Kato, Yujun Wang, J. Kobayashi, P. S. Julienne, and S. Inouye. Isotopic shift of atomdimer efimov resonances in K-Rb mixtures: Critical effect of multichannel Feshbach physics.
Phys. Rev. Lett., 118:163401, Apr 2017.
[76] Private communications with Dmitry Petrov.
[77] Benno Rem, Andrew Grier, Igor Ferrier-Barbut, U Eismann, Tim Langen, N Navon, Lev
Khaykovich, F Werner, D Petrov, F Chevy, and C. Salomon. Lifetime of the bose gas withresonant interactions. Phys. Rev. Lett., 110:163202, Apr 2013.
[78] Eric Braaten and H.-W. Hammer. Enhanced dimer relaxation in an atomic and molecular
Bose-Einstein condensate. Phys. Rev. A, 70:042706, Oct 2004.

150[79] Cheng Chin, Andrew Kerman, Vladan Vuletic, and Steven Chu. Sensitive detection of cold
cesium molecules formed on Feshbach resonances. Phys. Rev. Lett., 90:033201, Feb 2003.
[80] Cheng Chin and Rudolf Grimm. Thermal equilibrium and efficient evaporation of an ultracold
atom-molecule mixture. Physical Review A, 69:33612, Mar 2004.
[81] Jens Herbig, Tobias Kraemer, Michael Mark, Tino Weber, Cheng Chin, Hanns-ChristophNägerl, and Rudolf Grimm. Preparation of a pure molecular quantum gas. Science (NewYork, N.Y.), 301:1510–3, Oct 2003.
[82] José P. D’Incao, Javier von Stecher, and Chris Greene. Universal four-boson states in ultracold
molecular gases: Resonant Effects in dimer-dimer collisions. Phys. Rev. Lett., 103:033004,
Jul 2009.
[83] Sébastien Laurent, Xavier Leyronas, and Frédéric Chevy. Momentum distribution of a dilute
unitary bose gas with three-body losses. Phys. Rev. Lett., 113:220601, Nov 2014.
[84] P. O. Fedichev, M. W. Reynolds, and G. V. Shlyapnikov. Three-body recombination ofultracold atoms to a weakly bound s level. Phys. Rev. Lett., 77:2921–2924, Sep 1996.
[85] Benno Rem. The Road to the Unitary Bose Gas. PhD thesis, Paris, 2013-12.
[86] Eric Braaten, H.-W. Hammer, Daekyoung Kang, and Lucas Platter. Three-body recombination of identical bosons with a large positive scattering length at nonzero temperature. Phys.
Rev. A, 78:043605, Oct 2008.
[87] Christian Langmack, D. Hudson Smith, and Eric Braaten. Avalanche mechanism for theenhanced loss of ultracold atoms. Phys. Rev. A, 87:023620, Feb 2013.
[88] Private communications with Daekyoung Kang.
[89] J. P. D’Incao, H. Suno, and B. D. Esry. Limits on universality in ultracold three-bosonrecombination. Phys. Rev. Lett., 93:123201, Sep 2004.
[90] M. Zaccanti, B. Deissler, C. D’Errico, M. Fattori, M. Jona-Lasinio, S. Müller, G. Roati,
M. Inguscio, and G. Modugno. Observation of an Efimov spectrum in an atomic system.
Nature Physics, 5:586, Jul 2009.
[91] Scott E. Pollack, Daniel Dries, and Randall G. Hulet. Universality in three- and four-bodybound states of ultracold atoms. Science, 326(5960):1683–1685, 2009.
[92] Olga Machtey, Zav Shotan, Noam Gross, and Lev Khaykovich. Association of Efimov trimersfrom a three-atom continuum. Phys. Rev. Lett., 108:210406, May 2012.
[93] Olga Machtey, David A. Kessler, and Lev Khaykovich. Universal dimer in a collisionally opaque medium: Experimental observables and Efimov resonances. Phys. Rev. Lett.,
108:130403, Mar 2012.
[94] Christian Langmack, D. Hudson Smith, and Eric Braaten. Avalanche mechanism for atom
loss near an atom-dimer Efimov resonance. Phys. Rev. A, 86:022718, Aug 2012.
[95] Christian Langmack, D. Hudson Smith, and Eric Braaten. Atom loss resonances in a BoseEinstein condensate. Phys. Rev. Lett., 111:023003, Jul 2013.

151[96] Eric Braaten and H Hammer. Three-body recombination into deep bound states in a bosegas with large scattering length. Phys. Rev. Lett., 87:160407, Oct 2001.
[97] Paul M. A. Mestrom, Jia Wang, Chris H. Greene, and José P. D’Incao. Efimov–van der Waals
universality for ultracold atoms with positive scattering lengths. Phys. Rev. A, 95:032707,
Mar 2017.
[98] P. Dyke, S. E. Pollack, and R. G. Hulet. Finite-range corrections near a Feshbach resonanceand their role in the Efimov effect. Phys. Rev. A, 88:023625, Aug 2013.
[99] Noam Gross, Zav Shotan, Servaas Kokkelmans, and Lev Khaykovich. Nuclear-spinindependent short-range three-body physics in ultracold atoms. Phys. Rev. Lett., 105:103203,
Sep 2010.
[100] D. S. Petrov. Three-boson problem near a narrow Feshbach resonance. Phys. Rev. Lett.,
93:143201, Sep 2004.
[101] Bo Huang, Leonid A. Sidorenkov, Rudolf Grimm, and Jeremy M. Hutson. Observation of thesecond triatomic resonance in Efimov’s scenario. Phys. Rev. Lett., 112:190401, May 2014.
[102] T. Mukaiyama, J. R. Abo-Shaeer, K. Xu, J. K. Chin, and W. Ketterle. Dissociation anddecay of ultracold sodium molecules. Phys. Rev. Lett., 92:180402, May 2004.
[103] Peter Staanum, Stephan D. Kraft, Jörg Lange, Roland Wester, and Matthias Weidemüller.
Experimental investigation of ultracold atom-molecule collisions. Phys. Rev. Lett., 96:023201,
Jan 2006.
[104] N. Zahzam, T. Vogt, M. Mudrich, D. Comparat, and P. Pillet. Atom-molecule collisions inan optically trapped gas. Phys. Rev. Lett., 96:023202, Jan 2006.
[105] N. Syassen, T. Volz, S. Teichmann, S. Dürr, and G. Rempe. Collisional decay of 87 Rb
Feshbach molecules at 1005.8 G. Phys. Rev. A, 74:062706, Dec 2006.
[106] F. Ferlaino, S. Knoop, M. Mark, M. Berninger, H. Schöbel, H.-C. Nägerl, and R. Grimm.
Collisions between tunable halo dimers: Exploring an elementary four-body process withidentical bosons. Phys. Rev. Lett., 101:023201, Jul 2008.
[107] Shuta Nakajima, Munekazu Horikoshi, Takashi Mukaiyama, Pascal Naidon, and Masahito
Ueda. Nonuniversal Efimov atom-dimer resonances in a three-component mixture of 6 Li.
Phys. Rev. Lett., 105:023201, Jul 2010.
[108] T. Lompe, T. B. Ottenstein, F. Serwane, K. Viering, A. N. Wenz, G. Zürn, and S. Jochim.
Atom-dimer scattering in a three-component fermi gas. Phys. Rev. Lett., 105:103201, Sep
2010.
[109] Ruth S. Bloom, Ming-Guang Hu, Tyler D. Cumby, and Deborah S. Jin. Tests of universal
three-body physics in an ultracold bose-fermi mixture. Phys. Rev. Lett., 111:105301, Sep
2013.
[110] Ming-Guang Hu, Ruth S. Bloom, Deborah S. Jin, and Jonathan M. Goldwin. Avalanchemechanism loss at an atom-molecule efimov resonance. Phys. Rev. A, 90:013619, Jul 2014.

152[111] S Knoop, F Ferlaino, M Berninger, M Mark, H-C Nägerl, and R Grimm. Observation of an
Efimov resonance in an ultracold mixture of atoms and weakly bound dimers. Journal of
Physics: Conference Series, 194(1):012064, Nov 2009.
[112] Eric Braaten and H.-W. Hammer. Resonant dimer relaxation in cold atoms with a large
scattering length. Phys. Rev. A, 75:052710, May 2007.
[113] S Thompson, E Hodby, and C. Wieman. Spontaneous dissociation of 85 Rb feshbach molecules.
Phys. Rev. Lett., 94:020401, Feb 2005.
[114] K. Helfrich and H. Hammer. On the Efimov effect in higher partial waves. Journal of Physics
B: Atomic, Molecular and Optical Physics, 44, Jul 2011.
[115] Eric Braaten and H.-W. Hammer. Efimov physics in cold atoms. Annals of Physics, 322(1):120– 163, 2007. January Special Issue 2007.
[116] Tino Weber, Jens Herbig, Michael Mark, Hanns-Christoph Nägerl, and Rudolf Grimm. Threebody recombination at large scattering lengths in an ultracold atomic gas. Phys. Rev. Lett.,
91:123201, Sep 2003.
[117] Chris H. Greene, B.D. Esry, and H. Suno. A revised formula for 3-body recombination thatcannot exceed the unitarity limit. Nuclear Physics A, 737:119 – 124, 2004.
[118] Sanjukta Roy, Manuele Landini, Andreas Trenkwalder, Giulia Semeghini, Giacomo Spagnolli,
Andrea Simoni, Marco Fattori, Massimo Inguscio, and Giovanni Modugno. Test of the universality of the three-body Efimov parameter at narrow Feshbach resonances. Phys. Rev.
Lett., 111(5):053202, 2013.
[119] Noam Gross, Zav Shotan, Servaas Kokkelmans, and Lev Khaykovich. Observation of universality in ultracold 7 Li three-body recombination. Phys. Rev. Lett., 103:163202, Oct 2009.
[120] R. J. Wild, P. Makotyn, J. M. Pino, E. A. Cornell, and D. S. Jin. Measurements of Tan’scontact in an atomic Bose-Einstein condensate. Phys. Rev. Lett., 108:145305, Apr 2012.
[121] L. J. Wacker, N. B. Jørgensen, K. T. Skalmstang, M. G. Skou, A. G. Volosniev, and J. J.
Arlt. Temperature dependence of an Efimov resonance in 39 K. Phys. Rev. A, 98:052706, Nov
2018.
[122] G. Barontini, C. Weber, F. Rabatti, J. Catani, G. Thalhammer, M. Inguscio, and F. Minardi.
Observation of heteronuclear atomic Efimov resonances. Phys. Rev. Lett., 103:043201, Jul
2009.
[123] Jacob Johansen, BJ DeSalvo, Krutik Patel, and Cheng Chin. Testing universality of Efimov
physics across broad and narrow Feshbach resonances. Nature Physics, 13(8):731, 2017.
[124] Pascal Naidon, Shimpei Endo, and Masahito Ueda. Microscopic origin and universality classesof the Efimov three-body parameter. Phys. Rev. Lett., 112:105301, Mar 2014.
[125] P. Giannakeas and Chris Greene. Van der Waals universality in homonuclear atom-dimer
elastic collisions. Few-Body Systems, 58, Dec 2016.
[126] Mark D. Lee, Thorsten Köhler, and Paul S. Julienne. Excited Thomas-Efimov levels inultracold gases. Phys. Rev. A, 76:012720, Jul 2007.

153[127] Alexander O. Gogolin, Christophe Mora, and Reinhold Egger. Analytical solution of thebosonic three-body problem. Phys. Rev. Lett., 100:140404, Apr 2008.

Appendix A

Algorithms: Absorption Imaging Simulation

A.1

Absorption Imaging Simulation%% constants
h = 6.6260755e-34; hbar = h/(2*pi);c = 299792458; amu = 1.6605402e-27;
um = 1e-6; nm = 1e-9; MHz = 1e6; us = 1e-6;
Gamma = (2*pi)*6.035*MHz; lambda = 766.701*nm; %D2 line
m39 = 39*amu; sigma0 = 3*lambda^2/(2*pi);%% experiment parameters% intensities are in unit of Isat, Isat is in unit of count% delta is in unit of Gamma
I0_list = logspace( log10(1e-2),log10(3),15 );rho_pk = 2.3e17; size = 10*um; pulse = 40*us;%% stationary atom model
kr = 2*pi/lambda; vr = hbar*kr/m39;
N_photons = Gamma/(2*kr*vr); %=173 to shift by half a linewidthdz = size/1e3; dt = pulse/1e3;t_step = round(pulse/dt);

155
rho = rho_pk* exp( - [-size*5:dz:size*5].^2 / (2*size^2) );n_step = length(rho);
Int = zeros(n_step, t_step); delta = zeros(n_step, t_step);
OD_simu1_list = [];for i = 1:length(I0_list)
I0 = I0_list(i); Hf = 0;
Int(1,:) = I0; delta(:,1) = 0;
for t = 1:t_step
for n = 1:n_step
A = sigma0*rho(n)*dz;
B = kr*vr*dt/(2*sigma0*rho(n)*dz);
Int(n+1, t) = Int(n, t) - A* Int(n, t)./( 1 + 4*delta(n, t).^2 + Int(n, t) );delta(n, t+1) = delta(n, t) + B* ( Int(n, t) - Int(n+1, t) );end
Hf = Hf + Int(n_step+1, t)*dt;
end
OD_simu1 = -log( abs(Hf-I0*dt) / abs(I0*(t_step-1)*dt) )OD_simu1_list = [ OD_simu1_list OD_simu1 ];end%% moving atom modeli = 1; I0 = I0_list(i);kr = 2*pi/lambda; vr = hbar*kr/m39;
dz = size/1e2;
rho = rho_pk* exp( - [-size*5:dz:size*5].^2 / (2*size^2) );156
n_step = length(rho);z_frt = []; z_mid = []; z_bck = [];OD_simu2_list = [];for k = 1:length(pulse_list)
pulse = pulse_list(k)*us;
dt = pulse/1e3;t_step = round(pulse/dt);
Int = zeros(n_step, t_step); delta = zeros(n_step, t_step);z = zeros(n_step, 1);
z_cumu = [-size*2:dz:size*8]’;
Hf = 0; Int(1,:) = I0; delta(:,1) = 0;
for t = 1:t_step
for n = 1:n_step
A = sigma0*rho(n)*dz;
B = kr*vr*dt/(2*A);
Int(n+1, t) = Int(n, t) - A* Int(n, t)./( 1 + 4*delta(n, t).^2 + Int(n, t) );delta(n, t+1) = delta(n, t) + B* ( Int(n, t) - Int(n+1, t) );z(n) = dt*Gamma* delta(n, t+1)/kr/2;
end
z_cumu = z_cumu + z;
[z_cumu, z_order] = sort(z_cumu); rho = rho(z_order);Hf = Hf + Int(n_step+1, t)*dt;
end
OD_simu2 = -log( (Hf-I0*dt) / (I0*(t_step-1)*dt) );OD_simu2_list = [ OD_simu2_list OD_simu2 ];z_frt = [z_frt z_cumu(ceil(end/5))];

157
z_mid = [z_mid z_cumu(ceil(end/2))];z_bck = [z_bck z_cumu(ceil(end*4/5))];
end

A.2

OD Correction as Function of I/Isat%% experiment input parameters% intensities are in unit of Isat, Isat is in unit of count% delta is in unit of Gamma
I0_list = 0.1;
rho_pk = linspace( 2.3e17/32,2.3e17,50 );size = 10*um; pulse = 40*us;%% stationary atom model
kr = 2*pi/lambda; vr = hbar*kr/m39;
N_photons = Gamma/(2*kr*vr);dz = size/1e3; dt = pulse/1e3;t_step = round(pulse/dt);
OD_corr_list = []; OD_phys_list = [];for j = 1:length(rho_pk)
rho = rho_pk(j)* exp( - [-size*5:dz:size*5].^2 / (2*size^2) );n_step = length(rho);
Int = zeros(n_step, t_step); delta = zeros(n_step, t_step);
I0 = I0_list; Hf = 0;
Int(1,:) = I0; delta(:,1) = 0;
for t = 1:t_step
for n = 1:n_step

158
A = sigma0*rho(n)*dz;
B = kr*vr*dt/(2*sigma0*rho(n)*dz);
Int(n+1, t) = Int(n, t) - A* Int(n, t)./( 1 + 4*delta(n, t).^2 + Int(n, t) );delta(n, t+1) = delta(n, t) + B* ( Int(n, t) - Int(n+1, t) );end
Hf = Hf + Int(n_step+1, t)*dt;
end
OD_simu1 = -log( (Hf-I0*dt) / (I0*(t_step-1)*dt) )OD_phys = sigma0*sum(rho(:))*dz;
OD_correction = OD_simu1./OD_phys;
OD_corr_list = [OD_corr_list OD_correction];
OD_phys_list = [OD_phys_list OD_phys];
end

A.3

Fit Function for Isatx0 = [ 0.9 1237*2 ]; %initial guesses
pulse = 200;%% method I% modified Beer’s law[x1,resnorm,residual,exitflag,output,Lmda,Jmat] = ...
lsqcurvefit(@(x1,Int_list) OD_analytic(x1,Int_list), x0, Int_list, OD_list);
r21 = 1 - resnorm/sum((OD_list-mean(OD_list)).^2);
ci1 = nlparci(x1,residual,’jacobian’,Jmat); %confidence interval
fit_err1 = abs(ci1(:,2) - ci1(:,1))/1.96/2;

function OD_meas = OD_analytic(x, Int_list)%x(1) = OD_phys; x(2) = Isat159
OD_meas = x(1) - Int_list/x(2) + lambertw( exp(Int_list/x(2) - x(1)) ...
.*Int_list/x(2));
end%% method II% correct for recoil induced detuning[x2,resnorm,residual,exitflag,output] = ...
lsqcurvefit(@(x2,Int_list) OD_simu1(x2,Int_list,pulse), x0, Int_list, OD_list);
r22 = 1 - resnorm/sum((OD_list-mean(OD_list)).^2);ci2 = nlparci(x1,residual,’jacobian’,Jmat); %confidence interval
fit_err2 = abs(ci2(:,2) - ci2(:,1))/1.96/2;

Appendix B

Algorithms: Analysis of Optical Aberrations

B.1

PSF and MTF for an Ideal Pupil%% experiment parameter
diameter0 = 26.15e-3; radius0 = diameter0/2;
f1 = 31.03e-3; f2 = 1500e-3;
mag = f2/f1; px_img = 16e-6; px_obj = px_img/mag;
NA = 0.384; lambda = 767e-9; Abbe_f = 2*NA/lambda;%% put the pupil on a griddim = 5e3; pupil = zeros(dim,dim);[d1,d2] = size(pupil);x = -(d1-1)/2:(d1-1)/2; y = -(d2-1)/2:(d2-1)/2;[Y,X] = meshgrid(y,x);r = round(sqrt(Y.^2+X.^2));
theta = atan2(Y,X);
scal = 1e2; radius = dim/scal; pupil(r<=radius) = 1;%% get the units right
pupil_unit = radius0/radius;
PSF_unit = (2*pi/dim)*radius/(2*pi*NA/lambda);

161
MTF_unit = (2*pi/dim)/PSF_unit;%% PSF and MTF for ideal pupil
PSF = ifftshift(fft2(fftshift(pupil))); PSF = abs(PSF);
PSF = PSF.^2/max(max(PSF.^2));% modulus squared for fluorescence
PSF = ifftshift(fft2(fftshift(pupil))); PSF = real(PSF);
PSF = PSF/max(max(PSF));% real part for absorption
MTF = ifftshift(fft2(fftshift(PSF))); MTF = MTF/max(max(MTF));

B.2

Determination of Zernike Polynomial Coefficients%% fitting parameters
pk = 0.9575; xc = 0.96; yc = 1.18; NA = 0.384; offset = 0.024;
S0 = -0.0978; alpha1 = -0.0065; alpha2 = 0.0008; beta = 0.0027;
delta1 = 0.0053; delta2 = 0.0136; epsilon1 = -0.0021; epsilon2 = -0.0008;% these initial guesses are also the fitting results
x0 = [ pk, xc, yc, NA, offset, S0, alpha1, alpha2, beta, delta1, delta2, ...
epsilon1, epsilon2 ]; % initial guess parametersx0fix = [ 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0 ];%% cost function[x_free,resnorm,residual,exitflag,OptmInfo,Lmda,Jmat] = ...
lsqnonlin(@(x) abs( PSFZernike(interlace(x0,x,x0fix),xdata) - OD_trimmed ),...
x0(~x0fix),lb,ub);
x = interlace(x0, x_free, x0fix);

162%% fitting function
function F = PSFZernike(x, xdata)

pupil2 = pupil.*...
exp( 2*pi*1i* ( x(6)*sqrt(5)*( 6*(r/radius).^4 -6*(r/radius).^2 + 1 ) ...
% spherical+ x(7)*sqrt(6)*(r/radius).^2.*cos(2*theta) + x(8)*sqrt(6)*(r/radius).^2...
.*sin(2*theta) ...
% astigmatism+ x(9)*sqrt(3)*( 2*(r/radius).^2 - 1 ) ... % defocus+ gamma1*2*(r/radius).*sin(theta) + gamma2*2*(r/radius).*cos(theta) ... % tilts+ x(10)*sqrt(8)*(3*(r/radius).^3 - 2*(r/radius)).*sin(theta) ...
+ x(11)*sqrt(8)*(3*(r/radius).^3 - 2*(r/radius)).*cos(theta) ...
% coma+ x(12)*sqrt(8)*(r/radius).^3.*sin(3*theta) ...
+ x(13)*sqrt(8)*(r/radius).^3.*cos(3*theta) ));% trefoil

PSF2 = ifftshift(fft2(fftshift(pupil2))); PSF2 = abs(PSF2);
PSF2 = x(1)*PSF2.^2/max(max(PSF2.^2)) + x(5);F = interp2( Y1, X1, PSF2, ( xdata(:,:,1) + x(2) )*px/PSF_unit, ...
( xdata(:,:,2) + x(3) )*px/PSF_unit, ’spline’ ); % interpolationend%% auxiliary function

163
function a = interlace(a, x, fix)a(~fix) = x;endAppendix C

Algorithms: Finite-Temperature Effects in Atom-Dimer Relaxation

C.1

Parametrization of Atom-Dimer Scattering Amplitude%% Definitions
nK = 1e-9; cm = 0.01; s = 1; ms = 0.001;

amu = 1.660539040*10^-27; %NIST
m = 38.9637064864*amu; %NIST
h = 6.626070040*10^-34; %NIST
a0 = 0.52917721067*10^-10; %NIST
kB = 1.38064852*10^-23; %NIST
hbar = h/(2*pi); kHz = 1e3;s0 = 1.00624;

a = a1*a0; a_star = a_star1*a0;
kd = 2./(sqrt(3)*a);%% AD scattering amplitude f_AD(k,theta)
c1 = @(a,k) -0.22 + 0.39*k.^2.*a.^2 - 0.17*k.^4.*a.^4;
c2 = @(a,k) 0.32 + 0.82*k.^2.*a.^2 - 0.14*k.^4.*a.^4;
phi = @(a,k) 2.64 - 0.83*k.^2.*a.^2 + 0.23*k.^4.*a.^4;

165

f_AD_para = @(a,k,a_star,eta) ...
c1(a,k) + c2(a,k).* cot( s0*( log(0.19*a) - log(a_star) ) + 1i*eta + phi(a,k) );f_AD = @(a,k,a_star,eta) ...
( f_AD_para(a,k,a_star,eta) ./a - 1i.*k ).^(-1); %eqn(5)sigma_ela = @(a,k,a_star,eta) 4*pi*abs( f_AD(a,k,a_star,eta) ).^2; %eqn(2)sigma_tot = @(a,k,a_star,eta) 4*pi*imag( f_AD(a,k,a_star,eta) )./k; %eqn(3)g = @(a,k,a_star,eta) ..
3*hbar.*k/(2*m) .*( sigma_tot(a,k,a_star,eta) - sigma_ela(a,k,a_star,eta) );%eqn(11)k_rel_integrand = @(a,k,a_star,eta,T) ..
k.^2 .* exp( - 3*hbar^2*k.^2 / (4*m*kB*T) ) .* g(a,k,a_star,eta); %eqn(13)prefactor = @(ampl,T) ..
ampl* 4/sqrt(pi) *hbar^3 ./ (m*kB*T).^(3/2); %use Boltzmann distribution%% Integral
beta_out = [];for index = 1:length(a1)k_rel = logspace(log10(kd(index)/1e6),log10(kd(index)),1e3);
beta_out_scal(index) = ...
prefactor(ampl,T) .* trapz(k_rel, k_rel_integrand(a(index),k_rel,a_star,eta,T));
beta_out = [beta_out beta_out_scal(index)];
end

166

C.2

Atom-Dimer Chemical Equilibrium: a Toy Model%% Ebind function
rvdW = 64.49; %in unit of a0
abar = 0.955978*rvdW;
Eb = @(as) hbar^2 ./( ma*( (as - abar)*a0 ).^2 );lambda = @(temp) h* ( 2*pi*ma*kB *temp ).^(-1/2);%% AD chemical equilibrium (Saha-Boltzmann type, in a box)
PSD = @(na, temp) na.*lambda(temp).^3;
nd = @(na, temp, as) 2^(3/2)* na.^2 .*lambda(temp).^3 .*exp( Eb(as)./(kB*temp) );%% number and energy conservation
alist = logspace( log10(5e2), log10(8e3), 2 ); sol0 = [];%initial guesses
for j = 1:length(alist)
syms naf1
solg0 = vpasolve( nai1 == 2*nd( naf1, Ti1, alist(j)) + naf1 );j = j+1;
solg1(j) = solg0(2);endsolg = double(solg1);%solve coupled equations
for i = 1:length(alist)
syms naf1 Tf1
sol0 = vpasolve( [ nai1 == 2*nd( naf1, Tf1, alist(i)) + naf1, ...

167
nai1* (3/2)*kB*Ti1 == ...
nd( naf1, Tf1, alist(i) )*( (3/2)*kB*Tf1-Eb(alist(i)) )+naf1*(3/2)*kB*Tf1 ],...
[naf1, Tf1], [ solg(i), Ti1 ]);sol1(1,i) = sol0.Tf1;
sol1(2,i) = sol0.naf1;
sol1(3,i) = nd( sol0.naf1, sol0.Tf1, alist(i) );i = i+1;end

sol = double(sol1);
templist = sol(1,:);
nalist = sol(2,:);ndlist = sol(3,:);C.3

Fit Function for Dimer-Dimer Collision Rate Constant% all N’s are density% dNa/dt/Na =% - beta*Nd% A-D inelastic collisions to form deep dimers% - 2*delta’*(Na)(Nd)% ignore this for early decay time

% + 2*delta*(Nd)^2/Na;% D-D collisional dissociation, assume K2+K2=4K% dNd/dt/Nd =% - alpha% - gamma*Nd% - beta*Na% 1-body decay due to spin relaxation% D-D inelastic collisions to form deep dimers% ignore this for early decay time

168
% - delta*Nd% + delta’*(Na)^2;% D-D collisional dissociation, assume K2+K2=4K% ignore this for early decay timefunction y = AtomNumbODEmulti(A,tspan,B)% A(1) = beta, A(2) = delta;
y0 = A(3);[t,y] = ode45(@DifEq, tspan, y0);function dy = DifEq(t, y)ydot = zeros(1,1);
ydot = - A(1)*y.*double_expo(B,t) + 2*A(2)*double_expo(B,t).^2;
dy = ydot;
end
end

function y = DimerNumbODEmulti(B,tspan)% B(1) = alpha, B(2) = gamma, B(3) = delta;
y0 = B(4);[t,y] = ode45(@DifEq, tspan, y0);function dy = DifEq(t, y)ydot = zeros(1,1);
ydot = - y/B(1) - B(2)*y.^2 - B(3)*y.^2;
dy = ydot;
end
endAppendix D

Algorithms: Finite-Temperature Effects in Three-Body Recombination

D.1

L3 Fitting Protocol I%% sorting and statistics%sort by timing[t_list0, t_sort0] = sort(t_list0); %get sorted indices
N_list0 = N_list0(t_sort0); %sort N
T_list0 = T_list0(t_sort0); %sort T%do stats for each wait time
[t_list, ix, iy] = unique(t_list0); %get unique timings (non-repetitive)
N_list1 = [t_list, accumarray(iy, N_list0, [], @mean)]; %N mean
N_std1 = [t_list, accumarray(iy, N_list0, [], @std)]; %N std
T_list1 = [t_list, accumarray(iy, T_list0, [], @mean)]; %T mean
T_std1 = [t_list, accumarray(iy, T_list0, [], @std)]; %T std%final lists generated from stats
N_list = N_list1(:,2); %N mean list
N_std = N_std1(:,2);
T_list = T_list1(:,2); %T mean list, T = ( Tz + 2*Tx )/3
T_std = T_std1(:,2);%regulate
max_bound_inds = t_list(:) < t_max;

170
t_list = t_list(max_bound_inds);
N_list = N_list(max_bound_inds);
T_list = T_list(max_bound_inds);
min_bound_inds = t_list(:) > t_min;t_list = t_list(min_bound_inds);
N_list = N_list(min_bound_inds);
T_list = T_list(min_bound_inds);%% temperature fit function
x0 = [ 0.15 0.04 -1 ]; %initial fit parameter guesses[x,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat] = ...
lsqnonlin(@(x)(polyx(x,t_list)-T_list),x0);

function Temp = polyx(x, t_list)% 2nd-order/quadratic polynomial: (poly fit works better for short t hold range)% Temp = x(1) + x(2)*t_list + x(3)*t_list.^2;% or use exponential: (expo fit works better for long t hold range)
Temp = x(1) + x(2).*exp(t_list./x(3));
end%% number fit function[Nmax, Imax] = max(N_list); %get maximum N index
N0_guess = N_list(Imax); %guess initial N = maximum N
B0 = [ 0.0117 (2*K2)*cm^3/(T2V0*nK^(3/2)) ...
(0.1*1.5e-26*cm^6/(T2V*nK^3)) N0_guess ]; %initial guess values
fix = [ 1 1 0 0 ]; %fix fitting parameters[B_free,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat] = lsqcurvefit(@(B,t_list) ...

171
NumbODE3(interlace(B0,B,fix),t_list,T_para),B0(~fix),t_list,N_list,[],[]);
B = interlace(B0, B_free, fix); %defined in Appx A.2
r2 = 1 - Rsdnrm/sum((N_list-mean(N_list)).^2); %r-squared
ci = nlparci(B_free,Rsd,’jacobian’,Jmat); %confidence interval of the fitting
fit_err_free = abs(ci(:,2) - ci(:,1))/1.96/2;
fit_err = interlace([0 0 0 0], fit_err_free, fix);

function y = NumbODE3(B,tspan,para)% dN/dt = -K3 <n^2> N - K2 <n> N - alpha N% <n^2> = Int(n^3 dxdydz); <n^3> = Int(n^4 dxdydz)%dN/dt = - B(1)*N - B(3)*N^3/T^3 - B(2)*N^2/T^(3/2);
y0 = B(4); uK2nK = 1e3;[t,y] = ode45(@DifEq, tspan, y0);function dy = DifEq(t, y)ydot = zeros(1,1);
ydot = - B(1)*y - B(3)*y^3./((polyx(para,t)*uK2nK).^3) - B(2)*y^2./ ...
((polyx(para,t)*uK2nK).^(3/2));
dy = ydot;
end
end%% get the statistics% time-averaged density and temperature and collision ratet_grid = linspace( t_list(1), t_list(end), 1e3 );dt = ( t_list(end) - t_list(1) )/1e3;
N_grid = NumbODE3(B,t_grid,T_para); N_grid = N_grid’;
Temp_grid = polyx(T_para, t_grid);

172

T_sum = sum( Temp_grid *dt ); %average temperature from summation
T_std = std( T_list*1e3 );T_init = polyx(T_para, 0) *1e3; % extrapolate to 0
T_final = T_list(end) *1e3;
T_avg = T_sum /(t_list(end) - t_list(1)) *1e3; %average temperature in nK
T_output = [ T_init, T_avg, T_std ];den2_sum = sum( (N_grid.^2 ./ (Temp_grid*uK).^3 )*dt );den2_avg = den2_sum / ( t_list(end) - t_list(1) ) / T2V * (cm^6); % <n^2>den_pk_avg = sqrt( den2_avg * sqrt(27) ); %peak density
den_avg_avg = den_pk_avg / (2*sqrt(2)); %average density <n>n_avg_init = ...
(N_list(1).*w.^3.*(m39./(2.*pi.*kB.*T_list(1).*uK)).^(3/2))*cm^3/(2^(3/2));
den_combo = [ den_pk_avg, den_avg_avg, den2_avg ];n_avg_list = ...
(N_grid.*w.^3.*(m39./(2.*pi.*kB.*Temp_grid.*uK)).^(3/2))*cm^3/(2^(3/2));
col_rate_list = ...
n_avg_list./cm^3.*8.*pi.*a.^2.*sqrt((16.*kB.*Temp_grid.*uK)./(pi.*m39));
col_rate_init = ...
n_avg_init./cm^3.*8.*pi.*a.^2.*sqrt((16.*kB.*T_list(1).*uK)./(pi.*m39));
col_rate_avg = sum(col_rate_list *dt) /(t_list(end) - t_list(1));
col_rate_combo = [col_rate_init, col_rate_avg ];D.2

L3 Fitting Protocol II (Simultaneous Fit)%% number fit function

173%two conditions are put in two cells[Nmax1, Imax1] = max(N_list{1}); [Nmax2, Imax2] = max(N_list{2});N0_guess1 = N_list{1}(Imax1); N0_guess2 = N_list{2}(Imax2);
B0 = [ 0.5*(0.0117+0.0159) (2*K2)*cm^3/(T2V0*nK^(3/2)) ...
(0.1*1.5e-26*cm^6/(T2V*nK^3)) N0_guess1 N0_guess2 ]; %initial guess values
fix = [ 1 0 0 0 0 ]; %fix fitting parameters
fun_cell=@(B,t_list) [NumbODE3_simu1(interlace(B0,B,fix),t_list{1},T_para{1});...
NumbODE3_simu2(interlace(B0,B,fix),t_list{2},T_para{2}) ];NN = vertcat( N_list{:} );[B_free,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat] = ...
lsqcurvefit(fun_cell,B0(~fix),t_list,NN,[(2*K2)*cm^3/(T2V0*nK^(3/2)) 0 0 0],[]);
B = interlace(B0, B_free, fix); %defined in Appx A.2

D.3

S-Matrix Formalism%% fitting programB0 = [ 900, 0.20, 230*nK, 1 ]; % guess parameters, a+, eta, T, ampl
fix = [ 0, 1, 1, 0 ];options=optimset(’TolFun’,1e-60,’TolX’,1e-50,’MaxFunEvals’,2000,’MaxIter’,2000);
lb = [ 800, 0 ]; ub = [ 1100, 2 ];[B_free,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat] = ...
lsqnonlin(@(B) abs((1./L3_std)...
.*(L3_function_pos_a_for_fit(interlace(B0,B,fix),a_list)-L3_list)),B0(~fix),lb,ub,options); %cost function
B = interlace(B0, B_free, fix);
r2 = 1 - Rsdnrm/sum((L3_list-mean(L3_list)).^2); %r-squared
ci = nlparci(B_free,Rsd,’jacobian’,Jmat); %confidence interval of the fitting
fit_err_free = abs(ci(:,2) - ci(:,1))/1.96/2;

174
fit_err = interlace([0 0 0 0], fit_err_free, fix); %defined in Appx A.2%% fitting function
function [L3out_finiteT_all] = L3_function_pos_a_for_fit_final(B, a1)
a_plus = B(1); eta = B(2); T = B(3); ampl = B(4);complex_s11 = s11_abs .* exp( 1i*s11_arg );complex_s12 = s12_abs .* exp( 1i*s12_arg );complex_s22 = s22_abs .* exp( 1i*s22_arg );L3s_T0 = L3_function_pos_a_zero_T_shallow(a_plus, eta, a1, ampl);
L3s_out2 = L3_function_pos_a_shallow(a_plus, eta, a1, T, ampl, ka, ...
complex_s11, complex_s12, complex_s22);
L3d_T0 = L3_function_pos_a_zero_T_deep(a_plus, eta, a1, ampl);
L3d_out2 = L3_function_pos_a_deep(a_plus, eta, a1, T, ampl, ka, ...
complex_s11, complex_s12, complex_s22);

ka_cut = 0.4; %about 30 points removed
L3s_out2 = L3s_out2( ka(:)>=ka_cut, : ); L3d_out2 = L3d_out2( ka(:)>=ka_cut,: );ka = ka( ka(:)>=ka_cut );%concatenate
ka = [0; ka];
L3s_out2 = [L3s_T0; L3s_out2]; L3d_out2 = [L3d_T0; L3d_out2];%interpolate
ka_interp = logspace(log10(1e-2),log10(1e1),1e4)’;L3s_interp2 = interp2(a1, ka, L3s_out2, a1, ka_interp, ’spline’);L3d_interp2 = interp2(a1, ka, L3d_out2, a1, ka_interp, ’spline’);175%integrate
x = ka_interp; a = a1*a0;
expo = exp( - hbar^2.*(x./a).^2 /m/ (kB*T) );L3s_integrand = 2*hbar^7/m^4 .* (x.^5 ./ a.^2) .* expo .* ( L3s_interp2 );L3d_integrand = 2*hbar^7/m^4 .* (x.^5 ./ a.^2) .* expo .* ( L3d_interp2 );L3sout = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(x,L3s_integrand);L3dout = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(x,L3d_integrand);% J>0 terms[ka_out, KJ_out] = L3_function_pos_a_nonzero_J;
expoJ1 = exp( - hbar^2.*(ka_out(:,1)./a).^2 /m/ (kB*T) );L3J_integrand1 = 2*hbar^7/m^4 .* (ka_out(:,1).^5 ./ a.^2) .* expoJ1 .*KJ_out(:,1);
L3Jout1 = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(ka_out(:,1),L3J_integrand1);
expoJ2 = exp( - hbar^2.*(ka_out(:,2)./a).^2 /m/ (kB*T) );L3J_integrand2 = 2*hbar^7/m^4 .* (ka_out(:,2).^5 ./ a.^2) .* expoJ2 .*KJ_out(:,2);
L3Jout2 = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(ka_out(:,2),L3J_integrand2);

expoJ3 = exp( - hbar^2.*(ka_out(:,3)./a).^2 /m/ (kB*T) );L3J_integrand3 = 2*hbar^7/m^4 .* (ka_out(:,3).^5 ./ a.^2) .* expoJ3 .*KJ_out(:,3);L3Jout3 = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(ka_out(:,3),L3J_integrand3);
expoJ4 = exp( - hbar^2.*(ka_out(:,4)./a).^2 /m/ (kB*T) );L3J_integrand4 = 2*hbar^7/m^4 .* (ka_out(:,4).^5 ./ a.^2) .* expoJ4 .*KJ_out(:,4);L3Jout4 = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(ka_out(:,4),L3J_integrand4);

expoJ5 = exp( - hbar^2.*(ka_out(:,5)./a).^2 /m/ (kB*T) );L3J_integrand5 = 2*hbar^7/m^4 .* (ka_out(:,5).^5 ./ a.^2) .* expoJ5 .*KJ_out(:,5);176
L3Jout5 = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(ka_out(:,5),L3J_integrand5);
expoJ6 = exp( - hbar^2.*(ka_out(:,6)./a).^2 /m/ (kB*T) );L3J_integrand6 = 2*hbar^7/m^4 .* (ka_out(:,6).^5 ./ a.^2) .* expoJ6 .*KJ_out(:,6);
L3Jout6 = ampl /6/ ( 2*(kB*T)^3 ) .* trapz(ka_out(:,6),L3J_integrand6);

L3out_finiteT = L3sout + L3dout +...
(L3Jout1 + L3Jout2 + L3Jout3 + L3Jout4 + L3Jout5 + L3Jout6)*0;[Lmax,index] = max(L3out_finiteT);%interpolate to unitary value where s_ij lose accuracy
if index < length(a1)L3out_finiteT_truc = L3out_finiteT(1:index);
L3out_finiteT_cat = [L3out_finiteT_truc, L3max];
a_truc = a1(1:index);
a_cat = [a_truc,a1(end)];L3out_finiteT_all = interp1(a_cat,L3out_finiteT_cat,a1,’PCHIP’);else
L3out_finiteT_all = L3out_finiteT;end

L3out_finiteT_all = 3*(cm^6/s)^(-1)*L3out_finiteT_all;%3 is the number of atoms lostend%% sub-functions (zero-T)
function [ L3out ] = L3_function_pos_a_zero_T_shallow(a_plus, eta, a1, ampl)
L3_shallow = 768*pi^2*(4*pi-3*sqrt(3)) *( sin( theta ).^2 + sinh(eta).^2 ) ...

177./ ( sinh( pi*s0 + eta )^2 + cos( theta ).^2 );L3out = L3_shallow;
end

function [ L3out ] = L3_function_pos_a_zero_T_deep(a_plus, eta, a1, ampl)
L3_deep = (3/2)*128*pi^2*(4*pi-3*sqrt(3)) ./ sinh(pi*s0).^2 .* ...
sinh(2*pi*s0) .* sinh(2*eta) ./ ( sinh(pi*s0 + eta)^2 + cos( theta ).^2 );L3out = L3_deep;end%% sub-functions (finite-T)%eq3 of PRA 78,043605 (2008)
function [ L3out ] = L3_function_pos_a_shallow(a_plus, eta, a1, T, ampl, ka, ...
complex_s11, complex_s12, complex_s22)
x = ka; %change of variables. s11(ka) = s11(x)L3_shallow = ...
1 - abs( complex_s22 + complex_s12.^2 .* exp( 2*1i*theta - 2*eta ) ...
./ ( 1 - complex_s11 .* exp( 2*1i*theta - 2*eta )) ).^2 - ( 1 - exp(-4*eta) ) ...
.* abs(complex_s12).^2 ./ abs( 1 - complex_s11 .* exp( 2*1i*theta - 2*eta ) ).^2;
L3out = 144*sqrt(3)*pi^2 ./ (x.^4) .*L3_shallow;
end

function [ L3out ] = L3_function_pos_a_deep(a_plus, eta, a1, T, ampl, ka, ...
complex_s11, complex_s12, complex_s22)
theta = s0*log(a1./a_plus);%theta =x = ka;pi/2 + theta; %for fitting recombination maxima

178
L3_deep = ...
( 1 - exp(-4*eta) ) .* ( 1 - abs(complex_s11).^2 - abs(complex_s12).^2 ) ...
./ abs( 1 - complex_s11 .* exp( 2*1i*theta - 2*eta ) ).^2;
L3out = 144*sqrt(3)*pi^2 ./ (x.^4) .*L3_deep;endAppendix E

Data and Fits: Three-Body Recombination at Finite Temperatures

E.1

(0)

Extract a+ from 460 nK Data

Figure E.1: Finite-temperature fit for recombination minimum, 460 nK condition, J = 0. The
fitting results are consistent with zero-temperature fit within errors.

180
(0)

E.2

Extract ap from 410 nK and 230 nK Data

(a)

p

(b)

p

Figure E.2: Finite-temperature fit for recombination maximum, 410 nK and 230 nK conditions,
J = 0. We fix η to be 0.2 in this fit.

181
(1)

E.3

Exploration of a+

(a)

(b)

Figure E.3: Finite-temperature fit for recombination minimum, 230 nK and 210 nK condition,
J = 0.

182

(a)

(b)

Figure E.4: Finite-temperature fit for recombination minimum, 150 nK and 60 nK condition, J = 0.