# Derivation of density of bosons below Bose-Einstein condensation temperature

I am trying to understand the explanation of Bose-Einstein condensation for non-interacting bosons given in Piers Coleman's "Introduction to Many-Body Physics", pg. 85-86. Coleman first writes that the density of a gas of bosons at finite temperature is given by $$\rho = \int \frac{d^3 k}{(2\pi)^3} (e^{\beta (E_k - \mu)}-1)^{-1}, \tag{1}\label{1}$$ where we have $$E_k = \hbar^2 k^2/(2m)$$. He then points out that the chemical potential $$\mu$$ must be negative, presumably since this avoids a divergent integral yielding an infinite density, and that as the temperature is lowered ($$\beta$$ increases), $$\mu$$ must increase (become less negative) for a fixed particle density $$\rho$$. This is clear from a simple inspection of the integral---$$\mu(\beta)$$ increases monotonically (from a negative value) with $$\beta$$. As one increases $$\beta$$, at some point $$\mu(\beta) = 0$$. This is called the Bose-Einstein condensation temperature, and it can be computed in terms of $$\rho$$ by setting $$\mu = 0$$ in the integral and solving for $$\beta$$.

So far so good. Now I quote the part I don't understand:

Below this temperature, the number of bosons in the $$k = 0$$ state becomes macroscopic, i.e. $$n_{\epsilon = 0} = \frac{1}{e^{-\beta \mu} - 1} = N_0(T) \tag{2}\label{2}$$ becomes a finite fraction of the total particle number. Since $$N_0(T)$$ is macroscopic, it follows that $$\frac{\mu}{k_B T} = -\frac{1}{N_0(T)} \tag{3}\label{3}$$ is infinitesimally close to zero. For this reason, we must be careful to split off the $$k = 0$$ contribution to the particle density, writing $$N = N_0(T) + \sum_{k \neq 0} n_k \tag{4}\label{4}$$ and then taking the thermodynamic limit of the second term. For the density, this gives $$\rho = \frac{N}{V} = \rho_0(T) + \int \frac{d^3 k}{(2\pi)^3} \frac{1}{e^{\beta (E_k)} -1}. \tag{5}\label{5}$$ The second term is proportional to $$T^3$$. Since the first term vanishes at $$T = T_0$$, it follows that, below the Bose-Einstein condensation temperature, the density of bosons in the condensate is given by $$\rho_0(T) = \rho \left[1 - \left(\frac{T}{T_0}\right)^{3/2}\right].\tag{6}\label{6}$$

(a) I understand the use of the Bose-Einstein distribution to obtain the particle number in the $$k = 0$$ state. But at the condensation temperature $$T_0$$, $$\mu = 0$$, and this quantity becomes infinite. I know somehow this is going to be the point, but I don't understand precisely how: as we cross $$T = T_0$$ evidently something in our theory has broken down (for example, at $$T = T_0$$, $$\mu = 0$$, and we obtain an infinite particle density \eqref{1}). What actually limits the particle density from becoming infinite (what assumption in our theory breaks down here?)? And what does it mean to talk about what happens "below" $$T_0$$ within the equations we have written down? For example, within these equations, as $$T$$ goes below $$T_0$$, $$\mu$$ cannot cross zero and become positive, otherwise we get divergences. So in \eqref{2}, is $$\mu$$ positive or negative?

(b) How can one clean up the math leading to \eqref{5}, or at least give some reasonable explanation? I cannot make sense of this derivation. Can one take some integral, split it into an $$\epsilon$$-neighborhood of the origin and the remainder, evaluate, and take the limit $$\epsilon \to 0^+$$ to observe a finite contribution from the origin and the second term as written?

(c) In the last sentence, why does it say $$\rho_0(T)$$ vanishes at $$T = T_0$$? Isn't it exactly the opposite ($$\rho_0(T_0) = \infty$$, by dividing \eqref{2} by $$V$$ and sending $$\mu \to 0^-$$)?

In general, I think that your issue is in identifying how the different quantities scale with $$N$$ in the thermodynamic limit.

1. At finite $$N$$, there is no phase transition. In fact $$\mu$$ is a decreasing function in $$T$$, and always satisfies $$\mu<0$$. The subtlety is that it scales differently in $$N$$ below $$T_c$$ in the thermodynamic limit. Below $$T_c$$, the finite quantity is rather $$N\mu$$. In the thermodynamic limit:
• $$T_c\leq T<+\infty$$: $$\mu$$ is finite negative and increases from $$-\infty$$ to $$0$$ as $$T$$ decreases. Accordingly, $$N\mu=-\infty$$ in this range.
• $$0: $$\mu$$ is constant equal to $$0$$. However, according to (3) and (6), $$N\mu$$ increases from $$-\infty$$ to $$0$$ as $$T$$ increases.

You therefore never cross $$\mu=0$$. This is consistent with (2) as since $$N_0$$ becomes macroscopic, in the thermodynamic limit, $$N_0=\infty$$. Only the intensive quantity: $$\rho_0 = \frac{N_0}{V}$$ is finite for $$T\leq T_c$$, which is consistent with the fact that $$N\mu$$ has a finite limit.

1. The issue is that the Thomas-Fermi approximation breaks down. A unified way would be to use the discretised orbitals, in accordance to the first quantised particle in a box. In $$D$$ dimensions, using $$\hbar=m=1$$, $$L^D=N$$ (extensive volume): \begin{align} \rho N &= \sum_{k\in\mathbb Z^D}\frac{1}{e^{\beta (E_k-\mu)}-1} \\ E_k &= \pi\frac{2\pi}{N^{2/D}}\sum_{i=1}^D k_i^2 \end{align} You can check that the Thomas-Fermi approximation breaks down, and the term $$k=0$$ is not negligible anymore. Alternatively, you can also directly relate it to the usual formula with the polylogarithm $$\text{Li}_{D/2}$$ by resumming the expression in terms of number of particles: $$\rho N = \sum_{n=1}^\infty z^n\theta\left(\frac{2\pi\beta n}{N^{2/D}}\right)^D \\ z = e^{\beta\mu}$$ with $$\theta$$ being related to the Jacobi theta function: $$\theta(t) = \sum_{k\in\mathbb Z}e^{-\pi n^2t}$$ You can check that you recover your expression in the thermodynamic limit $$N\to\infty$$ using: $$\theta(t\to+\infty) = 1 \\ \theta(t^{-1}) = \sqrt t\theta(t) \\ \theta(t\to0^+) \sim t^{-1/2}$$ so by replacing each term by its limit in the sum: $$z^n\theta\left(\frac{2\pi\beta n}{N^{2/D}}\right)^D \sim N\frac{z^n}{(2\pi\beta n)^{D/2}} \\ \rho = \frac{1}{(2\pi\beta)^{D/2}}Li_{D/2}(z)$$ which is valid for $$\beta<\beta_c=\frac{\zeta(D/2)}{2\pi\rho^{2/D}}$$. Intuitively, this is not valid anymore because you also have the tail: $$z^n\theta\left(\frac{2\pi\beta n}{N^{2/D}}\right)^D \sim z^n$$ which will give the $$\frac{1}{1-z}$$. This contribution is negligible in the normal phase $$\beta<\beta_c$$ since $$z$$ has a finite limit and the polylogarithm contribution suffices. However, in the BEC phase, it is also important as $$z\to 1$$, so you get this additional term.

2. No, $$\rho_0$$ vanishes at $$T=T_c$$. You need to be careful with extensive and intensive quantities. $$n_0$$ diverges at the critical temperature in the thermodynamic limit. However, rescaling it by the volume, $$\rho_0$$ is intensive and has a finite value below $$T_c$$. You should view $$\rho_0$$ as the order parameter:

• $$T_c\leq T<+\infty$$: $$\rho_0=0$$, the normal phase
• $$0: $$\rho>0$$, the BEC phase. It increases according to (6) as $$T$$ decreases.

Hope this helps.

• I'm a bit confused by this answer, starting with (1). Don't you mean that at $T=T_c$, $\mu = 0$? And that $\mu$ increases monotonically with $T$, and varies from $-\infty$ to $0$? For the $T < T_c$ analysis, I simply do not understand which equation is my starting point. If $\mu = 0$ for $T < T_c$, then Eq. (2) is nonsensical; this equation is said to hold for $T < T_c$, so plugging in $\mu=0$ gives $\infty$. So I don't understand what comes after. Do I start with Eq. (1) plugging in $\mu = 0$? Then what? Jan 4 at 1:14
• getting $+\infty$ from (2) is normal in the thermodynamic limit, $N_0$ is macroscopic.
– LPZ
Jan 4 at 9:12