statistical mechanics - How to obtain the Bose-Einstein distribution from the canonical ensemble?

In the wikipedia page

https://en.wikipedia.org/wiki/Fermi%E2%80%93Dirac_statistics#Canonical_ensemble

the Fermi-Dirac distribution is obtained from the canonical ensemble in the following manner:

The average occupation number is given by

$$\begin{equation} \bar{n}_{i}=\frac{\displaystyle \sum_{n_1,n_2,...}n_{i}\,e^{-\beta(n_1E_1+n_2E_2+...)}}{\displaystyle \sum_{n_1,n_2,...}e^{-\beta(n_1E_1+n_2E_2+...)}} \end{equation}$$

where the summation is over the sets of $n_r$ that satisfy $\sum\,n_r=N$, being $N$ the total number of particles in the system. The derivation followed in Wikipedia then separates the state of energy $E_i$ and calls

$$\begin{equation} \bar{n}_{i}=\frac{\displaystyle \sum_{n_i} n_{i}\,e^{-\beta\, n_iE_i} \, Z_i(N-n_i)}{\displaystyle \sum_{n_i}\,e^{-\beta\, n_iE_i} \, Z_i(N-n_i)} \end{equation}$$

where $Z_{i}(N-n_i) = \sum_{n_1,n_n,...} \,e^{-\beta\,(n_1E_1+n_2E_2+...)}$. The summation here does not consider the state $n_i$. Since for the Fermi-Dirac case $n_i=0,1$, it is possible to write the summation easily (it has only two terms) to obtain

$$\begin{equation} \bar{n}_{i}=\frac{\displaystyle 1}{\displaystyle \bigg[Z_i(N)/Z_i(N-1)\bigg]\,e^{-\beta\,E_i} +1} \end{equation}$$

which will, in turn, result in the Dirac-Fermi distribution once the ratio between the partition functions is related to the chemical potential.

Now, I'm trying to follow the same procedure to obtain the Bose-Einstein distribution but in this case $n_i$ can vary between $0$ and $\infty$ and the summation cannot be done in any direct way I could think of. How can I approach this?

Thank you very much.

Answer

The expression:

$$\begin{equation} \bar{n}_{i}=\frac{\displaystyle \sum_{n_i} n_{i}\,e^{-\beta\, n_iE_i} \, Z_i(N-n_i)}{\displaystyle \sum_{n_i}\,e^{-\beta\, n_iE_i} \, Z_i(N-n_i)} \end{equation}$$ can be rewrite in a usual for statphysics manner: $$\bar{n}_{i}=-\frac{1}{\beta}\frac{\partial\ Ln {\ Z(N)}}{\partial E_i}$$ We denoted denominator as $Z(N)$ -canonical partition function of the system (do not confuse it with $Z_i(N)$ - partition function with "punctured" i-th one particle state. Then, denoting $z_i=e^{ \beta \mu}$ -"punctured" fugacity, we can write: $$Z(N)=\sum_{n_i=0}^ {n_i=N}\ e^{-\beta\, n_iE_i} \, Z_i(N-n_i) = \sum_{n_i=0}^ {n_i=N}\ (e^{-\beta\,E_i})^{n_i} (\frac{1}{z_i})^{N-n_i} =$$ $$=z_i^{-N}\sum_{n_i=0}^ {n_i=N}\ (z_ie^{-\beta\,E_i})^{n_i}$$ Then, sequentially we write down the sum formulae for a finite geometric series, take the logarithm, tend $N$ to infinity, and differentiate by $E_i$. Finally, acting in the wiki spirit, we replace $\mu_i$ =$\mu$

PS. After I wrote this answer, I looked at the link in Wikipedia, it turns out the derivation of the B-E distribution by a similar method is in the same link ((Reif 1965) on p. 342