Class: BoseEinstein Condensate
There are two types of particles in the world: fermions and bosons. Whereas fermions can never occupy the same quantum state due to the Pauli Exclusion principle, the same is not true for bosons. A BoseEinstein condensate (BEC) is a state of matter consisting of ultracold bosons that undergo a phase transition at a low critical temperature in which most bosons occupy the ground state of the system. It was theorized by Einstein and Bose in the 1920s as a new state of matter and produced for the first time in 1995 by Eric Cornell and Carl Wieman^{1}.
In this class, we simulate a BoseEinstein condensate in 1, 2 and 3 dimensions using the GrossPitaevski equation (GPE).
Variables and parameters
The primary field is the complex wave function \(\psi\)
The BoseEinstein condensate takes the same keyword arguments as the BaseSystem in addition toKeyword  Definition  Default value 

gamma 
dissipative factor  \(0.01\) 
Model
The BEC is in the mean field regime described by the GPE^{2} ^{3}. This is a nonlinear Schrödinger equation which reads
Here, \(\mu\) is the chemical potential, \(m\) is the mass of the bosons, \(g\) is an interaction parameter and \(\hbar\) is the Planc constant. \(\psi\) is the wave function describing the condensate phase and \(V_{ext}\) is an external potential. The GPE can be obtained from the variational principle ^{3}^{4}
with the Hamiltonian
We introduce dimensionless units for length \(\xi = \hbar/\sqrt{m\mu}\), time \(\tau = \xi/c\) and energy \(E=\eta\) and rescaling the wave function to \(\psi \rightarrow \sqrt{\frac{g}{\mu}}\psi\), in addition we include a dissipative factor \(\gamma\). This results in the dGPE on dimensionless form as ^{5} ^{6} ^{7} ^{8}
The danmped GrossPitaevski equation (evolve_dGPE
)
So we see that we can split into its linear and nonlinear parts as
The Hamiltonian and its density in dimensionless units are calculated by the functions
Approximation of Ground States
In simulations, it is often convenient to start in a configuration that is close to the ground state. We can estimate this ground state by noticing that the GPE dissipates energy when it is evolved in imaginary time \(t \rightarrow it\) ^{9} ^{3}. We refer to this as evolving the dGPE in imaginary time. Given an external potential \(V_{ext}\) we can therefore find an approximation to the ground state by starting with a guess and then removing energy from the guessed state by evolving the equations in imaginary time.
To get a guess of the ground state for a given potential \(V_{ext}\) we use the ThomasFermi approximation ^{2} ^{3}, where we assume that \(\psi\) is slowly varying so that we can neglect the Laplacian term. Looking for stationary solutions to the dGPE we obtain the equation
This has two solutions, \(\psi = 0\) and
In the case where \(V_{ext} > 1\) there is only one possibility, namely \(\psi = 0\). In the case of \(V_{ext} < 1\), both the stationary solutions exist, so we need to evaluate their energy. We can do this by considering the Hamiltonian
Inserting the anzats \(\psi =0\) we see that this vanishes, while if we insert the anzats \(\psi^2 = 1V_{ext}\), we get something negative. We therefore conclude that the ThomasFermi ground state is given as
This ground state can be initialized as
Once this guess is initialized we can propagate the wave function in imaginary time using the function
Note that the potential has to be a constant during this evolution.
Potentials
A lot of the interesting dynamics of a BEC come about when it is interacting with an external potential. This is included as the function
By default, this is assumed to be timeindependent and returns the field
The potential can be changed by the function
which can be used both to set it as a function or to set it as a constant potential depending on whether V_ext
is a function, a constant or an numpy array.
If additive =True
one add the constant V_ext
to the allready existing potential.
If V_ext
is a function it needs to be on the form
The evolver will then evaluate it using the bec.time
variable which is updated on the run.
An example using a Gaussian stirring potential is provided in the example folder.
To make life easier we have provided a couple of popular potentials. The harmonic potential is provided in the function
Here \(R_{tf}\) is the ThomasFermi radius ^{3}, and the harmonic potential takes the form \(\(V_H = \frac{r^2}{R_{tf}^2}\)\). The Gaussian potential is provided through the function
giving the potential
Here \(g\) is the strength, \(\sigma\) is the size and \(\mathbf{r}_p\) is the position.
Much of the interesting physics happens when the potential is timedependent. For this one can define a function as
and update the external potential to this function by calling
Hydrodynamics
The dGPE can be transformed to a hydrodynamic description of the BEC. The first step in doing this is to introduce the Madelung transformation \(\psi = \sqrt{\rho} e^{i\theta}\), where \(\rho = \psi^2\) is the superfluid density ^{3}. For \(\gamma = 0\) this density is conserved and satisfies the conservation equation
with the superfluid current given by
Here the superfluid velocity is introduced as \(\mathbf{v}_s = \nabla \theta\). The superfluid current can be found with the function
We can also put the Madelung transformation into the Hamiltonian to get ^{7} ^{10}
The first term here is the kinetic energy of the condensate. To calculate this it is convenient to introduce the density weighted velocity \(\mathbf{u} = \sqrt{\rho}\mathbf{v}_s\)^{7}. This have the advantage of not being singular at the centre of the topological defects. Using this we can write the kinetic energy as \(\(E_k =\int d \mathbf r \frac{1}{2}u^2.\)\) The density weighted velocity and the kinetic energy can be calculated by the functions
Further, if we insert the Madelung transformation into the dGPE and do some work we can map it into the NavierStokes equations ^{3} ^{7}
Notice that the condensate density is only conserved when \(\gamma = 0\).
Forces on external potential
To study how the BEC interacts with impurities, one can model the impurity as a Gaussian potential and measure the forces that are acting on it ^{11} ^{12} ^{13}. From the Ehrenfest theorem the forces on the condensate from the stirrer is given as
which means that the force acting on the stirrer is \(\mathbf{F} = \langle \nabla V_{ext}\rangle\). Written explicitly, this is
and is calculated by the function
Note that this calculates the average force on all of the external potential.
Quantized vortices
The topological defects in the BEC take the form of quantized vortices. This is because the velocity is given by the gradient of the phase \(\mathbf{v}_s = \nabla \theta\). This has the consequence that the circulation \(\(\Gamma = \oint d\mathbf{l} \cdot \mathbf{v} = \oint d\mathbf{l} \cdot \nabla \theta = \int d\theta\)\) is quantized in units of \(2\pi\), since the field \(\psi\) is singlevalued. In two dimensions this vortex is a point defect, while in three dimensions it is a line defect. For the BEC this method is given by the function
which finds the vortex nodes in two or three dimensions.
If \(\partial_t \psi\) is put in as the field dt_psi
, then the velocity of the defects are also found.
The defects can be created in multiple ways. In addition to putting them in by hand, one can create them by stirring the condensate with a potential or by relaxing a disordered initial state.
Comoving frame
When studying a BEC that is stirred by a potential it is in some cases convenient to consider the potential as stationary with the BEC flowing around. This can be done by boosting the dGPE so that it reads
where \(\mathbf{V_p}\) is the velocity of the boost. Note that a Galilean boost of the GPE is often accompanied by a phase shift of the wave function \(\psi \rightarrow \psi \exp{(\mathfrak i\mathbf V_p \cdot \mathbf r + \frac i 2 V_p^2 t)}\), which transforms the superfluid velocity to the new reference frame^{14}, leaving the GPE unchanged after the Galilean transformation. However, the equation with \(\gamma \neq 0\) is not Galilean invariant and we therefore do not include the phase factor in the transformation. This has the consequence that the superfluid velocity obtained from \(\psi\) is measured in the lab frame and not the comoving frame^{11}.
To reduce the recycling of excitations into the incoming flow, we introduce a buffer region around the computational domain where \(\gamma\) is large, similar to the one used in Ref. ^{15}. The dissipative factor becomes a function of space and is given by
in three dimensions and
in two dimensions. Here
and similar for \(\gamma_y(y)\) and \(\gamma_z(z)\). The constant \(\gamma_0\) is the value in the bulk, \(d\) is the size of the interface between the bulk and the buffer and \(w_x\) is the distance from the centre of the domain to the buffer in the \(x\) direction. Note that when \(\gamma\) is a function of space we can no longer put it into the nonlinear differential operator \(\omega\), and we have to move it into the nonlinear part. This is taken care of in the evolver
Here it is assumed that the boost is in the \(x\)direction, and that the dissipative factor is spatially dependent.

Anderson, M. H., Ensher, J. R., Matthews, M. R., Wieman, C. E., & Cornell, E. A. (1995). Observation of BoseEinstein Condensation in a Dilute Atomic Vapor. Science, 269(5221), 198–201. https://doi.org/10.1126/science.269.5221.198 ↩

Dalfovo, F., Giorgini, S., Pitaevskii, L. P. and Stringari, S. (1999). Theory of BoseEinstein condensation in trapped gases. Reviews of Modern Physics. 71, 3, 463. https://doi.org/10.1103/RevModPhys.71.463 ↩↩

Kevrekidis, P. G., Frantzeskakis, D. J. and CarreteroGonzález, R. (2008). Emergent nonlinear phenomena in BoseEinstein condensates: theory and experiment. Springer Science & Business Media. Berlin. ↩↩↩↩↩↩↩

Pitaevskii, L. and Stringari, S. (2016). BoseEinstein Condensation and Superfluidity. Oxford University Press. https://doi.org/10.1093/acprof:oso/9780198758884.001.0001 ↩

Gardiner, C. W. and Davis, M. J. (2003). The stochastic GrossPitaevskii equation: II. Journal of Physics B: Atomic, Molecular and Optical Physics. 36, 23, 4731. https://doi.org/10.1088/09534075/36/23/010 ↩

Rooney, S. J., Blakie, P. B. and Bradley, A. S. (2012). Stochastic projected GrossPitaevskii equation. Physical Review A. 86, 5, 053634. https://doi.org/10.1103/PhysRevA.86.053634 ↩

Bradley, A. S. and Anderson, B. P. (2012). Energy spectra of vortex distributions in twodimensional quantum turbulence. Physical Review X. 2, 4, 041001 https://doi.org/10.1103/PhysRevX.2.041001 ↩↩↩↩

Skaugen, A. (2018). A Unified Perspective on TwoDimensional Quantum Turbulence and Plasticity. PhD Thesis, University of Oslo. http://urn.nb.no/URN:NBN:no69394 ↩

Minguzzi, A., Succi, S., Toschi, F., Tosi, M. P. and Vignolo, P. (2004). Numerical methods for atomic quantum gases with applications to BoseEinstein condensates and to ultracold fermions. Physics reports. 395, 45, 223355. https://doi.org/10.1016/j.physrep.2004.02.001 ↩

Nore, C., Abid, M. and Brachet, M. E. (1997). Kolmogorov turbulence in lowtemperature superflows. Physical review letters. 78, 20, 3896. https://doi.org/10.1103/PhysRevLett.78.3896 ↩

Rønning, J., Skaugen, A., HernándezGarcía, E., López, C. and Angheluta, L. (2020). Classical analogies for the force acting on an impurity in a BoseEinstein condensate. New Journal of Physics. 22, 7, 073018. https://doi.org/10.1088/13672630/ab95de ↩↩

Astrakharchik, G. E. and Pitaevskii, L. P. (2004). Motion of a heavy impurity through a {BoseEinstein} condensate. Physical Review A. 70, 1, 013608. https://doi.org/10.1103/PhysRevA.70.013608 ↩

Pinsker, F. (2017). Gaussian impurity moving through a {BoseEinstein} superfluid. Physica B: Condensed Matter. 521, 3642. https://doi.org/10.1016/j.physb.2017.06.038 ↩

Pismen, L.M. (1999). Vortices in Nonlinear Fields: From Liquid Crystals to Superfluids, From NonEquilibrium Patterns to Cosmic Strings. Oxford university press. Oxford ↩

Reeves, M. T., Billam, T. P., Anderson, B. P. and Bradley, A. S. (2015). Identifying a superfluid Reynolds number via dynamical similarity. Physical Review Letters. 114, 15, 155302. https://doi.org/10.1103/PhysRevLett.114.155302 ↩