Confusion regarding the use of partition function Suppose we have a system filled with $N$ particles. There are $k$ energy levels in this system, labeled by $\epsilon_i$, each with a degeneracy of $g_i$. Let us imagine $n_j$ particles out of these $N$ particles occupy the energy level $\epsilon_j$.
The single particle partition function $Z_{sp}$ is given by : $Z_{sp}=\sum g_ie^{-\beta\epsilon_i}$
Let $Z_N$ be the multi-particle partition function. Its expression is a little more complicated, and can be given by: $Z_N=\sum G_ie^{-\beta E_i}$.
Here $E_i$ is the total energy of all the $N$ particles and we are checking over all the combinations. Similarly, $G_i$ is the degeneracy of these total energies.
Here is my question :
What is the difference between asking how many particles are in the energy level $\epsilon_j$, and the probability of the system to have a total energy $E_m$?
The second question is simple, and the answer is given by : $$P(E_m)=\frac{G_me^{-\beta E_m}}{Z_N}$$
However, how do I find out the answer to the first question, i.e. find $n_j$, the number of particles in the $\epsilon_j$ energy level ?
According to Wikipedia, the number of particles in a particular energy level is nothing but the probability of a single particle to occupy that level, multiplied by the total number of particles. Since the probability of a single particle being in energy level $\epsilon_j$ is given using the single-particle partition function, our final answer comes to be :
$$n_j=\frac{g_j e^{-\beta \epsilon_j}}{Z_{sp}}N$$
Is this correct?
Like for a system with only one particle, the probability of that one particle to be in a particular energy level is the same as the probability of the entire system having that energy. However, for a system with $N$ particles, it seems that the probability of a particle being in a particular state is very different from the probability of the entire system having a certain total energy.
Are the statements all correct ?
 A: It sounds like the main thing that would convince you is the check that I recommended in the last question :). I hedged my bets there by saying that the formula I gave for $N_j$ would agree with Wikipedia at large $N$. But it actually holds for general $N$.
The expected number of particles at level $j$ is found by summing $n_j$ over all microstates with their Boltzmann factor probabilities as weights. So again, that means
\begin{equation}
N_j = Z_N^{-1} \sum_{n_1 + \dots + n_k = N} n_j G(n_1, \dots, n_k) e^{-\beta (\epsilon_1 n_1 + \dots \epsilon_k n_k)}.
\end{equation}
But notice that this sum looks exactly like the the full partition function (which we know to be $Z_{\mathrm{sp}}^N$) differentiated with respect to $\epsilon_j$! Therefore,
\begin{align}
N_j &= -\beta^{-1} Z_N^{-1} \frac{\partial}{\partial \epsilon_j} Z_N \\
&= -\beta^{-1} Z_{\mathrm{sp}}^{-N} \frac{\partial}{\partial \epsilon_j} Z_{\mathrm{sp}}^N \\
&= -N\beta^{-1} Z_{\mathrm{sp}}^{-1} \frac{\partial}{\partial \epsilon_j} Z_{\mathrm{sp}} \\
&= N Z_{\mathrm{sp}}^{-1} g_j e^{-\beta \epsilon_j}
\end{align}
and everything checks out.
