Quantum pressure of a Bose gas in a harmonic trap: Why is the result divergent?

Question

The quantum pressure $E_{kin}$ of a Bose gas in the Thomas Fermi limit with contact interactions in a symmetric harmonic trap is determined by (Pitaevskij & Stringari, 2016):

$E_{kin}=\int d\mathbf{r}\frac{\hbar^2}{2m}\left|\nabla\sqrt{n}\right|^2$

In that same harmonic trap, the density profile is

$n_{TF}^0(\mathbf{r})=\frac{1}{g}\left(\mu_{TF}^0-V_{ext}(\mathbf{r})\right)$,

where $V_{{ext}}(\mathbf{r})=\frac{m\omega}{2}\left(x^2+y^2+z^2\right)$ is the trapping potential.

The following are derived constants;

• $g=\frac{4\pi\hbar^2a_s}{m}$ is the interaction parameter,
• $\mu_{TF}^0=\frac{\hbar\omega}{2}\left(\frac{15Na_s}{a_{ho}}\right)^{2/5}$ is the chemical potential,
• $a_{ho}=\sqrt{\hbar/m\omega}$ is the harmonic oscillator length,
• $R=a_{h o}\left(\frac{15 N a}{a_{h o}}\right)^{1 / 5}$ is the radius of the condensate.

Other terms in the equations are independent experimental constants.

---

To solve the integral, we can convert to spherical coordinates with $d \mathbf{r}=r^2 \sin (\varphi) d r d \theta d \varphi$, where $r=\sqrt{x^2+y^2+z^2}$. The limits of integration range from $0$ to $2\pi$ for $\theta$, $0$ to $\pi$ for $\varphi$, and $0$ to $R$ for $r$, where R is the radius of the condensate.

We get an integral that looks like this:

$E_{kin}=\int_0^{2\pi}d\theta\int_0^{\pi}sin(\varphi)d \varphi\int_0^{R} r^2\frac{\hbar^2}{2m}\left|\nabla\sqrt{n}\right|^2 dr$

The angular components can be separated from the rest of the integral since the density profile does not depend on the angle.

Let's expand that term near the end, $\left|\nabla\sqrt{n}\right|^2$. We get

$E_{kin}=\int_0^{2\pi}d\theta\int_0^{\pi}sin(\varphi)d \varphi\int_0^{R} r^2\frac{\hbar^2}{2m}\left(\frac{1}{2\pi R^2a_{ho}^4a_s}\frac{r^2}{1-r^2/R^2}\right) dr$.

The problem is that the $\frac{1}{1-r^2/R^2}$ part will blow up near $R$, which is where the limit of integration is. Putting the thing into wolfram alpha returns the message that the result does not converge, and solving it as an indefinite integral just returns an infinite result when you substitute $r=R$ into it.

How can this problem be resolved? It should be a textbook case but it is returning a divergent result?

---

Reference:

Pitaevskij, L.P. and Stringari, S. (2016) Bose-Einstein condensation and superfluidity. Oxford: Oxford University Press.

Answer 1

On page 170 of the book there is an explicit formulation that in the Thomas-Fermi limit $n^0_{TF}({\bf r})$ is obtained from the Gross-Pitaevsky equation (11.1) $$ \left(-\frac{\hbar^2}{2m}\nabla^2 + V_{ext}({\bf r}) + g|\Psi_0({\bf r})|^2 \right)\Psi_0({\bf r}) = \mu_0\Psi_0({\bf r}) $$ by ignoring the “kinetic energy term” and taking $\Psi_0({\bf r})= \sqrt{n^0_{TF}(\bf{r})}$. Therefore, I think that directly substituting $n^0_{TF}(\bf{r})$ into the integral expression for the quantum pressure $E_K$ makes little sense. Rather, you should think that $E_K\approx 0$. Indeed, as I see it, most of Chapter 11 of the book is devoted to considering the smallness of $E_K$. The smallness of the kinetic term in the GP equation is obviously violated near the boundary, where $n^0_{TF} \approx 0$, so additional consideration is required.

Answer 2

Mathematically integral of $\frac{1}{r^\alpha}$ diverges near zero for $\alpha > 1$ and diverges in infinity for $\alpha <$. For $\alpha=1$ it diverges in both cases.

The integral in question is of this type: $$\int\frac{r^4}{1-\left(r/R\right)^2}dr= \int\frac{r^4}{1-r/R}dr+ \int\frac{r^4}{1+r/R}dr $$ (for $r\rightarrow R$ the numerator is effectively constant $r^4\rightarrow R^4$.

Thus, either the derivation leading to this integral is incorrect or the physical model is not good.