# The Stochastic Gross-Pitaevskii Equation and some Applications

###### Abstract

The stochastic Gross-Pitaevskii equation represents a versatile approach for studying the dynamics of trapped degenerate ultracold Bose gases in the presence of large phase and density fluctuations. Following a brief review of the original formulation of Stoof, which highlights the benefits of this approach and its relation to alternative theories, we present selected applications for the dynamics of effectively one-dimensional systems, and briefly discuss the generalization to two-dimensional systems, highlighting certain potential pitfalls in their numerical implementations.

## I Introduction

The ultracold atomic samples now routinely produced in many laboratories worldwide, provide excellent means by which to investigate many-body quantum systems on a macroscopic scale. From a theoretical perspective, these systems are appealing as they support an array of rich dynamics, in many cases readily amenable to experimental observation. With this close relation between theory and experiment Stringari_Review ; General_Review_2 , many such experimental features have been subject to analytical and numerical investigations: much attention has been paid to numerical analyses based upon the Gross-Pitaevskii equation (GPE) in particular Frantzeskakis_Book , which has been shown to successfully describe condensate dynamics in the limit when the many-body wavefunction represents the condensate alone, thus corresponding to the low-temperature regime.

Most experiments on trapped ultracold gases feature, to a greater or lesser extent, a depletion of the condensate due to the effect of finite temperature. While in three dimensional systems, finite temperature effects play an increasingly important role as the critical temperature is approached, for low dimensional systems, thermal effects harbour greater importance still, due to the existence of enhanced phase fluctuations within systems under tight confinement in one or two directions Low_D_Experiment_1 ; Low_D_Experiment_2 ; Low_D_Experiment_3 . Moreover, the dynamics of the thermal component can play a significant role, for example in the case of an ultracold gas under strong external perturbations Jackson_Zaremba_1 ; Jackson_Zaremba_4 ; Excitations_ZNG ; Jackson_LasPhys .

A number of finite-temperature techniques have been developed to incorporate the effects of temperature into a description of trapped Bose gases, as recently summarised in MyReview ; MyTutorial ; NZ_OZ_tutorial ; Yukalov_SB , which feature a detailed comparison of the relative weaknesses and merits of such approaches.
In brief, various theories have been put forward and applied to the experiments, as discussed below: In mean-field type approaches, the condensate is treated as a classical object, separated from the effect of the thermal cloud; various such perturbative approaches have been formulated to second order in the effective interaction strength
ZNG ; Bijlsma_Zaremba_Stoof ; JILA_Kinetic_1 ; Proukakis_JPhysB , with the common feature of obtaining a set of coupled dynamical equations for the condensate and the non-condensate, which interact both via mean fields and by particle exchange.
The condensate obeys an appropriately generalised GPE, which,
in the most commonly adapted version following from early work of Kirkpatrick and Dorfman Kirkpatrick_1 ; Kirkpatrick_2 ; Kirkpatrick_3 , is coupled to a quantum Boltzmann equation for the thermal cloud; this model has been extensively used by Zaremba, Nikuni, Griffin and co-workers, who introduced the terminology ‘ZNG’ ZNG .
Although such approaches work remarkably well in predicting a large range of experimentally-relevant phenomena, they do not fully describe regimes where fluctuations in the condensate phase are large, as is appropriate for example in low-dimensional systems
Low_D_Experiment_1 ; Low_D_Experiment_2 ; Low_D_Experiment_3 .
To incorporate such effects, generalised zero temperature Low_D_Petrov and finite temperature Low_D_Stoof modified methods were developed, to describe the static properties of such systems in an *ab initio* manner.

However, in order to go beyond mean field effects in a dynamical manner, various alternative approaches have been developed and implemented. In first instance, the time-dependent mean field approach mentioned above was reformulated in a somewhat more justified number-conserving manner Morgan_Gardiner , following from early work Girardeau_Arnowitt ; Gardiner_NC ; Castin_Dum ; although this has led to good agreement with experiments in at least one case Excitations_Burnett_2003 , a full dynamical model of such theories is still lacking.

Classical field approaches are based upon the observation that Svistunov_3 ; Polish_Review , as an equation for a classical field, the GPE should be able to describe the classical aspects of an ultracold trapped Bose gas. This includes all highly occupied modes, for which it is possible to approximate the field operators for the creation/annihilation of particles, by c-number amplitudes PGPE_Trap . Consistency requires only the highly occupied modes to be propagated this way, which can be achieved by introducing a projector to the GPE PGPE . Temperature effects can be incorporated through the initial conditions, which are typically set to some non-equilibrium distribution, and subsequently allowed to equilibrate at some temperaturePGPE_Trap , via the (projected-)GPE, under the restriction of fixed particle number and energy.

A variant of this technique is known as the Truncated Wigner method, originally developed in the context of quantum optics and first used for Bose gases in Steel_TWA ; a key element of this method is the approximate incorporation of quantum effects, through the addition of a prescribed amount of quantum noise to the wavefunction initial conditions. The numerical procedure which results, is to propagate a set of stochastic matter fields, which serve to sample the Wigner distribution under this approximation, via an ‘appropriate’ equation of motion for the Wigner distribution. The terminology ‘truncated’ Wigner arises from the fact that the propagation is carried out in an approximate manner, via the (projected) GPE. A limitation of this approach is its inability to accurately describe thermal effects, although certain remedies of limited applicability have been proposed Castin_TWA . In the latter two approaches, a temperature can only be assigned at the end of the simulations, with the thermalisation arising respectively due to ergodicity PGPE_Trap , while simultaneously exhibiting an unphysical heating.

To gain better control over the temperature to which the system equilibrates, requires detailed consideration of the coupling between such modes, and the higher lying, thermal region in our system. Within such schemes, an appropriate generalisation has been formulated in PGPE_T . Earlier, rather distinct, yet analogous in spirit, means of achieving this was provided independently by Stoof Stoof_PRL ; Stoof_JLTP ; Stoof_Langevin ; Stoof_Duine ; Stoof_LesHouches , and by Gardiner, Davis and co-authors QK_V ; QK_VII ; SGPE_I ; SGPE_II ; SGPE_III , using very different theoretical approaches.

This led to the generalisation of the GPE to a Langevin equation Stoof_Langevin , often referred to as the Stochastic-Gross-Pitaevskii equation (SGPE) SGPE_II , a terminology that we also adopt here; these two independent formulations have many formal similarities and are closely related, up to some subtle difference discussed in MyTutorial ; SGPE_II ; QK_VII . The SGPE describes the stochastic evolution of the condensate plus highly occupied, low-lying thermal modes, under the influence of coherent and incoherent interactions with the particles occupying higher energy modes.

Given the diverse range of approaches currently available, our aim in this manuscript is to highlight the key elements of this latter approach, paying particular attention on some of its appealing features via somewhat idealised experimental scenarios. We start our discussion with a brief overview of one formulation of this approach, which has already attracted a relatively large literature Stoof_Langevin ; Stoof_Duine ; Low_D_PRA ; SGPE_AtomLaser ; Stoof_Vortex ; Proukakis_AtomChip ; Proukakis_Equilibrium , with related work discussed in SGPE_II ; SGPE_III ; SGPE_IV ; Davis_SGPE_New ; Davis_SGPE_Nature .

## Ii Methodology

### ii.1 The Stochastic Gross-Pitaevskii equation

The stochastic Gross-Pitaevskii equation is a non-equilibrium microscopic theory to describe the evolution of an ensemble of ultracold atoms in contact with a thermal cloud. Using the Keldysh non-equilibrium formalism, Stoof derived a Fokker-Planck equation for the the system in Stoof_JLTP ; such an equation is equivalent to the following Langevin equation Stoof_Langevin ,

(1) |

where here is the order parameter for the lowest, coherent system modes and the higher modes are assumed to obey their own coupled dynamics, being in general governed by a quantum Boltzmann equation as in ZNG ; Bijlsma_Zaremba_Stoof ; MyTutorial . In addition to the usual Gross-Pitaevskii terms for the kinetic energy, trapping potential () and nonlinear interaction term (where , parameterised by the s-wave scattering length ), this equation features two additional contributions: firstly, denotes a dissipative term which represents the coupling between the condensate and a thermal reservoir at a temperature , and accounts for the transfer of particles between the high and low energy modes of the system; the contribution is a noise term which accounts for the random nature of incoherent collisions within the system and is somewhat analogous to the random particle jitter in Brownian-type motion (see, e.g. Stoof_Duine ). Finally, denotes an effective chemical potential. The above noise contribution is characterized by Gaussian correlations of the form

(2) |

where denotes averaging over different realisations of the noise (and it is generally understood that the results of simulations will undergo appropriate averaging, although single-runs can also offer qualitative results as discussed later). The amplitude of the noise is in general position and time dependent, and is governed by the so-called Keldysh self-energy . For a thermal cloud which is sufficiently close to equilibrium, this takes the form Stoof_Langevin

(3) | |||||

Although the occupation numbers, , and therefore also energies are in general time-dependent (their values governed by a quantum Boltzmann equation and self-consistent Hartree-Fock energies Stoof_Langevin ; MyTutorial ), current numerical implementations only focus on near-equilibrium situations, for which the thermal cloud is treated as static; in this case, can be approximated by the usual Bose-Einstein distributions , where and labels the momentum of a particle in the th single particle energy level.

As one might expect, the strength of the dissipative effects arising from dynamical particle exchange between the two subsystems is ‘balanced’ by the magnitude of fluctuations present, a relation encapsulated by the fluctuation-dissipation relation for the system, upon noting explicitly the dependence of these quantities on the system energy , namely In order to make the solution of (1) numerically tractable, in first instance one can allow the system to equilibrate to a classical distribution by noting that for large occupation numbers, . Since the condensate energy is still an operator of the form , the Langevin equation takes the somewhat simplified form

The very first numerical implementation of this SGPE in the context of ultracold Bose gases was undertaken in 2001 by Stoof and Bijlsma Stoof_Langevin to model the reversible growth of a condensate, as seen in the MIT experiment of Stamper-Kurn et al. MIT_Dimple : in this experiment, a dimple potential was introduced in a reversible manner to the centre of a harmonic trap, thereby inducing localised phase space compression, and leading to a (periodic) crossing through the phase transition. Good qualitative agreement was demonstrated between the numerical results Stoof_Langevin and the observed lagging behind of the condensate growth, relative to the addition of the dimple trap, found in the experiment. The method, which is also amenable to analytic variational calculations of stochastic dynamics Stoof_Duine ; Stoof_Vortex , has since been applied (in collaboration with one of us, NPP) to numerical studies of coherence Proukakis_Equilibrium ; Low_D_PRA , atom laser SGPE_AtomLaser , and atom chip dynamics Proukakis_AtomChip . It should be noted that related work has been recently undertaken by Gardiner, Davis and co-workers SGPE_I ; SGPE_II ; SGPE_III ; SGPE_IV ; Davis_SGPE_New ; Davis_SGPE_Nature .

The aim of this paper is to briefly review the main concepts and illustrate the versatility of the SGPE formalism for describing systems in which fluctuations are of great importance; for numerical convenience, most of our discussion is restricted to effectively one-dimensional systems, for which we also consider the reduction of the SGPE to simpler theories, which incorporate finite temperature effects, to differing levels of approximation. The extension of this approach to two-dimensional systems is also briefly reviewed. Although, by construction, the results of such simulations are best interpreted after numerical averaging over the different noisy initial conditions, we shall see that important information is also contained in single runs, as also noted in Davis_SGPE_New ; Davis_SGPE_Nature .

## Iii One-Dimensional Applications

### iii.1 Spatiotemporal Evolution of Coherence in a One-dimensional Bose Gas

#### iii.1.1 Non-equilibrium Coherence

Low dimensional systems exhibit a richer phase diagram than their three-dimensional counterparts: in the case of a purely 1d homogeneous system, long wavelength fluctuations in the phase Popov , prevent the onset of Bose-Einstein condensation, even down to zero temperature. For a trapped, one-dimensional gas, however, the long-wavelength physics is altered due to the low energy cut-off imposed by the harmonic trapping potential, and the quantum degenerate regime is instead described by two characteristic temperatures, and , below which density and phase fluctuations, respectively, become suppressed Low_D_Petrov . The intermediate regime Low_D_PhaseDiagram , is found to be one in which density fluctuations no longer dominate, yet fluctuations in the phase remain; this is the regime in which the system is said to contain a quasi-condensate Popov .

Under tight transverse confinement, the description of the dynamics of highly elongated three-dimensional systems effectively reduces to consideration of a one-dimensional equation, as dynamics in the transverse direction becomes restricted to the lowest mode. As this is the case in many current experiments with atom chips Low_D_Experiment_2 , it is crucial to fully understand and characterize such fluctuations. The SGPE is well suited to investigations of phase-fluctuating, or quasi-condensates, as fluctuations about the mean field are implicitly retained within the order parameter, thus giving access to information on coherence, available through the evaluation of various temporal and spatial numerical auto-correlation functions.

A study of the lowest order correlation functions of the gas reveals crucial information about the system. In particular, information about the phase of the system is contained in the first order correlation function, with the next higher order correlation containing crucial details about density fluctuations. Such quantities can be directly measured in interferometric ‘juggling’ experiments Orsay ; Hannover1 ; Hannover2 , and can also be easily obtained numerically within this formalism. More specifically, the renormalised first order correlation function at time can be obtained in ‘symmetric form’ about a point , via the relation

(4) |

We start our discussion by the basic demonstration of the stochastic method, investigating the nonequilibrium evolution of the density and first-order correlation function during the growth of a one-dimensional quasi-condensate from a static thermal cloud. Fig. 1 shows snapshots of the time evolution of both the density (top row) and coherence (bottom row) towards their equilibrium values, for a gas in contact with a heat bath at a fixed temperature. As the atomic system is cooled and heads towards the new equilibrium, the spatial extent of the correlation function increases as phase coherence is established, with the final shape attained dictated by the temperature, for a given particle species, number and trap geometry. In this and all subsequent figures on coherence, distances are scaled to the effective quasi-condensate spatial extent at equilibrium. This quantity takes account of the varying system size at different temperatures due to quasi-condensate depletion, and reduces at to the ‘usual’ Thomas-Fermi solution . There are two ways to obtain , as discussed in Proukakis_Equilibrium : firstly, by considering the modified low-dimensional theory of Andersen et al. Low_D_Stoof , or equivalently by using the curvature of the density profiles to obtain a corresponding quasi-condensate size. The particular example considered here corresponds to the case where the coherence length of the system is comparable to the system size, with . The observed growth in the densities is a direct consequence of the transfer of particles from the heat bath to the system, and the noise is crucial in order to seed the growth process Stoof_JLTP ; Stoof_Langevin .

For temperatures such that , the first-order correlation function has been found experimentally to exhibit exponential decay, whereas for the corresponding spatial dependence becomes Gaussian Orsay ; Hannover1 ; Hannover2 . Thus, one possible way to quantify the degree of coherence and extract a uniformly-varying coherence length is based on a simple fitting function of the form Proukakis_Equilibrium , , in which is a parameter characterizing the crossover from exponential to Gaussian behaviour.

A comparison between the spatial coherence at the trap centre (i.e. for ) for early and late times, together with the results of the fitting procedure, is shown in Fig. 2 for two temperatures, corresponding roughly to the case of a ’pure’ condensate (left images) and a gas at the crossover from ‘pure’ to ‘quasi’-condensation (right images).

As is evident from this figure (right images), for the higher temperature results, corresponding to the parameters of Fig. 1, the equilibrium correlation function decays more rapidly than for the lower temperature limit. Although the spatial extent of the system differs between these two cases ( vs. ), the appropriately scaled correlation function plotted in the top row of Fig. 2 exhibit very similar behaviour, unlike the case at their respective equilibria (bottom row). Nonetheless, the chosen fitting function seems to appropriately capture the main features; importantly, this enables us to study in a universal manner the growth of coherence as the system is equilibrating, being driven by the heat bath.

Numerous experiments MIT_Formation ; QK_PRL_III ; Growth_Amsterdam ; Growth_Zurich have been performed to study the relation between condensate growth and evolution of coherence, and a detailed comparison to such experiments is pending (but see also QK_JPhysB ; Growth_Orsay and references therein). For our present illustrative purposes, we simply highlight that the growth of coherence appears to exhibit similar behaviour to typical quasi-condensate number growth curves, although these generally follow such curves with an additional small lagging time; most notably, once the norm has reached its steady-state value, the coherence is still growing, albeit at a very slow rate.

With this means of quantifying the degree of coherence across a sample, it is then interesting to assess how varies with time during the cooling process. The analysis of these results for the two different temperatures in the context of this somewhat idealised scenario is plotted in Fig. 3. As these simulations are based on a fixed value of , the final atom numbers are slightly increased in the high temperature case (top right image); note also the significantly increased relaxation time and observation of the initial spontaneous initiation process prior to the domination of the bosonic stimulation effect (top left).

#### iii.1.2 Dependence of Equilibrium Coherence on Position

We also investigate the position dependence of the correlation function at equilibrium, by comparing the evolution of at different positions , for a fixed temperature corresponding to the high temperature case presented in the above results. The results of this investigation are shown in Fig. 4 for and for various times during the growth to equilibrium. The final profile attained at this temperature, together with the Thomas-Fermi profile for the same number of atoms, is also shown, with the positions considered indicated. The coherence within the Thomas-Fermi radius decreases uniformly but at a relatively slow rate as soon as one probes off-centred regions, with the effective coherence length decreasing rapidly as soon as one moves beyond the Thomas-Fermi radius (right images).

#### iii.1.3 Potential ’Identification’ of a Quasi-condensate

An appealing feature of the stochastic approach that is often not fully appreciated, is that it actually generates total density profiles without needing to resort to a special mean-field treatment of the condensate mode. This means that, as in the experiments, one does not have an automatic ’visualisation’ of the condensate in the system. So, in order to analyse these results, one may consider applying bimodal fits precisely as done in the experiments. However, an alternative definition has been proposed Svistunov_QC , in which the quasi-condensate density is identified as
This approach has been implemented in Proukakis_Equilibrium ,
with the density profile contribution for the thermal cloud over low-lying modes determined self-consistently via the relation
within the region , and given simply by outside of the Thomas-Fermi radius.
The identification of these quantities enables spatial plots of the temperature dependent quasi-condensate density profiles to be obtained in an *ab initio* manner, as shown for different temperatures in Fig. 5.
Such a definition was used more recently to distinguish between a condensate and a quasi-condensate Davis_SGPE_New ; Bisset_Rot_cond

Having discussed some basic elements of finite temperature trapped, quasi-one-dimensional Bose gases via the SGPE, we now discuss certain limiting cases of the above theory which are commonly used in the literature, and which shed more light on the benefits and versatility of the SGPE.

### iii.2 Reduction to Alternative Theories

We have already argued that the two important features of the SGPE are that it includes the correct damping due to its dynamical coupling to the heat bath, and that it extends beyond mean field effects, by including a stochastic noise term. Removing either (or both) of these elements leads to simpler theories MyTutorial , which can nonetheless be useful in certain regimes, as discussed below.

#### iii.2.1 Removing Contact with the Heat Bath

In first instance, we can consider what happens when the coupling to the heat bath is abruptly (and perhaps somewhat artificially) removed. We recall that upon beginning a simulation with the SGPE, the effect of the term in (1) is to pump particles from the heat bath, and populate the low lying modes to their equilibrium occupations, in accordance with the classical equilibrium for a given set of trap frequencies, chemical potential and temperature. On reaching a dynamical equilibrium, the scattering rate becomes on average equal to zero, and there is no net particle exchange with the heat bath, although this continuous exchange can cause additional damping. However, any externally imposed perturbations of the system about this set of system parameters, would typically lead to the introduction or loss of particles, upon further evolution with (1).

In general, the SGPE is coupled to a quantum Boltzmann equation for the population of the higher-lying modes, and these two equations should be solved self-consistently. If this is not done, then modifying for example the trap geometry will lead to a change in the norm of the SGPE. If we are only considering the evolution of the low-lying modes, in order to keep the atom number within our system fixed, one could consider ‘fixing’ the norm in the numerics by removing the coupling to the heat bath. A better alternative might perhaps be to change the effective chemical potential in the SGPE to the final value corresponding to the same atom number in the new trap, although the observed dynamics will typically depend on how quickly the chemical potential is adjusted to its final value.

One procedure is thus to propagate initially the full stochastic solution in order to obtain the desired equilibrium state, before switching off the noise and damping terms of (1), and subsequently propagating the noisy initial condition via the usual GPE which enforces particle number conservation. This corresponds physically to the removal of contact with the heat bath upon reaching equilibrium.

As a method, this procedure has strong parallels with the Truncated Wigner method, and a comparative study of the two is under way currently SGPE_vs_TWA . The main difference arises in the fact that in the usual truncated Wigner implementations (see, e.g. Steel_TWA ; TWA_Norrie_PRL ; Scott_Hutchinson_Gardiner or MyTutorial for a more extensive list), only quantum (as opposed to thermal) noise is typically added in the initial conditions (but see also Castin_TWA ; TWA_Ruostekoski_2 ), whereas the SGPE in general contains both quantum and thermal fluctuations Stoof_JLTP . Although existing numerical implementations of the SGPE focus on the thermal effects (due to the ‘classical’ approximation mentioned earlier), the important advantage of this approach is that there is an explicit fluctuation-dissipation relation which guarantees convergence to the correct classical equilibrium at any pre-determined temperature; in other words, the initial thermal state for the unperturbed system can be correctly assigned. The validity of conventional truncated Wigner implementations is on the other hand limited to very low temperatures Castin_TWA and occupation numbers for which quantum fluctuations have a dominant role - see SGPE_vs_TWA for a more detailed discussion.

This finite temperature approach based on the removal of the heat bath once the correct equilibrium has been obtained via the SGPE was used to probe the dynamics of an initially equilibrated thermal stationary state following the addition of a deep Gaussian dimple trap Proukakis_AtomChip . A large and sudden perturbation is found to lead to the appearance of shock waves in the system, and in the averaged profiles typically studied, these manifest themselves in terms of an initial increase in the central density, followed by periodic damped decrease-increase cycles, as shown in Fig.6; this can be interpreted as direct evidence of the temperature-induced suppression of the initial shockwave amplitude induced by the addition of the new trap. It should be noted here that the introduction of the correct thermal noise in the initial conditions is crucial to the generation of the damping which would not be present in simulations of the ordinary GPE starting from a stationary condensed state.

#### iii.2.2 A Microscopically-based Dissipative GPE

An even simpler limiting case of the SGPE is based on removing only the noise term of (1); then, the SGPE reduces to a ‘dissipative GPE’, as noted for example in SGPE_III ; QK_PRL_IV . The form of the damping is given by the self energy, and as such is spatially dependent, consistent with the spatially inhomogeneous background of the harmonic trap. Such a model has been used to study vortex and soliton dynamics leading to a number of interesting predictions. However, the origin of such models can be traced to rather crude approximations in the more advanced microscopic models, such as the SGPE, with most such publications justifying this equation phenomenologically, and actually using a constant value for the damping Choi_Phenomenology ; Vortex_Lattice_Ueda_1 ; Vortex_Lattice_Ueda_2 ; Proukakis_Parametric .

The dissipative equation which results from the SGPE in general takes, within the ‘classical approximation’, the form

(5) |

where the damping term (for a system made up of a fixed total number of atoms) is given by

(6) |

although our discussion so far has additionally ignored the time-dependence for numerical simplicity.

It is very easy to see that the resulting evolution is no longer norm conserving, and it therefore becomes nessessary to renormalize the wavefunction as each time step to maintain a constant particle number.

For illustrative purposes, we apply this approach here to model the oscillations of an ‘ideally generated’ dark soliton within a one-dimensional condensate. The dynamics for three identical initial soliton solutions are shown at different temperature in Fig. 7 (top images), with each image showing the position of the dark soliton for three different temperatures, corresponding to a different evolution time. To best illustrate the different temperature-dependent decay rates, we give corresponding space-time plots of the density along the z-axis with such damping already observed in experiments Hannover_Soliton_Exp .

Although consideration of the 1D case enables us to bring out the main physical features of the SGPE model, comparison to experiments requires us to consider its properties in higher-dimensional systems, and we conclude our presentation by a simple illustration in two-dimensional geometries, and some related comments.

## Iv Extension to Two Dimensions

As in one dimension, for a uniform two dimensional Bose gas, thermal fluctuations remain even at very low temperatures, preventing the onset of Bose-Einstein condensation Mermin ; Hohenberg . For interacting trapped two dimensional gases, however, quasi-condensation is possible for low temperatures Petrov2d , with a Berezinskii-Kosterlitz-Thouless phase transition to a superfluid state associated with the the pairing of vortices of opposite circulation beneath some critical temperature. Recently, this has spawned much interest in the physics of two-dimensional trapped ultracold systems, particularly in assessing the nature of the phase transition at the critical point Low_D_Experiment_3 ; DalibardNJPhys ; CornellPRL1 ; PGPE_BKT_PRL .

### iv.1 Spontaneous Processes - Vortex Nucleation

Phase transitions are strongly associated with such topological defects, the appearance of which are predicted by the Kibble-Zurek mechanism Spontaneous_Vortex_Kibble ; Spontaneous_Vortex_Zurek . An appealing feature of the stochastic approach is that in including fluctuations about the mean field, excitations in the form of spontaneous vortices can be seen during the growth process, as investigated recently using an alternative formulation of the SGPE Davis_SGPE_Nature . Our discussion so far has been restricted to analysing predictions of the SGPE based on numerical averaging over different initial noisy conditions. However, the monitoring of the evolution within a single run can also provide important information on the system parameters, which may otherwise be lost through averaging - e.g. information related to spontaneous events, such as topological defect formation Spontaneous_Vortex_BEC .

Within the numerical SGPE applied here, condensate formation from the static thermal cloud is initiated for a positive value of the chemical potential in (II.1). A rapid quench can be simulated by an instantaneous change in the chemical potential, leading to a rapid condensate growth. Within simulations, the appearance of a (random) number of vortices is common, although their details vary from run to run. The appearance of these actually becomes ‘washed out’ following the averaging process, as, due to the random nature of the nucleation of vortices, the positions at which they appear varies between different realisations. Such studies should be contrasted to studies based on the GPE, in which a vortex is actually ‘artificially’ imprinted at a predetermined position; when modelling the ensuing dynamics, one often introduces damping by resorting to the dissipative GPE , with the damping term, , which is actually proportional to the two-dimensional self-energy at a particular temperature, typically treated as a constant.

In density plots, vortices appear as dark regions, corresponding to positions at which singularities in the phase of the system can be seen. An example of the results from a typical simulation based on the SGPE is shown in Fig. 8. In the density and phase snapshots shown, two spontaneously generated vortices (top left) can be seen to come together (middle) and then annihilate (right), illustrating the spontaneous dynamics captured within simulations. The averaged density profiles included show the result of averaging over runs, and the smooth density plots which result.

### iv.2 Numerical Issues - Dependence on Momentum Cutoff

The numerical means of evolving the system discussed so far make use of the macroscopic occupation of the low lying modes, assuming these states to be well described by a classical field. In two dimensions, known divergences arise, associated with large momenta, or equivalently small grid spacings, as the continuum limit is approached on the lattice. This is due to the high momentum cutoff, , related to the numerical grid spacing, , through the relation .

We present here preliminary results of an investigation into this dependence for the SGPE as a ‘caution’ to the reader. Here we make no adjustments to parameters to compensate for grid dependence and find the expected logarithmic divergence in the momentum cut-off, with the system norm (or number of particles in the classical region) increasing with decreasing spatial discretization , as shown in Fig. 9.

In principle such divergences can be removed from the problem through the appropriate renormalization of physical quantities from their original (‘bare’) values. For example in homogeneous systems, the lattice effects on temperature or density can be taken into account in an *a posteriori* manner, as described in Prokov'ev_renorm . An alternative method is to introduce counter-terms to the system Hamiltonian, which has been successfully applied to Langevin equations previously in, for example, GleiserPRE2000 ; GleiserNucB ; WojtasThesis ; HabibBettencourt .

## V Conclusions

The stochastic Gross-Pitaevskii is emerging as a key approach for the equilibrium and dynamical properties of degenerate Bose gases, and is particularly beneficial for describing regimes of large phase and density fluctuations, such as near (or at) the phase transition, or for weakly-interacting low-dimensional Bose gases. Single run results contain important (at least qualitative) information on spontaneous processes, such as defect formation characteristics, whereas averaged profiles are most suitable for determination of correlation functions, quasi-condensate fractions and smooth density profiles. Although these stochastic simulations have obvious advantages, their existing numerical implementation have not yet explored such theories to their full potential, with current models typically assuming a static (near-equilibrium) thermal cloud, and additionally discarding quantum fluctuations in simulations. Finally, the important issue of the cut-off dependence, already highlighted by other researchers should not be overlooked, particularly for dimensions higher than one, and various techniques are available for renormalising the observed profiles to their actual values, thereby eliminating (at least to leading order) such dependences. We believe that use of the stochastic Gross-Pitaevskii equation, which is also coupled to a quantum Boltzmann equation for the self-consistent dynamical treatment of the higher-lying thermal modes, should be able to describe essentially all main features of weakly-interacting Bose gases currently pursued experimentally, and therefore look forward to the exciting times which lie ahead.

## Acknowledgments

We are grateful to Henk Stoof for early collaboration and discussion on these ideas and to the UK EPSRC for funding.

## References

- (1) Dalfovo F, Giorgini S, Pitaevskii L P, and Stringari S 1999 Rev. Mod. Phys. 71 463
- (2) Bloch I, Dalibard J and Zwerger W 2008 Rev. Mod. Phys. 80 885
- (3) Emergent Nonlinear Phenomena in Bose-Einstein Condensates: Theory and Experiment, Ed. Kevrekidis P G, Frantzeskakis D J and R. Carretero-Gonzalez (Springer Verlag, 2008)
- (4) Görlitz A, Vogels J M, Leanhardt A E, Raman C, Gustavson T L, Abo-Shaeer J R, Chikkatur A P, Gupta S, Inouye S, Rosenband T and Ketterle W 2001 Phys. Rev. Lett. 87 130420
- (5) Ott H, Fortagh J, Schlotterbeck G, Grossmann A and Zimmermann C 2001 Phys. Rev. Lett. 87 230401
- (6) Hadzibabic Z, Krüger P, Cheneau M, Battelier B and Dalibard J 2006 Nature (London) 441 1118
- (7) Jackson B and Zaremba E 2001 Phys. Rev. Lett. 87 100404
- (8) Jackson B and Zaremba E 2002 Phys. Rev. Lett. 89 150402
- (9) Jackson B and Zaremba E 2002 Phys. Rev. Lett. 88 180402
- (10) Jackson B and Zaremba E 2002 Laser Physics 12 93
- (11) Proukakis N P Beyond the Gross-Pitaevskii Mean Field in Ref. [3] [ arXiv:0706.3541v1 ]
- (12) Proukakis N P and Jackson B 2008 J. Phys. B 41 203002
- (13) Blakie P B, Bradley A S, Davis M J, Ballagh R J and Gardiner C W 2008 Advances in Physics 57 363
- (14) Yukalov V I 2007 Las. Phys. Lett. 4 632
- (15) Zaremba E, Nikuni T and Griffin A 1999 Low J Temp. Phys. 116 277
- (16) Bijlsma M J, Zaremba E and Stoof H T C 2000 Phys. Rev. A 62 063609
- (17) Walser R, Williams J, Cooper J and Holland M 1999 Phys. Rev. A 59 3878
- (18) Proukakis N P 2001 J. Phys. B 34 4737
- (19) Kirkpatrick T R and Dorfman J R 1983 Phys. Rev. A 28 2576
- (20) Kirkpatrick T R and Dorfman J R 1985 J Low Temp. Phys. 58 301
- (21) Kirkpatrick T R and Dorfman J R 1985 J. Low Temp. Phys. ibid. 58 399
- (22) Richard S, Gerbier F, Thywissen J H, Hugbart M, Bouyer P and Aspect A 2003 Phys. Rev. Lett. 91 010405
- (23) Hellweg D, Cacciapuoti L, Kottke M, Schulte T, Sengstock K, Ertmer W, and Arlt J J 2003 Phys. Rev. Lett. 91 010406
- (24) Cacciapuoti L, Hellweg D, Kottke M, Schulte T, Ertmer W, Arlt J J, Sengstock K and Santos L 2003 Phys. Rev. A 68 053612
- (25) Petrov D S, Shlyapnikov G V, and Walraven J T M 2000 Phys. Rev. Lett. 85 3745
- (26) Andersen J O, Al Khawaja U and Stoof H T C 2002 Phys. Rev. Lett. 88, 070407
- (27) Gardiner S A and Morgan S A 2007 Phys. Rev. A 75 043621
- (28) Girardeau M D and Arnowitt R 1959 Phys. Rev. 113 755
- (29) Gardiner C W 1997 Phys. Rev. A 56 1414
- (30) Castin Y and Dum R 1998 Phys. Rev. A 57 3008
- (31) Morgan S A, Rusch M, Hutchinson D A W and Burnett K 2003 Phys. Rev. Lett. 91 250403
- (32) Kagan Y and Svistunov B V 1994 Zh. Eksp. Teor. Fiz. 105 353 [1994 Sov. Phys. JETP 78 187]
- (33) Brewczyk M, Gajda M and Rzazewski K 2007 J. Phys. B 40 R1
- (34) Blakie P B and Davis M J 2005 Phys. Rev. A 72 063608
- (35) Davis M J, Morgan S A and Burnett K 2001 Phys. Rev. Lett. 87 160402
- (36) Steel M J, Olsen M K, Plimak L I, Drummond P D, Tan S M, Collett M J, Walls D F and Graham R 1998 Phys. Rev. A 58 4824
- (37) Sinatra A, Lobo C and Castin Y 2002 J. Phys. B 35 3599
- (38) Davis M J, Ballagh R J and Burnett K 2001 J. Phys. B 34 4487
- (39) Stoof H T C 1997 Phys. Rev. Lett. 78 768
- (40) Stoof H T C 1999 Low J Temp. Phys. 114 11
- (41) Stoof H T C and Bijlsma M J 2001 J. Low Temp. Phys. 124 431
- (42) Duine R A and Stoof H T C 2001 Phys. Rev. A 65 013603
- (43) Stoof H T C in Coherent Atomic Matter Waves, Proceedings of the Les Houches Summer School Session 72, 1999, ed. Kaiser R et al. (Springer-Verlag, Berlin 2001).
- (44) Gardiner C W and Zoller P 2000 Phys. Rev. A 61 033601
- (45) Davis M J, Gardiner C W and Ballagh R J 2000 Phys. Rev. A 62 063608
- (46) Gardiner C W Anglin J R, and Fudge T I A 2002 J. Phys. B 35 1555
- (47) Gardiner C W and Davis M J 2003 J. Phys B 36 4731
- (48) Bisset R N, Davis M J, Simula T P, and Blakie P B 2009 Phys. Rev. A 79 033626
- (49) Al Khawaja U, Andersen J O, N. Proukakis P, and H.T.C. Stoof 2002 Phys. Rev. A 66, 013615 ; ibid. 66, 059902(E) (2002)
- (50) Proukakis N P 2003 Las. Phys. 13 527
- (51) Duine R A, Leurs B W A and Stoof H T C 2004 Phys. Rev. A 69 053623
- (52) Proukakis N P, Schmiedmayer J and Stoof H T C 2006 Phys. Rev. A 73 053603
- (53) Proukakis N P 2006 Phys. Rev. A 74 053617
- (54) Bradley A S, Blakie P B and Gardiner C W 2005 J. Phys. B 38 4259
- (55) Bradley A S, and Gardiner C W, cond-mat/0602162.
- (56) Bradley A S, Gardiner C W and Davis M J 2008 Phys. Rev. A 77 033616
- (57) Weiler C N, Neely T W , Scherer D R, Bradley A S, Davis M J, Anderson B P 2008 Nature 455 948
- (58) Stamper-Kurn D M, Miesner H-J, Chikkatur A P, Inouye S, Stenger J and Ketterle W 1998 Phys. Rev. Lett. 81 2194
- (59) Popov V N 1965 Sov. Phys. JETP 20 1185; Popov V N Functional Integrals and Collective Excitations (Cambridge University Press, 1987)
- (60) Al Khawaja U, Proukakis N P, Andersen J O, Romans M W J and Stoof H T C 2003 Phys. Rev. A 68 043603
- (61) Miesner H J, D.M. Stamper-Kurn, Andrews M R, Durfee D S, Inouye S and Ketterle W 1998 Science 279 1005
- (62) Köhl M, Davis M J, Gardiner C W, Hänsch T W and Esslinger T 2002 Phys. Rev. Lett. 88 080402
- (63) Shvarchuck I, Buggle Ch., Petrov D S, Dieckmann K, Zielonkowski M, Kemmann M, Tiecke T G, von Klitzing W, Shlyapnikov G V and Walraven J T M 2002 Phys. Rev. Lett. 89 270404
- (64) Ritter S, Öttl A, Donner T, Bourdel T, Köhl M and Esslinger T 2007 Phys. Rev. Lett. 98 090402
- (65) Hugbart M, Retter J A, Varon A F, Bouyer P, Aspect A and Davis M J 2007 Phys. Rev. A 75 011602(R)
- (66) Davis M J and Gardiner C W 2002 J. Phys. B 35 733
- (67) Prokof’ev N and Svistunov B 2002 Phys. Rev. A 66 043608
- (68) Cockburn S P, Negretti A, Henkel C and Proukakis N P 2009 (in preparation)
- (69) Norrie A A, Ballagh R J and Gardiner C W 2005 Phys. Rev. Lett. 94 040401
- (70) Scott R G, Hutchinson D A W and Gardiner C W 2006 Phys. Rev. A 74 053605
- (71) Ruostekoski J and Isella L 2005 Phys. Rev. Lett. 95 110403
- (72) Penckwitt A A, Ballagh R J and Gardiner C W 2002 Phys. Rev. Lett. 89 260402
- (73) Choi S, Morgan S A and Burnett K 1998 Phys. Rev. A 57 4057
- (74) Tsubota M, Kasamatsu K and Ueda M 2002 Phys. Rev. A 65 023603
- (75) Kasamatsu K, Tsubota M and Ueda M 2003 Phys. Rev. A 67 033610
- (76) Proukakis N P, Parker N G, Barenghi C F and Adams C S 2004 Phys. Rev. Lett. 93 130408
- (77) Burger S, Bongs K, Dettmer S, Ertmer W, Sengstock K, Sanpera A, Shlyapnikov G V and Lewenstein M 1999 Phys. Rev. Lett. 83, 5198
- (78) Mermin N D and Wagner H 1966 Phys. Rev. Lett. 22 1133
- (79) Hohenberg P C 1967 Phys. Rev. 158 383
- (80) Petrov D S, Holzmann M, and Shlyapnikov G V 2000 Phys. Rev. Lett. 84 2551
- (81) Hadzibabic Z, Krüger P, Cheneau M, Rath S P and Dalibard J 2008 New J. Phys. 10 045006
- (82) Schweikhard V, Tung S, and Cornell E A 2007 Phys. Rev. Lett. 99 030401
- (83) Simula T P and Blakie P B 2006 Phys. Rev. Lett. 96 020404
- (84) Kibble T W B 1976 J. Phys. A 9 1387
- (85) Zurek W H 1985 Nature (London) 317 505
- (86) Anglin J R and Zurek W H 1999 Phys. Rev. Lett. 83 1707
- (87) Prokof’ev N, Ruebenacker O, Svistunov B 2001 Phys. Rev. Lett. 87 270402
- (88) Borrill J and Gleiser M 1996 Nucl. Phys. B 483 416
- (89) Gagne C and Gleiser M 2000 Phys. Rev. E 61 3483
- (90) Wojtas D, 2006 MSc. Thesis, University of Canterbury, Christchurch
- (91) Lythe G, Bettencourt L and Habib S 1999 Physical Review D 60 105039