Gross–Pitaevskii equation

The simplest exact solution is the free particle solution, with ${\displaystyle V(\mathbf {r} )=0}$,

${\displaystyle \Psi (\mathbf {r} )={\sqrt {\frac {N}{V}}}e^{i\mathbf {k} \cdot \mathbf {r} }.}$

This solution is often called the Hartree solution. Although it does satisfy the Gross–Pitaevskii equation, it leaves a gap in the energy spectrum due to the interaction:

${\displaystyle E(\mathbf {k} )=N\left[{\frac {\hbar ^{2}k^{2}}{2m}}+g{\frac {N}{2V}}\right].}$

According to the Hugenholtz–Pines theorem, [4] an interacting bose gas does not exhibit an energy gap (in the case of repulsive interactions).

A one-dimensional soliton can form in a Bose–Einstein condensate, and depending upon whether the interaction is attractive or repulsive, there is either a bright or dark soliton. Both solitons are local disturbances in a condensate with a uniform background density.

If the BEC is repulsive, so that ${\displaystyle g>0}$, then a possible solution of the Gross–Pitaevskii equation is,

${\displaystyle \psi (x)=\psi _{0}\tanh \left({\frac {x}{{\sqrt {2}}\xi }}\right)}$,

where ${\displaystyle \psi _{0}}$ is the value of the condensate wavefunction at ${\displaystyle \infty }$, and ${\displaystyle \xi =\hbar /{\sqrt {2mn_{0}g}}=1/{\sqrt {8\pi a_{s}n_{0}}}}$ is the coherence length (a.k.a. the healing length[3], see below). This solution represents the dark soliton, since there is a deficit of condensate in a space of nonzero density. The dark soliton is also a type of topological defect, since ${\displaystyle \psi }$ flips between positive and negative values across the origin, corresponding to a ${\displaystyle \pi }$ phase shift.

For ${\displaystyle g<0}$

${\displaystyle \psi (x,t)=\psi (0)e^{-i\mu t/\hbar }{\frac {1}{\cosh \left[{\sqrt {2m\vert \mu \vert /\hbar ^{2}}}x\right]}},}$

where the chemical potential is ${\displaystyle \mu =g\vert \psi (0)\vert ^{2}/2}$. This solution represents the bright soliton, since there is a concentration of condensate in a space of zero density.

The healing length can be understood as the length scale where the kinetic energy of the boson equals the chemical potential[3]:

${\displaystyle {\frac {\hbar ^{2}}{2m\xi ^{2}}}=\mu =gn_{0}\,.}$

The healing length gives the shortest distance over which the wavefuction can change; It must be much smaller than any length scale in the solution of the single-particle wavefunction. The healing length also determines the size of vortices that can form in a superfluid; It's the distance over which the wavefunction recovers from zero in the center of the vortex to the value in the bulk of the superfluid (hence the name "healing" length).

In systems where an exact analytical solution may not be feasible, one can make a variational approximation. The basic idea is to make a variational ansatz for the wavefunction with free parameters, plug it into the free energy, and minimize the energy with respect to the free parameters.

Several numerical methods, such as the split-step Crank–Nicolson[5] and Fourier spectral[6] methods, have been used for solving GPE. There are also different Fortran and C programs for its solution for the contact interaction[7][8] and long-range dipolar interaction.[9]

If the number of particles in a gas is very large, the interatomic interaction becomes large so that the kinetic energy term can be neglected from the Gross–Pitaevskii equation. This is called the Thomas–Fermi approximation.

${\displaystyle \psi (x,t)={\sqrt {\frac {\mu -V(x)}{NU_{0}}}}}$

In a harmonic trap (where the potential energy is quadratic with respect to displacement from the center), this gives a density profile commonly referred to as the "inverted parabola" density profile[3].

Bogoliubov treatment of the Gross–Pitaevskii equation is a method that finds the elementary excitations of a Bose–Einstein condensate. To that purpose, the condensate wavefunction is approximated by a sum of the equilibrium wavefunction ${\displaystyle \psi _{0}={\sqrt {n}}e^{-i\mu t}}$ and a small perturbation ${\displaystyle \delta \psi }$,

${\displaystyle \psi =\psi _{0}+\delta \psi }$.

Then this form is inserted in the time dependent Gross–Pitaevskii equation and its complex conjugate, and linearized to first order in ${\displaystyle \delta \psi }$

${\displaystyle i\hbar {\frac {\partial \delta \psi }{\partial t}}=-{\frac {\hbar ^{2}}{2m}}\nabla ^{2}\delta \psi +V\delta \psi +g(2|\psi _{0}|^{2}\delta \psi +\psi _{0}^{2}\delta \psi ^{*})}$

${\displaystyle -i\hbar {\frac {\partial \delta \psi ^{*}}{\partial t}}=-{\frac {\hbar ^{2}}{2m}}\nabla ^{2}\delta \psi ^{*}+V\delta \psi ^{*}+g(2|\psi _{0}|^{2}\delta \psi ^{*}+(\psi _{0}^{*})^{2}\delta \psi )}$

Assuming the following for ${\displaystyle \delta \psi }$

${\displaystyle \delta \psi =e^{-i\mu t}(u({\boldsymbol {r}})e^{-i\omega t}-v^{*}({\boldsymbol {r}})e^{i\omega t})}$

one finds the following coupled differential equations for ${\displaystyle u}$ and ${\displaystyle v}$ by taking the ${\displaystyle e^{\pm i\omega t}}$ parts as independent components

${\displaystyle \left(-{\frac {\hbar ^{2}}{2m}}\nabla ^{2}+V+2gn-\hbar \mu -\hbar \omega \right)u-gnv=0}$

${\displaystyle \left(-{\frac {\hbar ^{2}}{2m}}\nabla ^{2}+V+2gn-\hbar \mu +\hbar \omega \right)v-gnu=0}$

For a homogeneous system, i.e. for ${\displaystyle V({\boldsymbol {r}})=const.}$, one can get ${\displaystyle V=\hbar \mu -gn}$ from the zeroth order equation. Then we assume ${\displaystyle u}$ and ${\displaystyle v}$ to be plane waves of momentum ${\displaystyle {\boldsymbol {q}}}$, which leads to the energy spectrum

${\displaystyle \hbar \omega =\epsilon _{\boldsymbol {q}}={\sqrt {{\frac {\hbar ^{2}|{\boldsymbol {q}}|^{2}}{2m}}\left({\frac {\hbar ^{2}|{\boldsymbol {q}}|^{2}}{2m}}+2gn\right)}}}$

For large ${\displaystyle {\boldsymbol {q}}}$, the dispersion relation is quadratic in ${\displaystyle {\boldsymbol {q}}}$ as one would expect for usual non interacting single particle excitations. For small ${\displaystyle {\boldsymbol {q}}}$, the dispersion relation is linear

${\displaystyle \epsilon _{\boldsymbol {q}}=s\hbar q}$

with ${\displaystyle s={\sqrt {ng/m}}}$ being the speed of sound in the condensate, also known as second sound. The fact that ${\displaystyle \epsilon _{\boldsymbol {q}}/(\hbar q)>s}$ shows, according to Landau's criterion, that the condensate is a superfluid, meaning that if an object is moved in the condensate at a velocity inferior to s, it will not be energetically favorable to produce excitations and the object will move without dissipation, which is a characteristic of a superfluid. Experiments have been done to prove this superfluidity of the condensate, using a tightly focused blue-detuned laser. [10] The same dispersion relation is found when the condensate is described from a microscopical approach using the formalism of second quantization.

The optical potential well ${\displaystyle V_{\rm {twist}}(\mathbf {r} ,t)=V_{\rm {twist}}(z,r,\theta ,t)}$ might be formed by two counter propagating optical vortices with wavelengths ${\displaystyle \lambda _{\pm }=2\pi c/\omega _{\pm }}$, effective width ${\displaystyle D}$ and topological charge ${\displaystyle \ell }$ :

${\displaystyle {E_{\pm }}(\mathbf {r} ,t)\sim \exp \left(-{\frac {r^{2}}{2D^{2}}}\right)r^{|\ell |}\exp(-i\omega _{\pm }t\pm ik_{\pm }z+i\ell \theta ),}$

where ${\displaystyle \delta \omega =(\omega _{+}-\omega _{-})}$.In cylindrical coordinate system ${\displaystyle (z,r,\theta )}$ the potential well have a remarkable double helix geometry: [11]

${\displaystyle V_{\rm {twist}}(\mathbf {r} ,t)\sim V_{0}\exp \left(-{\frac {r^{2}}{D^{2}}}\right)r^{2|\ell |}\left(1+\cos[\delta \omega t+(k_{+}+k_{-})z+2\ell \theta ]\right),}$

In a reference frame rotating with angular velocity ${\displaystyle \Omega =\delta \omega /2\ell }$, time-dependent Gross–Pitaevskii equation with helical potential is as follows: [12]

${\displaystyle i\hbar {\frac {\partial \Psi (\mathbf {r} ,t)}{\partial t}}=\left(-{\frac {\hbar ^{2}}{2m}}\nabla ^{2}+V_{\rm {twist}}(\mathbf {r} )+g\vert \Psi (\mathbf {r} ,t)\vert ^{2}-\Omega {\hat {L}}\right)\Psi (\mathbf {r} ,t),}$

where ${\displaystyle {\hat {L}}=-i\hbar {\frac {\partial }{\partial \theta }}}$ is the angular momentum operator. The solution for condensate wavefunction ${\displaystyle \Psi (\mathbf {r} ,t)}$ is a superposition of two phase-conjugated matter-wave vortices:

${\displaystyle \Psi (\mathbf {r} ,t)\sim \exp \left(-{\frac {r^{2}}{2D^{2}}}\right)r^{|\ell |}\times }$

${\displaystyle \left(\exp(-i\omega _{+}t+ik_{+}z+i\ell \theta )+\exp(-i\omega _{-}t-ik_{-}z-i\ell \theta )\right).}$

The macroscopically observable momentum of condensate is :

${\displaystyle \langle \Psi \vert {\hat {P}}\vert \Psi \rangle =N_{\rm {at}}\hbar (k_{+}-k_{-}),}$

where ${\displaystyle N_{\rm {at}}}$ is number of atoms in condensate. This means that atomic ensemble moves coherently along ${\displaystyle z-}$ axis with group velocity whose direction is defined by signs of topological charge ${\displaystyle \ell }$ and angular velocity ${\displaystyle \Omega }$: [13]

${\displaystyle V_{z}={\frac {2\Omega \ell }{(k_{+}+k_{-})}}}$

The angular momentum of helically trapped condensate is exactly zero:[12]

${\displaystyle \langle \Psi \vert {\hat {L}}\vert \Psi \rangle =N_{\rm {at}}[\ell \hbar -\ell \hbar ]=0.}$

Numerical modeling of cold atomic ensemble in spiral potential have shown the confinement of individual atomic trajectories within helical potential well.[14]

Vortex dipole trap with topological charge ${\displaystyle \ell =2}$ loaded by ultracold ensemble.