It's not clear to me, what's meant by this claim. The issue with the ideal Bose gas is as follows. To make it clear, you have to start with a finite volume. The most conventient one is to take a cube of length ##l## with periodic boundary conditions. Then the one-particle basis can be chosen as the momentum-spin basis. From your description I also assume that you deal with non-relativistic bosons.
The single-particle basis is ##|\vec{p},\sigma \rangle## with ##\vec{p}\in \frac{2 \pi}{l} \mathbb{Z}^3## and ##\sigma \in \{-s,\ldots,s \}##. The single-particle dispersion relation reads ##E=\frac{\vec{p}^2}/(2m)##.
Now we calculate the grand-canonical potential
$$Z=\mathrm{Tr} \exp \left [-\frac{\hat{H}-\mu \hat{N}}{T} \right ].$$
The trace is over all boson-Fock states, i.e., after some algebra you get
$$\ln Z=-\sum_{\vec{p}} \sum_{\sigma} \ln[1-\exp(-(\vec{p}^2/(2m)-\mu)/T)].$$
Now as long as you keep the volume finite, it's clear that ##T > 0## and ##\mu \leq 0##, because otherwise ##\ln Z## would get complex.
Now the temperature and chemical potential are used to adjust the mean total energy and particle number. To get these quantities it's more convenient to introduce the new independent variables ##\beta=1/T## and ##\alpha=\beta \mu=\mu/T##. Then you have
$$\ln Z=-\sum_{\vec{p}} \sum_{\sigma} \ln[1-\exp(-\beta \vec{p}^2/(2m)+\alpha)]$$
and
$$\langle E \rangle=-\partial_{\beta} \ln Z=(2s+1) \sum_{\vec{p}} \frac{\vec{p}^2}{2m} f_{\text{B}}(E-\mu), \quad \langle N \rangle = \partial_{\alpha} \ln Z=\sum_{\vec{p}} f_{\text{B}}(E-\mu),$$
where
$$f_{\text{B}}(E-\mu)=\frac{1}{\exp[(E-\mu)/T]-1]}, \quad E=\frac{\vec{p}^2}{2m}.$$
Now you can get any value ##\langle{E} \rangle \geq 0## and any ##\langle N \rangle \geq 0## if you have ##T>0## and ##\mu<0##. For ##T \rightarrow 0##, what happens is that all particles must occupy the ground state with ##E=0##, and indeed the smaller you make ##T## the larger you must make ##\mu## to keep the average particle number at the fixed given value.
Now make the box size ##l## very large. Then you can put the sum to an integral by noting that in any momentum-volume element ##\mathrm{d}^3 \vec{p}## there are ##l^3 \mathrm{\mathrm{d}^3 \vec{p}}{(2 \pi)^3}## states and thus the mean particle number
$$\langle N \rangle=V \int_{\mathbb{R}^3} \frac{\mathrm{d}^3 \vec{p}}{(2 \pi)^3} f_{\text{B}}(E-\mu).$$
Now, because ##\mathrm{d}^3 \vec{p}=4 \pi p^2 \mathrm{d} p## the integral stays finite for ##\mu \rightarrow 0## and becomes also arbitrarily small when taking ##T \rightarrow 0##, i.e., for a given fixed particle number at one point you cannot choose ##\mu \leq 0## anymore such that you get the correct average particle number, but from our finite-box discussion, it's clear what happens: there will be a finite number of particles occupying the ground state, and in the limit of infinite volume, this means you must take ##\mu=0## and write the distribution function as
$$f(\vec{p})=(2 \pi)^3 N_{\text{BEC}} \delta^3(\vec{p}) + f_{\text{B}}(E) \quad \text{if} \quad T< T_c.$$
On the other hand, if you can fix the particle number at the given temperature by choosing a chemical potential ##\mu \geq 0## you have to set
$$f(\vec{p})=f_{\text{B}}(E-\mu) \quad \text{if} \quad T>T_c.$$
The former case is known as Bose-Einstein condensation and was predicted to happen by Einstein in the early 1920ies. The experimental discovery was in 1995 (Nobel prize 2001).
