How Is the Fermi–Dirac Distribution Derived?

The Fermi–Dirac distribution gives the average occupation of a one-particle state for a system of fermions in thermal equilibrium. It applies to particles such as electrons, which obey the Pauli exclusion principle: each available one-particle state can contain either zero or one particle.

The derivation below counts the number of possible arrangements at each energy, then finds the most probable distribution while holding the average energy and particle number fixed. This is equivalent to the equilibrium result from the grand canonical ensemble.

Count the fermion arrangements

Suppose energy level $i$ has energy $\epsilon_i$ and contains $g_i$ distinct one-particle states. If $N_i$ fermions occupy those states, the number of ways to choose which states are occupied is

$$W_i=\binom{g_i}{N_i}.$$

The total number of arrangements across all energy levels is the product

$$W=\prod_i \binom{g_i}{N_i}.$$

This counting is where the exclusion principle enters. Because a state can be occupied at most once, the choices are combinations of $N_i$ states from $g_i$, not unrestricted assignments of particles to states.

For large numbers of states and particles, use Stirling’s approximation, $\ln n!\approx n\ln n-n$. Taking the logarithm of the multiplicity gives

$$\ln W \approx \sum_i \left[g_i\ln g_i-N_i\ln N_i-(g_i-N_i)\ln(g_i-N_i)\right].$$

Maximize subject to fixed energy and particle number

The distribution must satisfy the constraints

$$\sum_i N_i=N,\qquad \sum_i \epsilon_iN_i=E,$$

where $N$ and $E$ are the total particle number and energy. Maximize $\ln W$ subject to these constraints using Lagrange multipliers $\alpha$ and $\beta$. For each level, differentiation with respect to $N_i$ gives

$$\frac{\partial}{\partial N_i}\left(\ln W-\alpha\sum_jN_j-\beta\sum_j\epsilon_jN_j\right) =\ln\left(\frac{g_i-N_i}{N_i}\right)-\alpha-\beta\epsilon_i=0.$$

Rearranging,

$$\frac{N_i}{g_i}=\frac{1}{\exp(\alpha+\beta\epsilon_i)+1}.$$

The ratio $N_i/g_i$ is the average occupation per one-particle state at energy $\epsilon_i$. In thermal equilibrium, $\beta=1/(kT)$ and $\alpha=-\mu/(kT)$, where $T$ is absolute temperature, $k$ is Boltzmann’s constant, and $\mu$ is the chemical potential. Substitution yields the Fermi–Dirac distribution:

$$f(\epsilon_i)=\frac{N_i}{g_i}=\frac{1}{\exp\left((\epsilon_i-\mu)/(kT)\right)+1}.$$

The chemical potential adjusts to give the required average number of particles. For a system with a fixed particle number, it is generally determined from the condition $\sum_i g_i f(\epsilon_i)=N$.

Worked example

Consider a one-particle state with energy $\epsilon=\mu+kT\ln 3$. Its energy is above the chemical potential, so

$$f(\epsilon)=\frac{1}{\exp(\ln 3)+1}=\frac{1}{3+1}=\frac14.$$

Thus, in equilibrium, that state has a 25% average occupation probability. At $\epsilon=\mu$, the exponent is zero and $f=1/2$. A state well below $\mu$ is more likely to be occupied, while one well above $\mu$ is less likely to be occupied.

Why the same formula appears in the grand canonical ensemble

In that ensemble, a single fermion state of energy $\epsilon$ has two possible occupations: $n=0$ or $n=1$. Its grand partition function is

$$Z=1+\exp[-(\epsilon-\mu)/(kT)].$$

The probability that the state is occupied is its occupied-state statistical weight divided by $Z$:

$$P(n=1)=\frac{\exp[-(\epsilon-\mu)/(kT)]}{1+\exp[-(\epsilon-\mu)/(kT)]} =\frac{1}{\exp[(\epsilon-\mu)/(kT)]+1}.$$

Since $n$ can only be zero or one, its expected occupation equals $P(n=1)$. This gives the same formula directly. The counting derivation describes the most probable distribution for a large system; the grand canonical derivation gives the mean occupation of an individual state.

Common mistakes and limits

• Confusing a level with a state: If an energy level is degenerate, $g_i$ counts its distinct states. The formula gives occupation per state; the expected number of particles in the whole level is $g_i f(\epsilon_i)$.
• Dropping the plus one: The denominator’s $+1$ distinguishes fermions from the Bose–Einstein distribution, whose denominator has a $-1$.
• Using the wrong sign: The exponent is $(\epsilon-\mu)/(kT)$. Occupation must decrease as energy rises above the chemical potential.
• Treating it as a universal interacting-particle formula: The basic derivation assumes independent states, or a setting where particles can be treated as independent quasiparticles. Strong interactions can require a more detailed model.

As temperature approaches zero, the distribution approaches a sharp step: states below the chemical potential are occupied and those above it are empty. At nonzero temperature, the transition is smooth.