Fermionic System Trouble I'm  trying to solve a question from a test and i don't understand the reasoning behind the solution
The question:
what is $\frac{P(occupied)}{P(unoccupied)}$ ?
Given a fermionic system with a chemical potential $\mu$. whats the relation between the probability of finding a state occupied by a particle with energy $\mu + \delta$ and the probability of an unoccupied state with energy $\mu - \delta$?
My try:
using the Grand canonical ensemble -
$$ P(state) \propto e^{-\beta(\epsilon_s -\mu)}$$
then: $$\frac{P(occupied)}{P(unoccupied)} = \frac{e^{-\beta(\mu + \delta -\mu)}}{e^{-\beta(\mu - \delta -\mu)}} = \frac{1}{e^{2\beta\delta}}$$
The solution
Using the fermi-dirac distribution:
$$f_{F D}(\mu+\delta ; \tau, \mu)=\frac{1}{e^{\beta \delta}+1}=\frac{e^{-\beta \delta}}{1+e^{-\beta \delta}}=1-\frac{1}{e^{-\beta \delta}+1}=1-f_{F D}(\mu-\delta ; \tau, \mu)
$$
and so:
$$
\frac{p_{\text {occupied }}(\mu+\delta)}{p_{\text {unoccupied }}(\mu-\delta)}=\frac{f_{F D}(\mu+\delta ; \tau, \mu)}{1-f_{F D}(\mu-\delta ; \tau, \mu)}=\frac{\frac{1}{e^{\beta \delta}+1}}{1-\frac{1}{e^{-\beta \delta}+1}}=1
$$
Why are we using fermi-dirac for the solution when it itself is based on the Grand canonical ensemble? and why the solutions don't match?
 A: While the probability of being occupied compared
to unoccupied is proportional to $e^{-\beta(\epsilon_s-\mu)}$, the
proportionality constants are not the same for your two states.
That is writing the proportionality constant as $C$, the probability of
the state with energy $\epsilon_s$ being unoccupied is
$C$ and the probability of that state being occupied is
$Ce^{-\beta(\epsilon_s-\mu)}$. The total probability is
one, $C(1+e^{-\beta(\epsilon_s-\mu)}) = 1$. So the probability of being
occupied is
\begin{equation}
P(state) = \frac{e^{-\beta(\epsilon_s-\mu)}}{1+e^{-\beta(\epsilon_s-\mu)}}
= \frac{1}{1+e^{\beta(\epsilon_s-\mu)}} \,.
\end{equation}
Properly normalizing your grand canonical probabilities gives
the correct Fermi-Dirac result.
Added to address the question in the comments:
The relative probability for the many-particle microstate $S$
of the total system
with energy
$E_S$ and particle number $N_S$
is given by the grand-canonical Boltzmann factor $e^{-\beta (E_S-\mu N_S)}$.
This has a common normalization. So if you are calculating the ratio of
the probability of finding single-particle state $s_1$ occupied you must
sum over all the occupations of all the other single-particle
states, since you do not
measure them. Similarly to calculate the probability of finding $s_2$
unoccupied you must sum over all the other single-particle states.
So you can either calculate the normalized probability for each of those
as above, or you can calculate the ratio with the common normalization.
Since all the state sums except $s_1$ and $s_2$ will cancel, you will
find
\begin{equation}
\frac{P_{\rm occupied}(s_1)}{P_{\rm unoccupied}(s_2)} =
\frac{(1+e^{-\beta(\epsilon_{s_2}-\mu)})e^{-\beta(\epsilon_{s_1}-\mu)}}
{(1+e^{-\beta(\epsilon_{s_1}-\mu)})}
\end{equation}
where the numerator is  $e^{-\beta(E_s-N_s\mu)}$ summed over all the
occupations of all the states except $s_1$ which is occupied,
and the denominator
is  $e^{-\beta(E_s-N_s\mu)}$ summed over all the occupations of all the states except $s_2$ which
is unoccupied.
The terms not containing $s_1$ or $s_2$ cancel and I did not write them.
The result is the
same as the Fermi-Dirac result, since we derived the normalized Fermi-Dirac
distribution to not have to repeat doing this same many-occupation number
sum.
