Bose condensation in a three-dimensional harmonic trap (source code)

= Bose condensation in a three-dimensional harmonic trap
{c}
{title2=$k_BT_c=\hbar\bar\omega[N/\zeta(3)]^{1/3}$}

Consider $N$ noninteracting <bosons> in one internal <quantum state>, at <thermal equilibrium> in the three-dimensional <harmonic oscillator> potential
$$
V(\mathbf r)=\frac m2(\omega_x^2x^2+\omega_y^2y^2+\omega_z^2z^2),\qquad
\bar\omega=(\omega_x\omega_y\omega_z)^{1/3}.
$$
Here $m$ is the <mass>, the positive $\omega_i$ are <angular frequencies>, and $\bar\omega$ is their <geometric mean>. Measure the one-particle <energy> $\epsilon$ and <chemical potential> $\mu$ from the <ground-state energy>. The <quantum harmonic oscillator> levels then have $\epsilon=\hbar(\omega_x n_x+\omega_y n_y+\omega_z n_z)$ with nonnegative integer $n_i$ and reduced <Planck constant> $\hbar=h/(2\pi)$.

When $\epsilon\gg\hbar\omega_i$, counting <lattice points> can be approximated by the volume of a <tetrahedron> in the positive octant of $(n_x,n_y,n_z)$:
$$
\mathcal N(\epsilon)\simeq\frac{\epsilon^3}{6\hbar^3\omega_x\omega_y\omega_z},\qquad
g(\epsilon)=\frac{d\mathcal N}{d\epsilon}\simeq\frac{\epsilon^2}{2(\hbar\bar\omega)^3}.
$$
This is a <density of states>, with units of inverse <energy>. The factor $1/6$ is the unit tetrahedron volume; differentiating produces the factor $1/2$. Treat the single <ground state> separately: its macroscopic population would be lost by replacing every level with a smooth <density of states>.

For $k_BT\gg\hbar\omega_i$, integrate the <Bose-Einstein distribution> over excited <energy eigenstates>. As $\mu$ approaches zero from below, their maximum population is
$$
N_{\rm ex}^{\rm sat}(T)=\int_0^\infty\frac{g(\epsilon)\,d\epsilon}{e^{\epsilon/(k_BT)}-1}
=\zeta(3)\left(\frac{k_BT}{\hbar\bar\omega}\right)^3.
$$
Here $k_B$ is the <Boltzmann constant> and $\zeta$ the <Riemann zeta function>. A short proof uses $1/(e^x-1)=\sum_{\ell\ge1}e^{-\ell x}$ and $\int_0^\infty x^2e^{-\ell x}dx=2/\ell^3$. The nonnegative terms permit exchanging <summation> and <integration> by <Tonelli's theorem>, giving $\int_0^\infty x^2/(e^x-1)dx=2\zeta(3)$. Near zero the integrand is proportional to $x$, so its integral converges. That finite excited-state capacity is the mechanism of <Bose-Einstein condensation>.

Set $N_{\rm ex}^{\rm sat}(T_c)=N$. Particle conservation gives, to this continuum accuracy,
$$
k_BT_c=\hbar\bar\omega\left(\frac{N}{\zeta(3)}\right)^{1/3},\qquad
\frac{N_0}{N}=1-\left(\frac{T}{T_c}\right)^3\quad(0\le T\le T_c).
$$
Above $T_c$, the macroscopic <condensate fraction> is zero and the <chemical potential> adjusts below zero. This cubic power belongs to a three-dimensional harmonic trap; it is not the homogeneous three-dimensional <ideal gas> exponent. For $N=10^5$, $\bar\omega/(2\pi)=100\,\mathrm{Hz}$ and $\zeta(3)\simeq1.202056903$, the model gives $T_c\simeq209.5\,\mathrm{nK}$. At $T=0.8T_c\simeq167.6\,\mathrm{nK}$, $N_0/N=0.488$, or approximately $48\,800$ condensed <bosons>. These are chosen parameters and calculated values.

Finite $N$ gives a smooth crossover rather than a sharp <phase transition>. The exact <ideal gas> sum over excited <quantum harmonic oscillator> levels differs from the continuum integral, especially when $k_BT$ is not large relative to every $\hbar\omega_i$. A sharp trapped-gas limit takes $N\to\infty$, $\bar\omega\to0$ with $N\bar\omega^3$ fixed and fixed frequency ratios. Interactions also alter the <critical temperature>, <condensate fraction> and release profile; a dilute gas can have appreciable interaction <energy> compared with the trap level spacing. The cubic law is a useful baseline, not a fit guaranteed for either 1995 experiment.