BoseEinstein Condensation from a Rotating Thermal Cloud:
Vortex Nucleation and Lattice Formation
Abstract
We develop a stochastic GrossPitaveskii theory suitable for the study of BoseEinstein condensation in a rotating dilute Bose gas. The theory is used to model the dynamical and equilibrium properties of a rapidly rotating Bose gas quenched through the critical point for condensation, as in the experiment of Haljan et al. [Phys. Rev. Lett., 87, 21043 (2001)]. In contrast to stirring a vortexfree condensate, where topological constraints require that vortices enter from the edge of the condensate, we find that phase defects in the initial noncondensed cloud are trapped en masse in the emerging condensate. Bosestimulated condensate growth proceeds into a disordered vortex configuration. At sufficiently low temperature the vortices then order into a regular Abrikosov lattice in thermal equilibrium with the rotating cloud. We calculate the effect of thermal fluctuations on vortex ordering in the final gas at different temperatures, and find that the BEC transition is accompanied by lattice melting associated with diminishing long range correlations between vortices across the system.
pacs:
03.75.Kk,03.75.Lm,05.70.Fh,05.10.Gg,I Introduction
Rotating BoseEinstein condensation was first realized by Haljan et al.at JILA Haljan2001 (), who grew large vortex lattices by evaporatively cooling a rotating thermal cloud. The evaporation is performed in a prolate trap, along the axis of rotation so that escaping atoms make large axial excursions but do not stray far from the rotation axis. Escaping atoms carry away significant energy but little angular momentum, increasing the angular momentum per trapped atom. Using this cooling mechanism, which requires typically s of evaporative cooling and an extremely cylindrical trap, very large vortex lattices containing more than 300 vortices have been created Engels2004 (). This nucleation mechanism is distinct from more studied vortex formation mechanisms Matthews1999 (); Leonhardt2002 (); Madison2000 (); Hodby2002 (); AboShaeer2001 (); Inouye2001 (); Anderson2001a (); Dutton2001a (); Scherer2007 () because, rather than manipulating a condensate to generate vortices, vortices can be formed at the onset of condensation as the rotating gas reaches the critical temperature () for BoseEinstein condensation. This mechanism thus has connections with the KibbleZurek Kibble1976 (); Zurek1985 () and BerezinskiiKosterlitzThouless Berezinskii1971 (); Kosterlitz1973 (); Hadzibabic2006a () mechanisms, as it is fundamentally caused by spontaneous symmetry breaking and thermal fluctuations near a phase transition.
In this work we study the dynamics of vortex lattice formation during rotating evaporative cooling using a stochastic GrossPitaevskii equation (SGPE) Stoof1999 (); SGPEI (); SGPEII (). Our central aim is to include rotation, thermal fluctuations, and twobody collisions, by separating the system into a condensate band coupled to a (possibly rotating) thermal reservoir. The condensate band is described by a generalized mean field theory known as the truncated Wigner method Graham1973 (); Steel1998 (); Sinatra2001 () that includes the projected classical field method Davis2001b (); Davis2006a (); Simula2006a () as a special case in the regime of high temperatures, while the noncondensate band formalism is based on quantum kinetic theory QKV ().
We also provide expressions for the dissipation rates of the theory under the condition that the highenergy component is approximately BoseEinstein distributed, and details of our numerical algorithm that is central to the utility of the SGPE approach. These are the first SGPE simulations of vortex formation in a trapped Bose gas.
Ii System and equations of motion
The trapped system is divided into a condensate band of states with energy beneath the cutoff , and the remaining noncondensate band of highenergy states. The formal separation of the Hamiltonian into bands according to energy allows different methods to applied in the two bands. The condensate band is treated using the truncated Wigner phase space method Graham1973 (); Steel1998 () which is valid when the modes are significantly occupied Blakie2007a (). Within this approach the condensate band includes the condensate and the most highly occupied excitations. The noncondensate band is comprised of many weakly occupied modes and is treated as approximately thermalised; the effect of this band is to generate dissipative terms in the equation of motion for the condensate band. The highenergy states play the role of a grand canonical thermal reservoir for the low energy condensate band. In this section we provide a minimal overview of the formalism and simply state the final equations of motion to be used in this work. Refs. SGPEI (); SGPEII () contain the essential formalism, and Ref. SGPEIII () gives necessary modifications for the rotation system, along with a number of corrections to the derivation of the equations of motion.
ii.1 Hamiltonian and separation of the rotating system
The coldcollision regime for bosonic atoms of mass is characterized by the interaction parameter , where is the swave scattering length Dalfovo1999 (). The secondquantized manybody Hamiltonian for a dilute Bose gas confined by a trapping potential in the rotating frame defined by angular frequency is , where
(1) 
and
(2) 
The orthogonal projectors that separate the system are defined with respect to the single particle basis of the trap. This provides a convenient basis because (a) The many body spectrum becomes well approximated at highenergy by a single particle spectrum, allowing the cutoff to be imposed in this basis Davis2001b (), and (b) the numerical method used to solve the final equations of motion uses the single particle basis Blakie2005a ().
The condensate band projector takes the form
(3) 
where the bar denotes the highenergy cutoff and represents all quantum numbers required to index the eigenstates of the single particle Hamiltonian. To minimize notation we will use both for the operator and for its action on wave functions; the meaning will always be clear from the context. The field operator is decomposed as
(4) 
where the noncondensate field describes the highenergy thermal modes, and the condensate band is described by the field . In the spatial representation the eigensolutions of , denoted , generate the action of the projector on an arbitrary wave function
(5) 
The condensate band field theory generated by such a decomposition has a nonlocal commutator
(6) 
where . In the harmonic oscillator representation of a one dimensional system, the approximate delta function has a position dependent width, being narrowest in the center of the trap. This reflects the fact that trap wavefunctions become smoother near their semiclassical turning points. Note that since , ; put another way, this means that is a true Diracdelta function for any condensate band field.
There are two physical processes arising from the interactions between condensate and noncondensate bands, referred to as growth and scattering SGPEI (); SGPEII (); SGPEIII (), and shown schematically in Figure 2. The growth terms correspond to twobody collisions between condensate and noncondensate atoms whereby particle transfer can take place. The scattering describes interband collisions that conserve the populations in each band but involve a transfer of energy and momentum between bands.
The noncondensate band is described semiclassically in terms of the Wigner function
(7) 
In this work the trapping potential assumes the cylindrical form
(8) 
where and are the radial and axial trapping frequencies respectively. The noncondensate band Wigner function for a Bose gas confined in a cylindrically symmetric trap in equilibrium in the rotating frame at frequency , with chemical potential , is
(9) 
where
(10) 
defines the effective potential as , and is the distance from the symmetry axis of the trap. Since is the rotating frame momentum, we can work in terms of this variable, defining
(11) 
Hereafter we drop the subscript for brevity since we always describe the system in this frame. In the semiclassical approximation is only significant in the noncondensate region , which is the region of phase space where .
ii.2 Stochastic GrossPitaevskii equation
Following the theory in Refs. SGPEI (); SGPEII (); SGPEIII (), the equation of motion for the density operator describing the total system
(12) 
is reduced to a master equation for the condensate band density operator, which is the trace of over the noncondensate band degrees of freedom . The master equation is then mapped to a FokkerPlanck equation for the Wigner distribution describing the condensate band. In the truncated Wigner approximation, an equivalent stochastic equation of motion can be obtained that allows the dynamics to be efficiently simulated numerically QN (). Formally, there is a correspondence between symmetrically ordered moments of the field operator and ensemble averages over many trajectories of the projected field that evolves according to the SGPE to be defined below SGPEIII (). For example, for the density
(13) 
where denotes stochastic averaging over trajectories of the SGPE.
To state the SGPE we require some expressions from mean field theory, in particular the GPE obeyed by in the absence of any thermal atoms in the noncondensate band. The mean field Hamiltonian for corresponding to is
(14) 
with atom number given by .
The growth terms involve the GrossPitaevskii equation found in the usual way as
(15) 
that defines the evolution operator . To obtain this equation of motion we have used the projected functional calculus SGPEII () that plays the same role here as ordinary functional calculus in nonprojected field theory.
The scattering process couples to the divergence of the condensate band current
(16) 
where the second line is a rigid body rotation term arising from the transformation to the rotating frame. The limit gives the laboratory frame velocity field , corresponding to the irrotational system mimicking rigid body rotation. Persistent currents around the vortex cores are responsible for maintaining irrotationality.
ii.2.1 Full nonlocal SGPE
We finally arrive at the Stratonovich stochastic differential equation (SDE) SGPEIII ()
(17) 
where the growth noise is complex, and the scattering noise is real. The noises are Gaussian, independent of each other, and have the nonzero correlations
(18)  
(19) 
The growth and scattering amplitudes, and , are calculated in Appendix A for the case where the noncondensate band is described by a thermal BoseEinstein distribution.
The two validity conditions for this equation are (i) high temperature: where is the characteristic system frequency, and (ii) significant occupation of modes in the condensate band, necessitated by the truncated Wigner method. The first line of Eq. 17 is the PGPE developed by Davis et al Davis2001b () that, in essence, is the GPE formally confined to a low energy subspace of the condensate band. The condensate growth terms of the second line are closely connected to phenomenological dissipative theories of BECs based on the GPE Choi1998 (); Penckwitt2002 (); Tsubota2002 () but here have a specific form given in Appendix A. The remaining terms describe scattering processes that transfer energy and momentum between the two bands; these terms where first derived in Ref. SGPEII (), and have not appeared in previous SGPE theories.
Irrespective of the form of and , the equilibrium state for the Wigner distribution of the system has the grand canonical form SGPEII ()
(20) 
that is generated by the stochastic evolution of Eq. (17). Solutions of the SGPE also known to satisfy a set of generalized Ehrenfest relations which in the steady state reduce to a statement of the fluctuationdissipation theorem for the grand canonical Bose gas Bradley2005b ().
ii.2.2 Simple growth SGPE
A simpler equation that reaches the same equilibrium state can be obtained by noting that the scattering terms in Eq. (17) do not directly change the condensate band population. Such terms have been included in a kinetic theory of BoseEinstein condensation QKIV (), where they were found to play a role in the early stages of condensate growth by causing atoms in many modestly occupied states to merge more rapidly into a single highly occupied state which then dominates the growth process. Within the SGPE theory, the scattering terms generate an effective potential from density fluctuations. This can be seen by noting that upon neglecting the growth and scattering terms, that are relatively small compared to the mean field evolution, the continuity equation for the condensate band gives . Consequently, since is purely real, the deterministic part of the scattering terms generate an effective potential that will tend to smooth out density fluctuations.
Upon neglecting these terms the resulting stochastic equation is much easier to integrate numerically as it does not require the nonlocal noise correlations of the full equations. The simulations in this work are performed by numerically integrating the simple growth SDE
(21) 
using the mixed quadrature pseudospectral method detailed in Appendix B. The rate is given by Eq. (44), and the noise is complex Gaussian with nonvanishing correlator
(22) 
Iii Simulations of Rotating BoseEinstein condensation
The central question of interest for this work regards the formation of vortices during condensation. In particular, does a vortexfree condensate form and then vortices enter it from the edge, or does condensation proceed into a state with vortices already present? That this is nontrivial to answer is easily seen from the single particle energy spectrum of the cylindrically symmetric harmonic trap in three dimensions
(23) 
where are the radial, angular and axial quantum numbers. In the absence of interactions BEC will form in the ground state mode with energy . Vortices arise from occupation of states with nonzero angular momentum and the positive angular momentum part of the spectrum behaves as leading to neardegeneracy of positive angular momentum modes as . States with angular momentum are increasingly easy to excite as the rate of rotation increases and in the interacting gas there is a tradeoff between the interaction and rotational contributions to the energy due to the presence of a vortex that is known to favor the creation of vortices in the condensate once the system is rotating faster than the critical angular frequency, , for energetic stability of a vortex Fetter2001 (). For a rapidly rotating BEC transition (), we can expect that vortices will play a dominant role at all stages in the BEC growth process. It is then possible that atoms condense into a vortexfilled state. In general the vortex configuration can be disordered. It is not required, for example, that the vortices minimize the system energy by forming a regular Abrikosov lattice at all stages of condensate growth.
iii.1 Properties of the rotating ideal Bose gas and the noncondensate band
An important effect of rotation in BoseEinstein condensation is the reduction of Stringari1999 (); Coddington2004a (), that, along with the slow onedimensional evaporative cooling, accounts for the long evaporation time in the JILA experiments. The transition temperature for the ideal Bose gas is modified by the effective radial frequency in the rotating frame, , so that
(24) 
where is the nonrotating transition temperature for trapped atoms.
In this work we model the noncondensate region assuming the highenergy atoms are quickly thermalised into a BoseEinstein distribution, that we approximate as noninteracting since the density is usually small at high energies. We can find analytical expressions for the properties of the noninteracting gas, including the effect of the cutoff, by defining an incomplete BoseEinstein function
(25) 
where is the incomplete Gamma function. In analogy with the reduction to an ordinary Gamma function , we have , reducing to the ordinary BoseEinstein function. We then find for the noncondensate band density, for example
(26) 
where . Changing integration variable to gives
(27) 
where is the thermal de Broglie wave length. Setting , we recover the standard form for the semiclassical particle density of the ideal gas Dalfovo1999 (). The total number in the noncondensate band is given by
(28) 
where is the geometric mean frequency in the rotating fame. In this way the usual semiclassical expressions can be generalized to include a cutoff in terms of the incomplete BoseEinstein function.
The angular momentum and its fluctuations can also be calculated. We now write , and the variance per particle as . In rotating frame coordinates the angular momentum is , which in equilibrium gives the rigid body relation between angular momentum and rotational inertia . From this expression we find the angular momentum is
(29) 
where . Similarly, the mean squared angular momentum is
(30) 
where . In Figure 3 we give an example of the the variation of , and the angular momentum per particle together with its fluctuations as a function of the angular frequency of rotation. The suppression of by the rotation is seen to be associated with diverging angular momentum and angular momentum fluctuations near the critical frequency for stable rotation in the harmonic trap.
iii.2 Numerical procedure
The long cooling time and the large changes in length scales in the axial and radial dimensions during cooling makes modeling the full condensation dynamics challenging. However, as the gas cools and spins up, the system becomes progressively more twodimensional. We model the system by considering the two dimensional projected subspace of the quasitwo dimensional system. We thus use the SGPE of Eq. 21 with rescaled interaction parameter , where in this work we take the thickness through the oblate system as , obtained by projecting onto the ground state of the dimension. Note that a consequence of the approximation is that the growth rate is unchanged by reduction to the two dimensional effective equation.
We take as our starting point the gas after it has been evaporatively cooled to a rapidly rotating state above the transition temperature. We then simulate the dynamics of a rapid quench below threshold. The procedure is summarized as follows:

Initial states are generated by sampling the grand canonical ensemble for the ideal Bose gas for initial parameters , and .

To approximate a rapid quench, the temperature and chemical potential of the noncondensate band are altered to the new values and , while the rotation of the cloud is held fixed at , consistent with axial evaporation. The condensate band field is then evolved according to the SGPE.
Specifically, we aim to model a possible run of the JILA experiment consisting of a system of atoms held in a magnetic trap with frequencies , initially at temperature with chemical potential , and rotating with angular frequency . We examine a range of final quench temperatures with constant chemical potential . We note the ratio , where a value approaching unity indicates the onset of the mean field quantum Hall regime. The ratio , indicates the system is oblate but not highly twodimensional. We have found that simulations for parameters that are highly oblate are very challenging numerically, so we consider as our first application of the method these more modest parameters. For the results presented here we use a cutoff value of , giving a condensate band containing single particle states. To check the cutoff independence of our results we have also carried out simulations for (not shown here). We find qualitative agreement for single trajectory dynamics and quantitative agreement (at the level of a few percent) for the ensemble averaged condensate fraction shown in Figure 10. We are thus confident that our results are cutoff independent. For computational reasons we choose constant dimensionless damping for all simulations. Although we have derived a form for the damping rate given in Appendix A, the values of arising from Eq. (44) for our parameters are orders of magnitude smaller then our chosen value, so we have chosen a fixed damping such that the growth rate is comparable to the JILA experimental time of 50s. This affects the specific timescale for growth, but not the qualitative features of the simulations which merely require that the damping is much slower than any other timescale of the system. We emphasize that the equilibrium properties are independent of .
Since our numerical simulation method employs a projected eigenbasis of the associated linear Schrödinger problem, the initial BoseEinstein distribution is easily sampled. The grand canonical Wigner distribution for a noninteracting trapped Bose gas distributed over a set of single particle states with frequencies is
(31) 
where , the modes are defined in Eq. 55, and
(32) 
The energy cutoff is denoted by the upper limit of the product notation. The single particle spectrum of the condensate band is restricted by the cutoff to
(33) 
The distribution is sampled as , where are complex Gaussian random variates with moments , and . We simulate the time dynamics of the quenched system using the numerical algorithm described in Appendix B.
We note here that application of the PenroseOnsager criterion for obtaining the BEC Penrose1956 (); Blakie2005a () is hindered in this system because each run of the SGPE theory leads to vortices in different locations, and ensemble averaging the trajectories generates spurious results for . Specifically, the condensate number is significantly reduced below the condensate number that we find using a different technique of shorttime averaging single trajectories. Generalizing the PenroseOnsager criterion for BEC to multimode mixed states is not pursued in this work. Here we have taken the view that a given trajectory samples a subensemble of the full range of quantum evolution, consistent with spontaneously broken symmetry in the initial conditions. We present single trajectory results of our simulations of the SGPE that give a qualitative indication of what could be seen in individual experimental measurements of the atomic density time series.
iii.3 Initial states and dynamics
In Figure 4 we show the density of the condensate band for the thermal gas above the BEC transition. Individual samples of the Wigner function are seen to contain phase singularities. We also construct and diagonalize the one body density matrix, according the the PenroseOnsager criterion for BEC. By averaging samples from the Wigner function we find that the largest eigenvalue is and the associated eigenvector contains several phase singularities. Using samples washes out the singularities and gives the noninteracting ground state; the point here is that the high degree of degeneracy means that we must average over very many samples to remove the singularities. Put another way, this means that individual trajectories always contain many vortices, and the motion of these vortices destroys long range coherence associated with BEC.
In Figure 5 we plot the condensate number as a function of time after the quench for a range of different quench temperatures. The condensate is calculated by constructing the one body density matrix for the subensemble via short time averaging individual trajectories Blakie2005a (). We have found that by averaging over samples taken uniformly over an interval of radial trap periods, we obtain sufficiently good statistics to determine the condensate number , given by the largest eigenvalue of . The main point to note from this result is that at higher temperatures the absolute condensate number is slightly suppressed. However, as we will see in the next section, the condensate fraction is strongly suppressed, due to the increasing atom number in the noncondensate band for higher temperatures.
A typical time series is shown in Fig. 6, for the lowest temperature considered in Fig. 5, . The initial sample contains many phase singularities, including a few that are circulating counter to the cloud rotation. Countercirculating vortices only tend to arise in the initial state, are thermodynamically very energetic, and are quickly annihilated by positive vortices during postquench dynamics. The condensate, shown in Fig. 7, forms with many vortices already present, but in a highly disordered configuration. The vortices then assemble into a regular triangular lattice.
At low temperatures the equilibrium state of the system consists of a an Abrikosov vortex lattice with small thermal fluctuations of the vortices about their average positions. The time dynamics for the highest temperature of Fig. 5, , are shown in Fig. 8. We see that growth proceeds into a disordered state for which periodic order remains frustrated by significant thermal fluctuations.
iii.4 Equilibrium states
We now consider equilibrium properties of the rotating system, in the grand canonical ensemble generated by evolution according to the SGPE. In Fig. 10 we plot the condensate fraction of the final state for the temperatures used in Fig. 5. The total atom number in the system is calculated as , where is given by Eq. (28) evaluated for our chosen energy cutoff, and is the mean total number of atoms in the condensate band.
We note that the ideal gas relationship between condensate fraction and temperature, , is unaltered by rotation, provided the rotating frame transition temperature, given by Eq. 24, is used in place of the nonrotating transition temperature. The rotating ideal gas curve is plotted for comparison, and we see that the condensate fraction predicted by the SGPE simulations is suppressed relative to the ideal gas result. Although we are unable to simulate the nonrotating system with our twodimensional algorithm, since the ideal gas behavior is universal with respect to rotation we also compare with the experimental data for the nonrotating Bose gas Ensher1996a (). Our simulations show a condensate fraction significantly suppressed relative to the experimental data for the nonrotating system.
This difference presumably comes about from two effects. The most obvious is the loss of long range order caused by thermal motion of vortices that are absent from the nonrotating system. A less obvious and possibly less significant effect is that the reduction of condensate number by interactions is more pronounced for a rotating gas, due to the exclusion of atoms from vortex cores. Our highest temperature simulation, at , is right on the edge of BoseEinstein condensation, and this is reflected in the lack of long range order seen in Fig. 8.
The loss of long range order can be seen more clearly by looking at the relative positions of vortices in the condensate wavefunction at different temperatures Engels2003a (). In Figure 11 we plot the condensate wave function, along with the histogram of the separation of each pair of vortices in the system. We only consider vortices with radius . For the lowest temperature () the histogram is indistinguishable from the ground state histogram (not shown). The presence of many peaks for large separations reflects the long range spatial order in the system. As temperature increases, order is reduced starting at large vortex separation. At the highest temperature () we see that only two histogram peaks clearly persist, corresponding to nearest and nextto nearest neighbor order being preserved.
Iv Conclusions
We have developed and applied the SGPE theory to the rotating Bose gas. The formalism presented here can be equally applied to the nonrotating case and in that limit extends the prior SGPE theory SGPEI (); SGPEII (); SGPEIII () by providing expressions for the dissipation rates. We have developed an efficient algorithm for numerically integrating the SGPE in a rotating reference frame and applied it to the case of quenching a rotating Bose gas below threshold, modeling the rapidquench limit of the experiment at JILA Haljan2001 ().
Our simulations indicate that many vortices corotating with the thermal cloud are trapped in the emerging condensate, reflecting the large angular momentum of the system. In our calculations these primordial vortices are nucleated from vortices in the individual samples of the initial Wigner distribution and from fluctuations arising from the condensate growth terms in the SGPE. The process can be regarded as stimulated vortex production, in contrast with the nonrotating case where vortices can appear spontaneously but the total vorticity in the BEC must be zero on the average Drummond1999 ().
For low temperatures the vortices order into a regular Abrikosov lattice as the system approaches thermal equilibrium with the quenched thermal cloud. At sufficiently high temperatures lattice crystallization is frustrated by thermal fluctuations of the positions of vortex cores. This disorder is reflected in a loss of long range correlations in the relative positions of vortices that are most evident as the system approaches .
We note that the concept of rotation has not been investigated in the context of the KibbleZurek mechanism, which explains vortex nucleation during the phase transition in terms of spontaneous symmetry breaking and critical slowing down Anglin1999 (). Future quantitative studies of dynamics at the transition will allow for a deeper understanding of the condensation process, in particular the role of spontaneous symmetry breaking as a vortex nucleation mechanism in rotating BoseEinstein condensation.
We gratefully acknowledge B. P. Anderson and M. K. Olsen for their critical reading of this manuscript, and P. Jain, P. B. Blakie and M. D. Lee for helpful discussions about this work.
This research was supported by the Australian Research Council Center of Excellence for Quantum Atom Optics, by the New Zealand Foundation for Research Science and Technology under the contract NERFUOOX0703: Quantum Technologies, and by the Marsden Fund under contracts PVT202 and UOO509.
Appendix A Calculation of dissipation rates
Reservoir interaction rates closely connected to those required here have been previously computed within the context of the quantum kinetic theory (QKT) of BEC QKV (). We now derive the form of the growth and scattering rates including the effects of the cutoff, the external trapping potential and the rotation of the thermal cloud.
a.1 Growth
a.1.1 General form
From Eq. (84) of Ref. SGPEII (), the growth amplitude is
(34) 
where
(35) 
(36) 
and conserves energy and momentum. Note that effect of working in the rotating frame is only to introduce the effective potential term into the cutoff and Wigner functions in , which is otherwise unaltered. We will take the reservoir to be in equilibrium with BoseEinstein distribution at chemical potential and temperature . Integration over the momentum function reduces to
(37) 
where Eq. (11) gives the expression for the momentum cutoff
(38) 
Changing to variables , , and taking account of the deltafunction, we obtain
(39) 
where , and . We now show that the rate is composed of two regions: an inner region where the rate is spatially invariant, and an outer limit where there is a weak spatial dependence.
a.1.2 Inner region
From the from of the lower limits we find that wherever , . Proceeding similarly for we find the same condition to hold, so that in the inner region the growth rates are spatially invariant, taking the final form
(40)  
(41)  
where the Lerch transcendent is defined as . These terms are superficially similar to the condensate growth rates previously obtained using the quantum kinetic theory QKPRLII (), although the form of the cutoff dependence is different in the SGPE theory. Note that this form is independent of any collisional contribution to the effective potential.
a.1.3 Outer region
In the region where the rates depend on position, albeit quite weakly. At the edge of the condensate band where , we find
(42)  
(43)  
where is the modified Bessel function of the second kind.
In Fig. 12 we compute the variation of with radius for a spherical parabolic potential. For illustration we have also neglected the collisional contribution to the effective potential. To show the maximum variation of the total rate we also plot (inset) the ratio versus . For all data we use which is typically a good approximation for the energy at which the spectrum of the many body system becomes single particlelike SGPEII (). The rates are spatially invariant except for a small region near the edge of the condensate band where the scattering phase space in the harmonic trap is modified by the energy cutoff. Since the variation with radius is usually not significant, for our simulations we will use the approximate growth rate
(44) 
which is exact in the inner region and a good approximation elsewhere.
a.2 Scattering
Our starting point of the scattering amplitude calculation, that is related to the actual scattering via Fourier transform, is Eq. (86) of Ref. SGPEII ()
(45) 
where . Upon evaluating the momentum deltafunction we find
(46) 
Using the notation , where
(47) 
in the variable we have
(48)  
where . Writing to distinguish the term with one or two BEdistributions as before, and choosing the component of along , we find, for the term involving one BEdistribution
(49) 
The theta function generates the effective cutoff . Evaluating the integral gives
(50) 
which can be further simplified as follows. We want the condition for , which can be expressed as . The restriction of is determined by the requirement that , where