Theoretical Investigation of Cold Alkali Atom
Collisions
byJames Patrick Burke, Jr.
B.S., Georgia Institute of Technology, 1989
M.S., Georgia Institute of Technology, 1991

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
1999

This thesis entitled:
Theoretical Investigation of Cold Alkali Atom Collisions
written by James Patrick Burke, Jr.
has been approved for the Department of Physics

Chris H. Greene

John Cooper

Date

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

iii
Burke, Jr., James Patrick (Ph.D., Physics)
Theoretical Investigation of Cold Alkali Atom Collisions
Thesis directed by Prof. Chris H. Greene
Theoretical techniques are developed to solve the coupled Schrödinger equations
that describe cold magnetically trapped alkali atom dynamics. In particular, a multichannel quantum defect formulation is presented to separate cleanly the atomdependent from the atom-independent contributions to the final scattering matrices.
The resulting approach provides a more comprehensive understanding of cold alkali
collisions and improves substantially on the efficiency of calculations. Improved
Rb two-body interaction potentials are extracted from measured 87 Rb inelastic collision rates and from the measured position and width of a 85 Rb Feshbach resonance. Using these state-of-the-art potentials, specific atomic hyperfine states havebeen identified whose collisional properties suggest intriguing possibilities for interesting and novel degenerate gas studies. In addition, a detailed analysis of 39 K
photoassociation lineshapes has permitted the determination of accurate two-bodyinteraction potentials for this atom as well. Finally, three-body recombination of
doubly polarized trapped atoms is investigated. Two qualitatively different mechanisms which control the recombination rate are identified for positive and negative
two-body scattering lengths a. Surprisingly, the recombination rate is found to scale
approximately as a4 in both cases.

To Angela

Acknowledgements

This work would not have been possible without support from many people. First and foremost, I would like to thank my advisor Chris Greene. Chris has
shown uncommon patience and determination in his attempt to teach me the finer
points of atomic collision theory. I feel very fortunate to have had this opportunity
to study under a person of his stature.
I am also grateful that Chris guided me in the direction of cold atomic
collisions. In the course of this study I have crossed paths with a number of outstanding scientists. Notable among this group is John Bohn, who arrived at JILAat about the time this project started. John has proved to be an extremely patient“sounding-board” during my ongoing struggles with this subject and in the process,he has managed to teach me a great deal of physics. John has been involved tosome degree in almost all of the work presented here, and he has provided many
constructive comments on this dissertation.
The Greene group has had a number of postdocs and graduate students
from whom I have benefited. Of this group, I would like to single out Brett Esry
who has helped me appreciate the fine art of numerical computing and taught me
a great deal about many-body physics. Chapter 8 is the outgrowth of a project
Brett and Chris originally started. In addition, I would like to thank Hugo van derHart, Dörte Blume, William Clark, Kurt Meyer, and Brian Granger for many useful
discussions.
An advantage of working at JILA is the active experimental program thatis being pursued. In particular, I would like to thank the Cornell-Wieman-Jin experimental group studying BEC and degenerate Fermi gases for providing us withthe opportunity to model some really marvelous experiments. I would like to thankChris Myatt, Eric Burt, Jake Roberts, Neil Claussen, and Rich Ghrist for providing
a glimpse into the experimental world. In addition, I want to recognize the members of my thesis committee, John Cooper, Eric Cornell, Carl Wieman, and David
Nesbitt for constructive comments on this manuscript.
The photoassociation work presented in Chapters 6 and 7 resulted froma collaboration with the experimental group of Bill Stwalley, Phil Gould, and He
Wang. I would like to thank them for allowing us to analyze their data and for being
so cooperative with what probably seemed like an endless barrage of questions.

vi
Also, I would like to thank Paul Julienne, Carl Williams, and Eite Tiesinga for many
useful discussions regarding both photoassociation theory and cold collision theory.
On a personal note, I am indebted to my parents for their support andencouragement over all these years. My wife has probably worked harder for mydegree than I have. Without her sacrifice, I would have had to end my involvementwith this project long ago. Angela, thank you for hanging in there. And last but
certainly not least, I want to recognize my children, Julia and Joel, just for putting asmile on my face each day.
This work would not have been possible without financial support fromthe National Science Foundation.

vii

Contents

Chapter
1 Introduction

1

2 Coupled Schrödinger equations for two-body collisions

6

2.1

2.2

2.3

Scattering theory . . . . . . . . . . . . . . . . . . . . . . . . . . .

6

2.1.1

Derivation of the cross section and event rate formulas . . .

9

2.1.2

Relating the scattering amplitude to the scattering matrix . .

11Interaction potentials . . . . . . . . . . . . . . . . . . . . . . . . .

17

2.2.1

Constructing the Born-Oppenheimer potentials . . . . . . .

18

2.2.2

Atomic hyperfine Hamiltonian . . . . . . . . . . . . . . . .

21

2.2.3

Magnetic dipole interaction . . . . . . . . . . . . . . . . .

22

2.2.4

External magnetic field and Feshbach resonances . . . . . .

23Angular momentum representations . . . . . . . . . . . . . . . . .

24

2.3.1

Symmetries . . . . . . . . . . . . . . . . . . . . . . . . . .

29

2.3.2

Transformations . . . . . . . . . . . . . . . . . . . . . . .

30

3 Numerical techniques for solving the coupled Schrödinger equations
3.1

32

Finite element R-matrix approach . . . . . . . . . . . . . . . . . .

33

3.1.1

R-matrix equations . . . . . . . . . . . . . . . . . . . . . .

33

3.1.2

Finite element basis . . . . . . . . . . . . . . . . . . . . .

34

viii
3.1.3

Matching R-matrix boxes . . . . . . . . . . . . . . . . . .

36

3.2

Constructing the scattering matrix . . . . . . . . . . . . . . . . . .

39

3.3

Adiabatic representation . . . . . . . . . . . . . . . . . . . . . . .

40

4 Multichannel quantum defect treatment of cold collisions
4.1

43

Reference wave functions . . . . . . . . . . . . . . . . . . . . . . .

44

4.1.1

Milne equation . . . . . . . . . . . . . . . . . . . . . . . .

46

4.1.2

Calculating the MQDT parameters . . . . . . . . . . . . . .

48

4.1.3

Low energy behavior of the long-range parameters . . . . .

51

4.2

Constructing scattering matrices and wave functions . . . . . . . . .

55

4.3

Frame transformation approximation . . . . . . . . . . . . . . . . .

61

4.4

Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

64

5 A compendium of rubidium scattering properties

70

5.1

Suppressed 87 Rb inelastic cross sections . . . . . . . . . . . . . . .

71

5.2

Observation of a Feshbach resonance in 85 Rb collisions . . . . . . .

78

5.3

Predictions for 85 Rb-87 Rb collisions . . . . . . . . . . . . . . . . .

84

5.4

Scattering properties of the short-lived Rb isotopes . . . . . . . . .

91

6 Theory of photoassociation lineshapes

101

6.1

One color lineshape formula . . . . . . . . . . . . . . . . . . . . . 102

6.2

Symmetry rules for diatomic Hund’s case (a) wave functions . . . . 104

6.3

Excited and ground state Hamiltonians . . . . . . . . . . . . . . . . 108
6.3.1
6.3.2

6.4

S + P Hamiltonian . . . . . . . . . . . . . . . . . . . . . . 108
S + S Hamiltonian . . . . . . . . . . . . . . . . . . . . . . 113

Dipole matrix element . . . . . . . . . . . . . . . . . . . . . . . . 114

ix
7 Determination of potassium scattering lengths

121

7.1

Photoassociation experiment . . . . . . . . . . . . . . . . . . . . . 121

7.2

Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126

8 Three-body recombination

136

8.1

Adiabatic hyperspherical representation . . . . . . . . . . . . . . . 137

8.2

Recombination rates and mechanisms . . . . . . . . . . . . . . . . 141

9 Summary

155

Bibliography

157

Appendix
A FEM basis functions

167

B Adiabatic FEM R-matrix

168

C Thermal averaging

170

D Related publications

173

x

FiguresFigure2.1

Rb2 Born-Oppenheimer potentials . . . . . . . . . . . . . . . . . .

20

2.2

Magnetic field dependence of 85 Rb thresholds . . . . . . . . . . . .

25

2.3

Schematic of a Feshbach resonance . . . . . . . . . . . . . . . . .

26

2.4

Feshbach resonances in 85 Rb and 23 Na collisions . . . . . . . . . .

27

3.1

R-matrix box schematic . . . . . . . . . . . . . . . . . . . . . . . .

38

3.2

87 Rb f =2 adiabatic potentials and couplings . . . . . . . . . . . . .

42

4.1

Energy dependence of cold collision wave functions . . . . . . . . .

45

4.2

Rb2 long-range potential . . . . . . . . . . . . . . . . . . . . . . .

47

4.3

Reference wave functions and Milne parameters . . . . . . . . . . .

49

4.4

Standardized long-range parameters . . . . . . . . . . . . . . . . .

52

4.5

Threshold energy dependence of elastic cross section . . . . . . . .

54

4.6

Energy and field dependence of the short-range eigenphase shifts . .

57

4.7

85 Rb elastic s-wave partial cross section . . . . . . . . . . . . . . .

59

4.8

Rb and Na quantum defects . . . . . . . . . . . . . . . . . . . . . .

63

4.9

Comparison of calculations with an approximate reaction matrix . .

66

4.10 Scattering length comparison using standardized long-range MQDTparameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

67

xi
4.11 Inelastic rates derived from an approximate reaction matrix . . . . .

68

4.12 Comparison of FT approximation over entire s ; t plane . . . . .

69

5.1

87 Rb rate constants versus assumed singlet scattering length . . . . .

72

5.2

Block diagonal structure of the 87 Rb collision Hamiltonian . . . . .

74

5.3

Suppression of spin-exchange collisions in 87 Rb . . . . . . . . . . .

76

5.4

Singlet and triplet wave functions . . . . . . . . . . . . . . . . . .

77

5.5

Topology of 87 Rb f = 2(12) inelastic collisions . . . . . . . . . . .

79

5.6

Constraints on 85 Rb singlet-triplet parameter space . . . . . . . . .

81

5.7

Rb elastic cross sections versus energy . . . . . . . . . . . . . . . .

83

5.8

First minimum and maximum in s-wave partial cross sections . . . .

85

5.9

Field dependence of 85 Rb dipole loss rate . . . . . . . . . . . . . .

86

5.10 Feshbach resonances in 85 Rb-87 Rb collisions . . . . . . . . . . . . .

88

5.11 Energy-field evolution of two Feshbach resonances . . . . . . . . .

89

5.12 Block diagonal structure of 85 Rb-87 Rb Hamiltonian . . . . . . . . .

90

5.13 Magnetic-field suppression of 85 Rb-87 Rb loss rates . . . . . . . . .

92

5.14 Magnetic-field dependence of the 85 Rb-87 Rb short range reaction
matrix . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

93

5.15 Mixing angle dependence of 84 Rb-84 Rb and 86 Rb-86 Rb bound states

96

5.16 84 Rb-84 Rb and 86 Rb-86 Rb elastic cross sections . . . . . . . . . . .

98

5.17 84 Rb-84 Rb and 86 Rb-86 Rb Feshbach resonances . . . . . . . . . . .

99

5.18 83 Rb-87 Rb Feshbach resonances . . . . . . . . . . . . . . . . . . . 100
6.1

Block structure of the S + P Hamiltonian including the 0;
g=u electronic state . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110

6.2

Block structure of the S + P Hamiltonian including the 0+
g=u electronic state . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111

xii
6.3

Comparison of Franck-Condon factors . . . . . . . . . . . . . . . . 115

7.1

Schematic of the 39 K photoassociation 0;g experiment . . . . . . . . 123

7.2

Final state dependence of the 0;
g spectra . . . . . . . . . . . . . . . 124

7.3

Second color dependence of the 0;
g spectra . . . . . . . . . . . . . . 125

7.4

Schematic of the 39 K ground and excited state potentials and wave
functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127

7.5

Comparison of theoretical and measured 0;
g spectra . . . . . . . . . 129

7.6

Dependence of J 0 = 2 linewidth on position of d-wave shape resonance132

7.7

Allowed K singlet-triplet parameter space . . . . . . . . . . . . . . 133

7.8

Comparison of theoretical and experimental J 0 = 2 linewidths . . . . 134

8.1

Adiabatic hyperspherical potential curves for 7 Li +7 Li+7 Li . . . . . 142

8.2

Recombination length as a function of two-body scattering length . 144

8.3

Adiabatic hyperspherical potential curves and couplings for positive
and negative two-body scattering lengths . . . . . . . . . . . . . . . 145

8.4

Hyperradial position of the long-range avoided crossing for positive
two-body scattering lengths . . . . . . . . . . . . . . . . . . . . . . 146

8.5

4 He recombination length . . . . . . . . . . . . . . . . . . . . . . . 148

8.6

Comparison of analytic recombination lengths with numerical calculations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 150

8.7

Bound state dependence of rate limiting crossing . . . . . . . . . . 152

8.8

Bound state dependence of the three-body entrance channel barrier . 153

8.9

Energy dependence of three-body recombination rate . . . . . . . . 154

xiii

TablesTable
2.1

Born-Oppenheimer potential parameters . . . . . . . . . . . . . . .

19

4.1

Zero energy values of the long-range parameters . . . . . . . . . . .

53

5.1

Hamiltonian structure for unpolarized bosonic gases . . . . . . . . .

73

5.2

Rubidium scattering lengths and quantum defects . . . . . . . . . .

82

5.3

Rubidium atomic properties . . . . . . . . . . . . . . . . . . . . .

94

7.1

Potassium scattering lengths and quantum defects . . . . . . . . . . 130

7.2

Potassium atomic properties . . . . . . . . . . . . . . . . . . . . . 130

8.1

Recombination event rate constants . . . . . . . . . . . . . . . . . 150

Chapter 1
Introduction

The study of particle collisions has a long and storied history in the annalsof scientific exploration. Much of our understanding of the fundamental nature ofparticles and their underlying microscopic interactions has been derived from scattering experiments. A famous example is the work of Ernest Rutherford, who in
1911 scattered alpha particles off thin gold foil. He quickly realized that the only
explanation for the observed distribution of backscattered particles was the presenceof a charged massive core concentrated within a tiny spatial region. He thus discovered the nucleus and proposed a celestial model for the atom, with electrons orbitinga massive central nucleus. Soon thereafter, Neils Bohr combined his quantum theory of light with Rutherford’s model to build a theory of the hydrogen atom and withit modern atomic physics was born.
The existence of atoms (or at least something called atoms) had been accepted in some scientific circles hundreds of years earlier. One of the first gas dynamics experiment was performed by Robert Boyle who discovered in 1662 the inverse relationship between the volume and pressure of a gas. This relationship was
independently verified by Mariotte[1] five years later. However, it was probably
Bernoulli (1739) who first recognized that the pressure of a gas is a result of collisions between the particles and the container. The development of the kinetic theoryof gases by Maxwell, Boltzmann, Krönig, and Clausius, among others, was begun in
the mid-1800’s to combat the hopeless task of evaluating each individual particle’s
contribution to the bulk properties of the gas. Kinetic theory replaces the individual particle motion with average quantities based on the statistics of the sample.
From these considerations, Maxwell developed his theory of velocity distributions(1860). A key concept in kinetic theory is that any disturbance in one part of thegas is transmitted throughout by particle collisions. In this manner, energy (thermal
conductivity), momentum (viscosity), and mass (diffusion) can be moved aroundwithin the gas. Boltzmann (1872) (improving on an earlier result of Maxwell’s)
derived a nonlinear integrodifferential equation that describes the evolution of an
initially non-uniform gas and, from which, these transport properties can be calculated. Treating the gas particles as tiny classical “billiard balls” bouncing off one
another, Maxwell found that the viscosity of a gas is independent of its density. His

2
subsequent experiments validated this surprising result, which proved to be one of
the first successes for the classical kinetic theory of gases and thus for the underlying mechanism of particle collisions. (Later measurements showed that this result
is not exact for real gases but it is a good approximation at moderate densities).
With the development of quantum mechanics and a better understandingof the structure of atoms, it became clear that atomic collisions are a much more
complicated process than envisioned by Maxwell and his fellow 19th century scientists. However, classical kinetic theory continues to provide a reliable methodfor determining the macroscopic properties of a gas, provided quantum mechanicalcollision cross sections are used in place of the corresponding classical scatteringcross sections. In the particular case of dilute monatomic gases, Chapman (1916)
and Enskog (1917) were able to show that classically, the bulk transport propertiesare related to integrals over the interatomic potential energy surfaces. This provided
the first direct link between macroscopic properties of the gas and the underlying
microscopic forces that atoms exert on one another. With it came the obvious question: can a study of the bulk properties of a gas be used to understand the underlying
atomic interactions? Given the complexity of the relationship, the answer was far
from obvious, but over the last half-century great strides have been made (see forexample Ref.[2]). (Transport measurements are not the only avenue to obtain thisinformation but they are the most relevant for this particular work). In particular,
the temperature dependence of the viscosity of a gas has yielded useful information
about the interatomic potential. This is due in large part to the fact that viscosity canbe measured much more accurately than other transport properties. Measurements
over the temperature ranges from tens of degrees Kelvin to several thousand degrees
Kelvin have been made for both pure and binary mixtures of dilute gases, primarily
by the groups of Smith[3], Guevara[4], and Kestin[5].
Further departures from the classical physics of Maxwell’s day arose withthe development of laser cooling and trapping of neutral atoms (1985, see Ref.[6,7]). The translational motion of the atoms can now be slowed to the point where theeffects of quantum statistics begin to affect the overall behavior of the gas. Probably
the biggest single achievement in this regard has been the realization of a BoseEinstein condensate (BEC), first seen here at JILA (1995) in a dilute gas of 87 Rb[8]atoms and quickly followed by observations in gases of 7 Li[9] and 23 Na[10] atoms.
It has since been observed by nearly twenty groups (see http://amo.phy.gasou.edu/bec.html
for the latest BEC developments). Bose-Einstein condensation is a macroscopic
quantum mechanical phase transition in which a large fraction of bosons spontaneously occupies the quantum mechanical ground state as its temperature approaches absolute zero. A condensate thus exhibits purely quantum mechanicalproperties on a macroscopic scale, which makes it unique.
Collisions play an important role both in the road to condensation and in
the behavior of the condensate itself. Cooling atoms with laser light has thus far
proven insufficient to produce BEC. Instead, it has proven necessary to use evapo-

3
rative cooling techniques to maintain the phase space density required for condensation. With evaporative cooling, the atoms are loaded into a magnetic trap thatis approximately harmonic. A radio-frequency field is then applied, lowering thetrapping potential so that only the “hottest” (i.e., most energetic) atoms are removed
from the trap. The remaining atoms are left to re-thermalize to a colder temperature
via elastic collisions. These elastic collisions are loosely termed “good” collisions,since they are an essential part of the cooling process. On the other hand, magnetic
traps are extremely weak, with depths on the order of few milliKelvin. If a pair ofatoms should collide and “flip” a spin, they may end up in untrapped spin statesor the collision may release enough energy to eject the pair from the trap. Either
outcome constitutes a “bad” inelastic collision. A general rule of thumb devised byexperimentalists is that about 100 good elastic collisions are needed for every bad
inelastic collision if the cloud is to re-thermalize sufficiently quickly to condense
the gas before too many atoms are lost.
Once a condensate is produced, its properties are almost wholly determined by the s-wave two body scattering length a (which is related to the elastic
scattering cross section in the T ! 0 limit). In fact, the sign (positive or negative)of the scattering length determines whether condensates comprised of large numbers of atoms can even be produced[11]. A variety of mean-field and Hartree-Fockapproaches have been introduced for evaluating condensate properties, all of whichadequately approximate the strength of atom-atom interactions using a mean-fieldinteraction term that is proportional to a. (A good review of current BEC theories isprovided in Ref[12].) In addition, modern kinetic theories that describe condensate
formation[13] require atomic scattering input.
Another intriguing aspect of collisions in this cold regime is the prospectof using external fields to influence the outcome. One method proposed by Stwalley[14]is to induce a Feshbach resonance using a magnetic bias field. The resonance would
in principle allow control over the s-wave scattering length. This was recently realized in 23 Na condensates[15] and in thermal Rb[16, 17] and Cs[18] clouds. In
terms of condensates, the appeal of this idea is two-fold: 1) it provides a mechanism
for enhancing the elastic collision rate, which could allow condensation of isotopes
that are hard to cool otherwise (this is somewhat of a fine line since generally theinelastic collision rate is also increased near the resonance), and more importantly
2) it provides for “tunable” condensate properties. (Tunability was not achieved in
the Na experiment because of a large increase in the three-body recombination rate.
Whether this an inherent problem or instead is isolated to this particular Na experiment is still unresolved.) Magnetic-field tuning of the scattering length is possible
only because the average collision energy of the atoms is much less than the atomic
hyperfine splitting.
Extremely low collision energies also imply that the atoms interact with
one another for a long time during a collision event. This makes the study of cold
collisions an ideal method for understanding the weak long-range forces that atoms

4exert on one anther. Photoassociation spectroscopy is a tool that has been successfully applied in this regard. The idea was first proposed by Thorsheim, et al.[19].
In photoassociation, an excited molecular bound state is resonantly formed throughphotoexcitation during the collision of two atoms, (e.g. with a CW laser). The
low collision energies lead to well-resolved rotational (and in some cases hyperfine)
spectra with linewidths approaching the natural linewidth. In particular, this technique has been extremely useful for measuring nodes and antinodes of the groundstate colliding wave function (via the Franck-Condon overlap). In this way, accurate
values for the s-wave scattering lengths have been derived for most of the alkalis. In
addition, retardation effects[20] have been observed in excited Na2 bound states[21],and accurate values of atomic lifetimes for 23 Na[21] and 39 K[22] have been derived
from these bound state measurements.
This dissertation has developed and extended the theoretical tools necessary to calculate and understand two-body and three-body collisions of trappedalkali atoms. In brief, ground state two-body collisions comprise Chapters 2-5,
while photoassociation spectroscopy is the subject of Chapters 6-7, and three-bodycollisions of ground state atoms are treated in Chapter 8.
The chapters are broken down as follows: the coupled Schrödinger equations for two-body ground state alkali collisions are derived in Chapter 2. Here, definitions for cross sections, scattering event rates, and scattering lengths are provided
in terms of the scattering matrix. The interactions and the corresponding Hamiltonian for the colliding atom pair are also discussed in this chapter. Numerical techniques for obtaining the scattering matrix by solving the coupled Schrödinger equations are presented in Chapter 3. Chapter 4 develops a powerful alternative approachfor evaluating the scattering matrix. Here, the ideas of multichannel spectroscopy[23]are applied to the two-body collision problem. The outcome is a short-range reaction
matrix that incorporates the bulk of the scattering physics, yet remains essentiallyenergy-independent over nearly a degree Kelvin. Simple algebraic manipulations allow the scattering matrix to be evaluated in terms of this reaction matrix and a set oflong-range parameters. Moreover, these long-range parameters are “standardized”
in a form that makes them nearly atom-independent. In addition it is shown that thereaction matrix can be approximated, with reasonable accuracy sufficient for survey
calculation in most systems, entirely in terms of single channel parameters and a
recoupling matrix.
The theory developed in the previous chapters is then applied to a set ofcollision experiments conducted here at JILA. These results are discussed in Chapter
5. The new experimental information permitted us to refine the Rb interatomic interaction potentials. Based on these state-of-the art potentials, predictions of scatteringobservables are then provided for other Rb isotopes. The theory of photoassociationlineshapes is presented in Chapter 6. This theory is then applied in Chapter 7 toanalyze measured rovibrational spectra of the 39 K 0;
g state. This analysis extracted
39
the first accurate values of the K scattering lengths. Finally, Chapter 8 presents a5theory of three-body recombination, a process in which two atoms recombine into a
molecule along with the conversion of the excess binding energy into kinetic energy.
Three-body recombination is an extremely difficult problem and Chapter 8 will only
scratch its surface. However, from our investigation of doubly-spin-polarized atoms
using model potentials, we have been able to uncover universal systematics of therecombination process.

Chapter 2
Coupled Schrödinger equations for two-body collisions

This chapter presents the two-body scattering theory of cold ground statealkali atoms. Section 2.1 begins with the coupled Schrödinger equations that describe the two-body collision dynamics. In addition, the appropriate scatteringboundary conditions applied to these coupled differential equations are discussed.
Physical observables of primary interest include scattering cross sections, event
rates, and scattering lengths which are derived in subsections 2.1.1 and 2.1.2. In
particular, the identical particle aspects of scattering by spin-dependent potentials
are handled in detail. Section 2.2 discusses the physical interactions included in the
coupled Schrödinger equations. These are the Born-Oppenheimer potentials (subsection 2.2.1), the atomic hyperfine Hamiltonian (subsection 2.2.2), the magnetic
dipole-dipole interaction (subsection 2.2.3), and the effects of an external magnetic
field (subsection 2.2.4). The coupled Schrödinger equations can be solved in any ofseveral different angular momentum representations. Those of greatest importanceare discussed in section 2.3. Subsection 2.3.1 develops the symmetrization of each
representation while subsection 2.3.2 presents the unitary transformations between
representations.
2.1 Scattering theory
Quantum two-body scattering theory is a well established subject presented in numerous textbooks (see for example, Ref.[24, 25]). However, scatteringof identical particles with spin is still a tricky subject which can often lead to factorof two errors if one is not careful. I begin by briefly deriving a very general formof the coupled Schrödinger equations. The most general form of boundary conditions applied here is given in terms of a scattering amplitude f (; ). Subsection
2.1.1 relates the symmetrized scattering amplitude f (; ) + f ( ; ;  + ) to
the scattering cross section and and to the event rate. While these relationships arederived within the Born approximation, they should have applicability to nonperturbative processes. At the end of this subsection, the physical observables are still
given in terms of the Born integrals over the potential. These integrals are related to7symmetrized scattering S and transition T matrices in subsection 2.1.2. In particular, the four spin permutations of identical particle scattering for both isotropic andanisotropic potentials are examined carefully, tracking all factors of two. The final
relationships among the cross sections, event rates, and scattering lengths in termsof the symmetrized S - and T -matrix elements are given in this subsection.
The forces on an atom in a trap depend, to a good approximation, only
on the atom’s position in the trap and not explicitly on time. The time independent
Schrödinger equation

H =E

(2.1)

is therefore a proper starting point to describe the collision physics. A standard
approach for solving two-body collision problems is to separate the center-of-mass
motion from the relative motion of the atoms. This step is valid provided the forceson the two colliding atoms (labeled a and b) depend only on their relative orientation
R~ = ~ra ; ~rb and not their absolute position. Separating coordinates sounds a little
counterintuitive since the goal of this theory is to describe trapped atoms. The“degree” of coordinate separability can be established by considering the relevant
length scales in our problem. First, the range of a cold scattering event is roughlyR  103 Bohr. Typical magnetic trapping fields can be described accurately byan anisotropic harmonic potential which in the center-of-mass coordinate system isgiven by V = 12 (!x x2 + !y y 2 + !z z 2 ). Here R2 = x2 + y 2 + z 2 ,  is the reduced
mass of the pair, and ! the angular frequency of the trapping field. A pure quadratic
field is still separable as is evident from the coordinate dependence in V . It is only
the presence of anharmonic terms (e.g., V / R3 ) which would couple the centerof-mass and relative motions. However, these anharmonicities are several orders ofmagnitude smaller than the quadratic potential. Considering the isotropic case andplugging in typical values for 87 Rb collisions (!  100 Hz), the harmonic trapping
potential adds an additional 1K of energy only when the internuclear separationis larger than R  3(106 ) a.u., well beyond the distance one needs to consider forcollisions. The trapping fields can therefore be safely neglected, which permits thecenter-of-mass motion to be removed from the scattering problem.
A multicomponent wave function representing the relative motion of theparticles can be written in spherical coordinates as

= R;1

X
lmk

Flmk (R)Ylm(; )k :(2.2)

Here, the angular degrees of freedom are represented by an expansion in sphericalharmonics Ylm . (This is the partial wave expansion treated in textbooks). The spin
degrees of freedom are labeled k and Flmk (R) represents a set of radial functions to
be determined. Insertion of this form for into Eq. 2.1 and projection onto Yl m k
0

0

0

8
leads to an infinite set of coupled differential equations of the form (i  Ylm k )X   1 @2ff + 1)~
fi ; 2 @R2 + lf (2lR
2 ; E + Vfi (R) Ffi = 0 :(2.3)~ ) will be discussed in section 2.2. If the interaction
The interaction potentials Vfi (Rpotentials are isotropic, a finite set of coupled Schrödinger equations are obtainedfor each partial wave l, of which, only l=0-4 typically need to be considered forcollision energies below a few mK. The coupled Schrödinger equations (Eq. 2.3)
are solved subject to the following boundary conditions:
Ffi (0)=
R!1

Ffi(R) ;!

0

r 2  e;i(kiR;li=2) ; ei(kiR;li=2)S 
~2ki

p

fi

fi

2i

(2.4)

where ki = 2Ei =~2 . The amplitude of the wave scattered from channel i intochannel f is given by the scattering matrix element Sfi . All experimentally observable collision quantities can be expressed in terms of the scattering matrix. In
practice, it is more convenient to work with a real-valued reaction matrix K insteadof the complex-valued S matrix. This will be discussed in Section 3.2. The asymptotic boundary condition on Ffi implies that the total wave function behaves in thelarge R limit as:

!1
+ R;!
i

X  ik~f R~
f

f e



ikf R
e
fi + R ffi(; )(2.5)

Here, the scattering amplitude ffi (; ) takes on the role of the scattering matrix andthe exact relationship between the two quantities is derived in Section 2.1.2. In thecase of identical particle scattering, the total wave function must be an eigenfunctionof the nuclear permutation operator P12 , i.e. the wave function must be symmetrized.
The boundary condition on the symmetrized wave function + is then given by

!1
+ R;!
i

X

ik~f R~ + e;ik~f R~
e
p
f
fi
2
f
 (f (; ) + f ( ; ;  + )  eikf R 
pfi
+ fi
R
2

(2.6)

and the physical observables are determined by the symmetrized scattering amplitude ffi (; ) + ffi ( ; ;  + ). The scattering amplitudes will be related to thephysical quantities of interest in cold collisions (i.e., cross sections, event rates, andscattering lengths) in the next subsection.

9
2.1.1 Derivation of the cross section and event rate formulas
The Born approximation can be used to derive cross section and event
rate formulas, in a way that brings out the structure of the scattering amplitudes andtherefore makes the identical particle aspects of the scattering theory more transparent. First, consider the scattering of two non-identical spinless particles. The initialwave function in the center-of-mass coordinate system can be described by a plane
wave eiki Z Here, particles a and b approach from  = 0 and  =  , respectively.
The time-independent wave function that obeys outgoing-wave boundary conditions
obeys the integral equation[25]+ = eikZ +Z

d3R0 G+(R; R0) V (R0 ) +(R0 )(2.7)

where ~ki = kZ^ and the free particle Green function is

eikjR~ ;R~ j :(2.8)
4~2 jR~ ; R~ 0j~ ; R~ 0 j ' R ; R0R^  R^0 and the Green function becomes
In the limit of large R, jR
G+(R~ ; R~ 0 ) = ; 202 eikf R
G+(R~ ; R~ 0) 
= ; 4~2 R e;i~kf R~

0

R!1 ^
where ~kf ;! kR. Using the first Born approximation for
expression for the scattered wave function(2.9)+ , I obtain a final

 Z d3R0 e;i~kf R~ V (R0)eikZ :(2.10)
R 2~2The non-identical particle scattering amplitude f (; ) can then be extracted byikR!1 ikZ e
+ R;!
e ;

0

0

comparing equation 2.5 with 2.10. The integrals in equation 2.10 will be expressedin terms of T -matrix elements in subsection 2.1.2
I turn next to the adaptation of this derivation to the case of identical particle scattering. The permutation operator P12 commutes with the Hamiltonian of oursystem and can be used to obtain new eigenfunctions with the proper symmetry. For
the moment, the spin degrees of freedom will be suppressed. These will be taken
into account when actually evaluating the scattering amplitude integrals (subsection
2.1.2). For the sake of clarity, I will only treat collisions between identical bosons.
The new wave functions are1 + P 
p 12
2

i=

;iki Z
eiki Z +pe2(2.11)

10and1 + P 
p 12
2

f

ikf R
;! e R 2~2 

R!1Z(2.12)

i~kf R~ + e;i~kf R~
e
3
0
p
dR
V (R0) eikiZ :
0

0

0

2

^ ).
Here, P12 acts only on the unprimed coordinates (recall ~kf = kR
Equation 2.12 still needs a little work. First, note that the integral in Eq.
2.12 is symmetric in the Z 0 coordinate. Adding an integral over e;iki Z and thendividing by two, gives the following properly symmetrized expression
0f

eikf R  
;! p
2~2

R!1

Z

2R

d3R0

ei~kf R~ p+ e;i~kf R~
2
0

0

!

V (R0 )

 eikiZ + e;ikiZ 
0

0

p

2

(2.13):

Comparing Eq. 2.13 with Eq. 2.6, we find the following relationship for the symmetrized scattering amplitude

f (; ) + f ( ; ;  + ) = ; 2~2 M fi
whereM fi =

Z

d3R0

ei~kf R~ p+ e;i~kf R~
2
0

0

!

V (R 0 )(2.14)

 eikiZ + e;ikiZ 
0

p

0

2

: (2.15)

Mott and Massey[24] define the identical particle differential cross section as theratio of the scattered radial flux to the incident flux in a single beam. The result isd = kf jf (; ) + f ( ; ;  + )j2 :d f ki(2.16)

This expression for the differential cross section will be rewritten in terms of symmetrized T -matrix elements in subsection 2.1.2. The numerical solution of equation
2.3 to obtain these symmetrized T -matrix elements is the subject of chapter 3.
We are interested not only in the scattering cross section but also in theevent rate for collisions. To that end, consider the scattering of two particles froma different perspective. Place the particles into a large box of dimension L on eachside. The Born approximation treats scattering of particles in a box as a transition
between discretized eigenstates of the box (see, for example Ref[26]). The totaltransition rate ;fi from eigenstate i ! f is an integral of the state-to-state transition

11
rate Wfi over all distinguishable final states. The “Fermi golden rule” can be used
to determine Wfi , resulting in

;fi =

Z

kf2 d f dkf

 L 3
2

Wfi

 

Z
3
L
2

2
2
=
~ d f kf dkf jMfij (Ef ; Ei) 2

(2.17); 3 represents a density-of-states factor and the matrix element M iswhere 2L
given byfi!

Mfi = h f jV j ii =Z!

i~kf R~;i~kf R~
i~ki R~;i~ki R~
d3R e p+ e3
V (R) e p+ e3:
2L
2L(2.18)The integral over dkf in Eq. 2.17 is rewritten as an energy integral using dEf =~2kf dkf =, which gives;fi = 42~4L3Z

d f kf jM fij2

(2.19)

 fi = L3 Mfi . Equation 2.19 represents the probability per second for thewhere M
two particles to make a transition from state i to all final states f . The event rate
constant K is defined as the transition probability per pair of particles to make thetransition from state i to state f per volume per second. It is obtained by dividing
Eq. 2.19 by the number of pairs per unit volume in our box 1=L3 and is given interms of the scattering amplitude by

K2 =

Z

d f kf jf (; ) + f ( ; )j2 = ki  :(2.20)

2.1.2 Relating the scattering amplitude to the scattering matrix
The symmetrized scattering amplitude f (; ) + f ( ; ;  +  ) can now
be expressed in terms of a corresponding symmetrized T -matrix. In general, it isnecessary to include the spin degrees of freedom. I will consider the four different
possibilities of scattering identical particles in distinguishable or indistinguishable
initial and final spin states. The terminology “distinguishable” and “indistinguishable” is reserved for the spin degrees of freedom. Although the scattering of identical particles is considered here, an experiment can still identify these particles
according to their spin state. Identical atoms in the same spin state are therefore
truly indistinguishable. To simplify the presentation, consider collisions between

12
two spin-1 particles. The symmetrized initial plane waves are eikiZ + e;ikiZ pj10; 10i
2
iki Z+ e;ikiZ j1 ; 1; 11i
D = e j11; 1 ; 1i p
i
2
I
i

=

and final state wave functions are(2.21)!

i~kf R~ + e;i~kf R~ep=j10; 10i2
i~k R~+ e;i~kf R~ j1 ; 1; 11i
D = e f j11; 1 ; 1i pf
2If

(2.22)

where the I or D superscripts indicate indistinguishable or distinguishable spinstates. The spin ket notation is jfa ma ; fb mb i. The plane waves can be expanded
using the following identities[27]

eikZ

=

ei~kR~

=

1
X
l
l=0

i (2l + 1)Pl (cos )jl (kR)

1 X
+l
X
l=0 m=;l4il Ylm (f ; f )Ylm(; )jl (kR) :(2.23)

Analogous expansions of the complex conjugates of the above plane waves introduce additional (;1)l factors. Inserting the plane wave expansions, the initial andfinal state unperturbed wave functions become

and

I
f

=

D
f

=

I
i

=

D
i

=

p X li2

Xlil=evenili

4ilf Ylfmf (f ; f )Ylf mf (; )jlf (kf R)I

f =even mf
XlX
lf 

lf mf

(2.24)i (2li + 1)Pli (cos )jli (kiR)D

p X X
2

i (2li + 1)Pli (cos )jli (kiR)I

4i Ylf mf (f ; f )Ylf mf (; )jlf (kf R)D



(2.25)

where I = j10; 10i and D = p12 j11; 1 ; 1i + (;1)l j1 ; 1; 11i represent symmetrized spin basis functions. Note that the expansions for the indistinguishable

13
spin states are summed over even values of the partial wave. In the case of identicalfermions scattering in indistinguishable spin states the sums would include only oddl’s.
The first example will be outlined in detail. The steps are the same in allcases, so I will only quote final results for the last three. For concreteness consider
a spin-dependent isotropic potential V S (R);i) indistinguishable spin state ! indistinguishable spin state

X lf +li 0M fi = h If jV S j Ii i = 8i (2li + 1)Ylf mf (f ; f )Z

li ;lf ;mf

d P (cos )Y lf mf (; )liUsing Pl (cos  ) =

q 4

Z

(2.26)

R2dRjlf (kf R)hI jV S jI ijli (kiR) :2l+1 Yl;0 and the orthogonality properties of the spherical harmonics, the angular integral ( d = 4 ) can be evaluated givingM fi = 8RXZli =evenli ;lf(2l + 1)Pli (cos f )R2 dRjlf (kf R)hI jV S jI ijli (kiR) :(2.27)q2 kRj (kR), the radial
Introducing energy normalized radial wave functions[27] klintegral can be related to a scattering S - or transition T -matrix element[27]fT gfi = fS gfi2i; fi
= ;

Z

s

r

(2.28)

2 k Rj (k R)hI jV S jI i 2 k Rj (k R) :
dR k
f l f
ki i l if

The brackets on fT gfi indicate that the matrix element is calculated in a symmetrized spin basis. The final result for the matrix element is then given byXM fi = ; 4 p 1(2l + 1)Pli (cos f )fT gfi
kikf li =even iii) indistinguishable ! distinguishable(2.29)

14pX(2l + 1)Pli (cos f )fT gfi
M fi = ; 2 2 p 1
kikf li =even i(2.30)

iii) distinguishable ! indistinguishablepX(2l + 1)Pli (cos f )fT gfi
M fi = ; 2 2 p 1
kikf li =even i(2.31)iv) distinguishable ! distinguishable

1
X
2

1

(i)2li (2li + 1)Pli (cos f )fT gfi
Mfi =  p
kikf li=0(2.32)

R

The last step required to obtain a total cross section or event rate constant is evaluation of the integral over final states d f in equations 2.16 and 2.20. In the caseof distinguishable final spin states this integral is taken over all space. However,for indistinguishable final spin states we must restrict the integral over half space( d f = 2 ) to avoid double counting. The results for the four cases are quoted
below. To derive these, one must keep in mind the caveat on the angular integration
limits.
i) indistinguishable ! indistinguishable
 fi from equation 2.29 into equations 2.16 and 2.20 and using the orthogPlugging M
onality properties of the Legendre polynomials givesR = 8
K2 = 8X

li =even(2li + 1) jfTkg2fij

li =eveniXi(2li + 1) jfTkgfij

22(2.33)

This cross section expression can be checked in the s-wave limit. The scattering
length is defined as limki !0 a = (; tan =k ) [24], where  is the elastic phase shift.
In the zero energy limit, fT gfi !  and we find  = 8a2 as expected. In practice, I
calculate the scattering lengths from the S -matrix. The diagonal S -matrix elements
are related to the phase shift[24] by fS gii = e2ii . The scattering length in channel
i is then given by

lnfS gii :
lim
a
=
;
i
k !0
2ik
i

i

(2.34)

15ii) indistinguishable ! distinguishable

X(2li + 1) jfT g2fij
ki
li =even
2X
K2 = 8(2li + 1) jfT gfij
ki
li =even = 8

2(2.35)

iii) distinguishable ! indistinguishable

X(2li + 1) jfT g2fij
ki
li =even
2X
K2 = 4(2li + 1) jfT gfij
ki
li =even = 4

2(2.36)iv) distinguishable ! distinguishable

 = 4

1
X

2
jf
T
g
fij
(2l + 1)

i

ki2
1
X
jf
T
g
fij2
K2 = 4 (2li + 1) k
i
li =0
li =0

(2.37)In the case of non-identical particle collisions, the cross section and event rate formulas are given by case iv) and the T -matrix element is calculated in an unsymmetrized spin basis.
There is one more identical particle detail that must be considered. The
event rate K2 describes the rate at which pairs of atoms collide, whereby K2 inelastic events will generally lead to 2K2 atoms lost in unit trap volume per unit time.
Generally, the most useful parameter is the atom-loss rate constant, the coefficientof the density squared in the rate equation. This quantity is L2 = K2 =(1 +  ). The
denominator (1 +  ) is needed because the number of pairs depends on the initial
spin state. If the initial spin states are the same  = 1, otherwise  = 0. For example, a gas of N identical atoms initially in indistinguishable spin states will have
N (N ; 1)=2 pairs as opposed to the distinguishable case which will have NaNb
pairs. However, these additional factors of two are easily accommodated since K2
in equations 2.33 and 2.35 is a factor of two larger than in equations 2.36 and 2.37.

16
The definition for the atom loss rate constant will therefore be taken as(L2 )fi  4 (2li + 1) jfTkgfij :2(2.38)i

The rate equation controlling the number of

dn = X X(1 +  )[L (  !
;
2
dt
f g

spin state atoms would thus read:)n n ; L2 (!  )n n ](2.39)

where n is the density of atoms in the spin state. There is some confusion in
the literature regarding this terminology. Generally, the published results that weretermed “event rate constants”, e.g. in Ref.[28], are understood to be the coefficientof the density squared in the rate equations, (i.e., what we call the atom-loss rate).
However, our definition of L2 used in the above rate equation agrees with the “standard” definition of K2 .
The last situation that needs to be examined is the case of an anisotropic spin-dependent potential V S (R; ). Returning to the first indistinguishable
spin isotropic example (case i), equation 2.26 can be generalized as follows

X li+lfM fi = 8i (2li + 1)Ylf mf (f ; f )(2.40)

li ;lf ;mf

r 4 Z2li + 1R2 dRjlf (kf R)hI Ylf mf ( )jV S jYli0 ( )I ijli (kiR) :The new expressions for the perturbative S - and T -matrix elements arefT glfif mf li=

fS glfif mf li ; fi

= ;

Z

(2.41)

2i s

2 k R j (k R)dR kf lf ffr 2

hI Ylf mf ( )jV S jYli0 ( ) I i k kiR jli (kiR)iwhich can be substituted back into the equation above to give

X XXp
M fi = ; 4 p4
2li + 1 Ylf mf ( f )fT glfif mf li :kikf lf =even mf li=evenp(2.42)

 fi for cases ii) and iii) will again be a 2 smaller than the
The matrix element M
one given in Eq. 2.42 and case iv) a factor of two smaller. The cross sections and17event rates are obtained in the same fashion, although one should keep in mind the
caveats on the d f integration limits. The final results for the anisotropic potential
can be summarized as follows

 = 4

Xp
lf mf li2li + 1

fT glf mf li
ki

2(1 + )(2.43)

and2

4 n X p2l + 1 fT glf mf li :L2 = kii lf mf li(2.44)An anisotropic potential couples the different partial waves and therefore these must
be added coherently (with the appropriate limits on the sum over the initial and final
partial waves). Finally, the  function in Eq. 2.43 is one if the two atoms are initiallyin the same spin state and zero if the two atoms are initially in distinguishable spinstates.
2.2 Interaction potentials
The complete description of two interacting multi-component atoms isextremely complex. Levy and co-workers[29, 30] have identified a fourteen term
Hamiltonian valid for diatomic molecules in the Born-Oppenheimer limit[31]. However, most of these interactions are vanishingly small for two ground state alkali
atoms. In fact, most of the interesting cold collision physics comes about simply
because the collision energies are much smaller than the atomic hyperfine splittings. Consequently the collision partners can in some cases resonate (i.e. at a
Feshbach resonance in asymptotically closed channels). Moreover, if the partners
change hyperfine states, enough energy can be liberated to eject the atoms fromthe trap. Or if the trap is magnetic, the collision can produce untrapped hyperfine
states and again the atoms are lost. These hyperfine changing collisions are generally termed spin-exchange collisions, although one should keep in mind that otherinteractions besides electronic spin-exchange can also change the hyperfine states.
True spin-exchange collisions can be viewed as resulting from a transformation between angular momentum coupling schemes. At small internuclear separations R,
the two electronic spins tend to lock and precess around the total angular momentum
vector. At large R the electrostatic exchange splitting decays exponentially and theelectron spin on each atom now prefers to couple to the corresponding nuclear spin,
precessing around the total atomic angular momentum vector of each atom. However as the atoms recede, the small R evolution of the total electronic spin states
will not in general allow the spins to recouple in exactly the same hyperfine states

18
that initially approached. Spin-exchange coupling is generally strong, leads to large
trap loss rates, and is something to be avoided. This can be arranged either throughthe choice of symmetry (e.g., by trapping doubly polarized spin states, also known
as “stretched states”) or else through energetic considerations. A second and much
weaker loss process frequently occurs via the magnetic dipole-dipole interaction between the two valence electrons. It can usually be neglected unless one needs theloss rate from a spin-exchange forbidden collision.
Complete descriptions of these interactions are provided in this section.
Subsection 2.2.1 describes the construction of the Born-Oppenheimer potentials.
The hyperfine Hamiltonian is presented in subsection 2.2.2 and the magnetic dipoledipole Hamiltonian in subsection 2.2.3. Finally, magnetic field induced Feshbachresonances will be discussed in subsection 2.2.4.
2.2.1 Constructing the Born-Oppenheimer potentials
The motion of the constituent charged particles of the diatom is controlled
by two very different time scales. The electrons are much lighter than the nuclei
and therefore move at greater velocities. The Hamiltonian describing the complete
system, H = HN + He , can be approximately separated into a nuclear and electronic component. This is known as the Born-Oppenheimer limit[31]. The BornOppenheimer potential curves are constructed by freezing the interatomic separation
at a given value R, and then evaluating the energies of the resulting electronic eigenstates of He . These energies vary with R and constitute a set of potential energy
curves. The physics underlying these potentials is qualitatively different in regionsof small and large internuclear separations. The division between these two regions
is determined by the overlap of the atomic electronic clouds, which begins to become exponentially small beyond R  20 a.u. in the alkalis. Information regardingthese two regions is generally collected from different sources and must then besmoothly connected through the transition region near R  20 a.u.
The small R region is characterized by large overlaps of the electronic
clouds; here one must use molecular concepts to correctly describe the interactionof the collision pair. In the case of alkali atom ground state collisions, the potential energy depends solely on the internuclear separation R and the orientation of
the two atom’s valence electronic spins (si =1/2) which couple into singlet (S =0),
~ = ~sa + ~sb) or triplet (S =1) configurations. The Hund’s case (a) desigwhere (S3 +nations for these states are X1 +
g and a u for the singlet and triplet respectively.
The spin dependence of the potential energy is a consequence of the Pauli exclusion
principle which requires the wave function to have an odd permutation symmetry
under exchange of electrons[32]. As one can imagine, ab initio calculations to determine the interactions of multi-electron alkali atoms require sophisticated numerical techniques[33]. Fortunately, molecular spectroscopicists have devoted a great
deal of attention to measuring ro-vibrational spectra. In turn, these spectra can be19inverted to determine a corresponding potential using standard techniques such as
RKR[34] or inverted perturbation analysis (IPA)[35]. References for the small R
Born-Oppenheimer potentials used in this work are given in Table 2.1.
The physics simplifies as the overlap of the electron clouds decreases. The
residual overlap leads to a long range exchange term of the form V ex = AR exp; R
where + is used for the triplet state and ; for the singlet. Analytical values forand are derived in Ref.[36]. These expressions are given in terms of the valence
electron binding energy Eb , which for identical atoms becomes = 2p72E ; 1 andb= 2 . The value of the overall scaling coefficient A is determined from the difference VS =1 ; VS =0 near the transition R  20 a.u., (see for example Ref[37]). The
remaining large R electrostatic interactions are a result of “instantaneous” multipole
moments. The motion of the electrons leads to small fluctuations in the charge density surrounding each atom. In turn, one atom can polarize the other setting up an
momentary dipole moment which then attracts the first. This is the van der Waals
or dispersion force. In the case of two alkali atoms colliding in their atomic groundstates, the dispersion potential is given by V d = ;C6 =R6 ; C8 =R8 ; C10 =R10 . Ab
initio calculations[38, 39, 40] for the dispersion coefficients Cn have proven reliable
to within a few %. Ref.[38] is used for the C8 and C10 coefficients. C6 coefficientsare taken from a variety of sources which are listed in Table 2.1. The dispersion potentials will also be modified by retardation effects. However, current experiments
have not been able to separate ground state retardation effects from uncertainties
associated with other interactions. Retardation effects are therefore neglected. (Retardation effects have been observed in pseudo-bound excited state measurements,
see for example Ref.[21]). The Rb Born-Oppenheimer potentials used in this work
are shown in Fig. 2.1.
Table 2.1: Sources of information for constructing the Born-Oppenheimer potentials. The singlet and triplet columns list the method used in the correspondingreference to obtain the small R part of the potential. The C6 coefficients are thenominal values used in this work and are cited along with the corresponding reference.
Atom Singlet
TripletC6
Li
RKR[41]
RKR[41]1393.4[42]Na
RKR[43]
RKR[43]1538.9[44]K
IPA[45]
RKR[46]3897.0[39]Rb
IPA[47] ab initio[33] 4700.0[17]Cs
IPA[48] ab initio[33] 6851.0[39]The potentials as given are still not accurate enough for ultracold collision calculations. Typically, the highest bound state measured with conventional
molecular spectroscopy is more than 2 cm;1 below the dissociation threshold. Ex-

20

Rb2 Potentials
50.0
3

+

Σu

25.0
0.0

Thresholds

E (THz)

1

−25.0

+

a)

Σg

0.0
−50.0

2+2
∆87

2+1

∆87

1+1

b)

−75.0
0.0

∆85

∆87

∆85

−100.01000−125.0
0.0

10.0

3+2
2+2
3+1
2+1

1050
20.0

30.0

R (a.u.)Figure 2.1: Singlet and triplet Rb2 Born-Oppenheimer potentials. Hyperfine split
thresholds are indicated in the inset for a)two 87 Rb atoms and b) 87 Rb + 85 Rb collision. The atomic hyperfine quantum numbers fa + fb are indicated at each threshold.
The values of the energy splittings are 87 = 6.835 GHz and 85 = 3.036GHz

21
trapolating the Born-Oppenheimer potentials beyond this point has proven unreliable. In order to obtain collision parameters which agree with experiment we add acorrection to the inner wall of the singlet and triplet potentials of the form[49]V (R) =(CS tan;1h; R;Re 2i
R0; R > Re; R < Re

(2.45)

where Re denotes the minimum of the well and R is a width parameter generally
taken to be 2 Bohr. The energy dependence of the phase shifts has been tested
using other functional forms for the correction, such as quadratic, cubic, and quartic
polynomials in R, all of which give nearly identical results.
2.2.2 Atomic hyperfine Hamiltonian
The atomic hyperfine Hamiltonian represents the coupling of an electronwith the electromagnetic field created by the charged nucleus. In general, for nuclei
with spin i > 1/2 the Hamiltonian is quite complicated. For example, electric andmagnetic multipole terms of order k  2i are allowed[50]. However, expectation
values of these multipole terms average to zero for spherically symmetric atoms(L=0) which will be the case for alkali atoms in their ground state. The lone remaining term, commonly referred to as Fermi’s “contact term”, is a result of the
non-zero spatial extent of the nucleus and can be expressed as[51]H ahf = ;38 j e(0)j2 h~e  ~ni

(2.46)

H ahf = C (f~a2 ; ~sa2 ; ~ia2 )(2.47)

where ~e and ~n are the magnetic moments of the electron and nucleus, respectively. Evaluation of H in this representation is complicated by the existence of theelectronic wave function term j e(0)j2 . We can avoid this problem by introducing
measured hyperfine splitting energies. (The review paper by Arimondo, et al.[52]is an excellent resource for this information.) H can be rewritten in the following
form

where f~a = ~sa + ~ia and C is a constant which can be related to the hyperfine
splitting . H ahf is diagonal in the basis jfa (sa ia )mfa i. Requiring the barycenterof the alkali hyperfine Hamiltonian to be E =0, C is given by

C =  2i+ 1 :(2.48)

Here,  = +1 if the nuclear magnetic moment ~n is positive or  = ;1 if it isnegative. The sign determines whether the hyperfine structure is “standard” (i.e.,22

E (f = i + s) > E (f = i ; s)) or “inverted” (E (f = i + s) < E (f = i ; s)).
The total hyperfine Hamiltonian H hf for the collision complex taken to be simply
the sum of the two atomic hyperfine Hamiltonians H ahf(2.49)
H hf = a 2i+a 1 (f~a2 ; ~sa2 ; ~ia2 ) + b 2i+b 1 (f~b2 ; ~sb2 ; ~ib2):ab
The hyperfine split thresholds for both homonuclear and heteronuclear Rb collisions
are shown in the inset of Fig. 2.1. The intra-atomic hyperfine Hamiltonian as wellas any R-dependence in H hf are neglected.
2.2.3 Magnetic dipole interaction
The magnetic dipole Hamiltonian can be constructed by treating the interaction of a magnetic dipole (i.e., valence electron a) with the field created by asecond dipole (valence electron b). Starting with the classical expression for thisinteraction[51] and substituting quantum mechanical expressions for the electronic
moments leads to the following relationship (in a.u.)

H ss = ;~  B (~r) = ; 2

 3(^n  ~s )(^n  ~s ) ; ~s  ~s 
a

b
3
Ree

a

b

(2.50)

where is the fine structure constant. The vector ~r is drawn from dipole a todipole b and is defined as Ree n^. For our purposes, Ree can be approximated as
the internuclear separation R. We do not need to include the s-wave contact term
(which is important for the hyperfine interaction) because the two electrons repel
each other. Therefore, the electron pair wave function is exponentially small in thelimit Ree ! 0. H ss is most conveniently evaluated in the molecular frame with n^
defining the spin quantization axis. In the spin representation jS i, the Hamiltonianis diagonal, with matrix elements given by[53]

2

;



VS = ; 2R3 32 ; S (S + 1) :(2.51)

Here,  = a + b represents the total projection of the electronic spin onto theinternuclear axis.
The effect of the magnetic dipole term in the molecular frame is to lift
the spin projection degeneracy of the triplet state. This behavior is mimicked by anadditional effect known as the second-order spin-orbit interaction[53, 29, 54]. The
second-order spin-orbit interaction becomes important when the electronic clouds
overlap such that indirect perturbative coupling to excited electronic molecular states
can modify the direct spin-spin interaction. For ground state alkali collisions, only;
; +the i (L+i si + Li si )=2 part of the single electron spin-orbit Hamiltonian[53, 54]
contributes to the interaction. Here, L+ S ; are the standard raising and lowering op-P23
erators for the ith electrons orbital angular momentum and spin, respectively. Mies,et al.[54] have shown that the dominant interactions come from coupling the triplet
Born-Oppenheimer 3 u state with nearby excited (2S +1) u states. The interaction isshort-ranged and its strength scales with Z . The inclusion of this term is important
only for the heavier atoms such as Rb and Cs. Analytical expressions for Rb and Cs
have been derived in Ref[54].
Before moving on, it is instructive to examine the operator structure ofss
H in the lab frame. Derivations are given in Refs.[29, 50] and accordingly only
the final result is given here

p

"

#

2 X
H ss = ; R63(;1)q Cq2(s1a s1b )2;q :q(2.52)(s1a s1b )2;q represents the contraction of two first-rank spin tensors into a second-qrank tensor. The important result for cold collisions is the introduction of an angulardependence through the term Cq2 = 45 Y2q (; ). The second-rank nature of thespherical harmonic implies that incident s-waves will couple only with d-wave exit
channels. However, the H ss operator is an overall scalar and thus is rotationallyinvariant. Therefore, the squared total angular momentum F~ 2 (F~ = f~ + ~l) stillcommutes with the zero-field Hamiltonian.
2.2.4 External magnetic field and Feshbach resonances
The cold collisions community has energetically studied the possibility of
influencing cold collisions with external fields. In principle this would allow condensate properties to be tailored. Although several methods have been proposed[55,56, 14] probably the simplest and most robust approach is to induce a Feshbachresonance in the collision of two atoms by applying a static magnetic field. The
scattering length can be varied between 1 as the external magnetic field is swept
through the resonance field value, as has been demonstrated experimentally now, inNa condensates[15] and in thermal Rb[16, 17] and Cs[18] clouds.
The physics is rather straightforward. First, the magnetic field Hamiltonian is

H B = ;B ge(saz + sbz ) + N(gNa iaz + gNb ibz ) :(2.53)~ = B z^ is assumed to be constant and in the direction z^ that
The magnetic field B
defines the spin quantization axis. The constants ge and gN N include the gyromagnetic ratio g and magneton e = ;e~=2me for the electron and n = e~=2mP
for the nucleus. The individual atom spins are again labeled by a and b. The effect of Eq. 2.53 is two-fold: first the spin projection degeneracy is lifted whichshifts and splits the final state thresholds; second, the magnetic field couples the24entrance channel to additional spin configurations. H B conserves the total angularmomentum projection mF but not its squared magnitude F~ 2 . 85 Rb atomic hyperfine thresholds and the thresholds for a j2; ;2i85 + j2; ;2i85 collision complex are
shown as a function of magnetic field in Fig. 2.2. (The notation jf; mf iA will be
used throughout to define the isotope (A) and hyperfine spin state of each colliding atom.) A Feshbach resonance occurs when the collision energy becomes sufficiently close to the energy of a quasi-bound state associated with a closed channel
that is coupled to the entrance channel. As a first approximation, one can consider
these quasi-bound states as being “attached” to a particular threshold which are thenshifted by the magnetic field. A schematic of a magnetic field induced Feshbachresonance is shown in Fig. 2.3. Two situations can occur. A bound state can either be added or removed from the potentials; either case results in a change of thephase shift by  radians as a function of B . A scattering length that approaches;1 as B increases toward the resonance value BR (from low-field) indicates thatthe resonance came “down” out of the continuum (i.e., an additional bound state isadded for B > BR ). This is the general scenario for magnetically trapped atoms as
observed in Ref.[16, 17, 15]. The entrance channel threshold rises with increasing
field, which is why they are called “weak-field-seeking” (see Fig. 2.3). Strong-fieldseeking states can be trapped with optical methods; their threshold energies decrease
with increasing field strength, and the most likely scenario is therefore a bound statethat pops above the entrance threshold at B > BR . In this instance, the scattering
length a ! +1 as B approaches BR from below. This was observed for Na[15]and Cs[18] collisions. The scattering length behavior near the resonance for bothsituations is shown in Fig. 2.4. The scattering length near a Feshbach resonance isoften parameterized with the following expression[15]a = a0 1 ; B ;B(2.54)
R
where  is the width of the resonance and a0 is the zero-field value of the scattering
length.
2.3 Angular momentum representations
Different angular momentum coupling schemes play a key role in different aspects of the collision physics. In this section, I will define the different
coupling schemes and outline the advantages of each. Subsection 2.3.1 discusses
the symmetrization of each representation and presents the corresponding identicalparticle selection rules. Subsection 2.3.2 presents the unitary transformation matrix elements to transform back and forth from the short range representation to thevarious asymptotic representations.
Any collision event begins with atoms trapped in their individual atomic

25

0.2
a)

mf=3

E (K)

0.1
0
mf=−2

−0.1
−0.2

0

200

400 600
B (G)

800

1000

0.2
b)
|3,−2>+|3,−2>

E (K)

0.1

|3,−3>+|3,−1>0|3,−2>+|2,−2>|3,−3>+|2,−1>

−0.1
−0.2

|2,−2>+|2,−2>

0200

400 600
B (G)

800

1000

Figure 2.2: a) 85 Rb atomic thresholds versus magnetic field. The projection of fa islabeled for the highest energy state in each hyperfine manifold.
b) Thresholds versus magnetic field for a collision of two j2; ;2i85 atoms. The channel labels are given in terms of the atomic hyperfine quantum numbers jfa mfa i +jfbmfb i.

26

0
3

E (K)

−0.05

2

−0.1
1
−0.15
BR
−0.2

0

100

200

300
B (G)

400

500

600

Figure 2.3: Schematic of a Feshbach resonance. The resonance state (dashed lines)
is attached to threshold 2. The collision energy is relative to threshold 1 and becomesdegenerate with the resonance state energy at the field value labeled BR . In thisexample, the scattering length would approach ;1 as B is increased to BR .

27

3000.0a2,−2 (a.u.)

2000.01000.085

Rb

0.0−1000.0−2000.0−3000.0
135.0

155.0
B (G)

175.0

500
23

a1,+1 (a.u.)

250

Na

0
−250
−500
895

900

905
B (G)

910

915

Figure 2.4: Scattering length versus magnetic field in the presence of a Feshbachresonance.
a) Feshbach resonance observed[17] in the collision of two 85 Rb j2; ;2i atoms.
The scattering length behavior on resonance indicates a bound state is added to thepotentials.
b) Feshbach resonance observed[15] in the collision of two 23 Na j1; +1i atoms. The
scattering length behavior on resonance indicates a bound state is removed from thepotentials.

28
hyperfine state f~i = ~si + ~ii (i  fa; bg). This provides the first natural choice for anangular momentum representation and I will refer to it as theuncoupled hyper ne representation  jfamfa ; fbmfb ; lml i :Here, l represents the nuclear rotation (or partial wave) and m the projection of
the angular momentum vector on a space-fixed quantization axis. The collision
dynamics are simpler in a coupled representation. In the absence of a magnetic field(and when dipolar relaxation can be neglected), both l2 and f 2 (where f~ = f~a + f~b )are conserved quantities. The coupled Schrödinger equations therefore have a block
diagonal structure in the

coupled hyper ne representation  j(fafb )fmf ; lml i :This representation is still useful in the presence of a magnetic field that is a sufficiently small perturbation. If the dipole interaction is included in the Hamiltonian,thetotal angular momentum representation  j(fl)FmF i
is a better basis, because in the absence of a magnetic field the Hamiltonian is again
block diagonal.
The Hamiltonian expressed in any of the three representations discussed
so far becomes diagonal as R ! 1 (again neglecting magnetic field). In the presence of a magnetic field, the asymptotic Hamiltonian H (R ! 1) = H hf + H B
is diagonalized. The resulting eigenvectors U are used to transform the complete
potential matrix U T V (R)U into this eigenchannel representation. The eigenvalues
provide the channel thresholds. For the field strengths I will consider, the eigenchannel representation is very nearly the uncoupled hyperfine representation and will bereferred to as such. These representations are the most suitable for solving the coupled Schrödinger scattering equations since the asymptotic boundary conditions canbe easily applied. However, at small R the dominant interaction is provided by the
Born-Oppenheimer molecular potentials. These are diagonal in a total electronic~ = ~sa + ~sb) which I will refer to as thespin basis (Sshort range representation  jSmS ; ImI ; lml i :To help identify nuclear permutation symmetry, the nuclear spins are also coupled
to form a resultant I~ = ~ia + ~ib . The short range representation can also be written~ + I~, but generally what controls most of the coupledin terms of total spin f~ = S
channel dynamics is the projection of the hyperfine states onto the singlet and triplet
Born-Oppenheimer states.

29
2.3.1 Symmetries
Collisions of identical particles are of paramount interest and as such wemust consider the symmetrization of the total wave function in accord with the Pauliprinciple. The radial piece of the total wave function (Eq. 2.2) is automaticallysymmetric since it depends only on the magnitude of R. We therefore need only toconsider the effect of the symmetrization operator (1 + (;1)p P12 ) on the angularmomentum representations presented in the previous section. The (;1)p phase isincluded to handle both bosonic (p = 0) and fermionic (p = 1) symmetries. In
the center-of-mass representation, the spatial part of P12 is equivalent to the parity
operator, which applied to the nuclear rotation part of the wave function (jlml i =Ylm(; )) yields

P12jlml i = (;1)l jlml i :(2.55)The symmetrization of the uncoupled hyperfine representation is therefore straightforward(1 + (;1)pP12)jfa mfa ; fbmfb ; lml i = jfamfa ; fbmfb ; lml i +(2.56)p+l(;1) jfbmfb ; fa mfa ; lml i :Note that this relationship immediately implies that collisions of identical atoms inthe same hyperfine state (i.e., indistinguishable spin states) are allowed only for evenpartial waves for bosons and odd partial waves for fermions. The symmetry of thecoupled hyperfine representation can be obtained by first decomposing it into theuncoupled representation and applying P12(1 + (;1)pP12 )j(fafb )fmf ; lml i = (1 + (;1)pP12)

X

mfa ;mfb(2.57)

hfa mfa ; fbmfb jfmf ijfamfa ; fbmfb ; lml iUsing the Clebsch-Gordan identity[57]hfa mfa ; fbmfb jfmf i = (;1)fa+fb;f hfbmfb ; fa mfa jfmf i(2.58)

and recoupling the hyperfine momenta gives our final expression(1 + (;1)pP12)j(fa fb)fmf ; lml i = j(fa fb)fmf ; lml i +(2.59)
f
+
f
;
f
+
l
+
p
a
b
(;1)
j(fbfa )fmf ; lml i :30
In the case fa = fb definite selection rules are obtained in terms of f and l. A
similiar procedure can be used for the short range representation and one finds(1 + (;1)pP12 )j(sasb )Sms; (ia ib )ImI ; lml i =j(sasb )Sms; (ia ib )ImI ; lml i
s
+
s
;
S
+
i
+
i
;
I
+
l
+
p
a
a
b
+(;1) bj(sbsa )Sms; (ibia )ImI ; lml i :(2.60)P

The phase factor can be simplified by considering the case of identical alkali atoms(si = 1=2; ia = ib ) and noting that i (si + ii + p) will always be even. In thiscase, the short range representation is manifestly symmetric with states allowed
only when the phase factor (;1)S +I +l is equal to +1. Finally, the symmetry of thetotal angular momentum representation is given by the symmetry of its components(1 + (;1)pP12) j((fafb)fl)FmF i = j((fafb)fl)FmF i +(2.61)
f
+
f
;
f
+
l
+
p
a
(;1) bj((fbfa )fl)FmF i :2.3.2 Transformations
The four representations presented in the last section are connected byunitary transformations. Transformations between the different asymptotic hyperfine basis require one or at most two recoupling steps. However, connecting theshort range representation with the various asymptotic representations requires a little more effort. The derivations are not tricky, as they require only standard angularmomentum algebra techniques and consequently, only the results will be quoted
here.
i) uncoupled hyperfine representation $ short range representation

fhSmsImI l0 ml jfamfa fb mfb lml ig =0Xl;l ml ;ml
0

0

f;mf

(2.62)

hSms; ImI jfmf ihfamfa fb mfb jfmf ip(2fa + 1)(2fb + 1)(2S + 1)(2I + 1)

0s i f 1!
a a a
S +I +l
1
+
(1
;

)(
;
1)
f
;f
@ sb ib fb A
p ab
:
S I f

2 ; fa;fb

. . .

The fg brackets again indicate that this is a symmetrized matrix element and .. .. ..
represents a 9-J symbol.
ii) coupled hyperfine representation $ short range representation

31

fhSmsImI l0 ml j(fafb)fmf lml ig = l;l ml ;ml hSms; ImI jfmf i
0

p

0

0

(2fa + 1)(2fb + 1)(2S + 1)(2I + 1)

0s i f 1
!
a a a
S
+I +l
1
+
(1
;

)(
;
1)
@ sb ib fb A
p fa;fb:(2.63)

2 ; fa;fb

S I fiii) total angular momentum representation $ short range representation

fhSmsImI l0ml j(fl)FmF ig = l;l ml ;ml
0

0

X
0

(2.64)

mf ;mlhSms; ImI jfmf ihfmf lml jFmF ip(2fa + 1)(2fb + 1)(2S + 1)(2I + 1)

0s i f 1!
a a a
S +I +l
1
+
(1
;

)(
;
1)
@ sb ib fb A
p fa;fb:
S I fp

2 ; fa;fb

The unsymmetrized versions of these transformations are obtained by simply neglecting the (1 + (1 ; fa ;fb )(;1)S +I +l )= 2 ; fa ;fb term in each equation.

Chapter 3
Numerical techniques for solving the coupled Schrödinger equations

Solving coupled linear differential equations like Eq. 2.3 is a problem
found in all branches of physics and as such, extensive literature exists on the subject. In most cases one must resort to numerical techniques to solve these equations.
Many different numerical methods have been developed over the years and can bebroadly divided into two categories, explicit (sometimes called “propagators”) and
implicit approaches. Both approaches divide the domain of the independent variable
into grids or sectors. Explicit methods require the solution of the dependent variable
at the preceding grid point to determine its value at the next adjacent point. Implicit
methods typically require the solution of a large set of linear equations (or an eigensystem) and thereby determine the value at all grid points simultaneously. Explicit
approaches generally are easy to code, fast, and not memory limited since the solutions are calculated “on the fly”. The Numerov[58], Gordon[59], and log derivative
propagator[60] methods are examples commonly used in atomic collision physics;
a performance review of these algorithms is given in Ref.[61]. However, implicit
methods are inherently more stable[62] and often give better accuracy. Extremely
good linear algebra packages (e.g., LAPACK[63]) are now commonly available (andfree!). These packages are cpu-efficient and vastly reduce the programming complexity associated with the linear algebra steps. In addition, with recent advances incomputer technology, memory is not the prohibitive factor it once was.
In this chapter, I present an implicit method to solve the coupled-channel
Schrödinger equations. Section 3.1 develops the basic method which is based on a
finite-element (FEM) R-matrix approach. Section 3.2 connects the R-matrix to thephysical observables given in terms of a scattering S -matrix. Section 3.3 develops
an alternative adiabatic representation to describe the coupled equations. This ispresented mainly as a qualitative tool for understanding two-body collisions. However, the power of the adiabatic FEM R-matrix approach is brought to bear to solve
the three-body collision equations presented in Chapter 8.

33
3.1 Finite element R-matrix approach
The R-matrix method first introduced by Wigner and Eisenbud[64, 65, 66](1947) has been extended considerably over the years. The variant presented in Section 3.1.1 is the noniterative eigenchannel R-matrix[67, 68, 23] which is based on
the familiar Rayleigh-Ritz variational expression for energy eigenvalues. This is atechnique which provides a variational solution for the log derivative of a wave function on a given boundary. The R-matrix method has been used extensively in photoionization studies which are sometimes dubbed “half-collision” processes. Section 3.1.2 describes the finite element basis (FEM) used to solve the R-matrix equations. The use of finite elements to solve differential equations in atomic physics is
relatively new[69, 70]. This is one of the first implementations of a FEM R-matrix
approach (see also Ref.[71]), although a similiar approach using basis splines (Bsplines) has been presented elsewhere[72]. The last section 3.1.3 describes a technique to match R-matrix box solutions to the expected long-range form of the coupled channel solutions. This is useful if a prohibitively large number of channels or
sectors are needed to describe the collision complex, such that computer memory
becomes short. In this event, the problem can be side-stepped by dividing the radialdomain into several boxes and solving the equations in each box separately. An
appropriate transformation equates the asymptotic R-matrix with the individual box
solutions.
3.1.1 R-matrix equations
The eigenchannel R-matrix method solves the coupled Schrödinger equations (i.e., Eq. 2.3) within a finite reaction volume of configuration space, subject
to constant normal logarithmic derivative boundary conditions on the surface of . The collisional properties of the system, typically represented in terms ofan S -matrix, are easily obtained by a straightforward matching procedure once thenormal logarithmic derivative b = ;(@ =dn) ;1 is calculated. The present implementation restricts the discussion to a single coordinate R which represents theinternuclear separation. However, the techniques presented here are readily generalized to 2- or 3- dimensions[71].
Rearranging the Rayleigh-Ritz expression, one can obtain the following
variational expression for the normal logarithmic derivative b of the solutions of
Eq. 2.3 on the surface [67]:R 2 (E ; H ) dw ; R  @ dS
R  dS  @n :
b=


(3.1)

Here E represents the relative collision energy,  the reduced mass, and H the
Hamiltonian of the colliding pair. A typical approach for solving the above equation
is to approximate using a basis set expansion. A multi-component wave function

34P

can be represented by = m Fm (R)Ym ( )m , which upon introducing a basis
set expansion for the radial wave function becomes=Xici i (R)Yi ( )i :
0

0

0

0

(3.2)

0

Here, Y ( ) represents a spherical harmonic,  the spin degrees of freedom, m is anindex containing all quantum numbers needed to uniquely define a channel, and i0
incorporates the basis function index with m. Arbitrary basis functions are designated by i , and the coefficients ci are to be determined numerically. Insertion ofthis expansion (Eq. 3.2) into Eq. 3.1, leads to the following generalized eigenvalue
equation00;~c = b~c :(3.3)

The matrix elements of ; and  will be described in the following subsection. In
general, Eq. 3.3 can be solved directly obtaining the log derivatives and eigenvectors
needed to describe the solution on the boundary . However in practice, the matrices in Eq. 3.3 are quite large and direct diagonalization is inefficient. It has beenshown[67] that partitioning the matrices according to whether the basis functions
are non-zero (open  o) or zero (closed  c) on  ;cc ;co  cc 
;oc ;oo0= b 0 0oo

 cc coco
Eq. 3.3 can be reduced to a small (mP  mP ) eigensytemoo c~o = b oo c~o

(3.4)

(3.5)

where oo = ;oo ; ;oc(;cc );1 ;co . The main computational burden is thus shiftedto solving a set of linear equations;ccX co = ;co(3.6)

whose solution X co = (;cc );1 ;co provides the needed matrix inverse, after which
Eq. 3.5 can be solved efficiently.
3.1.2 Finite element basis
The R-matrix equations 3.5 and 3.6 were derived using an arbitrary radialbasis expansion. A particularly flexible choice is the FEM basis[73, 74, 69]. The
FEM method divides the radial domain into nmax sectors (or elements) and withineach sector defines a local basis. The six local basis functions uk (xn ) are fifth order
Hermite interpolating polynomials which are non-zero only in sector n. Here, xn

35
is a rescaled variable defined on the interval [-1,1] which is related to the physicalinternuclear separation R through the transformation R = an xn + dn with sector
coefficients an = (Rn+1 ; Rn )=2 and dn = (Rn+1 + Rn )=2. The FEM expansion
for the radial wave function becomes

F (R) =

Xifk;m;ng

ciui(xn ) ;(3.7)

where the set i contains the basis function index k , the channel index m, and thesector index n. The six basis functions are defined through the following boundary
conditions (exact expressions are given in Appendix A):

uk (;1) = 1k uk (0) = 3k uk (1) = 5k
duk (;1) =  duk (0) =  duk(1) =  :
2k
6k
4k
dx
dx
dx

(3.8)

The requirement that each channel component and its first derivative must be continuous across sector boundaries imposes the following constraints on the basis function coefficients ck;m;nn cc5;m;n = c1;m;n+1 c6;m;n = ana+12;m;n+1 :(3.9)

In addition, channel boundary conditions can be imposed quite simply by setting
the value of the appropriate coefficient. For example, the large R boundary condition in closed channels can be applied at the outset by simply setting c5;m;nmax = 0.
This eliminates problems with exponential growth in these channels and, if Rnmax ischosen well into the classically forbidden region, there is no loss of accuracy.
The matrix elements of ; and  in the finite element representation are
given by;ij = 2Z1

;1

ui(xn )(E ; H )uj (xn)andxn ; m;m n;nmax k;5k ;6=an
0

ij = m;m n;nmax k;5k ;50

00(3.10)(3.11)

where i  fk; m; ng and j  fk 0 ; m0 ; ng. H is the Hamiltonian of the system. Because the basis functions are non-zero only within a given sector the corresponding
matrices have a block diagonal structure. Each sector has an associated block whichis coupled only to its nearest neighbors through the continuity conditions (Eq. 3.9).
Note that the overlap matrix , whose elements are given by surface integrals overthe basis functions[67], is particularly simple in the FEM representation. It contains
only mP non-zero elements (all equal to one). The integrals representing the matrix
elements of ; are also particularly simple in the FEM representation. In fact, ex-36
cept for the integral over the interaction potential, all can be done once and for all
beforehand, significantly decreasing the cpu time necessary to construct the matrix.
In addition, the ; matrix is symmetric which reduces the effort required to construct
the matrix. The closed portion ;cc , which is by far the largest matrix involved, canbe written in banded format with dimensions (6mmax ; 1)  (4mmax nmax ), where thefirst term is the half-bandwidth. The actual storage required by LAPACK’s banded
format[63] is 3 times the half-bandwidth. The other large matrix required ;co hasdimensions (4mmax nmax )  mP . A general rule of thumb is to use at least four
sectors per shortest wavelength in the problem which has achieved roughly 6 digitsof accuracy in the phase shift, based on limited tests. The accuracy improves to8 digits using 8 sectors/wavelength. Increasing the number of sectors is mainly aconcern when memory becomes limited. For example, cpu time tests indicate thata single channel calculation using 1000 sectors requires less than 0.2 seconds on a
DEC Alpha 400 MHz workstation. This number of sectors is generally adequate todescribe most alkali two-body collisions. The cpu time scales roughly linearly withnmax and quadratically with mmax. Performance comparisons of this algorithm witha Numerov, Gordon, and integral equation[75] method are provided in Ref[76]. The
FEM R-matrix approach was found to be quite comparable to the integral equation
method (another implicit approach) in speed, memory requirements, and accuracy.
Both of these implicit approaches proved more accurate than the two explicit methods. No time comparisons were conducted with the Numerov and Gordon methods.
I use standard LAPACK[63] routines to solve equations 3.5 and 3.6. The
eigenvalues b and eigenvectors c~o completely specify the linearly independent solutions of the Schrödinger equation on the surface . This information is generally
packaged in terms of a R-matrix

Rmm = ;
0X

Zm b;1 (Z ;1) m0

where the columns of Z are given by the eigenvectors c~o . Heresolution index.

(3.12)
is the independent

3.1.3 Matching R-matrix boxes
The method outlined to this point assumes that equations 3.5 and 3.6 need
to be solved only once to obtain an R-matrix at sufficiently large R to apply theasymptotic boundary conditions. This is not always the case. Situations can arise
which require either an unusually large number of sectors or many tens of channels
to describe the collision process. This can result in a ;cc matrix which is too large
to reside in memory. However, it is not necessary to abandon the FEM R-matrix
approach. One need only divide the radial domain into “boxes” and solve the Rmatrix equations separately in each box. The resulting solutions can be combined

37
in such a way that a final solution is obtained which is equivalent to solving the
R-matrix equations in a single box encompassing the entire radial domain and thisis achieved with no loss of accuracy. This type of approach is sometimes called an“R-matrix propagator” method[77].
A description of the matching procedure is simplified by considering theexample shown in Fig. 3.1. The radial domain is divided into two boxes, A and B.
Boundary conditions are specified on the inner boundary of box A (R = 0) and left
unspecified (i.e., left open) at R = RI . In box B, the boundary conditions are left
open at R = RI and on the outer boundary R = RII , except for channels which arestrongly-closed at RII . In these strongly-closed channels, we must impose boundary
conditions to eliminate exponential growth. The resulting square solution matrix Z
A in A and (mB
has dimensions mmax  max
max + mP )  max in B, where mmaxis the total number of channels, mP is the number of channels open or weaklyclosed on RII . For our particular example, the number of independent solutions
A = m , and B = m
in A is max
max
max + mP in box B. The first step is tomax
require a linear combination of the independent solutions and their derivatives to becontinuous across the boundary (RI ) in each channel. This can be written as:
A
max
X
A

=1

Zm DA =0

X A A A
A
max

=1

B
max
X
B

=1

= ;

Zm b D

Zm DB
0

0

B
max
X
B B B

Zm b D
0

0

=1

0

0

(3.13)

where D represents the undetermined coefficients, the channel index m runs from
0 = Zm b . The
1-mmax and the wave function derivatives have been written as Zm
minus sign in the second condition comes from a reversal of the surface normal
direction. As it stands, Eq. 3.13 is not sufficient to fully specify the coefficients
D. A final condition is obtained by again requiring constant log derivatives on theouter boundary (RII ). This requirement selects the R-matrix eigenchannel solutions, which will be superimposed later to meet the physical boundary conditions
relevant for either a scattering or photoabsorption process. This statement can be
written as
B
max
X
B
0

=1

Zm bB DB = b f
0

0

0

0

B
max
X
B
0

=1

Zm D B :000(3.14)

where now the channel index m0 runs from mmax + 1 to mmax + mP . Equations
3.13 and 3.14 can be combined to form a generalized eigenvalue equation whose
non-infinite eigenvalues determine the new log derivatives b f on the boundary RII .
The corresponding eigenvectors D give the new solution through expressions such

38

B

A

0

Ι

ΙΙ

R
Figure 3.1: Division of the radial domain into two boxes A and B. The R-matrix
equations are solved separately from R = 0 ! RI in A and from RI ! RII in B.
Appropriate boundary conditions are applied at R = 0 and R = RII , leaving allchannels open across RI .

39
as:

Zmf f = N f
0

B
max
X
B
0

=1

Zm DB f000(3.15)

where N f is a normalization constant for each new independent solution. If it
becomes necessary to divide the radial domain into many sectors this procedurecan be iterated as many times as necessary to obtain the solution matrix and log
derivatives on the final boundary.
3.2 Constructing the scattering matrix
At this stage, the R-matrix represents mathematically correct solutionsof the coupled Schrödinger equations. The physically relevant solutions are obtained by applying asymptotic boundary conditions (Eq. 2.5). At large R, where“large” is typically of the order 103 Bohr (occasionally 105 Bohr for systems thatinvolve the long-range magnetic dipole interaction), V (R) ! constant and the twolinearly independent solutions are the standard spherical Bessel jl (km R) and Neumann nl (km R) functions in each channel. A set of linearly-independent solutions atenergy E can be written in this limit as[23]q!1
M (R) R;!
f ;g Kq(3.16)

where fmm = k2m km Rjl (km R)mm and gmm = k2m km Rnl (km R)mm arediagonal matrices of energy normalized spherical Bessel and Neumann functions.
The real-valued constant reaction matrix K determines the appropriate linear combination of f and g . It is evaluated by matching the log derivative of the wave
function integrated from small R (Eq. 3.12) with the asymptotic form given in Eq.
3.16. The result written in terms of the R-matrix is[23]
0

0

0





0



K = f ; f 0 R g ; g0 R ;1 :(3.17)The S -matrix, from which all physical observables can be evaluated, is related to
the K -matrix by the following[23]p

S = [1 + iK ] [1 ; iK ];1(3.18)

where i = ;1. Any scattering observable of interest can then be determined once
this S -matrix is constructed.

40
3.3 Adiabatic representation
Solving the diabatic two-body cold collision equations is the most straightforward approach for obtaining an S -matrix. However, it does not always provide
the most insight into the collision dynamics. Adiabatizing the equations in R canbe a useful step in this regard. A qualitative picture of the collision dynamics is
obtained immediately from the adiabatic potentials and couplings. From a numerical standpoint, an S -matrix obtained by solving the coupled adiabatic equations canbe just as accurate as solution determined in a diabatic representation. However, inthe two-body case it is actually more work to generate the adiabatic potentials andtherefore this approach is presented mainly as a qualitative tool. This is not the casefor the 3-body collisions presented in chapter 8, for which potentials and couplings
have been calculated using the adiabatic hyperspherical method. An adiabatic variant of the FEM R-matrix procedure has been developed to solve these equations andis outlined in Appendix B.
If we choose R as the adiabatic coordinate, the Hamiltonian can be separated into two terms H = TR + H ad , where the kinetic energy operator TR contains
all the derivative terms with respect to R. The adiabatic Hamiltonian H ad is givenbylm (lm + 1) + V
ad = 
Hmm
mm
mm
2R2
0

0

0

(3.19)

where Vmm includes the interaction potentials discussed in chapter 2.2 evaluated in
an appropriate spin basis. Adiabatic eigenfunctions  and eigenvalues (adiabatic
potentials) U (R) of H ad are calculated by parametric diagonalization of Eq. 3.19
as a function of R. The wave function can then be written in terms of the adiabatic
eigenfunctions as [78, 79]0=X

 ( ; ; R)M (R) ;(3.20)

where represents the angular degrees of freedom and  the spin degrees of freedom. Substituting this expansion back into the Hamiltonian gives the following
representation for the coupled radial equations:

"





#

@ + P 2 + U M~ = E M~ :; 21 1 @R

(3.21)All of the coupling (or nonadiabatic) effects are now contained in the real, antisym-

41
metric derivative coupling matrix P (R) given by





@
P =  j @R
00:

(3.22)The Born-Oppenheimer approximation neglects the P matrix. I emphasize that the
adiabatic formulation of the coupled Schrödinger equations is exact in principle, if
P is included and a complete angular momentum basis is used, which in the twobody case will always be possible in the Born-Oppenheimer limit[31], provided
external magnetic fields can be neglected. Using a standard identity, Eq. 3.22 canbe rewritten asP =0( h j@Had=@Rj i
U ;U
0

0;

0

;

6= 0 := 0(3.23)For two-body collisions, evaluating the P matrix is more straightforward with thisapproach since most of the terms in H ad are analytical functions of R. It is only
necessary to take numerical derivatives of the small R part of the Born-Oppenheimer
potentials.
As an example of the usefulness of this method, adiabatic potentials andderivative couplings are shown in Fig. 3.2 for the case of two 87 Rb atoms colliding in an f =2 total spin state. The two-body collision complex consists of three
coupled channels (one singlet and two triplet Born-Oppenheimer potentials) with
one channel converging to each of the three hyperfine split thresholds (see inset
Fig. 2.1). Magnetic dipole couplings have been neglected so this represents a purespin-exchange collision. The dynamics of spin-exchange collisions are extremely
simple. The derivative couplings shown in Fig. 3.2.b indicate that transitions between spin channels occur almost exclusively in the neighborhood of R  20 ; 25
a.u. This value of R is a fairly generic result for ground-state collisions in all thealkalis. The coupling region indicates the point at which exchange splitting in the
Born-Oppenheimer potentials becomes comparable to the hyperfine splittings. Outside this transition region the wave functions in each channel evolve freely in the
adiabatic potentials. At small R these potentials are basically the singlet or triplet
Born-Oppenheimer potentials. Because the couplings are so localized and all channel crossings occur near the same R, the opportunity arises for the wave functions
to interfere in such a way as to shut off spin-exchange losses. This is exactly thefortuitous situation with 87 Rb collisions, which will be discussed in more detail inchapter 5.

42

0.01

U (a.u.)

a)
0

−0.01
11.68

−0.02

0

10

20

30

40

50

0.5
(1,2)
(1,3)
(2,3)

Pij (a.u.)

0.4

b)

0.3
0.2
0.1
0

0

10

20
30
R (a.u.)40

50

Figure 3.2: Adiabatic potentials a) and P -matrix couplings b) for the s-wave collision of two 87 Rb atoms in a f =2 angular momentum state. The inset of a) shows
the splitting of the two triplet-like potentials which is approximately equal to one
hyperfine unit of energy (87 ' 10;6 a.u.). The adiabatic channel indices in b) aredefined in order of increasing energy i.e., channel 1 is a singlet-like state converging
to the 1+1 threshold. Channels 2 and 3 are predominately triplet states at small R,
which correlate adiabatically to the 1+2 and 2+2 thresholds, respectively at R ! 1.

Chapter 4
Multichannel quantum defect treatment of cold collisions

Chapter 2 presented the textbook approach for scattering theory of cold
atoms, deriving a set of coupled-channel Schrödinger equations which describe thecollision dynamics. Brute-force numerical techniques for solving these equations
were developed in Chapter 3. One job that falls to the cold collision practitioner isrefining the Born-Oppenheimer potentials such that theoretical calculations agreewith experiment. This usually requires the solution of the scattering equations
throughout a large parameter space that consists of the Born-Oppenheimer parameters: singlet as and triplet at scattering lengths, and the C6 dispersion coefficient.
In addition, the calculations often must be performed for a number of magnetic bias
fields B and collision energies E . Although the numerical algorithms presented in
Chapter 3 are reasonably fast, this still becomes a time consuming project whenmany thousands of calculations are needed to bound the parameter space, as is frequently the case.
To alleviate this difficulty, our understanding of cold collision physics canbe exploited to develop a more compact and efficient framework for solution of
the scattering equations. The method outlined here is a fully-developed theory forspin-exchange collisions. The magnetic dipole interaction is neglected throughout
this chapter. The techniques developed here should also allow dipolar effects to
be calculated conveniently within perturbation theory, though this has not yet been
implemented. Referring back to Fig. 2.1, it is obvious that the collision partners
experience two very different energy scales in the course of a scattering event. At
small R, the Born-Oppenheimer potentials are very deep compared to typical collision energies and the relative kinetic energy of the pair is large. At large R, thehyperfine energy sets the scale. The energy scale of interest for cold ground statecollisions ranges up to a few mK. Over this range, the small R wave functions do
not change appreciably (see Fig. 4.1) implying that, to a good approximation, thispart of the problem is nearly energy-independent. This in itself represents a huge
computational savings because the spin-exchange coupling occurs in this small R
region (see Fig. 3.2) and only needs to be calculated for a few values (or even justa single value) of E . The energy dependence of the final scattering matrices derives
from the large R evolution of the wave function where the channels are uncoupled;

44
the solutions in this region can therefore be obtained in terms of single-channel solutions to a simple radial Schrödinger equation. In addition, the large R solutions
can be parameterized such that the single channel equations need to be solved only
once as a function of energy. Moreover, we will see that the parameterization canbe made essentially atom-independent. The final scattering solutions can then be
obtained through the manipulation of a few small matrices.
This approach has its origins in Seaton’s[80] multichannel quantum defecttheory (MQDT) and in multichannel effective range theory[81]. I will refer to themethod as multichannel spectroscopy or sometimes simply MQDT. It has provento be a powerful yet simple theoretical tool for describing multichannel collisions
and half-collisions, in contexts such as Rydberg electron motion[23] in the field of astructured ionic core. The theory is well known and has been applied to cold alkali
collisions by other groups[82, 83, 84, 85, 86, 87, 88]. However, the formulation
presented here has two new key developments[89]; i) a new standardization of thelong-range field properties which allows a nearly atom-independent parameterization, and ii) a recoupling frame transformation (FT) approximation. The FT approximation combined with the long-range MQDT parameters reproduces the multichannel collision physics, with reasonable accuracy for most systems, completely
in terms of single channel parameters. In addition to increased numerical efficiency,
MQDT provides a great deal of insight. The techniques developed here are used inlater chapters to explain certain aspects of cold collision dynamics.
The structure of this chapter is as follows: In section 4.1, accurate solutions of the Schrödinger equation in the presence of a long-range potential aredescribed. The parameters that connect these solutions with free-particle solutions(jl ; nl ) at R ! 1 are defined and evaluated. In Section 4.2 the matrix equations
that allow one to connect the small and large R solutions and ultimately derive an
S -matrix are developed. In addition, it is shown how to construct energy-normalized
scattering wave functions from the reference wave functions and long-range parameters. A frame transformation approximation for the short-range reaction matrix ispresented in Section 4.3 and finally, example calculations illustrating the MQDT
approach are presented in Section 4.4.
4.1 Reference wave functions
In multichannel spectroscopy, the key idea is to express scattering observables in terms of a real, linearly independent base pair (f 0 ; g 0 ) of solutions to theradial Schrödinger equation in the appropriate long-range potential. The long-range
potential applicable for ground state spin-exchange collisions is

C6 ; C8 ; C10 :V lr (R) = ; R
6 R8 R10

(4.1)

45

2.0

E = 10 mKE = 1µK

a)

Ψ

1.0
0.0
−1.0
−2.0
0.0

100.0
200.0
R (a.u.)

300.0

0.20

Ψ

b)

0.00

−0.20
5.0

10.0

15.0
20.0
R (a.u.)25.0

Figure 4.1: Component of a multichannel wave function calculated at two collisionenergies E = 1K and E = 10mK. Top figure a) shows the strong energy dependence of the long range part of the wave function. Bottom figure b) is a blowup of
the small R part of the wave function. Except for an overall normalization constant
the two wave functions are almost identical.

46

V lr (R) is the exact form of the long range potentials for both singlet and triplet
Born-Oppenheimer states once the long-range exchange interaction becomes exponentially small. This occurs at roughly R  30 Bohr in all the alkalis (see Fig. 4.2).
An advantage to the numerical approach presented here is the fact that this methodology is not restricted to the form 4.1 of V lr (R). The theory can be easily extended
to incorporate retardation effects or other forms for the long range potentials.
In Section 4.1.1 the Milne phase amplitude method is used to generate
the energy-analytic base pair of radial solutions (f 0 ; g 0). This pair is then standardized making the MQDT parameters portable from atom to atom. Section 4.1.2
presents transformations that allow the base pair (f 0 ; g 0) to be written in terms ofan asymptotic energy-normalized base pair (f; g ) required by scattering theory. The
four long-range parameters needed for these transformations are defined and evaluated. Their low-energy properties are then detailed in section 4.1.3.
4.1.1 Milne equation
A convenient, exact representation of (f 0 ; g 0 ) can be obtained within theMilne[90] phase-amplitude method, details of which are provided in Ref.[91, 92]:

f 0("; l; R) =00g0("; l; R) = ;The Milne amplitudeq 2
R R dRq 2 (R) sin RRRx 2dR(R ) + bl  (R) cos

Rx 2 (R ) + bl :
0

(4.2)

0

(R) is a particular solution of the Milne nonlinear equation@ 2 + k2 (R) (R) = ;3(R) :
@R2

(4.3)

where the R-dependent wave number is defined ass+ 1) ; V lr ) :k(R) = 2(" ; l(2lR
2

(4.4)

Eq. 4.3 is solved with a predictor-corrector method[62] subject to WKB-like boundary conditions at some small fixed radius Rx (generally chosen to be Rx = 10 a.u.),
namely,(Rx) = k(Rx );1=20 (R ) = @ k(R );1=2  :
x
x
@R

(4.5)

WKB-type boundary conditions are not strictly necessary, however, Korsch andLaurent have shown[93] that this provides optimal smoothness in the radial andenergy variations of . The choice of Rx = 10 a.u. is somewhat arbitrary, all

47

−4

10

lr

V

−5

10

lr

|V | (a.u.)

C6
−6

10

C8
C10

−7

10

ex

V
−8

10

−9

10

20

25

30

35
R (a.u.)

40

45

50

Figure 4.2: Contributions to the the Rb2 long range potential V lr . These include thethree dispersion terms, each of which is labeled by its coefficient Cn . The long-range
exchange term V ex which is not included in V lr is also shown.

48
one requires is that V lr (Rx ) is deep enough that the semi-classical approximation is reasonably good where the boundary conditions are applied. An exampleof the base pair (f 0 ; g 0 ) derived from the “smooth” Milne amplitude and phase
(R) = ;2(R0)dR0 is shown in Fig. 4.3.
The energy independent phase shift bl in Eq. 4.2 standardizes the lowenergy asymptotic structure of (f 0 ; g 0). Its introduction allows us to develop a nearlyatom-independent parameterization of the long-range field properties (subsection
4.1.2). The phase shift bl is chosen such thatR!1
f 0(" = 0; l; R) R;!
const  Rl+1 :(4.6)

In the s-wave limit this is equivalent to demanding that the scattering length of f 0
must vanish. The standardization phase is determined by integrating the Milne equation from Rx to Rmax  200 Bohr at zero energy and then equating the logarithmic
derivative of f 0 (" = 0; l; Rmax ) with the logarithmic derivative of the zero energy
solutions of the Schrödinger equation in the presence of a van der Waals potential(i.e., V lr = ;C6 =R6 which is valid beyond R  200 Bohr). These solutions[24], R1=2 J (x), are given in terms of Bessel functions of fractional order  . One
independent solution, with  = ;(2l + 1)=4, has the proper asymptotic behavior
Rl+1. Here, x represents a rescaled variable x = C2 6 R;2. The final expression
for the standardization phase isq

bl = ;(Rmax ) +tan;1 2(R )[J (x) + 2R
max





(4.7)

2Rmax J (x)max J0 (x)] ; 2Rmax J (x) (Rmax ) 0 (Rmax )The derivative of the Bessel function J0 (x) is taken with respect to R. Setting large
R boundary conditions for f 0 is contrary to the “standard” MQDT approach[80,
82, 83, 84, 85, 86, 87, 88] which requires f 0 to vanish as R ! 0. However, thereference potential V lr is not physical inside R  30 a.u. and the base pair (f 0 ; g 0 )is never used there. This provides a very clean separation of the small and large R
physics, as well as a simple standardization of the long range parameters.
4.1.2 Calculating the MQDT parameters
We ultimately want a scattering solution normalized to  (E ; E 0 ). This isobtained by introducing a second energy normalized base pair (f; g ), which are nonanalytic functions of energy. The energy-normalized base pair is related to (f 0 ; g 0 )through the transformation

 f (R)   A 12
=
1
g(R)

01
;
A 2 G A; 2

 f 0(R) 

g0(R) :

(4.8)

49

250.0
0

a)

0.0

0

f,g

0

g

f

0

−250.010.0

20.0
R (a.u.)

30.0

10.0
b)

α

0.0
−10.0
−20.0
10.0

φ

20.0
R (a.u.)30.0

Figure 4.3: a) The base pair (f 0 ; g 0 ) are shown as a function of R. These wereconstructed from the Milne amplitude and phase  shown in b).

50
The MQDT parameters A("; l) and G ("; l) are independent of R and are essentiallynormalization factors. The new base pair are also mathematically correct solutionsof the Schrödinger equation in the presence of the long range-field V lr . Moreover,(f; g) can be expressed simply in terms of the energy-normalized spherical Bessel
and Neumann functions at large R,p!1
f (R) R;!
kR 2=k (jl (kR) cos  ; nl (kR) sin )p!1g(R) R;!
kR 2=k (jl (kR) sin  + nl (kR) cos ) :(4.9)

The non-analytic solutions (f; g ) incur a phase-shift  ("; l) (due to the presence ofV lr ) relative to the free particle solutions. These three parameters (A; G ; ) summarize all that is needed to characterize the positive energy solutions. They are
obtained by integrating Eq. 4.3 from Rx to R1  104 Bohr at energy ". It is only
necessary to construct the energy-analytic base pair (Eq. 4.2) and spherical Bessel
and Neumann solutions at R1 . Using Eq. 4.8 and 4.9, the MQDT parameters can be
related to the Wronskians W (: : :) of these two sets of linearly independent solutions
through the following relationships:

tan  = W (f 0; f s)=W (f 0; gs) W (g0; f s) ; tan  W (g0; gs) 

A = ; W (f 0; gs) + tan  W (f 0; f s)

 W (g0; gs) + tan  W (g0; f s) G = ; W (f 0; gs) + tan  W (f 0; f s) :(4.10)

Here, f s; g s are the energy normalized spherical Bessel and Neumann functions(defined in subsection 3.2). The last long-range parameter that will be needed is anegative energy phase ("; l), which represents the phase accumulated in V lr ,("; l) =Z1
Rx(R0);2dR0 + bl :(4.11)The standardization of f 0 allows us to tabulate these long-range field parameters once and for all as functions of the single parameter = (2)3 C6 "2 , ina.u. This scaling is strictly valid only for a pure van der Waals long-range potential,
but in practice holds to a good approximation for the above form of V lr . Examples of the standardized s- through f -wave parameters are shown in Fig. 4.4. The
standardization proposed here begins to break down for partial waves greater thanl = 2 for different atoms (i.e., when you vary the mass and V lr parameters) andbeyond l = 3 when one simply varies C6 . In addition, the current procedure fordetermining the nonzero l parameters at very low positive energies is numericallyunstable. Both of these problems are likely linked to the large exponential growth

51of the Milne amplitude in the classically forbidden region, which eventually leads
to a linearly dependent base pair. A possible remedy would be to incorporate theanalytical van der Waals solutions at finite energy derived in Ref.[88]. This would
circumvent most of the exponential growth because it would not be necessary tointegrate as far into the classically forbidden region. These issues will be the subjectof a future study.
4.1.3 Low energy behavior of the long-range parameters
In the limit of zero energy, the asymptotic solutions for both positive andnegative energy versions of the reference wave functions (f 0 ; g 0 ) should be equivalent. However, the long-range nature of the potential makes this difficult to achieve
numerically. Therefore, it is preferable to develop an analytical connection formula
to bridge this gap. The asymptotic form of the positive energy base pair (f 0 ; g 0 ) canbe written (using equations 4.8 and 4.9) as the following,r

2 A;1=2 sin(kR ; l=2 + )
(4.12)kr 2 A1=2 cos(kR ; l=2 + ) + A;1=2 G sin(kR ; l=2 + ) :R!1
0
g ;! ; k
In the zero energy limit, f 0 becomes simply!1
f 0 R;!

r2 A;1=2 [kl+1=2 Rl+1 cos  ; k;(l+1=2) R;l sin ] = R2l+1=C : (4.13)

Our standardization of f 0 forces its asymptotic behavior to be f 0 ! Rl+1 =C . Therefore in the low-energy limit, A ! C 2 k 2l+1 and  ! 0. In addition, we find empirically that G goes to a constant value, Gl ! Gl (0). The Wigner threshold propertiesof the S -matrix are ultimately obtained from this low-energy behavior of A. We can
now use these relationships to write down the low-energy limit for g 0 , which is
g0 ! RCl ; GC(0) Rl+1 :(4.14)The negative energy solutions are written in terms of the Milne phase ( =  + bl )and amplitude as originally given in Eq. 4.2. Using 2 2 = (f 0 )2 + (g 0 )2 and
0 R!1

f ;!

equations 4.13 and 4.14, we find

!1 R
2 R;!

2l+2
2 (0)][1+
G
2
C(4.15)

52

1.0

1.0
a)

1/4

A/γ

0.5
0.0

0.5

β/π

−1.0
−1.5
0.0

2.0

A/γ

1/4

b)

0.0

G

−0.5

G

−0.5

η/π

η/π

β/π

−1.0

4.0 6.0
1/4
γ

8.0

10.0

1.0

−1.5
0.0

2.0

4.0 6.0
1/4
γ

8.0

10.0

1.0
1/4

η/π

A/γ

G

c)

0.0

d)

0.5
β/π
G

−1.0

β/π

1/4

A/γ

0.0
η/π

−2.0
0.0

2.0

4.0 6.0
1/4
γ

8.0

10.0

−0.5
0.0

2.0

4.0 6.0
1/4
γ

8.0

10.0

Figure 4.4: The four long-range standardized MQDT parameters are shown for a)l = 0, b) l = 1, c) l = 2, and d) l = 3 partial waves as a function of the parameter
1=4 = (2C )1=4 jkj. On this scale, the parameters for l = 0 ; 2 are nearly indistin6
guishable for all the alkalis. As a point of reference, the energy range E = 1 K - 1
mK for a 87 Rb atom corresponds with the range 1=4 = 0.1-3.7.

53
and the negative energy reference function becomes

f0 !

r

2 Rl+1 1 + G 2(0)1=2 sin (0) : C(4.16)

Comparing f 0 in equations 4.13 and 4.16 provides the connection between the positive and negative energy parameters,

j sin (0; l)j = 1 + G 2 (0; l) ;1=2 :(4.17)

The correct threshold energy behavior of the elastic scattering cross section is obtained if satisfies this condition. This is illustrated in Fig. 4.5.
It is useful to have s-wave values for both A and G in the zero energy limit,
as these will be related to the scattering length in the next section. The constants
C 2 = c(2C6)1=4 , G (0) and (0) are tabulated for the nominal set of dispersion
parameters for each alkali atom (except Fr) in Table 4.1.
Table 4.1: Zero energy values of the s-wave long-range parameters. A is related tothe constant c through the expression A=k = c(2C6 )1=4 . Equation 4.17 is used todetermine (0).
Atom
c
G (0)
(0)
Li
0.96023 -1.0091 -0.78087
Na
0.95857 -1.0055 -0.78266
K0.95775 -1.0037 -0.78355
Rb
0.95722 -1.0026 -0.78410
Cs
0.95701 -1.0022 -0.78430

Analytical expressions for these constants can be obtained in a pure van der Waals
potential, (i.e., V lr = ;C6 =R6 ). Two linearly independent solutions of the single
channel Schrödinger equation, H = @ 2 =@R2 + 2C6 =R6 , are given byq

Y1 = R1=2 J;1=4 (x)Y2 = R1=2 N;1=4 (x)(4.18)

where, again x = C2 6 R;2 is a rescaled variable. Asymptotically expanding theBessel J and Neumann N functions and requiring W (f 0 ; g 0 ) = W (Y1 ; Y2 ) leads to54

8.10e−142

σ (cm )

8.05e−14

8.00e−14

7.95e−14

7.90e−14
−1.0

−0.5

0.0
E (µK)

0.5

1.0

Figure 4.5: Energy dependence of the 6 Li f = 0(1=2; 1=2) elastic cross section. The
collision complex is described by two channels which consist of both atoms being ineither their fa = fb = 1/2 or 3/2 hyperfine states. Zero energy indicates the positionof the 3/2 + 3/2 threshold. The (E < 0) solid line is a calculation using (0) from
Eq. 4.17. The (E < 0) dashed line uses the numerical (0) which results in a slight
artificial discontinuity in the energy dependence of the cross section.

55
the following expressions:s2 21=4;(5=4) R
Y1(R) ! k
(4.19)
1=
8(2C6);(3=4)ss!
1=8 ;(3=4)
1=4
(2

C
2
;(5
=
4)
2

6)
+
R :
Y2(R) ! ; k
21=4;(5=4) (2C6)1=8 ;(3=4)
Comparing equation 4.13 with 4.19 yields the following values for C 2 and G (0):C6)1=4 ;(3=4) = 0:95598(2C )1=4 = aC 2 = (2p62;(5=4)
G (0) = ;1:0 :(4.20)

This expression for C 2 is fairly well known, because it was derived previously inRef[94]. Gribakin and Flambaum call this parameter the characteristic scattering
length a that appears in their semi-classical expression for the physical scattering
length a. Also, the analytic constants derived from a pure van der Waals potential
are very nearly equivalent to the numerical values that include the effects of higher
order dispersion terms in V lr (see Table 4.1). This gives us some indication as towhy the standardization works as well as it does. Referring back to Fig. 4.2, theC8 and C10 terms contribute little to the potential beyond R  35 Bohr. Neglecting
the higher order dispersion terms, the Schrödinger equation can be rescaled in termsof the dimensionless parameter , which in the process removes all atom specificinformation from the equation and therefore from the long-range parameters.
4.2 Constructing scattering matrices and wave functions
The procedure outlined in chapter 3 requires the integration of a set of
N = No + Nc coupled Schrödinger equations outward from the origin to somelarge internuclear distance of order R  103 Bohr. However, we have now definedreference wave functions (f 0 ; g 0 ) that are mathematically correct solutions beyond
R0  35 Bohr. It is therefore only necessary to integrate the coupled Schrödinger
equations to R0 . The N  N solution matrix M (R) at a total energy E can be
written as a linear combination of (fi0 ; gi0 ) at all radii R  R0 . These N independent
solutions are standardized to have the following form:

M (R) = f 0 (R) ; g0(R)K sr ; R  R0 ;(4.21)

where Mij (R) is the i-th channel component of the j -th independent solution. Herefi0 and gi0 are diagonal matrices in the channel space, evaluated at the appropriate
channel energy "i = E ; Ei , with Ei the dissociation threshold energy, and li

56
the relevant orbital momentum in the i-th channel. The solution matrix still hascomponents in energetically closed channels, but this method usually works best if
all channels included at R  R0 are still “locally open”, i.e. with k 2 (R0 ) > 0: The
short range (atom-dependent) physics is completely encapsulated in the short range
reaction matrix K sr , which varies slowly with E or B , except near isolated poles.
The weak energy- and field-dependence of K sr is illustrated in Fig. 4.6.
The solution matrix (Eq. 4.21) does not yet satisfy the proper asymptotic scattering boundary conditions. In particular, M still contains closed channels.
The asymptotic boundary conditions for these channels must be imposed in a simple matrix manipulation at each desired final state energy, through the following
procedure[23]. First partition the K sr matrix into No open channels (those channels
i 2 P; for which "i = E ; Ei  0) and Nc closed channels (i 2 Q; for which"i = E ; Ei < 0). In this notation, K sr becomes:

K sr =

 K sr K sr 
PP
PQ
:
sr
sr
K QP K QQ(4.22)The exponentially-growing parts of the wave functions at R ! 1 are next “eliminated” by applying a N  No transformation matrix B (E ),

B=



1
;(tan + K srQQ);1 K srQP



(4.23)

where 1 represents an N0  N0 identity matrix and tan is a diagonal matrix ofthe negative energy phases defined in the preceding section (4.1.2), evaluated at theappropriate channel energy "i = E ; Ei . This transformation correctly assembles
the linear combination of independent solutions which forces the closed channel
components of the wave function to zero as R ! 1. The open channel subspace ofthe new “physical” solution matrix F 0 = M B is given by:

F 0 (E; R) = f 0 ; g0 K~ :(4.24)~ has dimensions No No and contains the potentially
The modified reaction matrix K
resonant influence of closed-channel pathways:;
K~ = K srPP ; K srPQ K srQQ + tan ;1 K srQP :(4.25)

Moreover, should all channels be closed, discrete bound states occur at roots of thefollowing determinantal equation:;

det K srQQ + tan

=0:(4.26)

Since K sr is only weakly-energy dependent and (") has been tabulated for the
long-range potential of interest, the search for bound states and resonances is easy.

57

0.300

µ = δ/π

a)
0.250
0.200
0.150
0.000 0.025 0.050 0.075 0.100
Energy (K)
0.300
b)

µ

0.250
0.200
0.150
0.0

250.0
B (G)

500.0

Figure 4.6: The eigenvalues of the short-range matrix K sr are plotted as functionsof a) energy and b) magnetic field. The eigenvalues of the reaction matrix, tan  ,are known as eigenphase shifts and these have been converted into quantum defects
through the relation  =  = . In both cases, the quantum defects are nearly
constant implying that K sr can be treated, to a very good approximation, as constantover this range of energy and magnetic field. This particular example is for thecollision of two j2; ;2i 85 Rb atoms.

58
We want the final scattering solutions to be energy-normalized. This isaccomplished by first relating the energy-analytic base pair (f 0 ; g 0) to the energynormalized pair (f; g ) (using Eq. 4.8):



F (E; R) = f ; g K :(4.27)(The solution matrices F 0 (E; R) and F (E; R) are related by an R-independent ma~ ), i.e. F 0(E; R) = F (E; R)C (E ), which can be ignoredtrix C (E ) = A1=2 (1 + G K
since it has no bearing on the final determination of the scattering matrix.) The
new No  No real, symmetric, “energy-normalized” reaction matrix in the (f; g )representation is obtained using[82]

K = A1=2 K~ ;1 + G

;1 1=2

A :(4.28)

This is still not quite the “physical” reaction matrix K phys of conventional scattering theory, owing to the phase shift  ("; l) of (f; g ) relative to (jl ; nl ) at R ! 1.
Nevertheless, K possesses other properties expected in K phys , such as the Wigner
threshold behavior which is reflected in the dependence A("; l) s "l+1=2 as " ! 0:
A1=2 and G are diagonal matrices in channel space, again evaluated at the appropriate channel energy "i = E ; Ei . At this stage, it is a simple matter to construct
the “physical” scattering matrix S . A scattering matrix in the (f; g ) representation
can be obtained by the relationship given in Eq. 3.18. The additional phase shift of(f; g ) relative to (jl ; nl ) is then accounted for when the physical scattering matrix isdetermined in the following manner,+ iK ei :S = ei 11 ;iK

(4.29)

Here,  is a diagonal matrix. The final S -matrix is given completely in terms of theshort-range reaction matrix K sr and the four long range parameters A("; l); G ("; l);  ("; l);and ("; l).
Fig. 4.7 compares the elastic s-wave partial cross section calculated withthe “essentially exact” FEM R-matrix approach and a quantum defect calculation
performed with an energy-independent K sr via equations 4.25, 4.28, and 4.29. The
MQDT results quantitatively reproduce all features of the spectrum over tens of mK
with a single K sr calculated at E =0, relative to the 2+2 threshold. This calculation
also predicts the energy of the zero-field Feshbach resonance[17, 16] (near 29 mK)
to an accuracy of  0.1 mK. The small shift in the resonance position indicates that
K sr is not strictly constant over this energy range. Effects associated with this slightenergy dependence can be incorporated by interpolation of K sr after it is tabulated
on a coarse energy grid.
Energy normalized scattering wave functions will be needed for the pho-

59

−11

10

FEM R−matrixsr
K (E=0)
Frame Transformation

−12

10

2

σl=0 (cm )

−13

10

−14

10

−15

10

−16

10

0.0

10.0

20.0
30.0
E (mK)

40.0

50.0

Figure 4.7: Elastic s-wave partial cross section for the collision of two 85 Rb atoms
in their j2; ;2i atomic hyperfine states. The feature near 29 mK is the zero-fieldposition of the Feshbach resonance observed in Refs.[17] and [16]. The solid curve
is a FEM R-matrix calculation. The dashed curve is an MQDT calculation using
a constant K sr calculated at E =0. The short dashed line is a calculation using theenergy-dependent frame transformation approximation for K sr .

60
toassociation study in Chapter 7, in the range beyond R  30 Bohr. This is achievedby right-multiplication of Eq. 4.27 using a normalization matrix N that tranforms
the wave function into scattering matrix states (i.e., incoming and outgoing waves)that satisfy the asymptotic boundary condition given in Eq. 2.4. The energynormalized standing-wave base pair (f; g ) are defined in terms of the energy-normalized
traveling-wave base pair (f + ; f ; ) through the following relationships[23]:fi = fi p; fi
2i+;gi = ;fip; fi ;2+;

(4.30)(4.31)and the base pair (f + ; f ; ) are defined as:fi+=

r 2 ei(kiR;li=2)p(4.32)2 ki~2r 2 e;i(kiR;li=2);p:
fi = ~2k(4.33)
2i
Substituting (fi+ ; fi;) for (fi ; gi ) into Eq. 4.27 and requiring (F N )ii to equal the0

asymptotic boundary condition in Eq. 2.4, we find the following expression for thenormalization matrix:

N = (1 ; iK );1 :(4.34)

An energy-normalized radial solution valid beyond R  35 Bohr is then represented
by(E; R) = F (E; R)N (E ) :(4.35)

Generation of these wave functions requires only K sr , the four long-range MQDTparameters, and the single channel reference wave functions (f 0 ; g 0). In the case ofa single open channel there is an alternative real representation for that is usefuland is given by:(E; R) = F 0 (E; R)N 0(E ) :(4.36)

Here, F 0 (E; R) is given by Eq. 4.24 and the normalization constant is equal to thefollowing expression:r
N 0 (E ) = sin
A1=2 K~(4.37)

61
where  r = tan;1 (K ).
Finally, we can derive a formula that connects the single channel scattering length a with the MQDT parameters. The physical phase shift is given in termsof the MQDT parameters by  =  r + i . In this context  r represents the phase shift
relative to (f; g ) which in the zero energy limit is equivalent to the physical phaseshift,  r =  . The scattering length a = ; tan =k can therefore be written in termsof the long-range parameters as follows

2~a = ;C K ~(4.38)
1 + G (0)K
where the low energy behavior of the parameter A ! C 2 k has been incorporated.
4.3 Frame transformation approximation
At this stage, the theory presented still requires a full solution of the coupled radial Schrödinger equations to obtain K sr . However, as discussed in section
3.3 the dominant spin exchange coupling is localized in the range R  25 Bohr(see Fig. 3.2). To a good approximation the atoms move at R . 20 a.u. in a setof uncoupled channels, labeled by their total electronic spin quantum number S .
We can therefore integrate the single channel Schrödinger equations with all hyperfine interactions omitted, after which each solution can be matched to the energyanalytic base pair (f o ; g o ) to obtain singlet (S =0) and triplet (S =1) quantum defects
S ("). (The quantum defects S are related to the singlet or triplet phase shift by
S = S =. The phase shift is given in terms of the single channel short-range reaction matrix by  = tan;1 K sr.) The short range reaction matrix is then approximatedby the frame transformation (FT) formulaKi;isr =X0

hiji tan (") hji0i:(4.39)

Here Xi = hiji represents the unitary transformation matrix that connects theshort range singlet-triplet basis ji  j(sa ; sb )S (ia ; ib )I FM i with the asymptotic
hyperfine basis jii  j(sa ; ia )fa (sb ; ib )fb FM i: Xi is understood to be the transformation matrix between symmetrized kets if the particles are identical. Thesetransformations are discussed in section 2.3.2. One important point to note is that
the FT approximation for K sr includes energetically closed channels as in previous
MQDT studies[23]. Therefore, our FT approximation includes the Feshbach resonance physics which distinguishes it from a related class of approximations (notably
the degenerate internal states method of cold collisions[95, 96], or equivalently the
adiabatic nuclei method[97, 98, 99] of electron-molecule scattering theory). We
have tried two different implementations of the frame transformation approximation. The first and simplest is an energy-independent approximation which uses

62
only the zero energy values of the quantum defects (i.e., " = 0 in Eq. 4.39). Ina different implementation of the energy-dependent FT than used previously, theenergy " at which  is calculated in Eq. 4.39 is chosen to be a weighted averageof the channel energies "i appropriate for eigenchannel ji,i.e. "  i "i jhijij2. In tests conducted to this point, the energy-independent FT approximation has
proven to be more accurate for the diagonal (elastic) S -matrix elements, while theenergy-dependent FT approximation has been more accurate for the off-diagonal(inelastic) S -matrix elements. However, the quantum defects have only a weakenergy-dependence (see Fig. 4.8) so neither approach has been remarkably better
than the other. Given the greater simplicity of the energy-independent FT approximation, I would recommend it above the other approach. In particular, only two
parameters are needed to determine K sr . The physical scattering matrix can thenbe obtained via equations 4.25, 4.28, and 4.29 as before, but now solely in termsof single channel quantities. Examples of the FT approximation are provided in thenext section.
As a last item, the quantum defect can be related to the scattering length.
In general, it is more useful to parameterize the singlet and triplet Born-Oppenheimer
potentials in terms of their quantum defect rather than their scattering length as iscommonly done. For one thing, the complete range is covered over 1=2 and accordingly they never experience a pole as scattering lengths can. But more importantly, given the quantum defects one can calculate an S -matrix that is reasonably
accurate for the entire entire range of energy and field relevant for cold collisions.
This is in contrast to the scattering length, which only allows one to predict an swave elastic cross section over a very limited energy range (E  0-few tens of K
in Rb). Returning to the scattering length equation (Eq. 4.38) and substituting in the~ gives the relationship between a and :
quantum defect for KP2
a = 1 +;CG (0)tantan :(4.40)

This equation is readily inverted to express  as a function of a: = 1 tan;1 C 2 +;aaG (0) :(4.41)

From this equation we find the following correspondence between the ranges of aand :

a=0 ; =0
 ;1 
1;1
a < 0 ; 0 <  <  tan G (0) ;1 
1;1
a = +1 ;  =  tan G (0) :(4.42)

63

0.50
singlettriplet

85

0.25

Rb

85

α

µ = δ /π

Rb

0.00

Na

α

23

87

−0.25
87

Rb

Rb
23

Na

−0.50
−1.00

−0.50

0.00
E (K)

0.50

1.00

Figure 4.8: The weak energy-dependence of the Rb and Na singlet and triplet quantum defects are shown. The minimum of the triplet potential is not as deep as the
singlet (see Fig. 2.1) and therefore the triplet quantum defects are a little moreenergy-dependent.

64
Values for G (0) are given in Table 4.1, from which we can see that a negative scattering length corresponds to 0 <  . 0:25. All other values of  (within its defined
range) correspond to a positive scattering length. Given no information about the
ground state Born-Oppenheimer potential, Gribakin and Flambaum’s semi-classicaltheory[94] predicts that there is a 75% chance that a van der Waals potential will
have a positive scattering length. We can find the deviation from this number causedby the inclusion of higher order dispersion terms in the potential by comparing therange of  for which a is positive with the total allowed range of , which is 
0.751/1.0. Therefore a potential which includes the higher order dispersion terms
has approximately a 75:1% chance of having a positive scattering length.
4.4 Examples
We have verified that no accuracy is lost in our S -matrix calculated whenthe multichannel solutions are represented beyond R = 35 Bohr by the semianalytic long-range reference functions (fi0 ; gi0 ), rather than the FEM R-matrix solution. However, if it required calculating K sr at each energy and field, the method
would not be that much of an improvement over the brute force numerical approach.
The real utility of the MQDT method comes from the fact that K sr can be treated,
to a good approximation, as constant in energy and field and therefore only needs
to be calculated once. Moreover, the frame transformation approximation allows
K sr to be calculated solely in terms of two single channel quantum defects and astandard recoupling matrix. The accuracy of this approach has already been illustrated in Fig. 4.7. Both approximations produce elastic s-wave partial cross sections
that agree within 2% of the FEM R-matrix calculations over a 50 mK energy range
except near the zero-field Feshbach resonance. The discrepancy in the energy position of the resonance is  0.1 mK with a constant K sr and by  1.4 mK in theenergy-dependent FT approximation. This will generally be the case. Both approximations for K sr produce accurate S -matrices provided no resonances or quantum
interferences are present. Using an approximate K sr will still predict these quantum
effects but at an energy or field value slightly shifted from its true location. Thisis examined in more detail in the examples below. (The possibility exists that our
method could predict shape resonances that are either spurious or at the wrong energy but we have seen no evidence for this in the partial waves (l=1-4) examined todate. This could occur if the reference potential V lr + l(l + 1)=(2R2) supported ashape resonance.)
Fig. 4.9 shows the magnetic field-induced Feshbach resonance measured
in Ref.[17]. Four methods have been used to calculate the scattering length a. The
FEM R-matrix is considered to be the “exact” solution with a resonance position(defined here to be simply where a ! 1) of Br = 155.4 G. The second method
uses a K sr calculated at E = 0; B = 0. The reaction matrix K sr is rotated into
the representation in which the asymptotic Hamiltonian is diagonal at each B (see65
section 2.3) and then equations 4.25, 4.28, and 4.29 are used to produce the final
S -matrix from which the scattering length is extracted (Eq. 2.34). Treating K sr as
independent of E and B predicts a resonance field value of Br = 156.4 G, which isshifted 1 Gauss higher than its true location. The other two methods are the energydependent and energy-independent FT approximations. These predict resonancefield values Br equal to 161.8 G and 156.1 G, respectively. The predicted resonance field value from the energy-independent FT is actually more accurate than the
one using a multichannel calculation for the constant K sr, although this is probably
just an accident. The coupled-channel version of K sr produces more accurate offresonance values for a. Finally the resonance is shown again in Fig. 4.10, which test
our standardization of the long-range parameters. Here, both calculations use thesame K sr (E = 0; B = 0) as in the previous graph. Two sets of standardized longrange parameters are used, one calculated with 39 K atomic properties (i.e., mass anddispersion coefficients) and the other with 85 Rb properties (these were used in theprevious figure). The off-resonance values of a agree to within  1%. However, thecalculation with the 39 K standardized parameters predicts a resonance field positionof Br = 136.2 G.
87 Rb inelastic loss rates are examined in Fig. 4.11. Again, four different
methods have been used, and we regard the FEM R-matrix calculation as providing
the “exact” benchmark rates. Here, we see that treating K sr as constant in energy
reproduces the inelastic rates to within 7% over a 10 mK energy range. However,in this example the FT approximations are disastrous. This is because the 87 Rb singlet as and triplet at scattering lengths are nearly equivalent producing a suppressedinelastic rate. This suppressed spin-exchange rate was observed in Ref.[100]. The
FT approximation will predict suppressed inelastic rates but only when as ' at (theequality holds exactly for the energy-independent FT). The 87 Rb inelastic scattering
rates are discussed in more detail in Section 5.1. The range over which the FT will
provide realistic inelastic rates is examined in Fig. 4.12. Here, the percentage difference of an off-diagonal mod-squared S -matrix element calculated with the energyindependent FT approximation and with the FEM R-matrix method are comparedover the entire 23 Na s ; t plane. There are no closed channels in this particularcase of scattering from an f = 2(2; 2) total spin state into an f = 2(1; 1) state. The
FT approximation agrees to within 10% of the exact value over most of the parameter space. The exception is along the line of suppressed rates s ' t , just as with
87 Rb collisions. The presence of resonances associated with closed channels will
modify these results just as it did with elastic scattering. The accuracy of the FT diminishes somewhat with increased hyperfine splitting and with increased mass. At
present, it has not been used for Cs calculations but the FT does in general, produce
excellent results for Li through Rb with the few exceptions noted in this section.

66

3000
2000

a (a.u.)1000
0
FEM R−matrixsr
K (B=0)
Energy Independent FTEnergy Dependent FT

−1000
−2000
−3000
135

145

155
165
B (G)

175

185

Figure 4.9: Scattering length a versus magnetic field B for the collision of two 85 Rb
atoms both of which are in their j2; ;2i hyperfine state. Three methods using anapproximate K sr are compared with an “exact” FEM R-matrix calculation. The
FEM R-matrix approach predicts the Feshbach resonance field position in accordance with Ref.[17]’s measurement. The collision energy was taken to be E = 1
K.

67

1000.0750.0

a (a.u.)

500.0

85
39

Rb (C6=4700)
K (C6=4000)

250.0
0.0
−250.0−500.0−750.0−1000.0
0.0

50.0

100.0 150.0 200.0 250.0 300.0
B (G)Figure 4.10: Scattering length a versus magnetic field B for the collision of two 85 Rb
atoms both of which are in their j2; ;2i hyperfine state. A constant K sr calculated
at E = 0, B = 0 is used for both calculations. One (solid line) uses long-range
standardized parameters calculated with 85 Rb atomic properties. The other (dashedline) uses long-range standardized parameters calculated with 39 K atomic properties.
The collision energy was taken to be E = 1 K.

68

−12

10

−13

3

L2 (cm /sec)

10

Energy−Independent FT

−14

10

Energy−Dependent FTsr

K (E=0)−15

10

FEM R−matrix
−16

10

0.0

2.5

5.0
E (mK)

7.5

10.0

Figure 4.11: The total inelastic loss rate L2 for a j2; 2i87 + j1; ;1i87 collision. Calculations using an approximate K sr are compared with the exact rates determinedfrom a FEM R-matrix solution.

69

Figure 4.12: Mod-squared off-diagonal S -matrix elements calculated with anenergy-independent FT approximation and with the FEM R-matrix approach arecompared in this contour plot. These matrix elements are plotted as the absolutevalue of the percentage difference. The 0:1 = 10% and 0:5 = 50% contour lines
are shown. The mod-squared S -matrix element for this particular case represents
the transition probability f = 2(fa = 2; fb = 2) ! 2(1; 1) in the collision of two
Na atoms. The physical values of the Na quantum defects are s = -0.087 and t =-0.400.

Chapter 5
A compendium of rubidium scattering properties

The low-energy scattering properties of Rb in a magnetic trap are now
well understood. This is a dramatic departure from the situation of just a few years
ago. It has taken a number of experiments and corresponding theoretical analyses
to provide a clear and consistent picture. Most notable among the early experiments were the observations of d-wave[101] and g -wave[102] shape resonances inscattering of doubly-polarized 87 Rb atoms and 85 Rb atoms, respectively. In bothexperiments, the observation of additional rotational lines in the photoassociation
spectra signaled the presence of enhanced l = 2 or 4 ground state collisions. Cloud
temperatures were in range of hundreds of K range, well below the centrifugal barrier whereby collisions in these higher partial waves should have been suppressed.
An analysis of the photoassociation lineshapes provided limits on both the triplet
scattering length at and the leading order dispersion coefficient C6 .
The singlet scattering length as was unknown until the experiment of Myatt, et al.[100], that obtained, surprisingly, a pair of Bose condensates coexistingin the same magnetic trap. Prior to this experiment it was generally expected thatthe two trapped spin-states (jfa ; mfa i  j2; 2i and j1; ;1i) would suffer huge spinexchange losses, leading to the destruction of the mixed cloud before it could condense. The anomalous suppression of the spin-exchange loss rate signaled the nearequivalence of the singlet and triplet scattering lengths, as was pointed out almost
simultaneously by three groups[103, 104, 105]. The other novel aspect of this experiment was that only atoms in one spin-state (j1; ;1i) was evaporatively cooled.
Atoms in the j2; 2i spin-state cooled sympathetically (i.e., through elastic collisions
with the j1; ;1i atoms). This provided a “proof-of-principle” for a new route tocondensation that might be applied to other hard to cool species. For example, thismight be useful for fermionic isotopes for which s-waves collisions between two
atoms in the same spin-polarized state is forbidden. An analysis of the Myatt, et al.
experiment will be the subject of Section 5.1.
The observation of a magnetically-induced Feshbach resonance by Roberts,
et al.[17] has provided the most definitive bounds on the Rb interaction parameters(as , at , and C6 ), to date. The same resonance had also been observed using photoassociation spectroscopy techniques[16] but uncertainties in the magnetic field

71
calibration and the difficulty of the corresponding theoretical analysis led to muchlarger uncertainties in the potential parameters extracted in Ref.[16]. An analysis of
the Roberts, et al. experiment will be detailed in Section 5.2.
These experiments, coupled with our theoretical analysis, have provided
a set of “state-of-the-art” potentials, from which reliable predictions can be madefor other Rb isotopes. In particular, we have identified selected hyperfine statesthat provide large elastic cross sections, small inelastic cross sections, and novel
quantum effects (i.e., interferences and resonances). These states offer prospects not
only for interesting collision studies but also for unique degenerate gas experiments,
involving both bosons and fermions. This is the subject of Sections 5.3 and 5.4.
5.1 Suppressed 87 Rb inelastic cross sections
Myatt, et al. measured[100] an inelastic spin exchange rate L2 of 2.2(9)10
cm3 /sec in collisions between the j2; 2i and j1; ;1i hyperfine states of 87 Rb. At the;14time of the experiment, the triplet scattering length was estimated[102] to be 110
10 a.u., whereas the singlet scattering length was completely unknown. Figure 5.1
shows the rate constant (using at = 110 a.u.) versus an assumed singlet scattering
length as . It is immediately apparent that the measured value of L2 places strict limits on as . After the published uncertainties in at were incorporated into our analysis,
we reported[103] bounds on the singlet scattering length of 72  as  102 (in a.u.).
We can interpret the strong suppression of spin-exchange collisions very
simply within the energy-independent frame transformation approximation. First,
the structure of the s-wave Hamiltonian in the coupled hyperfine representation(neglecting magnetic dipole and external magnetic field interactions) is as follows(see Fig. 5.2): the collision complex consists of six channels, three of which(f (fa ; fb ) = 1(1; 2), f = 3(1; 2), and f = 4(2; 2)), have purely triplet characterat small R. The remaining three channels constitute an f = 2 spin-exchange coupled block. The zero-field inelastic processes are controlled to a large extent by
this f = 2 Hamiltonian, which consists of one singlet and two triplet channels with
one channel converging asymptotically to each of the hyperfine split thresholds (seeFig. 2.1). (Most inelastic spin-exchange processes for the bosonic isotopes are controlled by a similar three channel Hamiltonian, see Table 5.1). The symmetrized
short-range reaction matrix for this f =2 Hamiltonian can be written in the FT limit(Eq. 4.39) as:

0 0:1875 + 0:8125 0:2652( ;  )
s
t
s
t
sr
@
0
:
2652(

;


)
0
:
3750

+
0
:
6250
t
K =
s
t
s
0:2864(s ; t )0:4050(s ; t )1

0:2864(s ; t )0:4050(s ; t ) A :(5.1)
0:4375s + 0:5625tand  = tan  . (The channel labels for K sr are 1  (fa ; fb f ) = (11)2, 2 (12)2, and 3  (22)2). The coefficients of the recoupling matrix (see Eq. 2.63)

72

10
10

3

L2 (cm /sec)

10
10
10
10
10

−9

−10

elastic
inelastic
Expt. Upper Bound

−11

−12

−13

−14

−15

−1000.0−500.00.0

500.0

1000.0as (a.u.)

Figure 5.1: Elastic and inelastic event rates L2 for the collision of j2; 2i87 + j1; ;1i87
atoms versus an assumed singlet scattering length as . The calculation is for an
incident energy of 1K above threshold with a fixed triplet scattering length of 110
a.u. The upper bound of the measured value L2 = 2.2(9) 10;14 cm3 /sec is also
displayed.

73
have been evaluated explicitly to show its structure. However, the transformation
between spin representations is unitary and therefore the off-diagonal structure (i.e.,
proportional to (s ; t )) is perfectly general for any coupled channel Hamiltonianin the FT limit. Evidently if the singlet (s ) and triplet quantum defects (t ) (orequivalently the two scattering lengths) are equal, then K sr and also the S matrix
will be diagonal. In such a case the inelastic rates vanish.
Table 5.1: The zero-field Hamiltonian describing an unpolarized gas of identicalbosons has a block diagonal structure in the coupled hyperfine representation (neglecting magnetic dipole effects). The total spin f labels and associated number ofspin-exchange coupled channels are identified for the naturally abundant alkali isotopes. The single channel blocks are pure triplet channels. The two coupled channel
block consists of a singlet and triplet channel converging to each of the two fa = fb
hyperfine thresholds. The three coupled channel block consists of one singlet and
two triplet channels with one channel converging to each of the three hyperfine split
thresholds (see Fig. 2.1).
Atom Nuclear Spin Single Channel Two Channels Three Channels
7 Li
3/2
f = 1,3,4
f =0
f =2
23 Na
3/2
f = 1,3,4
f =0
f =2
39 K
3/2
f = 1,3,4
f =0
f =2
41 K
3/2
f = 1,3,4
f =0
f =2
85 Rb
5/2
f = 1,3,5,6
f =0
f = 2,487 Rb
3/2
f = 1,3,4
f =0
f =2
133 Cs7/2
f = 1,3,5,7,8
f =0
f = 2,4,6
The physical picture is only slightly different from the FT approximation. Recall the discussion of spin-exchange in Section 3.3. There it was notedthat the localization of the spin-exchange coupling at an internuclear separation of
RI  20 ; 25 a.u. presented an opportunity for interference. Physically, the incident
wave splits when it propagates inward through the interaction region (RI ), sending
amplitude into both singlet and triplet channels at short range. These wave components evolve nearly independently, reflect from their inner turning points, thenrecombine in the interaction region to produce the final mixture of exit channels.
If the components meet with their original relative phase, the outgoing waves recombine constructively reproducing the original channel, in that case the scattering
is primarily elastic. Otherwise, the amplitudes recombine destructively requiring
scattering flux to exit in other channels. The suppression occurs when as  atsince singlet and triplet potentials quickly become degenerate for R > RI , meaningthat little additional phase difference accumulates. A similiar interpretation emerges
from the analysis of Dalgarno and Rudge[106]. This circumstance also explains thepoor quantitative performance of the FT approximation in Fig. 4.11. The FT pre-74

f=1(2,1)

f=2

f=3(2,1)

f=4(2,2)Figure 5.2: Block diagonal structure of the Hamiltonian describing the collision of
j2; 2i87 + j1; ;1i87 atoms. The 66 Hamiltonian is presented in a coupled hyperfine representation with the nonzero matrix elements indicated by the corresponding(fa ; fb)f hyperfine labels. The f =2 block consists of three coupled channels.

75
dicts the suppression of the spin-exchange rates, but the suppression is erroneously
predicted to be greatest when as = at (see Fig. 5.3).
This interference picture can be further solidified by examining the singlet
and triplet wave functions near RI . Figure 5.4 shows zero-energy singlet and triplet
wave functions. We see that the singlet wave function with as = 99 a.u. is very nearlyin phase with the triplet, while that with as = 553 a.u. is out of phase. Changing thesinglet scattering length from ;1 to +1 introduces one additional node into thesinglet wave function, pushing it through a  phase change in the interaction region.
The result is one region of inelastic suppression per as cycle, as shown in Fig. 5.1.
However, if the collision complex contains asymptotically closed channels
as is the case here, there is another interference mechanism that in principle can shutoff spin-exchange losses. As an example, consider the collision entrance channel tobe the middle threshold (see Fig. 2.1) with one atom in its fa = 1 hyperfine state andthe other in its fb =2 state (this is channel 2 in the above K sr , Eq. 5.1). In this case
we must follow the procedure outlined in Section 4.2 to eliminate the asymptotically~ is a 22 symmetric matrix
closed channel (channel 3). The new reaction matrix K
with off-diagonal elements given byK~ 12 = K12sr ; K13sr (tan ("3) + K33sr );1 K32sr :(5.2)~ 12 ! 0 the resulting S matrix is diagonal, again eliminating inelastic spinIf K
exchange processes. Feynman path language can be used to interpret the above
equation, where Kijsr is a measure of the flux transferred from channel j into channel
i. The inelastic loss of atoms in our example can occur along one of two pathways;
sr , or it iseither flux is lost from channel 2 directly into channel 1, described by K12
first transferred into the closed channel 3, evolves for awhile, and then hops fromchannel 3 into channel 1 and is lost, (described by the second term in Eq. 5.2). The
two pathways can interfere, shutting off the inelastic processes. Inserting the FT~ 12 ! 0 whenmatrix elements from Eq. 5.1, we find that K
tan ("3) = ; tan t :(5.3)

This relationship is general (within the FT approximation) for all three-channel
Hamiltonians listed in Table 5.1. This example of spin-exchange suppression depends primarily on the triplet quantum defect (i.e., scattering length) and the threshold splitting (through the relationship "3 = E ; E3 ). It is therefore conceivablethat a magnetic field could be employed to induce the interference, although a field
necessarily complicates this simple analysis. Evidence for this type of magneticfield induced interference will be provided in Section 5.3. Figure 5.5 shows the
87 Rb f = (12)2 inelastic rate over the entire  ;  plane. Here, we see thestline of suppression along s  t and also a horizontal line of suppression alongt = ; =  0:4. The topology shown in Fig. 5.5 also was discussed within an

76

−11

10

sr

K (E=0)FT Approximation−123

L2 (cm /sec)

10

−13

10

−14

10

−15

10

0

50

100
as (a.u.)

150

200

Figure 5.3: The inelastic rate constant L2 is plotted versus an assumed singlet scattering length as . The FT approximation predicts L2 ! 0 when as = at , here at isfixed at 106 a.u. However, the exact calculation predicts the minimum will occur atas ' 92 a.u. Our current best estimate of the 87Rb singlet scattering length is 90  1
a.u. The FT provides a poor approximation for 87 Rb inelastic rates, misrepresenting
the rate by orders of magnitude, because of this small shift in the minimum rate
position.

77Ψ

at = 110 a.u.
as = 99 a.u.
as = 553 a.u.

15.0

20.0
R (a.u.)25.0

Figure 5.4: Zero-energy single channel wave functions for the singlet and triplet(solid line) Rb Born-Oppenheimer potentials as functions of the internuclear distance R. The singlet wave function is calculated assuming two different scattering
lengths as = 99 and 553 a.u., while the triplet wave function corresponds to a scattering length of at = 110 a.u. The singlet (as = 99) and triplet wave functions are almost
locked in phase through the region of large spin-exchange coupling (RI =20-25 a.u.)that accounts for the suppressed inelastic collisions.

78
adiabatic framework in Ref.[104].
5.2 Observation of a Feshbach resonance in 85 Rb collisions
A Feshbach resonance in the collision of two j2; ;2i 85 Rb atoms, firstpredicted by Vogels, et al.[107], has now been observed independently by twogroups[16, 17]. Using spectroscopic techniques, Ref.[16] detected the resonancethrough the magnetic field dependence of the photoassociation spectrum. Basedon a theoretical analysis of the spectra, these authors reported the zero energy position of the resonance to be Bmax = 164  7 G with a width of 8.2  3.8 G. A
more precise measurement was reported in Ref.[17]. Here, the technique of “crossdimensional mixing”[108] was employed. That is, the rf fields used for evaporative cooling quickly “chop-off” part of the atom cloud in a single dimension. The
cloud equilibration time is directly related to the elastic collision rate. From thesemeasurements, Roberts, et al. found the maximum and minimum collision rates
occurred at field values of Bmax = 155.2  0.4 G and Bmin = 166.8  0.3 G.
The position of the resonance depends on both the singlet as and triplet atscattering lengths (or corresponding quantum defects) as well as the leading orderdispersion coefficient C6 . In addition, the width of the resonance reflects the coupling of bound and continuum channels and therefore provides an accurate measureof the scattering length difference as ; at . We performed a parameter search over the
singlet-triplet scattering length space for a number of dispersion coefficients. The
phaseshifts generated by the Born-Oppenheimer potentials were adjusted by adding
small corrections to their inner wall (Eq. 2.45) until the field values associated withthe theoretical maximum and minimum of the thermally averaged elastic cross sections agreed with the measured Bmax and Bmin within experimental uncertainties.
(The procedure for thermal averaging the rates and cross sections is presented inAppendix C.) An example of the allowed s ; t parameter space for a given C6 =4700 a.u. is provided in Fig. 5.6. As can be seen from this figure, the small measured
uncertainties of Bmax and Bmin tightly constrain the allowed quantum defect ranges.
In addition, we found that the area of this allowed parameter region was independentof the dispersion coefficients. Varying these coefficients simply shifted the regionin s ; t space. To bound the parameter space completely, we required that ourparameters also reproduced the measured g -wave shape resonance[102] within experimental uncertainties, which limits the van der Waals coefficient to C6 = 4700 
50. This value for C6 is slightly larger but with smaller uncertainty than the mostrecently published value[109] of 4550 100 a.u. We were unable to obtain useful
limits on C8 and C10 due to their small influence on the positions of both the shaperesonance and the Feshbach resonance. Our final scattering length bounds werefirst reported[17] in terms of a fractional bound state phase  to facilitate comparisons with Ref.[16]’s analysis. Here, the results are presented in terms of scattering
lengths and corresponding quantum defects (see table 5.2). The final bounds on the

79

Figure 5.5: Inelastic collision rate L2 for the process f = 2(12) ! f = 2(11) in
87 Rb. The rates are plotted versus singlet  and triplet  quantum defects. Thest
region of suppressed inelastic collisions are contained within the marked contour
lines. These regions are discussed in the text. The physical values of the 87 Rb
quantum defects are s = -0.294 and t = -0.354.

80
triplet quantities also include a 1 uncertainty[107] in the number of bound states
supported by the Born-Oppenheimer potential.
As a further check of the resulting potentials, we confirmed that they predict scattering quantities that are consistent with previous measurements. In particular, our triplet potential supports a 87 Rb d-wave shape resonance[101]. The predicted
87 Rb a =a2;1 1;;1 scattering length ratio[110] and thermally averaged j2; 2i + j1; ;1iinelastic rate constant[100] also agree with their measured values. In addition, we
find 10 of the 12 measured 85 Rb d-wave bound states[109] within the 2 experimental uncertainties. The new potentials also permit us to predict additional 85 Rb
Feshbach resonances: one at Bmax = 226.0  4.0 G with a width of  0:01 G, andanother at Bmax = 535.5  4.0 G with a width of 2.2  0.2 G.
The Feshbach resonance suggests a way to evaporatively cool 85 Rb, whichhas proven difficult at zero-field. Initially, the poor B = 0 evaporative cooling observed experimentally was surprising since it was known that the scattering length
was large (a  ;400), implying a large elastic cross section ( = 8a2 ) at threshold. However, the 85 Rb s-wave partial cross section suffers Ramsauer-Townsend[24]type minimum at low energies (E  375 K) (see Fig. 5.7). We can understandthe presence of this zero by studying the energy dependence of the elastic phaseshift l . Levinson’s theorem[24] states that the zero-energy phase shift l (0) = N
where N is the number of l-wave bound states supported by the potential. In the
E !1limit of infinite collision energy the phase shift l ;! 0. Each partial cross section (l / sin2 l ) will therefore suffer at least N -1 zeros as a function of energy.
However, at low energies the s-wave phase shift is related to the scattering lengthl=0 = ;ka. For negative scattering lengths, the phase shift rises initially and therefore incurs at least N zeros as a function of energy. Levinson’s theorem holds strictly
only for a single channel interaction. However, we find that the low-energy behaviorof the phase shift for multichannel collisions obeys this trend as well. Exceptions tothis rule could occur if resonances are present in the relevant energy range.
The energy position of the first s-wave zero in a single channel interaction
can be estimated using quantum defect theory. The physical phase shift is given interms of our standardized quantum defect parameters by

(E ) = tan;1

 A(E ) tan  
1 + G (E ) tan + (E ) :(5.4)

For a given value of  (which can be approximated as energy-independent in Rb up
to about E = 0.1 K), this equation can be inverted to find the energy at which thephysical phase shift obtains a specified value. In Fig. 5.8, the lines of constant  ==2 (i.e., unitarity limit) and  = 0 are shown as a function of the quantum defect and the standardization parameter 1=4 = (2C6 )1=4 k, where k = (2E )1=2. This
figure shows that the first zero in the s-wave partial cross-section will occur at lower
energies for a negative scattering length than for a positive one. This does not turn

81

Figure 5.6: Lines of constant maximum  (Bmax) and minimum  (Bmin ) cross section versus magnetic field are shown in singlet-triplet quantum defect space. The
magnetic field values are labeled on each line, (dashed=nominal, solid=bounds).
Dispersion coefficient values of C6 =4700, C8 =523070, and C10 = 7665(104 ), in a.u.,were used for this particular calculation.

82

Table 5.2: Singlet and triplet scattering lengths a (in a.u.) and quantum defects (dimensionless) for the different Rb isotopic pairs given a constant C6 = 4700
a.u. and C8 = 550600 a.u. The uncertainty in the singlet quantum defect S is0.0014 in all cases. Variations of the scattering lengths with C6 and C8 can be
determined by the following formulas:  =  + b6 (C6 ; 4700) ; b8 (C8 ; 550600)
and a = ;C 2 tan()=(1 + G (0) tan()). Here, C 2 = 0.957217(2mC6 )1=4 , m is
the reduced mass of the atom-pair (in a.u.), and G (0) = -1.0026. The constant b6
equals 7.510;5 for the singlet and 1.410;4 for the triplet. The constant b8 is thesame for both potentials and is given by 2.010;8 .
Pair
as
s
at
t
82-82

-381

0.0612+15110-0.4910.019

83-83

+661

-0.2021+813-0.2620.013

84-84+1421-0.4678+153-0.0330.007

85-85+2400+600;350

0.2603-369160.1940.001

86-86

+71

-0.0144+21170.4210.007

87-87

+901

-0.2939+1064-0.3540.003

82-87

+551

-0.1568-4040.0640.004

83-87-990+60;700.2263-29570.1840.001

84-87+1171-0.3952+550+45;350.3020.004

85-87

+111

-0.0233+21370.4180.007

86-87+33640.3443+1434-0.4670.010

83

10

−10
87

Rb

10

−11
85

87

2

σ (cm )

Rb− Rb

10

−12
85

Rb

10
10
10

−13
85

Rb s−wave

−14

−15

0.0

200.0

400.0 600.0
E (µK)

800.0

1000.0Figure 5.7: Total elastic cross sections for doubly polarized Rb collisions are plottedversus collision energy. The 85 Rb s-wave partial cross section suffers a minimum
at E = 375 K. In addition, a g -wave shape resonance (E = 650 K) is evident in
the 85 Rb total cross section. The 87 Rb cross section exhibits a broad d-wave shaperesonance near E = 400 K and a p-wave shape resonance (E  100 K) enhances
the mixed-isotope cross section.

84
out to be a problem for 7 Li (at = -27.6) simply because of its small mass. There,
the first s-wave zero does not occur until E > 1 mK, well beyond the temperature
range where evaporative cooling techniques are employed.
Evaporative cooling near a Feshbach resonance must tread a fine line,
balancing increased inelastic as well as elastic collision rates. Although zero-fieldmagnetic dipole loss rates for collisions occurring at the lowest hyperfine thresholdare exceedingly small ( 10;18 cm3 /sec in Rb), near a resonance these rates can increase dramatically. The field-dependence of the 85 Rb j2; ;2i dipole loss rate L2 isshown in Fig. 5.9. At the measured resonance field value (Bmax = 155 G) L2 peaks
at a loss rate of 10;13 cm3 /sec. However, the peak loss rate is extremely energydependent, growing prohibitively large at energies below 1K. The good news isthat the lineshape is asymmetric, with the minimum loss occurring on the positive
scattering length side of the resonance. (Negative scattering length ramifications fora 85 Rb condensate are discussed at the end of Section 5.3.) It therefore looks promising for evaporative cooling to work if the bias field is set in the range B  160-165
G. The dipole loss rates also show much more structure than the elastic rates. The
broad resonances at B >300 G are associated with d-wave spin-exchange coupledchannels. A measurement of this structure could provide even tighter constraints on
the Rb potentials than we have to date.
5.3 Predictions for 85 Rb-87 Rb collisions
A second possible way to Bose-condense 85 Rb is to cool it sympathetically in a gas of 87 Rb atoms. The 85 Rb-87 Rb triplet scattering length is large andpositive ( 200 a.u., see table 5.2) implying a large elastic collision rate between
the atoms (see Fig.5.7). Generally, collisions between upper and lower manifold
hyperfine states suffer large spin-exchange losses. (Note, however, the exception
presented in Section 5.1.) Therefore, the two best choices for mixed-isotope condensation require trapping j3; 3i85 atoms with j2; 2i87 or trapping j2; ;2i85 with
j1; ;1i87 atoms, as neither combination suffers spin-exchange losses. Of these twocandidates, the second is the more interesting because the scattering length between
isotopes can be “tuned” via Feshbach resonances. In principle, this should allow
control of the overlap between the two condensates and produce interpenetrating quantum fluids. In this section, I provide updated values for the mixed-isotope
scattering observables that were first reported in Ref.[111], reflecting our improved
knowledge of the Born-Oppenheimer potentials.
The scattering length a for j2; ;2i85 + j1; ;1i87 collisions in a magnetic
field is shown in Fig. 5.10. Two pronounced resonances are apparent at field valuesof Bmax = 288.5  15 and 420  35 G with widths of 4  1 and 47  10 G,
respectively. The energy-field evolution of these resonances is shown in Fig. 5.11.
Here, the time delay operator[27] Q = 2d=dE is introduced. The signature of aresonance is a rapid increase by  radians in the physical phase shift  . Therefore,85

First max/min s−wave σ
0.5
δ=0δ = π/2

0.4
0.3
0.2

µ

0.1
0.0
−0.1
−0.2
−0.3
−0.4
−0.5
0.0

1.0

2.0

3.0

4.0
1/4
γ

5.0

6.0

7.0

8.0

Figure 5.8: Constant lines of the physical s-wave phase shift  are plotted versus
quantum defect  and the standardization parameter 1=4 = (2C6 )1=4 k, in a.u..
Negative scattering lengths correlate with 0 <  . 0:25 (see Eq. 4.42). This
figure shows that atom-pairs with negative scattering lengths suffer a zero in thes-wave cross-section at lower energies than positive scattering length pairs. Also,
the unitarity limit ( = =2) is obtained at lower energies for positive scattering
lengths. The physical 85 Rb triplet quantum defect is 0.194. From this graph one
finds the first zero occurs at roughly 1=4 = 2.2 that is an energy of 366 K, in good
agreement with the full calculation shown in Fig. 5.7.

86

3

L2 (cm /sec)

10
10
10
10
10

−12

a)

E=1µK

−13

−14

−15

−16

0

100

200
B (G)

300

400

500

a2,−2 (a.u.)

b)
0

−500

−1000

0

100

200
B (G)

300

400

Figure 5.9: a) Magnetic field dependence of the inelastic dipole loss rate L2 forcollisions of two 85 Rb atoms in their j2; ;2i hyperfine state and b) the corresponding
scattering length dependence a2;;2 . The dipole loss rates exhibit a large increase at
the measured Feshbach resonance (B = 155 G) field. The sharp feature in a) at B 
225 G is the narrow Feshbach resonance predicted previously. The broad features ina) B > 300 G are resonances in d-wave spin-exchanged coupled channels.

87
near a resonance, Q can be fit to a Lorentzian profile

Q = (E ; E )2=+2 (=2)2 ;
r

(5.5)

from which the maximum value Qmax immediately identifies the width  = 2=Qmax
and resonance position Er . It is sometimes useful to define a similiar operator withrespect to the change in field strength QB = 2d=dB . (This analogous “time-delay”
like definition of resonance position and width is not generally used for magnetic
field resonances in the cold collision community. Instead, the resonance position
is defined as the field Bmax at which a(Bmax ) ! 1 and the width is given by = jBmax ; Bminj where a(Bmin)=0 (see Eq. 2.2.4). Positions and widths quoted
in this work (including the ones given in this section) will adhere to these commonly
accepted definitions unless otherwise stated.) The Feshbach resonances originate
from f =3 and f =4 quasi-bound states that are nearly degenerate at zero field. The
block diagonal structure of the 88 Hamiltonian that describes the s-wave collision complex is provided in Fig. 5.12. In zero-field, the f = 3 entrance channel
is uncoupled from the f =4 spin-exchange coupled channels and therefore only one
resonance appears in Fig. 5.11. As the bias field is turned on, the f =4 resonance immediately appears with an extremely narrow width “on top” of the f =3 resonance.
As the coupling (HB ) is increased, the two states repel one another, eventually becoming completely separated in energy before crossing the entrance channel threshold. In addition the f =4 resonance broadens considerably as the coupling strength
is increased.
Another interesting magnetic-field dependence is exhibited in the spinexchange loss rates of j3; 3i85 + j1; ;1i87 collisions. Ordinarily, these two hyperfinestates would quickly expel one another from the trap due to their large spin-exchange
rate, L2 = 3:4(10;12 ) cm3 /sec. However, this rate turns out to be nonresonantly suppressed in the presence of a magnetic field. The field dependence of L2 is shown in
Fig. 5.13. The minimum rate of  10;15 occurs at B = 297.5 G. A 1212 Hamiltonian, with five channels open at R ! 1, describes these s-wave collisions. It istherefore surprising (and nontrivial to explain) how an interference effect can shut
down all four inelastic channels at the same value of the magnetic field. MQDT
can be used to gain some insight into the dynamics. Figure 5.14 shows the shortrange reaction matrix elements between the inelastic channels (labeled j ) and the~ j5) eliminating closed
entrance channel (labeled 5), both before (Kjsr5 ) and after (Kchannels. The elements of Kjsr5 are nonzero and almost constant in field. However,when the closed channels are eliminated from the reaction matrix all off-diagonal~ j5 go to zero over the same small range of magnetic field. Therefore,
elements of K
the closed-pathways are interfering with the direct paths in much the same fashion
as discussed at the end of Section 5.1. The loss rate L2 does not exactly vanish since
the off-diagonal elements vanish at slightly different field strengths. This destruc-881000.0a (a.u.)

500.0

0.0

−500.0−1000.0

0

100

200
300
B (G)

400

500

Figure 5.10: The magnetic-field dependence of the j2; ;2i85 + j1; ;1i87 scattering
length. The zero-energy position of the Feshbach resonances are 287.5  15 and420  35 G with widths of 4  1 and 47  10 G, respectively.

89

Figure 5.11: The evolution of two Feshbach resonances in energy and field areshown for collisions of j2; ;2i85 and j1; ;1i87 atoms. The time delay operator isdefined as Q = 2d=dE . Its peak value identifies the position and width  =2=Qmax of each resonance.

90

f=3

HB

HB

f=4

HB

HB

f=5

Figure 5.12: Block diagonal structure of the j2; ;2i85 + j1; ;1i87 s-wave Hamiltonian. The f =3 spin-exchange coupled channels consist of 1 singlet and 3 triplet
channels with one channel converging to each of the four hyperfine-split thresholds(see Fig. 2.1). The f =4 block consists of three channels, 1 singlet and 2 triplets,with one channel converging to each of the three highest energy hyperfine thresholds. The last f =5 block is a single triplet channel converging to the 3+2 hyperfine
threshold. A magnetic field HB introduces additional coupling within each spinexchange block and couples the f = 1 channels.

91
tive interference raises the possibility that a triple condensate could be formed in amagnetic trap with the species j1; ;1i87 , j2; 2i87 , and j3; 3i85. Triple condensates of
23 Na atoms in each Zeeman sublevel of the f =1 hyperfine state have been observedain an all optical dipole trap[112].
A mixture of 85 Rb atoms and 87 Rb atoms exhibits interesting cold collision properties. However, this system suffers a drawback from the point of view ofcondensate mixtures. Namely, the large negative scattering length of 85 Rb implies
that its condensate is unstable at B = 0 beyond  100 atoms. The 85 Rb Feshbachresonance can in principle alter this situation by applying a bias field to tune its
scattering length to a positive value, which should then allow large number condensates to be formed. Unfortunately, the 85 Rb scattering length is negative at the fields
required to tune the 85-87 interaction. Therefore, experiments designed to exploit
the 85-87 interaction must account for the small 85 Rb condensates. Detecting andstudying such small condensates will be no small feat but may be possible[113].
5.4 Scattering properties of the short-lived Rb isotopes
Experimentalists have recently developed the ability to cool and trap shortlived radioisotopes. Experiments with francium[114, 115], sodium[116], as well as
82 Rb[117] have been reported. These experiments open up new territory, offering
the possibility that these short-lived species can be used in studies of degenerate
Fermi and Bose gases. In this section, I highlight some of the interesting scattering
properties of the Rb radioisotopes. In particular, mixed isotope collisions of 83 Rb87 Rb and mixed spin state collisions of the individual fermionic isotopes 84 Rb and86 Rb will be investigated. We first reported these results in Ref.[118]. That work
showed that the 82 Rb scattering properties are only marginally suitable, at best, fordegenerate gas studies and will not be discussed here.
First, consider the scattering properties of the fermionic isotopes. Oneof the most important considerations for degenerate gas studies is the magnitude ofthe elastic scattering cross section that enables evaporative cooling. The Pauli exclusion principle prohibits s-wave collisions between identical spin-polarized fermions.
Therefore, formation of a degenerate Fermi gas will require sympathetic cooling of
the atoms, by using either a different spin state or else a completely different atom.
Both 84 Rb and 86 Rb are ideal candidates for sympathetic cooling with 87 Rb, since
they have large positive triplet scattering lengths (see Table 5.2). In addition, theyboth have an “inverted” hyperfine structure, making f = 5=2 the lower energy
state. (The atomic properties of the Rb isotopes considered here are provided intable 5.3.) This means that spin exchange losses are forbidden in collisions between
j5=2; 5=2iA and j5=2; 3=2iA atoms, where A = 84 or 86. Evaporative cooling of
like-fermions in different spin states has been demonstrated[119] in 40 K, which has
a similarly inverted hyperfine structure. However, the scattering length is negative
for both 84-84 and 86-86 collisions, despite the positive sign of the singlet and triplet

92

10

10

−12

−13

3

L2 (cm /sec)

10

−11

10
10
10

−14

−15

−16

0

100

200
300
B (G)

400

500

Figure 5.13: Field-dependence of the j3; 3i85 + j1; ;1i87 spin-exchange loss rate L2
is plotted. A nonresonant suppression of the loss rate is observed at 297.5  40 G.

93100.000a)

10.000sr|K j5|

1.000
0.100
0.010
0.001
0.000
250

275

300

325

350

100.000b)

10.000

|Kj5|

1.000
0.100
0.010
0.001
250

275

300
B (G)

325

350

Figure 5.14: Field-dependence of the reaction matrix elements are plotted for collisions of j3; 3i85 + j1; ;1i87 atoms. In a) the elements of K sr that connect the entrance channel (labeled 5) to the inelastic channels are shown as solid lines and thesr is represented by a dashed line. In b) the closed channels havediagonal element K55
been eliminated from K sr (Eq. 4.25). The elements of the new reaction matrix K~j 5
are plotted here. (solid = off-diagonal elements, dashed = diagonal element). The
destructive interference between the direct and closed inelastic pathways is evident
in Fig. b)

94scattering lengths (see Table 5.2). Such a result, contrary to the “Degenerate InternalStates” (DIS) model of hyperfine state scattering [106], bears exploring here.
Table 5.3: Nuclear spin I , ground-state atomic hyperfine splittings , and halflife  of the Rb isotopes. A negative value of  indicates an “inverted” hyperfine
structure. The half-life units are m=minutes, d=days, y=years.
Isotope I  (GHz)82
1
1.54741.273 m83
5/2 6.3702
86.2 d84
2 -3.08316
32.9 d85
5/2 3.03573
stable86
2 -3.94688
18.63 d87
3/2 6.83468 4.88(1010 ) y

If the magnetic dipole coupled channels are neglected, this is a coupled
two channel problem with the j5=2; 5=2iA + j3=2; 3=2iA channel energetically closed
at R ! 1. The unitary frame transformation connecting the short range basis withthe uncoupled hyperfine basis implies that the entrance channel has 80 per centtriplet character. The DIS model, which neglects hyperfine energies, would thenpredict (a = 0:8at +0:2as ), i.e., that the 84-84 and 86-86 scattering lengths are bothreasonably large and positive.
In fact, the positions of bound states relative to the appropriate hyperfine
thresholds are crucial for determining the actual scattering length. We can explore
the present situation by introducing a mixing angle  to represent the coupling inthis system[24]. The model for our potential is thenV = U ()V S U T() + E(5.6)

where U ( ) is a standard 22 real orthogonal rotation matrix, V S is a diagonal matrix of singlet and triplet Born-Oppenheimer potentials, and E is a diagonal matrix
of hyperfine energies. The model in the uncoupled limit ( =0) gives simply a triplet
potential connected to the lowest hyperfine threshold (i.e., the entrance channel) anda singlet channel with one additional unit of hyperfine energy.
The -dependences of the 84 Rb-84 Rb and 86 Rb-86 Rb bound states are shownin Fig. 5.15. In 86 Rb (Fig. 5.15a), we find an extremely high-lying bound state in
the uncoupled limit (=0) that accounts for the large triplet scattering length. As thecoupling is turned on, a nearly degenerate pair of bound states, one singlet and one
triplet, begin to repel each other. Eventually the singlet state “pushes” the high-lying
triplet state above threshold, resulting in a negative scattering length. In 84 Rb (Fig.
5.15b) the highest-lying bound state is rather deep  100 MHz for  = 0, in accord
with this isotope’s small positive triplet scattering length. As the coupling is turned

95
on, this bound state interacts with a “singlet” state lying above threshold pushing it
still deeper into the potential, ultimately producing a negative scattering length. In
both cases, the position of the singlet states relative to the upper hyperfine threshold
are the determining factors in the physical scattering length.
Finally, we can compare the DIS scattering length result with that predicted by the frame transformation approximation. The short-range reaction matrix
in the frame transformation approximation is given by

K sr =

 0:2 + 0:8 0:4( ;  ) 
s
t
s
t
:
0:4(s ; t ) 0:8s + 0:2t(5.7)

Again, the notation  = tan  is used. Channel 2 is then eliminated (see Eq.
4.25) leaving a single elastic matrix element

sr sr
12 K21:K~ 11 = K11sr + K sr +Ktan
(;)
11

(5.8)

The scattering length (Eq. 4.38) is given in terms of this matrix element by2~
a = ;C K11~ :1 + G (0)K11(5.9)

Inserting values for the quantum defects (Table 5.2), the zero energy constants C 2 ,G (0) (Table 4.1), and the phase parameter (" = ;) (a table of values for the
long-range MQDT parameters is provided in Ref.[120]), the FT approximation predicts scattering lengths (in a.u.) of a = -103 and -286 for 84 Rb-84 Rb and 86 Rb86 Rb, respectively. The frame transformation approximation reproduces the scattering lengths to within 20% of their physical values of a = -115 (84 Rb) and a =-237 (86 Rb), including the correct sign. The difference in the two models can beseen immediately from Eq. 5.8. In MQDT language, the DIS model is equivalent to~ 11 = K11sr , which completely neglects the influence of the closed channel.
setting K
The closed-channel physics is crucial in this particular case and it is the inclusionof this physics that distinguishes the frame transformation approximation from the
DIS model.
The (surprisingly) negative scattering length implies a net attraction between the two spin states j5=2; 5=2iA + j5=2; 3=2iA, which could have importantconsequences for forming Cooper pairs in these fermionic systems [121, 122]. On
the other hand, the negative scattering lengths in these species again produce zeros
in their s-wave cross sections, as shown in Fig. 5.16. The good news is that 84 Rb
exhibits a Feshbach resonance in the presence of modest-sized magnetic fields (Fig.
5.17a). This resonance could then be used to move the position of the s-wave partial
cross section zero to higher collision energies (Fig. 5.17b), allowing the atoms to beevaporatively cooled into the degenerate regime. The extremely large width of this

96

250.0
86

a) Rb
E (MHz)0.0

T−250.0S−500.0T−750.0
0.0

10.0
20.0
Θ (deg)

30.0

750.0
84

b) Rb
E (MHz)

500.0
S
250.0
0.0
−250.0
0.0

T
10.0
20.0
Θ (deg)30.0

Figure 5.15: Bound and quasi-bound state positions versus the coupling parameter
. The physical coupling strength is given by =26.57. The labels S, T refer tosinglet or triplet states in the  = 0 limit. Zero energy defines the entrance channel
threshold.
a) j5=2; 5=2i86 + j5=2; 3=2i86.
b) j5=2; 5=2i84 + j5=2; 3=2i84 .

97resonance eliminates the need for accurate control of magnetic field strengths. The
resonance would also allow the experimentalist a means to study both repulsive andattractive effective interactions between atoms in a single system.
In the 86 Rb case, a magnetic field will also influence the scattering length,
but will change its sign only at very large fields  2800 gauss (Fig. 5.17a). The
position of the cross section zero can be moved to a higher collision energy (Fig.
5.17b), but to a lesser extent than in 84 Rb. However, in this case the s-wave zerois not a major problem, since the enhanced p-wave scattering should enable theexperimentalist to evaporatively cool the 86 Rb mixed spin states without a magnetic
field bias. However, introducing a bias field of a few hundred gauss would increase
the low energy total cross section by roughly a factor of 3.
Finally, consider the bosonic isotope 83 Rb. The singlet and triplet 83 Rb
scattering lengths are relatively large, positive and very nearly the same (see table
5.2). These are remarkably similiar to 87 Rb implying that 83 Rb by itself adds little to the field of BEC. However, mixtures of 83 Rb and 87 Rb atoms do offer someexciting possibilities. First, the 83 Rb-87 Rb triplet scattering length is large and negative (see Table 5.2), which is not surprising since the 83 Rb-87 Rb reduced mass isvery nearly the same as two 85 Rb atoms. Therefore, a double condensate with large
numbers of atoms could in principle be formed from the hyperfine states j3; 3i83 +j2; 2i87. Within the Thomas-Fermi approximation, a double condensate composedof equal mass atoms with scattering lengths a1 and a2 is unstable whenever thepmutual scattering length a12 satisfies ja12 j > a1 a2 [123, 124]. This relationshipis strictly true only for isotropic like-species condensates. However, a more general derivation[125] for anisotropic mixed-isotope double condensates shows thatthe instability remains. The nature of this instability has yet to be fully interpreted,particularly in the case of a12 negative. The 83-87 mixture would provide a goodsystem to explore this phenomenon. In addition, Feshbach resonances arise in collisions between atoms in the lower hyperfine manifold (see Fig. 5.18) , i.e. j2; ;2i83+ j1; ;1i87 . Thus, in contrast to the 85 Rb-87 Rb system, it is possible to envision
two large condensates with a tunable interspecies interaction. This capability could
enable detailed studies of double condensates, all the way from completely overlapping to utterly immiscible[126, 127], in particular near the stability limits where
ja12 j  pa1 a2.

98

84

a) Rb

−11

2

σ (cm )

10

−13

10

−15

10

0.0

250.0
E (µK)

500.0

250.0
E (µK)

500.0

−11

2

σ (cm )

10

−13

10

86

b) Rb
−15

10

0.0

Figure 5.16: Elastic cross sections versus collision energy using our nominal potentials. Solid lines represent the total cross section, dashed lines indicate the s-wave
contribution. The energy position of the first zero in the s-wave partial cross section
is 20030 K for a) and 29525 K for case b).
a) j5=2; 5=2i84 + j5=2; 3=2i84.
The feature near 250 K in the total cross section is a f -wave shape resonance.
Although we find the position of this resonance is uncertain to  150 K, it does
not contribute over a broad enough energy range to compensate for the zero in thes-wave cross section.
b) j5=2; 5=2i86 + j5=2; 3=2i86 .
We find a broad p-wave shape resonance near 60K. The height of the p-wave centrifugal barrier is roughly 80K. In this case, the strong p-wave scattering is probably enough to compensate for the zero in the s-wave cross section.

99

3000.0a)a (a.u.)1500.00.0−1500.0−3000.0
0.0

1000.02000.0B (G)3000.02000.0E0 (µK)

1500.0b)1000.0
500.0
0.0
0.0

250.0 500.0 750.0 1000.0
B (G)Figure 5.17: Solid lines represent j5=2; 5=2i84+j5=2; 3=2i84 collision, dashed lines
represent the same spin state collision for 86 Rb.
a) Scattering lengths versus applied magnetic field. The zero energy positions of thelow-field 84 Rb resonance is 1069 G.
b) Energy position E0 of the first s-wave partial cross section zero versus appliedmagnetic field.

100

400.0

a (a.u.)0.0−400.0−800.0
0.0

250.0
B (G)

500.0

Figure 5.18: Scattering length versus applied magnetic field for a collision of
j2; ;2i83 + j1; ;1i87 atoms. The zero energy positions of the resonance peaks are
13810 G, 1938 G, and 3713 G. The low-field resonance has a width of  =6.52.5 G. The two higher field resonances are extremely narrow (  1 G).

Chapter 6
Theory of photoassociation lineshapes

Photoassociation (PA) spectroscopy has played a significant role in understanding cold collisions. In Chapter 1, photoassociation was described as theprocess by which an excited molecular bound state is resonantly formed throughphotoexcitation during the collision of two atoms. Generally, experiments are conducted with the atoms initially in their atomic ground state so that an analysis of thephotoassociation lineshapes provides information, via the Franck-Condon overlap,
on the colliding ground state wave function. In addition, transitions to the so called“long-range” molecular states[128, 129] that typically have inner turning points in
the range  30 ; 70 Bohr, provide information on the ground state wave function
at large R. In this way, the s-wave scattering lengths can be determined accurately.
In fact, most of the alkali scattering lengths have been determined (at least initially)
from an analysis of PA spectroscopy.
Photoassociation spectroscopy was first proposed by Thorsheim,
Weiner and Julienne[19]. In that work, the authors developed the theoretical framework needed to analyze a PA spectrum. However, they were mostly concernedwith experiments at the comparatively “warm” cloud temperatures of tens of mK.
Photoassociation theory has since been extended[130] to include threshold effects
of atoms colliding at magneto-optical trap (MOT) temperatures of hundreds of K(see also Refs.[131, 132, 133] for a discussion of photoabsorption line shapes). In
addition, a further development[134, 135] of the theory now includes the effect ofa second laser. An extension of the one-color formulation of Ref.[130] to include
multichannel scattering along the lines of Ref.[136] is presented in this chapter. This
theoretical framework is then applied in Chapter 7 to analyze the measured rovibrational spectrum of the 39 K 0;
g state, and to extract the s-wave scattering lengths.
This chapter is organized as follows: in Section 6.1, the basic lineshapeformula is derived. The main component of this formula is a dipole matrix element that provides the radiative coupling between the ground state and the excitedstate. It is important to construct molecular wave functions that have the proper
symmetries in order to evaluate the radiative coupling matrix element. Symmetryrules for diatomic wave functions in a Hund’s case (a) representation[137] are presented in Section 6.2. The excited state Hamiltonian and ground state Hamiltonian

102
are constructed in Section 6.3. Here, the “ground state Hamiltonian” describes theinteraction of two ground state atoms (S ), while the “excited state Hamiltonian” describes the interaction of one ground state atom (S ) and one excited state atom (P ).
Finally, the dipole matrix element that connects these states is evaluated in Section
6.4.
6.1 One color lineshape formula
The centerpiece of photoassociation theory is an expression for the thermal averaged rate constant Kp that describes the rate at which molecular boundstates are resonantly formed during the collision of two atoms. The general approach for evaluating Kp is to treat the photoassociation problem as an inelastic
scattering process. The inelastic rate constant Kp for the formation of molecules(i.e., coefficient of squared density in the rate equation) is then proportional to

Kp =

*X

(2l + 1) 1+l v
l

+

(6.1)

where the brackets denote a thermal average over a distribution of incident velocitiesv, l is the total inelastic cross section for scattering from an incident continuum
channel of partial wave l into a molecular bound state channel, and  = 1 if theinitial pair of atoms are identical and in indistinguishable spin states, or  = 0
otherwise (see Section 2.1.2). Photoassociation experiments are usually conductedwith dilute atom clouds at temperatures of a few hundred K. Therefore, the thermal
averaging procedure can be treated by the classical method presented in Appendix
C. This expression for the rate constant assumes that the atom density n can befactored out of the thermal average as will be the case for most atom traps. The
result of thermal averaging, after writing the cross section in terms of a scattering2
matrix element =(1 +  ) = v
k2 jSpj , is [130]Z1
X
1
Kp = hQ(2l + 1)
jSpj2e;E=kB T dE ;
T l
0

(6.2)

where kB is the Boltzmann constant, T represents the cloud temperature, and QT =(2kB T=h2)3=2 is the translational partition function. Perturbative expressions forthe scattering matrix element jSp j2 , valid for the low laser intensities generally usedin these experiments, have been derived in Refs. [130, 135]. The general result canbe written as

jSp(E; l; !)j2 = (E + h 0;sE(E;)2l)+ ( =2)2
1

b

(6.3)

103
where 1 is the laser frequency, Eb is the bound state energy, is the total width (inenergy units) of the bound state, and 0 =~ represents the rate at which observableproducts are produced (i.e., atoms lost from the trap by spontaneous emission or
atoms (molecules) ionized by a second laser). The main component of equation 6.3
is the radiative width, which can be expressed in terms of a coupling matrix elementasE
2s(E; l) = 2 jh b (R)jVrad(R)j l (R)ij :(6.4)

Here, b (R) represents a unit normalized bound state wave function, E
l (R) is theenergy-normalized continuum wave function in the incident channel, and Vrad (R) =(2I=c)1=2d(R~ ) represents the radiative coupling that is proportional to the square~ ). (The molecularroot of the laser intensity I and the molecular dipole operator d(Rdipole operator is described in Section 6.4.)
The last step is to substitute expressions 6.3 and 6.4 into the rate constantequation (6.2) and generalize the result for multichannel scattering. The final resultfor Kp is given by[136]4 2 I
Kp(v 0; J 0 ; !; T ) = hcQ

T

X

p flp

nf (2l + 1)

(6.5)

Z 1 e;E=kB T 0jh p J (R)jdJflpp (R)j Eflp (R)ij2
:
dE
0

0

v0
0

0

0

0

0

0

0

(E + h1 ; Ev J )2 + ( =2)20

0This expression (6.5) is used for the analysis presented in Chapter 7. It includes asum over all unresolved final states, an average over all initial states and a thermal
average over incident energy. In this multichannel context, the bound molecularwave function b is labeled by its vibrational quantum number v 0 , rotational quantum number J 0 , parity p0 , and a set of additional quantum numbers 0 required todefine the wave function uniquely. For example 0 includes the nuclear spin degrees
of freedom. The continuum wave function E is energy-normalized and satisfies
outgoing-wave boundary conditions. It is labeled by the total spin f of the colliding atoms, the partial wave l, parity p, and includes all other quantum numbers
required to define the initial quantum state uniquely. Also present is a constant nf
representing a weighting factor for the population of atoms in the f entrance channel. For a randomly polarized gas of identical particles, this weight is given by
nf = N(2ffa+1)
where Nffa ;fb g = (2fa + 1)(2fb + 1 + fa ;fb )=(1 + fa ;fb ) is the total;fb
number of distinct symmetrized spin channels. It should be noted that averaging the
dipole matrix elements over initial states cancels both the (2f + 1) factor above andthe (2l + 1) term in Eq. 6.5.
f

g

104
6.2 Symmetry rules for diatomic Hund’s case (a) wave functions
The field-free diatomic Hamiltonian possesses a number of discrete symmetries beyond the translational and rotational dynamical symmetries. The Hamiltonian commutes with the electron permutation operator P^e , the parity operator ^ ,and for homonuclear molecules, the nuclear permutation operator P^n . In the BornOppenheimer limit[31], the electronic part of the factored wave function possesses
additional symmetries. These include reflection of the electronic wave function
through a plane containing the internuclear axis, whose operator is denoted ^v ; andfor nuclei of equal charge, inversion of the electronic wave function through thecenter of mass in the molecular frame, whose operator is denoted ^i.
We can construct an eigenfunction of the molecular Hamiltonian in the
Born-Oppenheimer limit from a set of Hund’s case (a) basis functions. These basis
functions have the following form:

(J; MJ ; ; ; S; ; I; MI ; R; w) =r

2J + 1 DJ  ( ; ; )(; w)
MJ ;4
j(sasb)S ij(ia ib)I MI i FJ;(R) : (6.6)

J
Here, DMJ ; is a matrix element of the (2J + 1)  (2J + 1) unitary rotation
matrix[57]. It is a function of the Euler angles ( ; ; ) which define the relationshipof the space-fixed lab frame with that of the body-fixed molecule frame. Zare’s[57]definition for the Euler angles, which agrees with that of other standard angular momentum texts (for example, Brink and Satchler[138]; see also Larsson[139] for a
discussion of the different conventions) is adopted here. The normalization factor
implies that = 0, which is a standard convention for linear molecules[57]. The~ + S~ + ~l, the total angularquantum numbers of the D -function consist of J~ = L
momentum (neglecting nuclear spin), MJ is the projection of J~ on the space-fixedz axis; and =  + , the projection of J on the body-frame z0 axis (taken to
be the internuclear axis). The electronic spatial wave function is represented by ,~ ontowhich depends on the projection of the total electronic angular momentum L
the internuclear axis  and on the the body-fixed inversion quantum number w that
distinguishes gerade (g) and ungerade (u) states. The complete expression for  isworked out later. The total electronic spin and its projection onto the internuclear
axis is represented by jS i. Since, the 0;
g experiment did not resolve hyperfine
structure our analysis is simplified by quantizing the nuclear spin jIMi i in the labframe. Finally, FJ;(R) represents the radial wave function.
First, the electronic orbitals that comprise  must be constructed such that
 is an eigenfunction of both P^e and ^i. We are primarily interested in the molecularwave function at large R, where there is negligible overlap of the electronic charge
distributions. There,  can be expressed in terms of properly symmetrized linear

105
combinations of atomic orbitals ULc [140]:(; w) = N [ULac (r~1 )ULb d (r~2 )] + (;1)n[ULac (r~2 )ULb d (r~1 )]+(6.7)o(;1)n [ULad (r~1)ULb c (r~2 )] + (;1)n [ULad (r~2 )ULb c (r~1)] :000

Here, N is a normalization constant and the phase factors n; n0 ; n00 are to be determined. The notation [ULac (r~1 )ULb d (r~2 )] indicates an atomic orbital with quantum
numbers Lc (Ld ) centered on atom a(b) that is a function of electron 1(2)’s coordi~ = L~ c + L~ d onto the internuclear axis.
nates and  denotes the total projection of L
Each atomic orbital has an associated spherical harmonic to represent its angulardegrees of freedom, YLi  (e ; e ), where the angles e ; e are defined relative to thebody-fixed coordinate system defined by the Euler angles ; ; , with the origin
translated to coincide with the atomic center. The inversion of the electronic wave
function through the center-of-mass of the body-fixed frame does two things: first it
switches the orbitals, then it inverts each orbital about the new atomic center. This isequivalent to letting e !  ; e and e !  + e , which amounts to adding an additional phase factor to the spherical harmonic, ^iYLi  (e ; e ) = (;1)Li YLi  (e ; e ).
Therefore, the effect of ^i on each component of  (Eq. 6.7) can be determined fromthe following:^i[ULac (r~1)ULb d (r~2)] = (;1)Lc+Ld [ULad (r~1)ULb c (r~2 )] :(6.8)

We can make  in Eq. 6.7 an eigenfunction of ^i by simply attaching (;1)Lc +Ld +w
phase factors to the third and fourth components. The eigenvalue of ^i is (;1)w ,where w = 0 defines a gerade (g) symmetry and w = 1 defines an ungerade (u)symmetry.
The electronic spin degrees of freedom must be included to antisymmetrize the electronic wave function with respect to electron exchange P^e . The
result when applied to a single component of  (Eq. 6.7) isP^e [ULac (r~1)U b Ld (r~2)] j(sasb)S i = (;1)sa+sb;S [ULac (r~2)ULb d (r~1)]
j(sasb )S i ;(6.9)
where the phase factor arises from the properties of Clebsch-Gordan coefficients.
The electronic wave function must be odd under the permutation of two electrons,i.e.

P^e (; w) j(sasb)S i = ;(; w) j(sasb)S i :(6.10)In the case of two alkali atoms sa = sb = 1=2, this equality requires an additional(;1)S phase factor to be included in the second and fourth components of Eq. 6.7.
The final result for the electronic orbital angular momentum wave function, which

106
is now an eigenfunction of both ^i and P^e , is(; w) = N [ULac (r~1)ULb d (r~2 )] + (;1)S [ULac (r~2 )ULb d (r~1 )]+(;1)Lc+Ld+w [ULad (r~1 )ULb c (r~2)] +(;1)Lc+Ld+w+S [ULad (r~2)ULb c (r~1 )] :(6.11)

With the complete basis function defined, we now turn to the effect of theparity operator. Hougen[141] has shown that the parity operator acting on the total^v , provided that ^v acts on all degrees ofwave function is equivalent to ^  = 
freedom of the wave function. There is more than one plane that can be chosen
to define the reflection operator ^v . We adopt Hougen’s convention of reflecting
through the x0 z 0 plane (see Ref.[139] for a discussion of the different approaches forthis symmetry operation). First, the symmetry rules for the action of ^x z on eachcomponent of  are provided and then a brief outline of their derivations is given.
(The separation of the electronic orbital and spin degrees of freedom, as done here,is strictly valid only in one- and two-electron molecules[142]). The symmetry rules
for the reflection operator are as follows:00^x ;z DMJ  ( ; ; ) = (;1)J ; DMJ ; ( ; ; )^x ;z (; w) = (;1)(;; w)
^x ;z jS i = (;1)S;jS ; i^x ;z jIMI i = jIMI i :
0

0

0

0

0

0

0

(6.12)

0

There is an additional (-1) phase factor in the relationship for , if  is a ; electronic state. However, these states are ignored here because they do not factor into
our analysis of the 0;
g spectrum. The parity of a single basis function (Eq. 6.6) isfound from these rules (Eq. 6.12) to be (;1)J +S for alkali atoms. A wave function
with definite parity (;1)p , which is also an eigenfunction of the symmetry operators^i and P^e, is constructed with the following linear combination of basis functions:

mol

=

s

DJ  ( ; ; ) (; w) jS i +2J + 14(2 ; ;0;0) M(1 ; ;0;0 )(;1)J +S+pDMJ ; ( ; ; ) (;; w)jS ; ijIMI i FJ;(R) ;(6.13)

provided  is constructed according to Eq. 6.11.
The effect of the reflection operator x ;z , which is a sense-reversing symmetry operation, on the rotational degrees of freedom is nontrivial. Hougen disJ =cusses the problem at length and eventually adopts the following scheme: x ;z DM000JC2(y )DM . Here, C2(y ) defines a  radians rotation about the y axis. The outcome of this convention is that ^x ;z  is equivalent to a  rotation of the molecule0

0

0

0

0

0

107
fixed coordinate system about an axis perpendicular to the internuclear axis (y 0 ),followed by the reflection of all other wave function coordinates through the plane
defined by (x0 ; z 0 ). The  rotation about y 0 amounts to changing the Euler angles by!  + , !  ; , and !  ; . Standard identities[57] for the D-functions
can be used to determine the symmetry rule for the rotational piece of the factored
wave function,

DMJ  ( + ;  ; ;  ; ) = (;1)(J ; ) DMJ ; ( ; ; ) :(6.14)To assess the effects of ^x ;z on (; w ), note that reflection of electronic orbital angular momentum through the x0 ; z 0 is equivalent to changing theangle e ! ;e . Each spherical harmonic then transforms as YL (e ; ;e ) =YL (e; e) = (;1) YL ;(e; e). The total wave function acquires the phase factor (;1)(+ ) = (;1) . Finally, the action of ^x ;z on the electronic spin would
have no effect if it was quantized in the lab frame (as is the case with jI; MI i).
Rotating the electronic spin wave function into the space-fixed frame,
0

0

0

0

jSMS i =0

X S

DMS ;( ; ; )jS i ;



(6.15)

the reflection ^x ;z of the spin wave function becomes00^x ;z jS; MS i =
0

0

=

X

[^x ;z DMS S ( ; ; )][^x ;z jS i]
0

X



0

0

0

(;1)S;DMS S ;( ; ; )[^x ;z jS i] :0(6.16)

0Thus the lab frame spin wave function will remain invariant only if 
^x ;z jS i =S;(;1) jS ; i.
The nuclear permutation operator acting on the total wave function provides nuclear selection rules. Hougen has shown[141] that this operation is equivalent to P^n = ^i^ P^I , where P^I only permutes the nuclear spins. We have alreadyevaluated ^ , as well as ^i, which acts only on the electronic orbital. Therefore,
applying P^n to the wave function in Eq. 6.13, gives:
0P^n mol = (;1)p+w+I ;ia;ib mol :0(6.17)For the bosonic alkali atoms, which have half-integer nuclear spins, the wave function must be odd under permutation of the nuclei and therefore the equality (;1)p+w+I =1 must be satisfied. It is necessary to include an additional [1 + (;1)p+w+I ]=2 termin the wave function (Eq. 6.13) to ensure that it is also an eigenfunction of P^n . The
final expression for the wave function, which has definite parity, and is an eigen-

108function of ^i, P^e , and P^n is then

s





2J + 1
1 + (;1)p+w+I jIM i F (R)(6.18)mol =I J;
4(2 ; ;0;0 )2
DJ  ( ; ; ) (; w) jS i +M
(1 ; ;0;0)(;1)J +S+pDMJ ; ( ; ; ) (;; w)jS ; i :We are primarily interested in the consequences of these symmetry relationships for the 39 K 0;
= 0. The rotational levels of ang electronic state, i.e.
= 0 state have a definite parity and therefore have definite nuclear spin selection
rules. This is simply because one way to construct an = 0 wave function is witha single basis function (note the ;0 ;0 function in Eq. 6.18). These individual
basis functions (Eq. 6.6) have parity (;1)J +S and the bosonic nuclear permutation
symmetry requires (;1)J +S +w+I = 1. In the next section, it will be shown that the0;g electronic state is a triplet (S = 1) spin state. Therefore, the odd J rotationalenergy levels have positive parity and the even J rotational energy levels have negative parity. In addition, the nuclear permutation symmetry requires (;1)J +I = ;1,
which limits the nuclear spin states that can contribute to each rotational level. These
selection rules play a large role in the analysis of the measured 0;
g spectra.
6.3 Excited and ground state Hamiltonians
This section is divided into two parts, the bulk of which is presented in
6.3.1. Here, the interactions and general matrix structure for the excited state Hamiltonian are presented. First, the general excited state molecular Hamiltonian is developed and the approximations employed in the analysis of the 0;
g state are detailed.
Next for completeness, the ground state Hamiltonian needed for the analysis presented in Chapter 7 is provided in Section 6.3.2.
6.3.1

S + P Hamiltonian

The interactions that comprise the S + P Hamiltonian include those listed
in Section 2.2 and the atomic spin-orbit interaction H so . The zero-field total Hamiltonian is thus composed of the following components:

H = H BO + H so + H l + H ss + H hf

(6.19)

where H BO includes the Born-Oppenheimer potentials, H l is the rotational Hamiltonian, H ss represents the magnetic dipole-dipole interaction, and H hf is the hyperfine
interaction. The Hamiltonian is constructed in a Hund’s case (a) molecular representation, which will be a convenient representation for evaluating the dipole matrix
element. The Born-Oppenheimer potentials that correlate asymptotically with the S

109~ , the projection of the electronic angular mo+ P threshold depend on the total spin S
mentum on the internuclear axis , and also the inversion symmetry of the electronic
wave function in the molecular frame, labeled g (even) or u (odd). There are thus
eight unique Born-Oppenheimer potentials for the alkali atoms, labeled 2S +1 g=u ,and H BO is an 8  8 diagonal matrix in the Hund’s case (a) representation.
The good quantum numbers in the Hund’s case (c)[137] angular momentum coupling scheme (neglecting nuclear spin) are the total angular momentum
J~ = L~ + S~ + ~l and its projection onto the internuclear axis =  + . The
nuclear rotation (or partial wave) ~l is oriented at right angles to the internuclear axis.
The Hund’s case (c) labels become appropriate when the strength of the spin-orbit
interaction H so becomes at least comparable to the energy splittings between the
Born-Oppenheimer potentials H BO . Spin-orbit matrix elements in a Hund’s case(a) representation have been evaluated in Ref.[143] and are reproduced here ( isthe atomic spin-orbit splitting):

h3 w ; jH soj3 0w ; 0 i = 3 w;w [; ; + ; ;1; +1
0

0

0

0

0

0

+ ; +1; ;1 ]h1 w ; jH soj3 0w ; 0 i = ; 3 w;w [; ;0 ; ; ;1 ;;1+ ; +1;+1]
1
so
1
0
0
h w ; jH j w ;  i = 0 :
0

0

0

0

0

0

0

0

(6.20)

The atomic spin-orbit interaction H so breaks the electronic spin projection degeneracy, resulting in 24 Hund’s case (c) electronic states. These 24 states can be broken
down into four different six-channel blocks which correspond to a definite parityand definite g/u symmetry. In addition, the Hamiltonian H = H BO + H so is block
diagonal with respect to the value of . Nonvanishing contributions of H BO + H so
are illustrated in figures 6.1 and 6.2 for opposite parity states.
1 ~l2 in a Hund’s case
Evaluation of the rotational Hamiltonian H l = 2R2(a) representation is simplified by first decomposing the nuclear rotation angular~ ; S~ , and then using the commutation properties ofmomentum vector ~l = J~ ; Leach Cartesian component of these operators, i.e. [Ji ; Si ] = 0. The final expressionfor H l , written in terms of the angular momentum operators J; L; S , is[53]1 (J 2 ; J 2 ) + (L2 ; L2 ) + (S 2 ; S 2 )H l = 2Rzzz2+(L+ S; + L;S+) ; (J+L; + J;L+ ) ; (J+S; + J;S+)] (6.21)
000where J = Jx  iJy represents standard raising and lowering operators definedin the body-fixed coordinate system, (L and S are defined similarly). There are
two difficulties with this expression. First, Hund’s case (a) wave functions are not
0

0

110

−

0 g/u parity
3

3

3

+

BO

so

l

so

3

Πg/u(Ω=0) H +H +H H +H3

Σ g/u(Ω=0)

+

3

3

Πg/u(Ω=1)

+

3

so

H +H
H

l

l

l

BO

l

H +H
H

l

l

l

H

l

H

l

3

so

1

Πg/u(Ω=1)00H

Πg/u(Ω=2)

0

0

H

H

H

BO

H

H +H

so

l

l

0

0

l

0

0

H

H +H

Σ g/u(Ω=1)

3

1

+

Πg/u(Ω=0) Σ g/u(Ω=0) Πg/u(Ω=1) Σ g/u(Ω=1) Πg/u(Ω=1) Πg/u(Ω=2)

l

l

so

H +H

BO

l

H +H

H

so

H

l

l

H

so

H

so

BO

H +H

0

H

H

l

l

l

0

BO

so

H +H +Hl

Figure 6.1: Block Structure of the 6 6 S + P Hamiltonian. The parity (;1)p =(;1)J +1 of the eigenfunctions composing the matrix elements is determined by thesymmetry rules for the 0;
g=u electronic state. The structure of this matrix is the samefor both g and u symmetries. The origin of the nonzero matrix elements are indicatedby H BO  Born-Oppenheimer potentials, H so  spin-orbit interaction, and H l 
the rotational Hamiltonian.

111

+

0 g/u parity
3

1

3

+

BO

so

l

l

3

Πg/u(Ω=0) H +H +H

1

Σ g/u(Ω=0)

H

l

H +H

3

H

l

H

3

+

Πg/u(Ω=1)

+

3

l

H

BO

l

l

l

H

l

H

l

3

so

1

Πg/u(Ω=1)00H

Πg/u(Ω=2)

0

0

H

H

H

BO

H

H +H

so

l

l

0

0

l

0

0

H

H +H

Σ g/u(Ω=1)

3

1

+

Πg/u(Ω=0) Σ g/u(Ω=0) Πg/u(Ω=1) Σ g/u(Ω=1) Πg/u(Ω=1) Πg/u(Ω=2)

l

l

so

H +H

BO

l

H +H

H

so

H

l

l

H

so

H

so

BO

H +H

0

H

H

l

l

l

0

BO

so

H +H +Hl

Figure 6.2: Opposite parity block as that shown in Fig. 6.1(with same notation fornonzero matrix elements). The eigenfunctions composing these matrix elements
have a (;1)p = (;1)J parity which can be derived from the symmetry rules for a0+
g=u electronic state. In the absence of rotation, the  parity eigenstates are degenerate for = 1 and 2

112~ 2. However, electronic orbitals constructed in the separated atom
eigenfunctions of L
~ 2 . The
limit, as in Section 6.2, can be treated as approximate eigenfunctions of L
second difficulty is that the total angular momentum J in the body-fixed frame obeys
anomalous commutation rules[57], i.e. [Jx ; Jy ] = ;iJz . A standard approach[57]for resolving the anomalous commutation problem is to redefine the raising andlowering operators such that J ! J . Matrix elements for H l in a Hund’s case (a)representation have been evaluated in Ref.[141].
The matrix structure of H l in the Hund’s case (a) representation can be
determined from the operators comprising the Hamiltonian (Eq. 6.21). The first
three terms in Eq. 6.21 are diagonal in the Hund’s case (a) representation. The fourthterm (/ L+ S; ) gives rise to rotational coupling within a constant block of the
Hamiltonian. The last two terms are responsible for electronic-rotational couplings
and are often denoted as Coriolis couplings. The rotational Hamiltonian preserves
the g/u symmetry of the electronic wave function and therefore preserves the block
diagonal matrix structure of H . This is illustrated in figures 6.1 and 6.2. Finally thematrix values of H ss , which are diagonal in a Hund’s case (a) representation, were
given in Section 2.2.3. These matrix elements are relatively small in the excitedstate (roughly 10;5 the value of the Born-Oppenheimer potentials at large R) andare neglected. In addition, H hf is neglected in the analysis presented in Chapter 7.
In Ref.[21], it was shown that the bound state energies of the lowest few vibrational
levels of the 0;
g state were essentially independent of the hyperfine interaction.
In practice, 0;g potentials that accurately describe the lowest few vibrational levels are obtained by considering only H BO + H so , along with the diagonal
elements of the rotational Hamiltonian H l (the =0 off-diagonal rotational couplingelements for 39 K are of the order  1= times smaller than the spin-orbit coupling).
This results in the following 22 matrix in the representation (3 g ;3 +
g ):
0

H=

V +pVl ; 3
2
3 

0

0

p2

3 

!

V + Vl(6.22)

where  = 57.706 cm;1 represents the atomic spin-orbit splitting. The Born-Oppenheimer
potentials at large R are well-approximated by a multipole expansion,

=
=
=
V = ;f = CR3 3 ; CR6 6 ; CR8 8 ;0(6.23)and the 39 K dispersion coefficients were taken from Ref.[144]. Retardation effects
are also incorporated in f = . These have the form[145]f  = cos(R=) + (R=) sin(R=) ; (R=)2 cos(R=)f  = cos(R=) + (R=) sin(R=) :(6.24)

113
Here,  = 4S !4P =2 is the transition wavelength. The rotational matrix elementsVl incorporate only the first three terms from Eq. 6.21 and are given explicitly by0~2 [J 0(J 0 + 1) + 2]Vl

= 2R2Vl= 2R2~2 [J 0(J 0 + 1) + 4] :(6.25)J
The J 0 -dependent 0;
g potentials, V , are determined by first diagonalizing
Eq. 6.22 and then incorporating non-adiabatic diagonal corrections[78]. Bound stateenergies Ev J and single channel wave functions v J (R) are now easily obtained
for the potentials V J . The radial part of the excited state wave function b is thenconveniently written as an R-dependent linear combination of the two Hund’s case(a) basis functions
0

0

0

0

0

0

Fb(R) =XbJ  (R)v J (R)
0

0

0

0

(6.26)

0where bJ (R) is the eigenvector obtained from the diagonalization of Eq. 6.22 that
correlates with the potential whose asymptote is S1=2 + P3=2 . This form (6.26)
facilitates the evaluation of the dipole matrix element.
0

6.3.2

S + S Hamiltonian

The ground state Hamiltonian has already been described in detail in Section 2.2. Here, I briefly outline the pertinent details needed for the analysis presented
in Chapter 7. The 39 K photoassociation experiment[146] was conducted with unpolarized fa = fb =1 ground state atoms trapped in a dark-spot MOT with cloud
temperatures of  400 K. The experiment resolved J 0 =0-4 rotational peaks for anumber of vibrational levels. Since no rotational lines were observed with J 0 > 4, it~ + S~ ,was only necessary to include partial waves up to l=2 in the analysis (J~ = ~l+ L
L=1, S =1). Neglecting the magnetic dipole interaction, the Hamiltonian for collisions among fa =1 atoms can be separated into five pieces (see Table 5.1). These
include f =0,2 s-wave and d-wave components, and an f = 1 p-wave component.
The analysis required thermal averaged radial wave functions for each Hamiltonian.
However, the nonzero radial overlap of the ground state wave function and the 0;g
excited state wave function occurs at R > 40 a.u.. Therefore, the MQDT method
outlined in Chapter 4 provides an extremely efficient method for generating thesewave functions. A single K sr matrix calculated at E = 0 for each Hamiltonian(which are functions of the singlet and triplet scattering length and the C6 coefficient) proved sufficient to obtain accurate radial overlaps over the entire energy
range E = 0-3 mK required for the thermal average. A comparison of the integrated
radial overlaps using MQDT-generated ground state wave functions and FEM R-

114
matrix generated wave functions is shown in Fig. 6.3.
6.4 Dipole matrix element
The molecular dipole matrix element needed to evaluate the rate constant

Kp (Eq. 6.5) is constructed in this section. This derivation provides “nuts and bolts”

intermediate steps. The reader more interested in the photoassociation analysis canskip ahead to the next chapter without loss of continuity. The approach detailed
here is to first evaluate the molecular dipole matrix in a Hund’s case (a) molecular
representation and then to project the ground state molecular wave function onto
coupled hyperfine states.
The excited 0mol and ground state mol molecular wave functions are assumed to have the form given in Eq. 6.18. The radial piece of the wave function will
be neglected in the intermediate steps since it has no bearing on the angular momentum algebra. The electronic orbital angular momentum  piece of the factored wave
function can be constructed from Eq. 6.11. The excited state orbital is given by:n(0; w0) = 21 [S a(r~1)P b(r~2)] + (;1)S [S a (r~2)P b(r~1)] +00(6.27)o0(;1)w +1[S b(r~1 )P a(r~2)] + (;1)S +w +1 [S b(r~2)P a(r~1 )]
0

0

0

0

0

;

and the ground state orbital is given by the expression:(0; w = S ) = p1 [S a(r~1 )S b(r~2)] + (;1)S [S a (r~2)S b(r~1)] : (6.28)
2Here, the orbitals are designated by the Li quantum numbers, i.e. S a represents
an atomic S -orbital on atom a. Primed and unprimed quantities are used to denote
excited state and ground state quantum numbers, respectively. The molecular dipole~ ) is approximated as the sum of two atomic dipole operators, which
operator d(Rcan be rotated into the body-fixed frame as follows:d(R~ ) =X 1
q

Dmq [dq (r~1 ) + dq (r~2 )] :(6.29)

Here, dq (r~i ) is an atomic dipole operator that acts only on electron i. This ensuresthat only one electron absorbs a photon at a time. The quantum numbers q and mrefer to the projection of the dipole operator onto the body-fixed and space-fixed
frames, respectively. The dipole matrix element for an S ! P transition in the

115| Radial Overlap |

1500
(a)

s−wave1000

f=0(1,1)0(1,1)

f=0(2,2)0(2,2)

500

0
0.0

1.0

E (mK)

2.0

3.0| Radial Overlap |

4000
3000

(b)

d−wavef=0(2,2)

2000

0(1,1)
0(2,2)
0(1,1)

1000
0
0.0

1.0

2.0

3.0

E (mK)
Figure 6.3: 39 K Franck-Condon overlaps of the v 0 =0, J 0 =2, 0;
g rovibrational levelwith a f = 0 s-wave state (a) and a f = 0 d-wave state (b) are compared using
two different methods to calculate the ground state wave functions. These were: 1)
MQDT approach (Eq. 4.35) using a constant K sr calculated at E =0, denoted by
the solid and dashed lines and a 2) FEM R-matrix approach. The FEM R-matrix
calculations are denoted by the closed circles. Each graph shows the absolute valueof the four components that comprise the total overlap. Each component depends onthe fa ; fb ground state quantum numbers and the excited state 0 quantum number.
The f (fa ; fb ) quantum numbers are labeled on the graph. A solid line represents the0 = 0 state, and a dashed line is used for the 0 = 1 state. The presence of a f = 0
d-wave shape resonance near 1.2 mK is evident from the radial overlaps in (b).

116
Hund’s case (a) representation is given by the following expression:s

h 0moljd(R~ )j moli = I ;I MI ;MI S ;S  ; (2 ;   1 )(2 ;  )  1 + (;1)p +I +w  1 + (; 1);0 p+I;0+w  ;0(6.30)

0

0

0

0

0

0

0

0

0

2
2
(2J 0 + 1)(2J + 1) 
162

Z
X
0
0
1  DJh( ; w )jd(R~ )j(0; w)i d DMJ J Dmq
MJ

r

0

0

q

0

+ (;1)p +p+J +J (1 ;  ;0 ;0)(1 ; ;0) 

Z
0
0
J
1

J
h(; ; w )jd(R~ )j(0; w)i d DM ; Dmq DMJ ;
0

0

0

0

0

J
0

0

First, the electronic dipole piece of this expression is evaluated:

h
i
d


;q
S
+1+
w
= p 1 + (;1):
2h(0; w0)jd(R~ )j(0; w)i00(6.31)

Evaluating the second dipole term in Eq. 6.30 leads to the same result after replacing ;q ! ; ;q . The constant d represents the atomic dipole matrix element. Next,
the integrals over the rotational functions are evaluated:00r(2J 0 + 1)(2J + 1)
162

s

Z

1  DJ
d DMJ J Dmq
MJ
0

0

0

=

(6.32)(2J + 1) hJM ; 1mjJ 0M 0 ihJ ; 1qjJ 0 0 i
JJ(2J 0 + 1)

The same result is obtained for the second integral in Eq. 6.30 except ! ; and0 ! ; 0 . The final expression for the Hund’s case (a) dipole matrix element is117
given by:





d  
1 + (;1)p +I +w 
~
p
j
d
(
R
)
j
i
=


mol
I
;I
S
;S

;

M
;M
mol
I I
2
2
 1 + (;1)p+I+w  S +1+w 
1
+
(
;
1)
2

h 0

0

0

0

0

0

0

0

0

h

i

1 + (;1)p +p+1(1 ;  ;0 ;0)(1 ; ;0) 
0

0

s

0

s

(2J + 1) 
(2 ;  ;0 ;0)(2 ; ;0 ) (2J 0 + 1)100

hJMJ ; 1mjJ 0MJ0 ihJ ; 10jJ 0 0 i :(6.33)

The first two factors in parentheses enforce the nuclear permutation symmetry rules
for the excited state and ground state, respectively. The third factor in parentheses
guarantees the g$u dipole selection rule is obeyed. This can be seen from the(;1)S+1+w phase factor and recalling that (;1)S determines the ground state g/u
symmetry. The fourth factor in parentheses ensures that only terms of opposite
parity are coupled. The first Clebsch-Gordan coefficient guarantees that the dipole
selection rules m + MJ = MJ0 and J 0 = J; J  1 (excluding J 0 = J = 0 transitions)
are satisfied.
Next, the transformation from the Hund’s case (a) basis into the coupledhyperfine representation is needed for the ground state wave function. Each component (fa ; fb ) of the multichannel wave function in the coupled hyperfine representation has the form
0

hf

=

 1 + (;1)f +l 2

p2 ;1 fa ;fbF(Efb fa)fl (R) j(f f )fm ; lm i + (;1)f +l;fa;fb (1 ; 
a b

f

l

(6.34)

;fa ;fb )j(fb fa )fmf ; lml iwhere F(Efb fa )fl (R) represents a component of the energy-normalized radial wave
function. The transformation matrix element that connects the ground state molecular wave function mol(J; MJ ; ; S; ; I; MI ; p; w) with the ground state coupled

118
hyperfine wave function hf ((fa fb )fmf ; lml ) is given by:

h molj hf i =s(2fa + 1)(2fb + 1)(2S + 1)(2I + 1)(2l + 1) 
(2J + 1)(2 ; ;0)(2 ; fa;fb )(6.35)

1 + (1 ;  ) + (;1)l+S+I (1 ; 

fa ;fb );0+ (;1)l+S+I (1 ; ;0)(1 ; fa;fb ) hS ; l0jJ i 

0s i f 1@ sab iab fab A XhSMS ; IMI jfmf ihSMS ; lmljJMJ i :
S I f

MS

The final result for the dipole matrix element that couples the excitedstate wave function 0mol (J 0 ; MJ0 ; 0 ; 0 ; S 0 0 ; p0 ; w0 ) with the ground state hyper-

119fine wave function hf is obtained from equations 6.33 and 6.35:

X

h 0moljd(R~ )j hf i =h 0moljd(R~ )j molih molj hf i

J;MJ ; ;S;

(6.36)

X 1 + (;1)p +I +w 1 + (;1)l+I +S
= pd
2
2
2 J
 1 + (;1)f +l  

S
+
w
+1
1 + (;1)

2
0

0

0

0

h

0

0



0

i

1 + (;1)J +p +J +l (1 ;  ;0 ;0)(1 ;  ;0) 
0

0

0

h

0

0

1 + (1 ;  ;0) + (;1)l+S +I (1 ; fa;fb )
0

0

0

i

+(;1)l+S +I (1 ;  ;0)(1 ; fa;fb ) 000

hJ MJ0 ; m; 1mjJ 0MJ0 ihJ 0 ; 0; 10jJ 0 0 i 
hS 0 Mf ; MI0 ; I 0MI0 jfmf ihS 0 Mf ; MI0 ; lMl jJ MJ0 ; mi 

0s i f 1@ sab iab fab A hS 00; l0jJ 0 ; 0i 
S I f
Z
dR bJ  (R)J (R)F(Efa fb)fl (R) :000Here, the excited state radial wave function has been written in terms of an adiabatic eigenvector, FJ0 ; (R) = bJ ; (R)J (R) (see Eq. 6.26). The complete modsquared dipole matrix element that is needed for the rate constant equation (Eq. 6.5)
is given by a coherent sum of the dipole matrix elements in Eq. 6.36 over the internal quantum numbers of the excited and ground state wave functions, an incoherent
average over all initial states, and an incoherent sum over all unresolved final states.
For the particular experiment presented in Chapter 7, this expression is:
0

X

0

p flp
0

1

0

0

0

0

nf (2l + 1)jh vp J (R)jdJflpp (R)j Eflp (R)ij2 =

X

0

0

0

0

X

0

X

0

0

0
~
Nffa;fbg M ;I ;M f;Mf ;l;Ml  ;fa;fbh moljd(R)j hf i
J
0

0

I
0

0

2

(6.37):

120
As a final test of this formulation, we have checked that all dipole selection rules areobeyed and that the lifetime of each rotational level (J 0 , MJ0 ) of the 0;
g state is equal
to one half the atomic lifetime[128]. In addition, to ensure that all pieces of thetheory were integrated properly, we have reproduced the Na 0;
g spectra illustrated
in Ref.[136].

Chapter 7
Determination of potassium scattering lengths

The theory of photoassociation (PA) lineshapes developed in Chapter 6 isherein applied to analyze measured 39 K 0;
g rovibrational spectra. In particular, synthetic spectra are generated that reproduce most features of the experimental spectra. Requiring that the theoretical relative intensities and rotational linewidths agreewith the experimentally measured counterparts, within experimental uncertainties,
allowed us to place limits on the singlet and triplet s-wave scattering lengths. These
results were first reported in Ref.[147]. Previous estimates of the K scattering
lengths had in fact already appeared in the literature[148, 149]. However, the older
works were based on conventional molecular spectroscopy. These data alone proved
insufficient to reliably extrapolate the potentials to the dissociation threshold. The
analysis presented here, which makes use of both the conventional spectroscopic
data as well as the new photoassociation data, should therefore provide a more reliable prediction for the scattering lengths. In addition, an independent analysis[150]of the 39 K 1u spectra confirms many of the results presented here. This chapter isbroken into two sections. First, an overview of the experiment is provided in Section
7.1. The analysis of the spectra and the resulting bounds on scattering lengths arepresented in Section 7.2.
7.1 Photoassociation experiment
The experiment discussed here was conducted by H. Wang, P. L. Gould,
and W. C. Stwalley at the University of Connecticut. They have published complete
details of the experimental setup in Ref.[144, 146]. This section provides a brief description of those experimental details that are necessary to understand the analysis
presented in Section 7.2. In particular, the experimental checks taken to ensure thequality of the data are emphasized.
A sample of 107 39 K atoms, at a density greater than 1011 cm;3 and atemperature of  400K, are prepared in a “dark spot” vapor cell magneto-opticaltrap (MOT). The atoms are trapped predominately in their fa = 1 hyperfine state(for 39 K, ia = 3=2). Separate single-mode tunable ring Ti:sapphire lasers provide
the trapping beams and induce the PA transitions. The total output power of the122trapping laser is  300 mW at 766.5 nm while the PA laser intensity is typically
50 mW/cm2 . The analysis presented in Section 7.2 concentrates on the v0 =0-6
vibrational levels of the 0;
g state. For these low lying states, trap loss spectra are noteasily observed because the local kinetic energy is not large enough to allow atoms
to escape from the trap when the excited vibrational state decays through boundfree spontaneous emission of a photon. The PA resonances can then be detectedby recording ion production following a double-resonance transition from the 0;grovibrational state to an autoionizing 1u state at an energy below the highly excited
4S + nD (n = 5 or 6) asymptote[151]. In particular, a tunable ring dye laser provides
the second photon, at frequency 2 and typical intensity 12 W/cm2 , for the doubleresonance process and a channeltron multiplier records the subsequent productionof ions. A schematic of the transition steps is provided in Fig. 7.1.
Unfortunately, there is essentially no information on the excited state potentials 4S + nD , which prevents us from applying the two-color theory presented
in Ref.[134, 135]. It was therefore necessary to ascertain the influence of the second laser on the PA spectrum of the 0;
g state experimentally. These tests included
recording spectra using final autoionizing states below both the 4S +5D and 4S +6D
asymptotes. Representative examples of these spectra are shown in Fig. 7.2. The
autoionizing final states are structureless and broad[151] with widths ranging from 1 GHz for levels near the 4S+5D threshold to  2 GHz for levels near the 4S+6D
asymptote. General trends in the PA spectra of the 0;
g state, such as weak or missing
odd J 0 rotational levels, anomalous J 0 = 2 linewidths, and large J 0 = 4 relative intensities were consistent in both cases. The only noticeable differences occur in theJ 0 = 4 relative intensity which was typically a factor of 1.5-2 bigger in the spectra
observed using final states below the 4S + 6D threshold. This difference can apparently be attributed to slightly hotter atoms in the MOT for the 4S + 6D experiment(based on linewidths), and to the broader final state that experiment.
PA spectra of the 0;
g state were measured by fixing the second laser fre;00quency 2 to a 0g (v J = 2) ! 1u (v) transition near the 4S+nD asymptote, thenscanning the PA laser frequency (1 ). Altering this setup so that 2 was tuned from
J 0 = 0 or J 0 = 4 to the final autoionizing state again produced no qualitative changes
in the measured 0;
g spectra. The only quantitative difference again appeared in the0
J = 4 relative intensities. These were roughly a factor of 1.5 larger when the second laser was tuned from J 0 = 4 instead of J 0 = 0 or 2 (see Fig. 7.3). Changing
the power of laser 2 only affected the total ion production rate. In addition, theexperimentalists varied the PA laser intensity by a factor of two and observed no
noticeable power broadening. Based on this near independence of the spectra on thesecond laser, we attribute the bulk of the relative peak intensities to the first transition. The single color theory presented in Chapter 6 should therefore provide a good
approximation to the measured lineshapes.

123

30190
30180

2

2

2

2

2

2

4 S1/2 + 5 D3/2

30170

−1

E (cm )

30160

hν2

13044

4 S1/2 + 4 P3/2

−

0g
13036

hν1

4 S1/2 + 4 S1/2

0
−0.5

25

50

75

100

125

150

R (a.u.)Figure 7.1: Schematic of the 39 K photoassociation transition sequence. Groundstate atoms are first excited (h1 ) to a rovibrational level of the 0;
g electronic state.
The free-bound transition is detected by exciting the molecule at frequency 2 toan autoionizing state below either the 4S + 5D or 4S + 6D threshold and thenmeasuring the resulting production of ions.

124

Ion Rate (arb units)4000.02
J’=02000.0

4S+5D

1
3

4

Ion Rate (arb units)

0.0
4S+6D
1000.0

0.0
0.0

500.0

1000.0Relative PA Laser Frequency (MHz)0
Figure 7.2: Representative examples of the 0;
g v =0 spectra obtained using autoionizing final states below the 4S + 5D and the 4S + 6D thresholds. The rotational
assignment J 0 for each peak is indicated in the top graph.

Relative Transition Amplitude

125

1.5

v’=3

1.0

J’=0

a)
2

4

0.5
0.0

Relative Transition Amplitude

−0.5
1.5

v’=3

b)

1.0
0.5
0.0
−0.5
0.0

500.0

1000.0Relative PA Laser Frequency (MHz)0
Figure 7.3: Representative examples of the 0;
g v =3 spectra obtained for different
0values of 2 . The peaks correspond to the J = 0, 2, 4 rotational levels, respectively.
0a) 2 fixed to a 0;
g (J =2) ! 1u(v) transition near the 4S + 6D asymptote.
0b) 2 fixed to a 0;
g (J =4) ! 1u(v) transition near the 4S + 6D asymptote.

126
7.2 Analysis
The analysis of the 0;
g spectra is simplified somewhat by its symmetryrules. These were discussed briefly in Section 6.2 and are elaborated on here. The
J +1 , while the ground state wave
parity of a 0;
g rotational wave function is (;1)
function parity is simply (;1)l . Therefore, dipole selection rules restrict transitions
from ground state collisions in even partial waves l to even rotational levels J 0 andlikewise odd partial waves to odd J 0 ’s. More specifically, s-wave collisions contribute to J 0 = 0,2; d-wave collisions contribute to J 0 = 0,2,4; and p-wave collisions
contribute only to J 0 = 1,3. In addition the 0;
g electronic state for like atoms has a
nuclear permutation symmetry which requires (;1)J +I to be odd.
The measured spectra show several distinct trends. First, the p-wave features are extremely weak or missing from all vibrational scans except v0 = 0 and
2. Second, the relative intensity of the J 0 = 4 (pure d-wave) peak is found to bequite large in all the scans. This is surprising since the height of the 39 K d-wave
centrifugal barrier is roughly 1.5 mK and at cloud temperatures around 400 K one
would expect, in the absence of resonances, the d-wave features to be suppressed by
the Wigner threshold law. Finally, the J 0 = 2 linewidths are quite large  35 MHz
compared to the other rotational lines  25 MHz. The analysis of the measured data
will focus on these particular features.
Suppressed p-wave features were also observed[136] in the Na 0;
g spectra. However, in sodium the Franck-Condon overlap was controlled essentially bya node in the open channel wave function. In potassium, the closed channel wave
functions can contribute substantially to the transition, as can be seen from Fig. 7.4.
The Franck-Condon contribution from each channel must be added coherently (Eq.
6.37), which can contribute to the suppressed p-wave features. The outer turning
points of the closed channel potentials occur at roughly the same internuclear separation as the inner turning point of the excited state potential. This means that theclosed channels, which are on their last half cycle at the Condon point, can have alarge effect on the radial overlap integral. Moreover, the shapes of the closed channel wave functions near the Condon point are essentially “fixed” in R, as they must
begin their exponential decay. Varying the singlet and triplet scattering lengths thus
only affects the relative amplitudes in the channels. For example, the amplitudes
will be quite small when as  at just as 87 Rb spin exchange inelastic rates[100] aresuppressed by this criterion. It is therefore not obvious (or unique) where the open
channel node needs to be placed in order to reproduce the observed suppression of
the p-wave features.
This complication prompted us to perform a complete search of singlettriplet scattering length space to ensure that our final fits to the spectra provided aunique set of parameters. The ground state potentials were systematically varied(Eq. 2.45) such that the singlet as and triplet at scattering lengths ranged between
1. This search was also conducted for a number of different C6 coefficients rang0

0

0

127

v’=6

v’=0

2+2

2+1

1+1

10.0

30.0

50.0
R (a.u.)

70.0

90.0

Figure 7.4: Schematic of the ground and excited state potentials (dashed lines) andwave functions (solid lines). The ground state thresholds are labeled by their fa + fb
quantum numbers. The ground state wave functions are computed at zero energy.
The excited state wave functions correspond to the lowest and highest vibrational
levels (v0 = 0 and 6, respectively) used in the analysis. Notice that the closed channel
components of the ground state wave function can contribute significantly to theFranck-Condon overlap with the excited state wave function. The potentials are notshown to scale.

128
ing from 3600-4000 a.u., to allow for the stated 5% uncertainty given in Ref.[38].
Ground state wave functions for each set (as ; at ,C6 ) of potentials were generatedusing the MQDT method outlined in Section 4.2 (see also Section 6.3.2 for a discussion of the ground state wave functions needed for this analysis). For each setof potentials, a thermally averaged synthetic rotational J 0 =0-4 series was generated
for vibrational levels v 0 =0-6 (equations 6.5 and 6.37). A “match” to the measured
data required that the theoretical relative intensities and linewidths of a rotational
progression within a given v0 agree within experimental bounds simultaneously foreach vibrational level v0 =0-6. A constant cloud temperature, estimated to be 400 
100 K from fits to the measured spectra, was assumed for each vibrational level.
Absolute intensities are quite uncertain, therefore each vibrational spectrum was
normalized to Kpmax (J 0 =2) = 1. In addition, the absolute PA laser frequency 1 wasadjusted such that the blue edge of the synthetic and measured J 0 =2 lines are aligned
for each vibrational spectrum. The absolute intensity and absolute frequency are theonly two adjustable fit parameters incorporated into the theory.
The data obtained from the 4S + 6D autoionizing final state was used exclusively in this analysis. The uncertainties associated with the second laser were
better characterized for this data and therefore more precise bounds could be placedon the scattering lengths. A comparison of spectra calculated with our “nominal”
set of potentials is shown in Fig. 7.5. The J 0 =4 relative intensity is the main discrepancy between experiment and theory. This can be attributed to optimizing the
double resonance experiment using the J 0 =2 intermediate level which then under
represents the J 0 =4 intensities (see Fig. 7.3). The final values for the singlet andtriplet scattering lengths, as well as quantum defects, are given in Table 7.1. The
table also translates these results into scattering lengths for other K isotopes. In
addition, a table of K atomic parameters is provided (Table 7.2).
In the following discussion, the contribution of each spectral feature to
the final limits placed on either the singlet or triplet scattering length is slightly
artificial. Bear in mind that it is only through the combination of these featuresthat unique limits on as and at are obtained. The large J 0 = 2 linewidths are aparticularly interesting feature. We find that these linewidths are a result of a broad
f =0 d-wave shape resonance, at an energy ER = 1:15 0.35 mK with a width  =0.5  0.2 mK (the time-delay definition of position and width has been used here,
see Eq. 5.5). The J 0 = 2 spectral width provides a fairly sensitive probe of theposition of the shape resonance. Figure 7.6 shows the effect on the J 0 = 2 linewidth
as the position of the shape resonance is varied. Recall that atoms scattering inboth s and d partial waves contribute to the intensity of the J 0 = 2 rotational level.
Moving ER to larger energies shifts the d-wave contribution to the red of the s-wave
contribution, thereby increasing the overall width. However, if ER becomes toolarge relative to the cloud temperature, the thermal averaging suppresses the d-wave
intensity contribution effectively reducing the linewidth. Similarly, if ER is moved
to smaller energies not only is the red shift reduced but the d-wave contribution

Relative Transition Amplitude

129

1.25

v’=3

v’=4

v’=5

v’=6

0.75
0.25

Relative Transition Amplitude

−0.25
1.25
0.75
0.25
−0.25
0.0

500.0

1000.0 0.0

Relative PA Laser Frequency (MHz)500.0

1000.0Relative PA Laser Frequency (MHz)Figure 7.5: Comparison of experimental (dashed line) and theoretical (solid line)
spectra. The synthetic spectra were calculated using our nominal set of potentials,
as = 140 a.u., at = -17 a.u. and C6 =3800 a.u., assuming a cloud temperature of
400 K. The prominent peaks in each vibrational spectra are the rotational levels
J 0 = 0, 2, 4, reading from left to right. It is believed that the red wing of the J 0= 0 experimental peaks correspond to additional flux obtained from fa =1 + fb =2
collisions. These collisions were not theoretically modeled since they contributelittle to the spectrum.

130

Table 7.1: Singlet and triplet scattering lengths a (in a.u.) and quantum defects (dimensionless) for collisions among different pairs of potassium isotopes assuming 6 =3800 a.u. The scattering length is related to the quantum defect by:a constant C2
a = ;C tan()=(1 + G (0) tan()) where C 2 = 0:95775(2mC6)1=4 , m is the
reduced mass of the atom-pair and G (0) = -1.0037. The singlet quantum defect:019
uncertainties are +0;0:011 for each collision pair, independent of the C6 value. Finally,
the uncertainties associated with each triplet quantum defect and its variation with:04 6) +0
C6 is given by: t = t + 8:0(10;5 )(C6 ; C;0:06 .
Isotopes
as
s
at
t

39 + 39

140+6
;9

0.460

;17  250.039

40 + 40

105+3;5-0.445

194+114;35

0.388

41 + 4185  3-0.36665+13;8-0.26839 + 40

;1+4
;8

0.002

-460+330;10.212

39 + 41 113  5-0.474205+140;40

0.379

40 + 41 ;50+9
;18

0.089

104+20;11-0.441

Table 7.2: Mass in 12 C amu, nuclear spin i, nuclear g-factor gn , and ground-stateatomic hyperfine splittings of the K isotopes. A negative value of  indicates an“inverted” hyperfine structure.
IsotopeMassign (MHz)39
38.963707 3/2 0.39146
461.7240
39.963999 4-1.298 -1285.79
41
40.961825 3/2 0.21487
254.01

131
with its extra (2l+1) phase space factor overwhelms the s-wave contribution, again
resulting in a smaller linewidth. The large J 0 = 4 peaks are also a by-product of thisshape resonance and their relative intensities are quite sensitive to the position of theresonance.
These constraints force ER to lie in the range 1.15 0.35 mK, which
creates an “allowed” strip in Fig. 7.7. In particular, the position of the resonancedepends sensitively on the singlet scattering length and thus imposes tight bounds
on this quantity. However, because the resonance is near the top of the centrifugal
barrier and therefore quite broad, the spectral linewidths are relatively insensitive
to the C6 value. In fact, we find the bounds on the singlet scattering length areunchanged when C6 is allowed to vary between 3600-4000 a.u. Allowing a quarter
linewidth uncertainty on each J 0 = 2 peak, places the following bounds on as =0
140+6
;9 a.u., which is in good agreement with Ref.[149]. Comparison of the J = 2
measured linewidths with our “best” fit potentials are shown in Fig. 7.8. Finally,it should be noted that this f = 0 d-wave shape resonance manifests itself in otherfeatures seen in the 1u spectra[150]. However, we do not agree on the exact positionof the resonance and therefore have different constraints on the singlet scattering
length. We believe this discrepancy indicates a limitation of the “naive” PA theory(which was also applied to the 1u spectra) that disregards the effect of the ionizing
laser. This point is discussed in more detail in Ref.[150].
There is, in fact, another bound on ER . Raising ER above  1.3 mK introduces an additional f = 2 g -wave shape resonance[150]. The existence of a g -wave
0
resonance should manifest itself in the 0;
g spectrum as an additional J = 6 peak (as0well as contributing to J = 2 and 4). However, there is no experimental evidence foran additional rotational peak. Also, the presence of a g -wave resonance would beinconsistent with the 1u spectra[150]. The exclusion of the g -wave resonance thus
reduces slightly our allowed parameter space (see Fig. 7.7).
Placing bounds on at is more complex. The suppression of the p-wave
features ultimately controls the allowed at range, but it is only after limiting theavailable parameter space by considering the d-wave shape resonance that we findan unambiguous result. However, these combination of features require that thewave function nodes are bounded as follows: f =1, p-wave node 68+9
;6 a.u. and the+4
f =2, s-wave node is 62;3 a.u. The positions of these nodes in our allowed as-at
parameter space are controlled almost exclusively by the triplet scattering length.
Unfortunately the suppression of the p-wave features is fairly robust in this region
which is largely responsible for the bigger uncertainty in the triplet scattering length.
Variations in C6 require a corresponding change in the triplet scattering length topreserve the nodal positions. We found changing the value of the C6 coefficient
simply shifted our allowed parameter region along the triplet axis of the as -at plane
while preserving the area. The triplet scattering length is therefore parameterized 6)  25 a.u., C 6 = 3800. The twoin the following manner at = -17 - 0.045(C6 - C
most important scattering lengths for BEC are the triplet scattering length at and theRelative Transition Amplitude132

ER = 1.1 mK
ER = 0.9 mK
no resonance

1.00

0.50

0.00
325.0

375.0

425.0

475.0

Relative PA Laser Frequency (MHz)Figure 7.6: Comparison of the J 0 = 2 linewidths versus the position ER of the f =0
d-wave shape resonance. The no resonance case indicates the s-wave contribution to
the linewidth. The position of the resonance controls the relative d-wave contribution to the peak and the amount it is red-shifted relative to the s-wave contribution.

133

Figure 7.7: The dashed lines indicate the energy position ER (in mK) of the dwave shape resonance in the as -at plane. The solid rectangle represents our allowed
individual singlet and triplet parameter ranges (C6 = 3800 a.u.) including the  100
K uncertainty in cloud temperature. The hatched region indicates the parameter
space excluded by the absence of the f = 2 g -wave shape resonance.

Relative Transition Amplitude

Relative Transition Amplitude

134

v’=0

1.0

1.0

0.5

0.5

0.0

0.0

v’=3

1.0

v’=2

1.0

0.5

v’=4

0.5

0.0
300.0

350.0

400.0

450.0

Relative PA Laser Frequency (MHz)

0.0
300.0

350.0

400.0

450.0

Relative PA Laser Frequency (MHz)Figure 7.8: Comparison of experimental (closed circles) and theoretical (solid line)
J 0 = 2 linewidths. The synthetic spectra was calculated using our nominal set of potentials, as = 140 a.u., at = -17 a.u. and C6 =3800 a.u., assuming a cloud temperatureof 400 K.

135

f = 2 s-wave scattering length a2s . From our final potentials, we find a2s = -20  6 ) +35
0.055(C6 - C;45 a.u. Both at and a2s are in agreement with the results from the1u analysis[150]. In addition, our 40 K triplet scattering length is consistent with thevalue extracted from a recent elastic scattering measurement[119]. Finally, we haveused these refined Born-Oppenheimer potentials to survey the threshold scattering
properties of a number of K hyperfine-selected states. These considerations werereported in Ref.[152].

Chapter 8
Three-body recombination

As greater densities are attained in ultracold atom experiments, their description in terms of two-body interactions is no longer sufficient. In fact, condensate densities limited by three-body inelastic collisions[153] have been achieved. It
is conceivable that, in some regimes, an accurate theoretical treatment of mean-fieldenergetics in a dense condensate will need to account for both two- and three-bodyelastic collisions. Probably the most dramatic three-body effect measured to datein a cold atom trap is the catastrophic loss of atoms from a Na condensate that wasmagnetically tuned to a Feshbach resonance[15]. It is believed that this loss was
driven by three-body recombination, in which three free atoms coalesce allowing
two to form a molecule while the third atom carries away most of the excess binding
energy. Qualitatively, one might reason that the probability of three atoms colliding
increases near a two-body resonance simply because the effective range of interaction between the atoms, characterized by the two-body scattering length a, becomesquite large. However, the Na experiment measured dramatic recombination rates
even near a  0, causing losses that severely limited the ability of that experiment
to study the condensate near the resonance field value. It is not known whether thiseffect is unique to the Na atom or if it is instead an inevitable consequence of thethree-body recombination physics.
Three-body recombination has been a largely neglected area of theoretical ultracold collision physics. A recent theoretical effort[154] predicted an a4
scaling for the field-free recombination rate. However, the approximations usedin that treatment are valid only in the limit of large positive a. Recombination of
spin-polarized hydrogen[155, 156] was studied in the ’80’s and more recently inother alkali atoms[157], mostly within the confines of the Jastrow approximation.
Those treatments based the recombination rates on a perturbative formulation thatapproximates the initial three-body wave function as a product of two-body wave
functions. Using the Jastrow approximation, Moerdijk, et al.[157] obtained an a2
scaling for the recombination rate which is at odds with the scaling proposed in
Ref.[154]. In short, the physical processes that control three-body recombination
have been poorly understood even at a qualitative level.
In this chapter, our first steps towards understanding the recombination

137
process are presented. In particular, we develop a quantitative method that can beused to calculate recombination rates that are stable as a function of both hyperradius and three-body continuum channels. Based on calculations using a wide rangeof two-body model potentials, we have been able to uncover universal systematics of the recombination physics. This work was first reported in Ref.[158]. The
chapter is divided into two sections. In Section 8.1 the numerical method is developed. In Section 8.2 the results of our calculations are given and mechanisms forthe recombination process are discussed.
8.1 Adiabatic hyperspherical representation
The formulation presented here has its origins in the work of Delves, whoshowed[159] that the three-body continuum can be represented by a countable set ofchannels in a hyperspherical representation. Most importantly the scattering matrix
S is symmetric and unitary in this representation. The hyperspherical coordinate
method has been applied successfully to a number of diverse quantum mechanicalthree-body problems. (A review of these treatments can be found in Ref.[160]).
The basic approach involves a coordinate transformation that first separates out theunimportant center-of-mass motion. The remaining six degrees of freedom are represented by three Euler angles and three internal coordinates, the hyperradius R and
two hyperangles,  and . It is unlikely that a resonance will occur over the rangeof collision energies available to cold magnetically trapped atoms. The Wigner
threshold law leads to K3 / E J near threshold therefore, it is sufficient to concentrate on collisions with total angular momentum J = 0. In this case the collisions
are isotropic and the three-particle dynamics are described by a three-dimensional
Schrödinger equation, which can be written in terms of hyperspherical coordinates.
There are numerous conventions for defining the hyperradial coordinates.
In the case of identical particle scattering, the Smith and Whitten “democratic” coordinate system[161, 162] is the most useful. The definitions for the hyperangles (; )and hyperradius (R) are best developed in two stages, as shown in Ref.[163]. After
separating out the center-of-mass motion, the remaining internal degrees of freedom are defined in terms of Jacobi coordinates (~1 ; ~2 ). The relationship between
the Jacobi coordinates and the space-fixed atom coordinates ~ri is given by[163]~1 = (~r2 ; ~r1 )=d(m
1~r1 + m2~r2 )~1 = d ~r3 ; (m + m )
1

(8.1)

2

where mi represents the atomic mass of atom i and the mass-weighting factor d is

138
given byd2 = (mm3=+)(mm1++mm2)
1

2

(8.2)

3

m2 m3 :
2 = m m+1m
+m
1

2

3

The hyperradius and hyperangles are defined in terms of the Jacobi coordinates asfollows[161]:

R2 = 21 + 22

(8.3)

and

(1 )x = R cos  cos  ; (1 )y = ;R sin  sin (2 )x = R cos  sin  ; (2 )y = R sin  cos  :
0

0

0

0

(8.4)

Here, (x0 ; y 0 ) are Cartesian coordinates referred to the principal axes of the bodyfixed system. Qualitatively, the hyperradius is a measure of the size of the threebody system. With these definitions, the hyperspherical coordinates have convenient
symmetry properties that simplify the application of boundary conditions for identical particle scattering. In fact, for identical particle scattering it is only necessary
to define the hyperangles over the range 0    =4 and 0    =6.
In general, the three-body interaction V (R; ; ) can be approximated accurately by a sum of two-body interactions v (r12 ) + v (r23 ) + v (r31 ). The first correction to this expression, the Axilrod-Teller term[32], has an R;9 dependence andshould be unimportant for our problem. The three-body potentials used in this work
are discussed in more detail in the next section. The two-body internuclear radii rij
are defined in terms of the hyperradial coordinates by[163]r12 = 3;1=4 R[1 + cos(2) cos(2)]r23 = 3;1=4 R[1 + cos(2) cos(2 + 2=3)]
r31 = 3;1=4 R[1 + cos(2) cos(2 ; 2=3)] :(8.5)

The Schrödinger equation that describes the three particle dynamics isgiven in the hyperspherical representation by the following expression[163]:

 ;~2  1 @

R5 @





2(; ) + V (R; ; );2 R5 @R @R
R2

= E :

(8.6)

The first derivative operator can be removed from this expression by introducing a

139rescaled wave function

 ;~2  @2

= R5=2 . The Schrödinger equation for2 + 15=4 + V (R; ; );2 @R2R2

is then

= E :

(8.7)

Here, the grand angular momentum operator 2 is given by the following expression:@ sin 4 @ ; 1 @ 2 :2 = sin;14 @@ cos2 2 @2(8.8)

A standard approach to solve this three-dimensional Schrödinger equation begins
by writing the three-body wave function as an adiabatic channel decomposition(R; ; ) =

1
X
 =0

F (R) (R; ; ) :(8.9)

Here represents the two hyperangles and  represents the spin degrees of freedom,
while  indexes the channel functions. The adiabatic channel functions  (R; ; ),which vary smoothly in R except near avoided crossings, are eigenfunctions of the
fixed-R hyperspherical Hamiltonian ~2(2 + 15=4)2+ V (R; )  (R; ; ) = U (R) (R; ; ) : (8.10)

The adiabatic potentials U (R) form a set of effective potentials that can be used togain qualitative insight into the particle dynamics. Further, insertion of the channel
decomposition of (R; ; ) into Eq. 8.7 leads a set of coupled one-dimensional
equations ;~2  d2

d +Q +U ;E F =0:
1+
2P2 dR2dR(8.11)d  i and Q = h j d2  i represent coupling matrices. In
Here, P = h j dR
 dR2 2 = ;h d  j d  i component of Q (see Eq. B.2) is neededpractice, only the P



dR
dR
to solve the coupled equations, which can be added to the potential matrix to form
W = U ; 2~2 P 2.
Once the symmetric matrix W and antisymmetric matrix P are determined,
the coupled hyperradial equations can be solved using essentially the same approachthat was presented in Chapter 3 for two-body scattering. The coupled equations (Eq.
8.11) are solved using the adiabatic FEM R-matrix code described in Appendix B
subject to F (0) = 0 boundary conditions. The scattering matrix is then derived
from the R-matrix solutions according to the procedure developed in Section 3.2.
However, we do find it necessary to integrate out to surprisingly large R  105 a.u.

140
before matching to energy-normalized channel functions. A code that calculates W
and P by first expanding the angular coordinates in a B-splines basis[164] followed
by a diagonalization of the Hamiltonian using banded iterative techniques[165] hasbeen developed by Dr. B. Esry. In practice, we use this code to calculate W and P
over the hyperradial range 0  R < 8a, where a is the two-body scattering length.
Beyond, R  8a the matrix elements of W and P are smooth and can be fitted to
an inverse polynomial series. We use the following series expansions to define theasymptotic values of these matrix elements[166]:

W; = [c1 + Rc22 + Rc33 ]W; = [ Rc12 + Rc23 + Rc34 ]W ; = [ Rc12 + Rc23 + Rc34 ]P; = [ Rc11 + Rc22 + Rc33 ]P; = [ Rc51=2 + Rc72=2 + Rc93=2 ]P ; = [ Rc12 + Rc23 + Rc34 ] ;(8.12)

0

0

0

0

0

where ; 0 label the recombination channels and ; 0 index the three-body continuum channels. The coefficients ci are determined by a least squares fit to matrix
elements calculated with the B-splines code at large R. Tests indicate that the coefficients ci are relatively insensitive to the range of R used in the fit, provided that itcorresponds to hyperradii greater than  6a.
The last step needed is to relate the scattering matrix to the recombination
rate. A generalization of the two-body cross section derivation given in Section 2.1
has been developed[167], which defines the three-body recombination cross section
for identical particle scattering in the following manner: jS23 = 1152
k5 A2+A A+A+Aj :2

p(8.13)

Here, k =2E=~2 is the hyperradial wavenumber in the incident three-bodycontinuum channel. The large coefficient in Eq. 8.13 is a result of Bose-Einstein
statistics and is a (3!)2 -fold enhancement of the cross section compared to one thatwould be obtained for nonidentical particle scattering. The generalized cross section
has dimensions (length)5 as is expected for radial flux scattered in six dimensions.
The cross section defined here follows the Mott and Massey[24] convention for
identical particle scattering. That is, the cross section is defined as the ratio of
scattered radial flux to the incident flux in one of the six permutations of the incident
plane wave. As with two-body scattering, the event rate constant K3 per atom triplet

141
is defined in terms of the cross section as

K3 = ~k 3 :(8.14)

Finally, the atom-loss rate is given by

L3 = 3K3=6 :(8.15)

Here, the factor of three denotes the loss of both the diatom and the final state atomfrom the trap. The division by six comes from the fact that there are N 3 =6 atom
triples in a fixed volume. These formulae (equations 8.13, 8.14, and 8.15) were
derived for collisions in a thermal gas and need to be divided by a factor of six todescribe collisions in a condensate[168].
8.2 Recombination rates and mechanisms
The three-body interaction potential used is a sum of triplet Born-Oppenheimer
two-body potentials. This choice is appropriate for collisions of doubly-polarizedatoms. The major numerical challenge of a fully quantum mechanical approach
to solving the recombination problem is the sheer multiplicity of two-body boundstates. This is evident from Fig. 8.1, which shows the adiabatic hyperspherical potential curves for the collision of three doubly-polarized 7 Li atoms. Each bound state
supported by the two-body interaction is represented by exactly one hyperspherical three-body potential that converges asymptotically to a dimer bound state. For
doubly-polarized 7 Li, this amounts to 118 recombination channels and the number
grows still larger for the heavier alkalis, with 1466 recombination channels available
for doubly-polarized 87 Rb collisions.
However, it is physically plausible that the dominant pathways for recombination will involve vibrational states closest to threshold, particular those with low
angular momentum l. In fact, previous studies[154, 157, 169] suggest that the twobody scattering length a ultimately controls ultracold recombination rates. Basedon this assumption, we rescale the two-body triplet potentials by an overall scaling
factor , where 0   1. This scheme allows us to reduce the total number of recombination channels, which alleviates many of the numerical complexities, whileallowing us to control a. This premise is tested by treating two-body potentials that
support different numbers of bound states but which generate the same scattering
length. In addition, a number of calculations have been performed using a potentialof the form v (r12 ) = D sech2 (r12 =r0 ), which was utilized in the three-body studyof Ref.[169]. If the recombination rate is controlled by a, then the exact form of theshort-range potential should not matter.
We have performed calculations with roughly 125 different two-body potentials, including a wide range of scattering lengths and different numbers of two-

142

Figure 8.1: The J =0 hyperspherical potential curves for 7 Li are shown as a functionof hyperradius R. These potential curves describe doubly-polarized collisions of
three 7 Li atoms, i.e. all in the same fa = 2; mfa = 2 spin state. All of the roughly
120 bound and continuum channels shown must be included in order to carry out afull calculation of the scattering matrix. Data courtesy of Dr. B. Esry.

143
body bound states. The potentials used include Na-Na, K-K, and Rb-Rb rescaled
Born-Oppenheimer triplet potentials, and sech2 potentials using r0 =5 and 15 a.u.
Fig. 8.2 shows the resulting scattering length dependence of the “recombination
length”. The recombination length is defined asp 1=4
3 = ~ K3 :(8.16)

Here,  = m= 3 is the three-body reduced mass. This definition is intended topresent a quantity largely independent of the atomic properties other than the scattering length, and to reduce to a number roughly comparable to the scattering length
in order of magnitude.
The most striking result of Fig. 8.2 is the manner in which the calculated
recombination lengths cluster along a single curve. This implies that the recombination rates are in fact predominantly controlled by the two-body scattering length.
In addition, Fig. 8.2 shows two additional features of the recombination processthat had not appeared before in the literature. Specifically, an interference minimum
occurs near a  290 a.u. and a shape resonance occurs near a  ;175 a.u. Thesefeatures can be explained qualitatively by the two following mechanisms.
Recombination for positive scattering lengths is controlled primarily by adistant avoided crossing between the three-body entrance channel and the highests-wave recombination channel that occurs near Rpeak  3a. Figure 8.3(a) shows
the hyperspherical potential curves for the sech2 two-body potential scaled to give a= 100 a.u. but with only a single s-wave bound state. In this case the hypersphericalpotentials are very simple and their main features are consistent with Refs.[169,170]. Figure 8.3(a) also shows the crossing near Rpeak  3a, as a peak in thedimensionless coupling strength parameter, P 2 =U . The approximate hyperradial
peak position of P 2 =U is plotted in Fig. 8.4 for a number of rescaled alkali
potentials. The position of the crossing is consistent for the Rb-Rb and K-K rescaled
potentials but approaches Rpeak  4a for the large positive Na scattering lengths.
It is not yet understood whether this deviation from the “nominal” behavior of thecrossing position is an isolated occurrence.
The suppression of the recombination rates near a  290 a.u. is analogous
to that of the 87 Rb spin exchange rates discussed in Section 5.1. Namely, two competing pathways interfere destructively, which effectively “turns off” the flux exitingalong the inelastic path. One minor difference between these two cases of so-called“Stückelberg oscillations” is that here the kinetic energy in the entrance channel isnegative at the crossing. This case has been treated by Nakamura[171] who calls
this process “nonadiabatic tunneling”. The scattering probability for recombination

144

6000
Rb−Rb
2
D sech (rij/5)
2
D sech (rij/15)
Eq. 8.19
Na−Na
K−K

5000

cs ρ3 (a.u.)

4000
3000
2000
1000
0
−750

−500

−250
0
cs a (a.u.)

250

500

Figure 8.2: Numerically calculated recombination lengths for a wide range of twobody potentials, which are discussed in the text. We find that each set of two-bodypotentials has a characteristic length scale, labeled r0 . In the case of the sech2 potentials, r0 is part of the definition of the potential, i.e. v (r12 ) = D sech2 (r12 =r0 ),r0= 5 and 15 used here . The value of r0 for the alkali potentials is determined“empirically” from the calculated position of the interference minimum and a zeroof the phase factor sin2 [ln(3a=2r0 )]. (The phase factor is discussed in the text.) We
find r0 (Rb) = 18.8, r0 (K) = 16.8, and r0 (Na) = 12.4 (in a.u.). The dimensionless
scale factor is set to cs = r0 (Rb)/r0 (Vi ), i.e. cs = 1 for the Rb potentials. Using this
rescaling, each set of potentials produces an interference minimum near cs a  290
a.u. and a shape resonance near cs a  175 a.u. The solid line is the analytical formula given in Eq. 8.19 with the sin2 phase factor set equal to 12 . These calculations
were performed at E = 1K, except for the range a  ;300 a.u., where an energyof 0.1K was used in order to reach the ultracold limit.

145

0.03

2.5
0.0
0.01
−2.5
−5.0

2

0.02

P /(µ∆U)−8

U (10 a.u.)

5.0(a) a=+100 a.u.
0

100

200

300

400

0.00
500

R (a.u.)
0.08
0.01

0.06
0.00
100 200 300 400

0

2

0.04

P /(µ∆U)−8

U (10 a.u.)50(b) a=−100 a.u.

−50
0

50

100

150

0.02
0.00
200

R (a.u.)Figure 8.3: Different qualitative mechanisms control ultracold recombination at a >0 (a) and a < 0 (b), as is clear from the lowest two adiabatic hyperspherical potential
curves (solid line) depicted along with their nonadiabatic coupling strength (dashedline). We define a dimensionless measure of the nonadiabatic coupling strength
2 (R), thein terms of the the squared first derivative coupling matrix element P12
difference in the adiabatic hyperspherical potential curves U (R), and the reduced
mass . These were calculated using the sech2 two-body potential (r0 = 15 a.u.)
discussed in the text, for the simplest situation of a single two-body s-wave boundstate. The basic nature of these potentials remains the same even as the number of
two-body bound states is increased. The inset of (b) shows the potential barrier inthe three-body entrance channel discussed in the text.

146

4.00

Rpeak/a

3.50

3.00
K−K
Na−Na
Rb−Rb

2.50

2.00

0

50

100 150 200 250 300 350 400
a (a.u.)Figure 8.4: Hyperradial position of the long-range avoided crossing Rpeak as a function of two-body scattering length a for a number of alkali potentials rescaled tosupport only one s-wave bound state.

147
into the highest s-wave dimer channel is given approximately by the expression:

jSA2 +A A+A+Aj2;! 0:067 (ka)4 sin2 ;
a ! +1(8.17)

where  is the phase difference between the amplitudes associated with the twocompeting pathways. An estimate of  can be made in the limit of large a byapproximating the hyperspherical potential curve converging to the highest s-wave
dimer potential as an “Efimov potential curve”[169], i.e. V ! ;0:631=R2 , forhyperradii less than the crossing point. In this limit, the phase factor at zero energy
is approximately

sin2 

;!
a ! +1sin2

  3a 
ln 2r
0

;

(8.18)

where r0 defines a characteristic length scale for the two-body potentials. (Thisis discussed in more detail in the caption of Fig. 8.2.) We expect that sin2  canbe approximated by 12 in most systems having numerous two-body bound states
because the amplitude of the interfering returning wave is likely to be dissipatedinto other recombination channels at small R. This point will require further study.
However there is one atom, 4 He, for which the physical two-body potentials do
support only one s-wave bound state and in fact, our calculations show that theinterference minimum occurs near the physical value of the scattering length. Thisis shown in Fig. 8.5.
Returning to Fig. 8.3, the (b) panel shows the hyperspherical potential
curves, for a sech2 two-body potential that generates a scattering length of a = ;100a.u., which supports a single s-wave bound state. Again, the hyperspherical potential
curves are simple, but for the negative a case we find that the rate-limiting crossing
occurs at small R. In fact, the hyperradius of the crossing and its coupling strength
are largely independent a. The scattering length dependence of the recombination
rates in the range a < 0 is the result of quantum mechanical tunneling under abarrier (previously identified in Ref.[169]), whose peak occurs at RB  2a with a:079 in a.u. (see the inset of Fig. 8.3(b)). In ordermaximum height of W (RB )  0a2
to reach the avoided crossing(s) at small R, the incoming wavefront in the threebody recombination channel must tunnel inside the barrier. If the incident particle
does tunnel successfully, then we expect that it has a high probability to undergo anonadiabatic Landau-Zener transition into one or more recombination channels. On
the other hand, the probability per collision to tunnel into the small R region is astrong function of the barrier height and the collision energy. In fact, the probability
varies as (ka)4 at energies well below the barrier peak. This leads to an approximate
a4 dependence of the recombination rates for both positive and negative scattering
lengths even though the recombination mechanisms are different. The existence of

148

800

ρ3 (a.u.)

600

400

200

0
100

125

150

175
200
a (a.u.)225

250

Figure 8.5: The recombination length versus scattering length a for collisions of
three 4 He atoms (closed circles). The physical scattering length is a=172 a.u. The
solid line is the analytical recombination length defined in Eq. 8.19 using r0 =9.4.

149
a potential barrier also provides the opportunity for a shape resonance in the threebody continuum to dramatically increase the recombination rates. This effect isdemonstrated in Fig. 8.2 near a  ;175 a.u. However, in the absence of a fullcalculation the positions of these resonances will be hard to predict for any real
alkali system.
We have developed the following expression for the recombination length
based on the roughly 125 calculations using a variety of two-body potentials:

9
88 ;
1
h

i

<
=
>
4
>
5:0 a ; 32 r0 2 sin2 ln 23ra0 ; a > 0 ;<;th3 (a) ' max >: ;7:7a; a < 0>:(8.19)

4r0;

where the slopes in this expression are accurate to about 10%. The lower limit
is set to 3 = 4r0 because we never see the recombination length decrease below
50-100 a.u. in the alkalis. This effect is presumably due to other recombination
channels becoming more important as the dominant s-wave channel weakens sufficiently. A comparison of our analytic expression for 3 with the numerically calculated +a recombination lengths is shown in Fig. 8.6.
Table 8.1 provides our best estimate of the recombination rates and compares these with the available experimental information and other theoretical predictions. For the estimates of alkali rates, we have replaced the sin2  phase factor
in Eq. 8.19 with its average value of 1/2. This was not necessary for 4 He since
the calculation of its recombination rate used a realistic two-body interaction potential. Table 8.1 presents the results in terms of event rates K3 . Where necessary
the quoted experimental loss rates L3 have been converted to event rates using Eq.
8.15. In addition, the experimental results quoted for atoms in a condensate havebeen multiplied by a factor of 6 to compare with our theoretical results for atoms
in a thermal cloud[168]. Our estimates are in reasonable agreement with most of
the measured rates. The lone experiment that deviates substantially from our predicted dependence of K3 on a is the magnetic field Feshbach resonance studied inNa[15, 172]. In particular, that experiment observed an increase in the recombination rate as a was decreased from its nominal value of a = 54:6 a.u. to a = 0. We
have not observed that type of behavior in our calculations, which might suggestthat the mechanism controlling the Na recombination differs from those discussed
here.
The majority of the calculations performed to date used two-body potentials scaled to include only a single s-wave bound state. In order to assess the generality of our results, it is important to study the dependence of K3 on the numberof s-wave bound states, nS , for a given a. Our preliminary results indicate that thedependence is not particularly strong. For example, the recombination event rate

150

2000
Rb−Rb
2
D sech (rij/5)
2
D sech (rij/15)
Eq. 8.19
Na−Na
K−K

cs ρ3 (a.u.)

1500

1000

500

0

0

100

200
300
cs a (a.u.)

400

500

Figure 8.6: Numerical recombination lengths are compared with the +a analytic
formula given in Eq. 8.19 (r0 = 18:8). The numerical data are discussed in thecaption of Fig. 8.2.

Table 8.1: Recombination event rate constants
Atom
4 He
7 Li
23 Na(1,1)
39 K
85 Rb
87 Rb(1,;1)87 Rb(2,2)
133 Cs(3,3)

presenttha (a.u.) K3 (cm6/s)
172
4:2  10;27;27.6 2:5  10;28
54.6
1:1  10;28;33 9:2  10;29;370 6:6  10;25
106
2:0  10;28
106
2:0  10;28;1250 5:6  10;23

experiment
K3exp (cm6/s)—–< 6  10;27[173]1  10;27[172]—–—–9  10;29[168]2:2  10;28[174]< 10;23[18]

other
th
K3 (cm6 /s)
2:5  10;27[154]5:2  10;28[157]
4:0  10;28[157]—–—–
8  10;30[157]
1:7  10;29[154]1  10;28[175]151

K3 was found to decrease by a factor of three when the potentials were rescaled
from nS = 1 to nS = 3 while maintaining a = +100. Examining the bound state

dependence of the position and strength of the rate-limiting crossing suggests that
this difference will not grow much larger (see Fig. 8.7) as nS is increased farther.
The negative scattering length recombination rates have an even weaker dependence
on nS . In fact, our preliminary tests indicate that, excluding resonances, the recombination rate is essentially converged with nS = 1, (see Fig. 8.8). However, theseresults are not conclusive and we plan to conduct more tests in the future.
It is also important to understand the energy dependence of K3 and theenergy range over which the threshold values of K3 are applicable. For positive
scattering lengths, we find that K3 is reasonably constant over collision energies
ranging from E = 0 to several hundred K. However, owing to the low energy
barrier in the entrance channel, K3 is far more energy-dependent for the negative
scattering lengths. In fact, the threshold limit for K3 is achieved at a < 0 only
for collision energies well below the barrier maximum. At higher energies, thescattering probability saturates and K3 decreases due to the k ;4 factor (see equations
8.13 and 8.14). The energy dependence of the recombination rates for a are shownin Fig. 8.9. This behavior should have a dramatic effect on the 133 Cs recombination
rates if a  ;1250 as quoted in Ref.[18]. In this case, the peak of the barrier in thethree-body entrance channel is roughly W (RB ) = 0.1K. Our experience suggeststhat K3 will decrease several orders of magnitude as the collision energy rises to
10K. The resulting rate would then be consistent with the experimental bound set
in Ref.[18] even though our threshold value is greater than the bound quoted in
Ref.[18](see Table 8.1).

152

0.03

1

P /(µ∆U)

0.02

2

2

3
6

0.01

0
150

200

250
300
R (a.u.)

350

4

400

Figure 8.7: Coupling between the three-body entrance channel and the highest swave recombination channel is shown as a function of hyperradius. Each curve
represents a Rb-Rb potential rescaled to hold nS s-wave bound states with a=+100.
The labels on each curve denote the number nS .

153

40.0
30.0

W(RB) (µK)

1
20.0

3

10.0
0.0
−10.0
100

300

500
R (a.u.)700

900

Figure 8.8: The centrifugal barrier in the three-body entrance channel is shown as afunction of hyperradius. Each curve represents a Rb-Rb potential rescaled to holdnS s-wave bound states with a = -100 a.u. Curves with nS = 1 and 3 are labeled,nS = 2 is the unlabeled curve. The height and position of the barrier are essentiallyindependent of nS . The event rates K3 calculated for these three potentials varied
by roughly 30%.

154

−25

10

a=−1006

K3 (cm /sec)

−26

10

a=+195
−27

10

−28

10

−7

10

−6

10

−5

10
E (K)

−4

10

−3

10

Figure 8.9: The energy dependence of K3 is shown for a. These data used RbRb potentials rescaled to support one s-wave bound state. The negative scattering
length recombination rate exhibits an E ;2 dependence for collision energies above
the barrier in the three-body entrance channel (W (RB )  30K for this example).

Chapter 9
SummaryThis dissertation has extended the theoretical methods needed to calculate
and understand two-body and three-body collisions of trapped alkali atoms. The“conventional” theory of two-body collisions in a magnetic trap was developed in
Chapter 2. Numerical techniques for solving the resulting coupled differential equations that describe the collision dynamics were presented in Chapter 3. Chapter 4
developed an alternative formalism to describe the two-body dynamics. Here, multichannel quantum defect methods were applied to the two-body collisions, yielding amore comprehensive understanding of the collision physics as well as substantially
improving the efficiency of the calculations. These theoretical techniques were applied in Chapter 5 to a set of Rb collision experiments. From an analysis of theseexperiments, we were able to derive state-of-the art Rb2 Born-Oppenheimer potentials. In addition, we used the refined potentials to investigate the scattering properties of several Rb isotopes. From these surveys, we have been able to identify
specific hyperfine states whose collisional properties could lead to novel degenerate
gas experiments. The theory of two-body collisions in a magnetic trap seems well inhand. However, with the recent demonstration of trapping atoms outside the alkalifamily[176, 177, 178] there are certainly more surprises to come.
Photoassociation spectroscopy has been an essential tool in cold collision
physics. A theory of photoassociation lineshapes was presented in Chapter 6. This
theory was applied in Chapter 7 to interpret measured 39 K rovibrational spectra,
which permitted us to derive improved ground state potentials. The analysis ofphotoassociation lineshapes provides vital information on the ground state scatteringwave functions. However, attempts to better quantify these results have sufferedfrom uncertainties associated with the second, ionizing, laser. Theories to include
multiple transitions[134, 135] have been devised but remain untested. These issues
will need to be addressed in the future.
A new theoretical description of three-body recombination in collisionsof doubly-polarized alkali atoms, was presented in Chapter 8. This theory is a fully
quantum mechanical treatment of the recombination process. From our numericalcalculations we have been able to identify two separate mechanisms that control therecombination physics depending on the sign of the two-body scattering length a.

156
Remarkably, we find that the recombination rate scales as a4 in both cases, although
the coefficient is much larger at a < 0. The treatment presented here used rescaled
two-body potentials to simplify the numerical calculations. In the future, it will benecessary to test our results by using more realistic two-body potentials. In addition,
it seems conceivable to introduce both hyperfine and magnetic field interactions into
our formalism. This would then allow a more realistic treatment of the Na Feshbachresonance experiment. Three-body recombination physics remains in its infancy,
but Chapter 8 develops a framework upon which future studies can build.

Bibliography[1] M. H. Shamos, Great Experiments in Physics, (Dover, New York, 1987).
[2] D. W. Gough, G. C. Maitland and E. B. Smith, Mol. Phys. 24 (1972) 1551.
[3] R. A. Dawe and E. B. Smith, J. Chem. Phys. 52 (1970) 693.
[4] F. A. Guevara, B. B. McInteer and W. E. Wageman, Phys. Fluids 12 (1969)
2493.
[5] J. Kestin, H. E. Khalifa, S. T. Ro and W. A. Wakeham, Physica 88A (1978)
242.
[6] W. D. Phillips, J. V. Prodan and H. J. Metcalf, J. Opt. Soc. Am. B 2 (1985)
1751.
[7] J. Dalibard and C. Cohen-Tannoudji, J. Opt. Soc. Am. B 2 (1985) 1707.
[8] M. H. Anderson, J. R. Ensher, M. R. Matthews, C. E. Wieman and E. A.
Cornell, Science 269 (1995) 198.
[9] C. C. Bradley, C. A. Sackett, J. J. Tollet and R. G. Hulet, Phys. Rev. Lett. 75(1995) 1687.
[10] K. B. Davis, M. O. Mewes, M. R. Andrews, N. J. van Druten, D. D. Durfee,
D. M. Kim and W. Ketterle, Phys. Rev. Lett. 75 (1995) 3969.
[11] P. A. Ruprecht, M. J. Holland, K. Burnett and M. Edwards, Phys. Rev. A 51
(1995) 4704.
[12] F. Dalfovo, S. Giorgini, L. P. Pitaevskii and S. Stringari, Rev. Mod. Phys. 71(1999) 463.
[13] C. W. Gardiner, M. D. Lee, R. J. Ballagh, M. J. Davis and P. Zoller, Phys.
Rev. Lett. 81 (1998) 5266.
[14] W. C. Stwalley, Phys. Rev. Lett. 37 (1976) 1628.

158[15] S. Inouye, M. R. Andrews, J. Stenger, H. J. Miesner, D. M. Stampur-Kurnand W. Ketterle, Nature (London) 392 (1998) 151.
[16] Ph. Courteille, R. S. Freeland, D. J. Heinzen, F. A. van Abeelen and B. J.
Verhaar, Phys. Rev. Lett. 81 (1998) 69.
[17] J. L. Roberts, N. R. Claussen, Jr. J. P. Burke, C. H. Greene, E. A. Cornell and
C. E. Wieman, Phys. Rev. Lett. 81 (1998) 5109.
[18] V. Vuletić, A. J. Kerman, C. Chin and S. Chu, Phys. Rev. Lett. 82 (1999)
1406.
[19] H. R. Thorsheim, J. Wiener and P. S. Julienne, Phys. Rev. Lett. 58 (1987)
2420.
[20] H. B. G. Casimir and D. Polder, Phys. Rev. 73 (1948) 360.
[21] K. M. Jones, P. S. Julienne, P. D. Lett, W. D. Phillips, E. Tiesinga and C. J.
Williams, Europhys. Lett. 35 (1996) 85.
[22] H. Wang, J. Li, X. T. Wang, C. J. Williams, P. L. Gould and W. C. Stwalley,
Phys. Rev. A 55 (1997) R1569.
[23] M. Aymar, C. H. Greene and E. Luc-Koenig, Rev. Mod. Phys. 68 (1996)
1016.
[24] N. F. Mott and H. S. W. Massey, The Theory of Atomic Collisions, (Clarendon Press, Oxford, 1965).
[25] L. S. Rodberg and R. M. Thaler, Introduction to the Quantum Theory of
Scattering, (Academic Press, New York, 1967).
[26] S. Gasiorowicz, Quantum Physics, (John Wiley and Sons, New York, 1974).
[27] U. Fano and A. R. P. Rau, Atomic Collisions and Spectra, (Academic Press,
New York, 1986).
[28] J. Wiener, V. S. Bagnato, S. Zilio and P. S. Julienne, Rev. Mod. Phys. 71(1999) 1.
[29] D. H. Levy, Adv. Mag. Res. 6 (1973) 1.
[30] A. Carrington, D. H. Levy and T. A. Miller, Adv. Chem. Phys. 18 (1970) 149.
[31] M. Born and J. R. Oppenheimer, Ann. Phys. (Leipzig) 84 (1927) 457.
[32] G. C. Maitland, M. Rigby, E. B. Smith and W. A. Wakeham, Intermolecular
Forces, (Clarendon Press, Oxford, 1981).

159[33] M. Krauss and W. J. Stevens, J. Chem. Phys. 93 (1990) 4236.
[34] A. L. G. Rees, Proc. Phys. Soc. Lond. A 59 (1947) 998.
[35] C. R. Vidal, Comments At. Mol. Phys. 17 (1986) 173.
[36] B. M. Smirnov and M. I. Chibisov, JETP 21 (1965) 624.
[37] R. Côté, A. Dalgarno, H. Wang and W. C. Stwalley, Phys. Rev. A 57 (1998)
R4118.
[38] M. Marinescu, H. R. Sadeghpour and A. Dalgarno, Phys. Rev. A 49 (1994)
982.
[39] A. Derevianko, W. R. Johnson, M. S. Safronova and J. F. Babb, Phys. Rev.
Lett. 82 (1999) 3589.
[40] S. H. Patil and K. T. Tang, J. Chem. Phys. 106 (1997) 2301.
[41] W. T. Zemke and W. C. Stwalley, J. Chem. Phys. 97 (1993) 2053.
[42] Z.-C. Yan, J. F. Babb, A. Dalgarno and G. W. F. Drake, Phys. Rev. A 54
(1996) 2824.
[43] W. T. Zemke and W. C. Stwalley, J. Chem. Phys. 100 (1994) 2661.
[44] M. Marinescu, (private communication).
[45] C. Amiot, J. Verges and C. E. Fellows, J. Chem. Phys. 103 (1995) 3350.
[46] G. Zhao, W. T. Zemke, J. T. Kim, B. Ji, H. Wnag, J. T. Bahns, W. C. Stwalley,
L. Li, A. M. Lyyra and C. Amiot, J. Chem. Phys. 105 (1996; and references
therein) 7976.
[47] C. Amiot, J. Chem. Phys. 93 (1990) 8591.
[48] W. Weickenmeier, U. Diemer, M. Wahl, M. Raab and W. Demtröder, J. Chem.
Phys. 82 (1985) 5354.
[49] E. Tiesinga, C. J. Williams, P. S. Julienne, K. M. Jones, P. D. Lett and W. D.
Phillips, J. Res. Natl. Inst. Stand. Technol. 101 (1996) 505.
[50] C. Cohen-Tannoudji, B. Diu and F. Laloe, Quantum Mechanics, (John Wileyand Sons, New York, 1977).
[51] J. D. Jackson, Classical Electrodynamics, (John Wiley and Sons, New York,
1999).

160[52] E. Arimondo, M. Inguscio and P. Violino, Rev. Mod. Phys. 49 (1977) 31.
[53] H. Lefebvre-Brion and R. W. Field, Perturbations in the Spectra of Diatomic
Molecules, (Academic Press, New York, 1986).
[54] F. H. Mies, C. J. Williams, P.S. Julienne and M. Krauss, J. Res. Natl. Inst.
Stand. Technol. 101 (1996) 521.
[55] P. O. Fedichev, Yu. Kagan, G. V. Shlyapnikov and J. T. M. Walraven, Phys.
Rev. Lett. 77 (1996) 2913.
[56] J. L. Bohn and P. S. Julienne, Phys. Rev. A 56 (1997) 1486.
[57] R. N. Zare, Angular Momentum, (John Wiley and Sons, New York, 1988).
[58] B. Numerov, Publ. Observatoire Central Astrophys. Russ. 2 (1933) 188.
[59] R. G. Gordon, J. Chem. Phys. 51 (1969) 14.
[60] B. R. Johnson, J. Comp. Phys. 13 (1973) 445.
[61] L. D. Thomas, M. H. Alexander, B. R. Johnson, Jr. W. A. Lester, J. C. Light,
K. D. McLenithan, G. A. Parker, M. J. Redmon, T. G. Schmalz, D. Secrestand R. B. Walker, J. Comp. Phys. 41 (1981) 407.
[62] W. H. Press, S. A. Teukolsky, W. T. Vetterling and B. P. Flannery, Numerical
Recipes in Fortran: The Art of Scientific Computing, (Cambridge University
Press, New York, 1992).
[63] E. Anderson, Z. Bai, C. Bischof, J. Demmel, J. Dongarra, J. Du Cros, A.
Greenbaum, S. Hammerling, A. Mckenney, S. Ostouchov and D. Sorensen,
LAPACK Users Guide, (SIAM, Philadelphia, 1992).
[64] E. P. Wigner, Phys. Rev. 70 (1946) 15.
[65] E. P. Wigner and L. Eisenbud, Phys. Rev. 72 (1947) 29.
[66] T. Teichman and E. P. Wigner, Phys. Rev. 87 (1952) 123.
[67] C. H. Greene, Phys. Rev. A 28 (1983) 2209.
[68] H. Le Rouzo and G. Raseev, Phys. Rev. A 29 (1984) 1214.
[69] L. R. Ram-Mohan, S. Saigal, D. Dossa and J. Shertzer, Comput. Phys. 4(1990) 50.
[70] J. Shertzer and J. Botero, Phys. Rev. A 49 (1994) 3673.

161[71] K. W. Meyer, C. H. Greene and B. D. Esry, Phys. Rev. Lett. 78 (1997) 4902.
[72] H. W. van der Hart, J. Phys. B 30 (1997) 453.
[73] K. J. Bathe, Finite Element Procedures in Engineering Analysis, (Prentice
Hall, Englewood Cliffs, NJ, 1976).
[74] K. J. Bathe and E. Wilson, Numerical Methods in Finite Element Analysis,(Prentice Hall, Englewood Cliffs, NJ, 1976).
[75] R. A. Gonzales, S. Y. Kang, I. Koltracht and G. Rawitscher, (submitted to J.
Comput. Phys.).
[76] G. H. Rawitscher, B. D. Esry, E. Tiesinga, Jr J. P. Burke and I. Koltracht,(submitted to J. Chem. Phys.).
[77] T. G. Schmalz, E. B. Stechel and J. C. Light, J. Chem. Phys. 70 (1979) 5660.
[78] H. Klar, Phys. Rev. A 15 (1976) 1452.
[79] H. Klar and U. Fano, Phys. Rev. Lett. 37 (1976) 1132.
[80] M. J. Seaton, Rep. Prog. Phys. 46 (1983) 167.
[81] C. M. Lee, Phys. Rev. A 11 (1975) 1692.
[82] F. H. Mies, J. Chem. Phys. 80 (1984) 2514.
[83] F. H. Mies and P. S. Julienne, J. Chem. Phys. 80 (1984) 2526.
[84] F. H. Mies and M. Raoult, (unpublished).
[85] M. J. Cavagnero, Phys. Rev. A 50 (1994) 2841.
[86] M. J. Cavagnero, (unpublished).
[87] B. Gao, Phys. Rev. A 54 (1996) 2022.
[88] B. Gao, Phys. Rev. A 58 (1998) 4222.
[89] J. P. Burke, C. H. Greene and J. L. Bohn, Phys. Rev. Lett. 81 (1998) 3355.
[90] W. E. Milne, Phys. Rev. 35 (1930) 863.
[91] B. Yoo and C. H. Greene, Phys. Rev. A 34 (1986) 1635.
[92] C. H. Greene, A. R. P. Rau and U. Fano, Phys. Rev. A 26 (1982) 2441.
[93] H. J. Korsch and H. Laurent, J. Phys. B 14 (1981) 4213.

162[94] G. F. Gribakin and V. V. Flambaum, Phys. Rev. A 48 (1993) 546.
[95] A. Dalgarno and M. R. H. Rudge, Proc. R. Soc. London, Ser. A 286 (1965)
519.
[96] H. T. C. Stoof, J. M. V. A. Koelman and B. J. Verhaar, Phys. Rev. B 38 (1988)
4688.
[97] D. M. Chase, Phys. Rev. 104 (1956) 838.
[98] R. K. Nesbet, Comments At. Mol. Phys. 11 (1981) 25.
[99] H. Gao and C. H. Greene, J. Chem. Phys. 91 (1989) 3988.
[100] C. J. Myatt, E. A. Burt, R. W. Ghrist, E. A. Cornell and C. E. Wieman, Phys.
Rev. Lett. 78 (1997) 586.
[101] H. M. J. M. Boesten, C. C. Tsai, J. R. Gardner, D. J. Heinzen and B. J.
Verhaar, Phys. Rev. A 55 (1997) 636.
[102] H. M. J. M. Boesten, C. C. Tsai, B. J. Verhaar and D. J. Heinzen, Phys. Rev.
Lett. 77 (1996) 5194.
[103] J. P. Burke, J. L. Bohn, B. D. Esry and C. H. Greene, Phys. Rev. A 55 (1997)
R2511.
[104] P. S. Julienne, F. H. Mies, E. Tsienga and C. J. Williams, Phys. Rev. Lett. 78(1997) 1880.
[105] S. J. J. M. F. Kokkelmans, H. M. J. M. Boesten and B. J. Verhaar, Phys. Rev.
A 55 (1997) R1589.
[106] A. Dalgarno and M. R. H. Rudge, Proc. R. Soc. London, Ser. A 286 (1965)
519.
[107] J. M. Vogels, C. C. Tsai, R. S. Freeland, S. J. J. M. F. Kokkelmans, B. J.
Verhaar and D. J. Heinzen, Phys. Rev. A 56 (1997) R1067.
[108] C. R. Monroe, E. A. Cornell, C. A. Sackett, C. J. Myatt and C. E. Wieman,
Phys. Rev. Lett. 70 (1993) 414.
[109] C. C. Tsai, R. S. Freeland, J. M. Vogels, H. M. J. M. Boesten, B. J. Verhaar
and D. J. Heinzen, Phys. Rev. Lett. 79 (1997) 1245.
[110] M. R. Matthews, D. S. Hall, D. S. Jin, J. R. Ensher, C. E. Wieman and E. A.
Cornell, Phys. Rev. Lett. 81 (1998) 243.

163[111] J. P. Burke, J. L. Bohn, B. D. Esry and C. H. Greene, Phys. Rev. Lett. 80(1998) 2097.
[112] H.-J. Miesner, D. M. Stamper-Kurn, J. Stenger, S. Inouye, A. P. Chikkatur
and W. Ketterle, Phys. Rev. Lett. 82 (1999) 2228.
[113] E. Cornell, (private communication).
[114] Z.-T. Lu, K. L. Corwin, K. R. Vogel, C. E. Wieman, T. P. Dinneen, J. Maddi
and H. Gould, Phys. Rev. Lett. 79 (1997) 994.
[115] J. E. Simsarian, L. A. Orozco, G. D. Sprouse and W. Z. Zhao, Phys. Rev. A
57 (1998) 2448.
[116] Z.-T. Lu, C. J. Bowers, S. J. Freedman, B. K. Fujikawa, J. L. Mortara, S.-Q.
Shang, K. P. Coulter and L. Young, Phys. Rev. Lett. 72 (1994) 3791.
[117] R. Guckert, X. Zhao, S. G. Crane, A. Hime, W. A. Taylor, D. Tupa, D. J.
Vieira and H. Wollnik, Phys. Rev. A 58 (1998) R1637.
[118] J. P. Burke and J. L. Bohn, Phys. Rev. A 59 (1999) 1303.
[119] B. DeMarco, J. L. Bohn, J. P. Burke, M. Holland and D. S. Jin, Phys. Rev.
Lett. 82 (1999) 4208.
[120] http:/condon.colorado.edu/ chg/Collisions.
[121] H. T. C. Stoof, Phys. Rev. Lett. 76 (1996) 10.
[122] A. G. K. Modowi and A. J. Legget, J. Low Temp. Phys. 104 (1997) 625.
[123] T. L. Ho and V. B. Shenoy, Phys. Rev. Lett. 77 (1996) 3276.
[124] B. D. Esry, C. H. Greene, J. P. Burke and J. L. Bohn, Phys. Rev. Lett. 78(1997) 3594.
[125] B. D. Esry, Many-Body Effects in Bose-Einstein Condensates of Dilute
Atomic Gases, PhD thesis University of Colorado, Boulder 1997.
[126] H. Pu and N. Bigelow, Phys. Rev. Lett. 80 (1998) 1130.
[127] B. D. Esry and C. H. Greene, Nature (London) 392 (1998) 434.
[128] M. Movre and G. Pichler, J. Phys. B 10 (1977) 2631.
[129] W. C. Stwalley, Y.-H. Uang and G. Pichler, Phys. Rev. Lett. 41 (1978) 1164.
[130] R. Napolitano, J. Weiner, C. J. Williams and P. S. Julienne, Phys. Rev. Lett.
73 (1994) 1352.

164[131] K. Burnett, J. Cooper, R. J. Ballagh and E. W. Smith, Phys. Rev. A 22 (1980)
2005.
[132] K. Burnett and J. Cooper, Phys. Rev. A 22 (1980) 2027.
[133] K. Burnett and J. Cooper, Phys. Rev. A 22 (1980) 2044.
[134] J. L. Bohn and P. S. Julienne, Phys. Rev. A 54 (1996) R4637.
[135] J. L. Bohn and P. S. Julienne, Phys. Rev. A 60 (1999).
[136] E. Tiesinga, C. J. Williams, P. S. Julienne, K. M. Jones, P. D. Lett and W. D.
Phillips, J. Res. Natl. Inst. Stand. Technol. 101 (1996) 505.
[137] G. Herzberg, Molecular Spectra and Molecular Structure, Vol. 1, (Krieger,
Malabar, Fl., 1989).
[138] D. M. Brink and G. R. Satchler, Angular Momentum, (Clarendon Press,
Oxford, 1968).
[139] M. Larsson, Physica Scripta 23 (1981) 835.
[140] R. W. Heather and P. S. Julienne, Phys. Rev. A 47 (1993) 1887.
[141] J. T. Hougen, N. B. S. Monograph 115 (1970) 1.
[142] I. Roeggen, Theor. Chim. Acta. 21 (1971) 398.
[143] M. Ya. Ovchinnikova, Theoret. Exptl. Chem. 1 (1965) 22.
[144] H. Wang, P. L. Gould and W. C. Stwalley, J. Chem. Phys. 106 (1997) 7899.
[145] W. J. Meath, J. Chem. Phys. 48 (1968) 227.
[146] H. Wang, P. L. Gould and W. C. Stwalley, Phys. Rev. A 53 (1996) R1216.
[147] J. P. Burke, C. H. Greene, J. L. Bohn, H. Wang, P. L. Gould and W. C. Stwalley, (submitted to Phys. Rev. A).
[148] R. Côté, A. Dalgarno, H. Wang and W. C. Stwalley, Phys. Rev. A 57 (1998)
R4118.
[149] H. M. J. M. Boesten, J. M. Vogels, J. G. C. Tempelaars and B. J. Verhaar,Phys. Rev. A 54 (1996) R3726.
[150] C. J. Williams, E. Tiesinga, P. S. Julienne, H. Wang, W. C. Stwalley and P. L.
Gould, (submitted to Phys. Rev. A).

165[151] H. Wang, X. T. Wang, P. L. Gould and W. C. Stwalley, Phys. Rev. Lett. 78(1997) 4173.
[152] J. L. Bohn, J. P. Burke, C. H. Greene, H. Wang, P. L. Gould and W. C. Stwalley, Phys. Rev. A 59 (1999) 3660.
[153] D. M. Stampur-Kurn, M. R. Andrews, A. P. Chikkatur, S. Inouye, H. J. Miesner, J. Stenger and W. Ketterle, Phys. Rev. Lett. 80 (1998) 2027.
[154] P. O. Fedichev, M. W. Reynolds and G. V. Shlyapnikov, Phys. Rev. Lett. 77(1996) 2921.
[155] Yu. Kagan, I. A. Vartan’yants and G. V. Shlyapnikov, Zh. Eksp. Teor. Fiz. 81(1981) 1113.
[156] L.P.H. de Goey, H.T.C. Stoof, B. J. Verhaar and W. Glöckle, Phys. Rev. B 38
(1988) 646.
[157] A. J. Moerdijk, H.M.J.M. Boesten and B.J. Verhaar, Phys. Rev. A 53 (1996)
916.
[158] B. D. Esry, C. H. Greene and J. P. Burke, (submitted to Phys. Rev. Lett.).
[159] L. M. Delves, Nuclear Physics 20 (1960) 275.
[160] C. D. Lin, Phys. Rep. 257 (1995) 1.
[161] R. C. Whitten and F. T. Smith, J. Math. Phys. 9 (1968) 1103.
[162] A. Kupperman, J. Phys. Chem. 100 (1996) 2621.
[163] J. G. Frey and B. J. Howard, Chem. Phys. 99 (1985) 415.
[164] C. de Boor, A Practical Guide to Splines, (Springer, New York, 1978).
[165] R.
B.
Lehoucq,
D.
C.
Sorensen
http://www.caam.rice.edu/software/ARPACK.

and

C.

Yang,

[166] J. Macek and E. Nielsen, (private communication).
[167] C. H. Greene, (notes).
[168] E. A. Burt, R. W. Ghrist, C. J. Myatt, M. J. Holland, E. A. Cornell and C. E.
Wieman, Phys. Rev. Lett. 79 (1997) 337.
[169] B. D. Esry, C. H. Greene, Y. Zhou and C. D. Lin, J. Phys. B 29 (1996) L51.
[170] B. D. Esry, C. D. Lin and C. H. Greene, Phys. Rev. A 54 (1996) 384.

166
[171] H. Nakamura, Int. Rev. Phys. Chem. 10 (1991) 123.
[172] J. Stenger, S. Inouye, M. R. Andrews, H. J. Miesner, D. M. Stampur-Kurnand W. Ketterle, Phys. Rev. Lett. 82 (1999) 2422.
[173] J. M. Gerton, C. A. Sackett, B. J. Frew and R. G. Hulet, Phys. Rev. A 59
(1999) 1514.
[174] J. Söding, D. Guéry, P. Debiolles, F. Chevy, H. Inamori and J. Dalibard,(preprint cond-mat/9811339).
[175] E. Tiesinga, A. J. Koerdijk, B. J. Verhaar and H. T. C. Stoof, Phys. Rev. A 46
(1992) R1167.
[176] T. P. Dinneen, K. R. Vogel, E. Arimondo, J. L. Hall and A. Galllagher, Phys.
Rev. A 59 (1999) 1216.
[177] J. D. Wienstein, R. deCarvalho, J. Kim, D. Patterson, B. Friedrich and J. M.
Doyle, Nature (London) 395 (1998) 148.
[178] M. Walhout, U. Sterr, C. Orzel, M. Hoogerland and S. L. Rolston, Phys. Rev.
Lett. 74 (1995) 506.
[179] M. Holland, (private communication).

167

Appendix A
FEM basis functions

The six Hermite interpolating polynomials used as a local FEM basis are given by

u1(x) = x2 ; 1:25x3 ; 0:5x4 + 0:75x5
u2(x) = 0:25(x2 ; x3 ; x4 + x5 )
u3(x) = 1:0 ; 2:0x2 + x4
u4(x) = x ; 2:0x3 + x5
u5(x) = x2 + 1:25x3 ; 0:5x4 ; 0:75x5
u6(x) = 0:25(;x2 ; x3 + x4 + x5 ) :(A.1)

Appendix B
Adiabatic FEM R-matrix

The equations for the adiabatic FEM R-matrix method can be derived by inserting
the adiabatic channel decomposition (Eq. 3.20) into Eq. 3.1. Using the relationships, H = TR + H ad , H ad  ( ;  ; R) = U  ( ;  ; R), Eq. 3.22 for the derivative coupling, and the orthonormality of the  ( ;  ; R)’s leads to the following
equation for b

b

X

XZ

M (R0)M (R0 ) =

;
;

X
X
Z;

2M (E ; U )M dR

M (R0 )M 0 (R0 )M (R0 )P ; (R0 )M (R0)000; 2 TR dw(B.1)

where R0 is the internuclear separation on the boundary. In order to evaluate the last
integral, I need to introduce the second derivative coupling matrix Q[78]@2 
Q =  j @R
2
0



(B.2)

0

which can be related to P using[78]Q = (P 2) + @P =@R :
0

0

0

XZ

@

The last integral in Eq. B.1 becomes

;

Z

2 T

R

dw =

;

0

M @R + P

(B.3)

2
0

M dR
0

(B.4)

The original diabatic equations are recovered exactly provided @=@R is replaced

169
with @=@R + P . The matrix elements of  remain the same, however the new ;matrix elements are given in the FEM representation by

Z1;ij = 2 ui(xn)(E ; H )uj (xn)andxn;1; ; n;nmax k;5k ;6=an;n;nmax k;5k ;5P ; 0 :
0

0

(B.5)

0

Again, the indices i and j include the basis function index k , the channel index, and the sector index n, (i  fk ng and j  fk 0 0 ng). This equation can be
improved upon by first squaring the kinetic energy operator TR and using Eq. B.3

2(TR ) ; =  @ 2 + 2P ; @ + @ P ;
@R
@R
@R
0

0

0



0

+ (P 2) ; :
0

(B.6)

@ + @ P term by parts
Plugging in the FEM basis expansion and integrating the P @R
@R
results in the following expression for the ; matrix elements

;ij =

 @2

2
uk (xn) @ (a x )2 + (P ) + 2(E ; U ) uk (xn )andxn;1
Z 1  @un n @u Z10

P uk @xk ; @xk uk dxnnn;1; ; n;nmax k;5k ;6=an :

+

0

0

0

0

0

0

0

(B.7)

Numerical derivatives of the P matrix are no longer required in this formulation
which generally improves the accuracy of the calculations. From this point on, theprocedure to solve the R-matrix equations is exactly as outlined in chapter 3.1.

Appendix C
Thermal averaging

The total number of two-body collisions N between atoms of type a andb per unit trap volume V per unit time t is given by na nb (vr)vr, where ni is thedensity of i atoms, vr the relative velocity of atom a with respect to atom b, and  isthe cross section. To get an average expression for the collision rate as a function of
the atom cloud temperature T , we must consider the distribution f (~xa ; ~xb ; ~pa ; p~b ) of
atoms in a trap volume element dV with momenta (measured in the lab frame) p~aand ~pb . Further, we require the distribution function f to be unit normalized, i.e.

Zf (~xa ; ~xb; p~a ; ~pb)d3xa d3xb d3pa d3pb = 1 :(C.1)

The thermally averaged value is then simplyZ

hna nb(vr)vriT = na (~xa )nb (~xb)(vr)vrf (~xa; ~xb ; ~pa; p~b)d3xa d3 xbd3pa d3pb
(C.2)
If the thermal de Broglie wavelengths of the atoms are small compared to the average spacing between atoms, then the gas can be treated with classical MaxwellBoltzmann statistics. This condition can be stated numerically as 1 1=3nh p3mkT
B(C.3)

where m is the atomic mass, kB is the Boltzmann constant, and T the cloud temperature. A 1 K cloud of Rb atoms with density of the order 108 atoms/cm3 easily falls
into the classical regime. Assuming the cloud is in thermal equilibrium (or at least
close to it) the canonical distribution function is given by the f = C exp(;H=kB T ),where H is the Hamiltonian of the system and C is a normalization constant. Moreover, in this classical regime we can neglect particle interactions. Accordingly, the171
Hamiltonian for two particles in an isotropic harmonic trap is given by

1 p~ 2 + 1 ~p 2 + 1 m !2~x 2 + 1 m !2~x 2
H = 2ma
2m b 2 a a 2 b b
a

b

(C.4)

where ! represents the trap frequency, and all coordinates and momenta are definedin the lab frame. (Trap anisotropies only effect the spatial density distribution of
atoms, which is an experimentally measured quantity.) Next, H must be converted
into the center-of-mass frame in the usual way by introducing center-of-mass and~ cm = ma~xa + mb~xb =(ma + mb ), and ~xr = ~xa ; ~xb ,relative coordinates: where, X
M = ma + mb ,  = ma mb =(ma + mb). The transformed Hamiltonian is given by
the expression:

1 P~ 2 + 1 p~2 + 1 M!2X~ 2 + 1 !2~x2 :H = 2Mcm 2 r 2
cm 2
r

(C.5)

Equation C.2 can now be evaluated by inserting the center-of-mass Hamiltonian (Eq.
C.5) into the distribution function, resulting in the following expression

R n (~x )n (~x ) e;V (Xcm;xr)=kB T d3X d3xa Ra b b
cm r
hna nb (vr)vri =
;
V
(
X
+
x
)
=k
T
3
3
cm
r
B
e
d Xcmd xr
R (v ) v e;T (pr)=kB T d3p
R re;Tr(pr)=kB T d3p r :
r

(C.6)

Here, the Hamiltonian is written as sum of three terms, H = T (Pcm ) + T (pr ) +V (Xcm; xr). The center-of-mass momenta drops out of Eq. C.6 since the densities,
cross section, and relative velocity have no explicit P~cm dependence. Eq. C.6 can beseparated into a term that depends only on the spatial coordinates and another termthat depends only on the relative momenta. The spatially dependent term is effectively the total number of atoms N times the density-weighted-density hDWD i,

R n (~x )n (~x ) e;V (Xcm;xr)=kB T d3X d3x
N hDWDi = a Ra e;b V (Xb cm ;xr)=kB T d3X d3x cm r
cm

r

(C.7)

which is an experimentally determined quantity. The remaining term is the thermal
averaged rate constant K2

R (v ) v e;T (pr)=kB T d3p
hK2 i = R re;Tr(pr)=kB T d3p r :
r

(C.8)

This expression can be simplified by evaluating the denominator and converting the
integral over momenta to an integral over energy, E = p2r =2. The final expression

172
for the thermal averaged rate constant is given byr2

hK2i = 2 

(kB T )3=2

Z10

(E )E e;E=kB T dE :(C.9)

Cold collision experimental observables are often the rate for some process to occur.
Therefore, the standard definition for the thermal averaged cross section is taken to
be

Z1
h
K
2i
;
2
hi = hv i = (kB T )(E )E e;E=kB T dE :r
0

(C.10)

Finally, in cross-dimensional mixing experiments[108] such as the one used to observe the Feshbach resonance discussed in Section 5.2, the experimental observable
is actually the rate of energy transfer[179]. In this case, the relevant thermal average
is hEK2 i, i.e. the integral in equation C.9 must include an additional factor of E .

Appendix D
Related publications

“Recombination of three atoms in the ultracold limit”, B. D. Esry, C. H. Greene,
and J. P. Burke, Jr., (submitted to Phys. Rev. Lett.).
“Determination of 39 K scattering lengths using photoassociation spectroscopy of
the 0;
g state”, J. P. Burke, Jr., C. H. Greene, J. L. Bohn, H. Wang, P. L. Gould, andW. C. Stwalley, (submitted to Phys. Rev. A).
“Comparison of numerical methods for the calculation of cold atom collisions”, G.
H. Rawitscher, B. D. Esry, E. Tiesinga, J. P. Burke, Jr., and I. Koltracht, (submitted
to J. Chem. Phys.).
“Observation of p-wave threshold law using evaporatively cooled Fermionic atoms”,
B. DeMarco, J. L. Bohn, J. P. Burke, Jr., M. Holland, and D. S. Jin, Phys. Rev. Lett.
82, 4208 (1999).
“Ultracold scattering properties of the short-lived Rb isotopes”, J. P. Burke, Jr. and
J. L. Bohn, Phys. Rev. A 59, 1303 (1999).
“Collisional properties of ultracold potassium: consequences for degenerate Boseand Fermi gases”, J. L. Bohn, J. P. Burke, Jr., C. H. Greene, H. Wang, P. L. Gould,
and W. C. Stwalley, Phys. Rev. A 59, 3660 (1999).
“Resonant magnetic field control of elastic scattering in cold 85 Rb”, J. L. Roberts,
N. R. Claussen, J. P. Burke, Jr., C. H. Greene, E. A. Cornell, and C. E. Wieman,
Phys. Rev. Lett. 81, 5109 (1998).
“Multichannel cold collisions: simple dependences on energy and magnetic field”,
J. P. Burke, Jr., C. H. Greene, and J. L. Bohn, Phys. Rev. Lett. 81, 3355 (1998).
“Prospects for mixed-isotope Bose-Einstein condensates in rubidium”, J. P. Burke,Jr., J. L. Bohn, B. D. Esry, and C. H. Greene, Phys. Rev. Lett. 80 2097 (1998).
“Hartree-Fock approximation for double condensates”, B.D. Esry, C.H. Greene,
J.P. Burke, Jr., and J.L. Bohn, Phys. Rev. Lett. 78, 3594 (1997).
“Impact of the 87 Rb singlet scattering length on suppressing inelastic collisions”,
J.P. Burke, Jr., J.L. Bohn, B.D. Esry, and C.H. Greene, Phys. Rev. A 55, R2511(1997).