The Grand Canonical Ensemble and Chemical Potential

Derivation

The Grand Canonical Ensemble and Chemical Potential

D-159 Home PU-203 Threads chance · matter Depends on The Partition Function and F = -kT ln Z, The Boltzmann Distribution, Equilibrium as Free-Energy Minimization

Statement

A small system in thermal and diffusive contact with a large reservoir at temperature T and chemical potential μ occupies microstate s (energy E_s, particle number N_s) with probability P_s = e^{−β(E_s − μN_s)} / Ξ, where the grand partition function Ξ(T,V,μ) = Σ_s e^{−β(E_s − μN_s)} and β = 1/kT. The chemical potential μ = (∂F/∂N)_T,V is the free-energy cost of adding one particle at fixed temperature and volume.

Why it matters

The canonical ensemble fixes particle number, but most real systems — adsorbed layers, electrons in a metal, gas above a liquid, reacting mixtures, semiconductor carriers — freely exchange particles with a surrounding bath. The grand canonical ensemble is the natural language for all of these, and it is the only tractable route to quantum statistics: Fermi–Dirac and Bose–Einstein distributions fall out of a per-mode grand partition function almost trivially, whereas the fixed-N sum is a combinatorial nightmare.

Just as temperature is the intensive variable that equalises when energy flows, the chemical potential is the intensive variable that equalises when particles flow. Phase coexistence, chemical equilibrium, and osmosis are all statements that μ is uniform. Making μ the free-energy price of one particle turns these into one-line conditions.

Assumptions

The reservoir is much larger than the system in both energy and particle capacity.If dropped, the first-order Taylor expansion of the reservoir entropy is not enough; second derivatives (finite reservoir heat capacity and finite particle number) survive, T and μ drift as the system samples states, and no single Gibbs factor exists.

System and reservoir are weakly coupled, so total energy and total number are additive: E_tot = E_s + E_R, N_tot = N_s + N_R.If dropped, interaction energy at the boundary cannot be assigned to either subsystem; the reservoir multiplicity no longer depends only on E_tot−E_s and N_tot−N_s, and the factorised probability breaks down. Relevant for long-range forces and small systems.

The combined isolated supersystem is in equilibrium and obeys equal a priori probabilities over its accessible microstates.If dropped, we cannot equate the state probability to a reservoir microstate count, and the entire ensemble construction loses its microcanonical foundation.

Entropy is smooth and extensive in E and N at the reservoir's operating point.If dropped — e.g. the reservoir sits at a first-order phase transition — the derivatives (∂S/∂E) and (∂S/∂N) are ill-defined or discontinuous, and T, μ are not sharp.

Derivation

P_s ∝ Ω_R(E_tot − E_s, N_tot − N_s)

The system-plus-reservoir is isolated (microcanonical), so every joint microstate is equally likely. Fixing the system in state s leaves the reservoir with Ω_R ways to carry the remaining energy and particles; the probability is proportional to that count. A

P_s ∝ exp[ S_R(E_tot − E_s, N_tot − N_s) / k ]

Boltzmann's definition S_R = k ln Ω_R inverts to Ω_R = e^S_R/k. Working with the entropy lets us expand a slowly varying quantity instead of a wildly varying exponential-large multiplicity. A

S_R ≈ S_R(E_tot, N_tot) − E_s (∂S_R/∂E)_V,N − N_s (∂S_R/∂N)_E,V

Because the reservoir is huge, E_s ≪ E_tot and N_s ≪ N_tot; Taylor-expand to first order about the full reservoir. Higher-order terms scale as (system size / reservoir size) and vanish in the reservoir limit. B

(∂S_R/∂E)_V,N = 1/T (∂S_R/∂N)_E,V = −μ/T

These are the thermodynamic definitions read off dU = T dS − P dV + μ dN, rearranged as dS = (1/T)dU + (P/T)dV − (μ/T)dN. The reservoir imposes its T and μ on the system; note the minus sign on the number derivative. B

P_s ∝ exp[ −E_s/kT + μN_s/kT ] = e^{−β(E_s − μN_s)}

Insert step 4 into step 3, divide by k. The constant term S_R(E_tot, N_tot)/k is the same for every system state, so it is an overall factor absorbed into the normalisation. A

P_s = e^{−β(E_s − μN_s)} / Ξ, Ξ ≡ Σ_s e^{−β(E_s − μN_s)}

Probabilities must sum to one; the normalising denominator is defined as the grand partition function Ξ. The sum now runs over all microstates of every allowed particle number, not just a fixed N. A

Ξ = Σ_N=0^∞ e^βμN Σ_i(N) e^−βE_i = Σ_N=0^∞ z^N Z(N,V,T)

Group the microstate sum first by particle number N, then over the states i at that N. The inner sum is exactly the canonical partition function Z(N,V,T) (prior result), and z ≡ e^βμ is the fugacity. The grand partition function is a fugacity-weighted generating function of the canonical ones. B

Φ(T,V,μ) ≡ −kT ln Ξ = F − μN̅ , μ = (∂F/∂N)_T,V

Define the grand potential Φ as the log of Ξ, the analogue of F = −kT ln Z. It is the Legendre transform of F with respect to N: from dF = −S dT − P dV + μ dN the coefficient of dN at fixed T,V is μ, so μ is the free-energy paid per added particle. C

Result

Ξ(T,V,μ) = Σ_s e^{−β(E_s − μN_s)} = Σ_N z^NZ(N,V,T), Φ = −kT ln Ξ, μ = (∂F/∂N)_T,V

Reading. Opening the walls to particle flow replaces the Boltzmann weight e^−βE with e^{−β(E−μN)}: each particle the system holds is “rebated” by μ against its energy cost, because holding it denies the reservoir that particle. The reservoir sets the price μ; the system decides how many particles to buy. Averages follow by differentiating ln Ξ: N̅ = kT(∂ ln Ξ/∂μ)_T,V = −(∂Φ/∂μ)_T,V and U̅ − μN̅ = −(∂ ln Ξ/∂β)_βμ.

Units check. The exponent β(E − μN) must be dimensionless: β is 1/energy, E is energy, and μN is (energy)×(pure number), so μ carries units of energy — joules per particle, or eV. Hence Ξ is a pure number, ln Ξ is dimensionless, and Φ = −kT ln Ξ has units of energy, matching F − μN.

Limiting cases

Fixed particle number. If only one value N₀ is allowed, Ξ = z^N₀Z(N₀,V,T) and Φ = F − μN₀; the canonical ensemble is recovered exactly.
Dilute / classical limit (z → 0). Then Ξ ≈ 1 + zZ₁, N̅ ≈ zZ₁ ≪ 1; multiple occupancy is negligible and Fermi–Dirac and Bose–Einstein both collapse to Maxwell–Boltzmann.
Thermodynamic limit. Relative number fluctuations σ_N/N̅ ∼ N̅^−1/2 → 0, so the grand canonical and canonical predictions for intensive quantities coincide (ensemble equivalence).
Ideal quantum gas per mode. A single-particle level of energy ε gives Ξ_ε = Σ_n e^{−β(ε−μ)n}: n∈{0,1} yields Fermi–Dirac, n∈{0,1,2,…} yields Bose–Einstein.

Breaks when

The reservoir is not large. For a nanoscale bath the first-order entropy expansion (step 3) is inadequate; finite reservoir heat capacity and particle capacity make T and μ fluctuate, and the simple Gibbs factor is only an approximation.
Long-range or strong interactions couple system and reservoir. When surface/interaction energy is not negligible against bulk energy, E_tot ≠ E_s + E_R; the reservoir multiplicity no longer factorises and the ensemble is ill-defined (gravitating systems, small clusters).
At a first-order phase transition or condensation point. N̅(μ) becomes non-analytic, (∂N̅/∂μ) → ∞, and fluctuations diverge; the grand canonical average can differ from the canonical result (ensemble inequivalence, e.g. ideal Bose gas condensate number fluctuations).

Failure modes

Wrong sign on the chemical-potential term: writing e^−β(E+μN). The physics minus sign comes from (∂S/∂N)_E,V = −μ/T; adding particles the reservoir loses lowers S_R only for μ>0.
Expecting μ>0 for a classical gas. For a dilute ideal gas μ = kT ln(nλ³) with nλ³ ≪ 1, so μ<0; a positive μ signals quantum degeneracy.
Forgetting Ξ sums over N too. Treating it as one more canonical sum drops the fugacity weighting and gives wrong occupancies.
Confusing potentials: using F = −kT ln Z where Φ = −kT ln Ξ is required; Φ = F − μN, not F.
Double-counting identical particles. The canonical Z(N,V,T) inside Ξ already carries the 1/N! (or full quantum symmetrisation); omitting it inflates Ξ and corrupts μ.

Discussion

The chemical potential plays for particle exchange exactly the role temperature plays for energy exchange. Two systems that can trade energy reach equilibrium when their temperatures match; two systems that can trade particles equilibrate when their chemical potentials match. This single statement, μ₁ = μ₂, governs phase coexistence (liquid–vapour lines), chemical equilibrium (Σ ν_iμ_i = 0 for a reaction), membrane and osmotic equilibrium, and the alignment of Fermi levels when two metals touch. The grand canonical ensemble is where that condition is derived rather than assumed.

The reading of μ as a free-energy cost is worth internalising. Because dF = −S dT − P dV + μ dN, adding one particle at fixed T,V changes the Helmholtz free energy by μ. Equivalently μ = (∂U/∂N)_S,V = (∂G/∂N)_T,P; the last form shows μ is just the Gibbs free energy per particle, which is why for a pure substance G = μN. The three partial derivatives agree because each holds fixed the natural variables of its own potential.

Number fluctuations are physically observable, not an artefact of the ensemble: σ_N² = kT(∂N̅/∂μ)_T,V, and this connects directly to the isothermal compressibility κ_T via σ_N²/N̅ = (N̅kT/V)κ_T. Critical opalescence — the milky scattering of light near a critical point — is exactly the divergence of these density fluctuations. The ensemble that lets N float is the natural home for such phenomena.

Structurally, Φ(T,V,μ) is the double Legendre transform of the energy that trades the extensive pair (S,N) for the intensive pair (T,μ), keeping V. Its natural variables are all but one intensive, so by Euler's theorem on the extensive energy U = TS − PV + μN one gets the compact identity Φ = −PV. Thus ln Ξ = PV/kT directly delivers the equation of state: the grand partition function is, quite literally, the pressure. This is why virial expansions and quantum gas thermodynamics are almost always launched from Ξ rather than Z.

Common misconceptions. Chemical potential is not restricted to chemistry and has nothing intrinsically to do with reactions; it is defined for any conserved particle number, including photons (where μ=0 because photon number is not conserved). It is also not the potential energy of a particle — it is a free energy that includes an entropic term, which is why it can be negative even when every single-particle energy is positive.

Worked examples

Example 1 — occupancy of a single adsorption site (Langmuir isotherm).

Ξ = e⁰ + e^{−β(−ε − μ)} = 1 + e^β(ε+μ)

The site is empty (E=0, N=0) or holds one molecule with binding energy −ε (E=−ε, N=1). Two terms only. A

N̅ = e^β(ε+μ) / (1 + e^β(ε+μ)) = 1 / (1 + e^−β(ε+μ))

Weight the occupied state by N=1 and divide by Ξ; the empty state contributes zero to N̅. Note the Fermi–Dirac form. A

kT = (8.617×10⁻⁵ eV/K)(300 K) = 0.02585 eV, ε+μ = 0.20 + (−0.30) = −0.10 eV

Insert numbers: binding ε = 0.20 eV, gas reservoir at μ = −0.30 eV (dilute gas, μ<0), T = 300 K. A

β(ε+μ) = −0.10 / 0.02585 = −3.868, e^−β(ε+μ) = e^3.868 = 47.85

Evaluate the exponent, then the exponential needed in the denominator. A

N̅ = 1 / (1 + 47.85) = 0.0205 ≈ 2.0% occupancy

Units check. Both ε and μ in eV, kT in eV, ratio dimensionless; N̅ is a pure occupation probability between 0 and 1. Raising μ (denser gas) toward −ε drives the site toward half-filling, as a Langmuir isotherm should.

Example 2 — chemical potential of a classical ideal gas (argon at STP-like conditions).

Ξ = exp( zV/λ³ ) ⇒ N̅ = z(∂ ln Ξ/∂z) = zV/λ³

For the ideal gas Z(N) = (V/λ³)^N/N!, so Ξ = Σ_N(zV/λ³)^N/N! = e^zV/λ³. Differentiate ln Ξ to get the mean number. B

z = e^βμ = N̅λ³/V = nλ³ ⇒ μ = kT ln(nλ³)

Solve for the fugacity, then take the log. λ = h/√(2πmkT) is the thermal de Broglie wavelength, n = N̅/V the number density. B

n = P/kT = 101325 / (1.381×10⁻²³ × 300) = 2.446×10²⁵ m⁻³

Number density from the ideal-gas law at P = 1 atm, T = 300 K. A

λ = h/√(2πmkT), m = 39.95×1.6605×10⁻²⁷ = 6.634×10⁻²⁶ kg

Argon molar mass 39.95 g/mol converted to per-atom mass. A

2πmkT = 1.727×10⁻⁴⁵ ⇒ λ = 6.626×10⁻³⁴/4.156×10⁻²³ = 1.594×10⁻¹¹ m

Evaluate the square root in the denominator, then divide Planck's constant by it. A

λ³ = 4.05×10⁻³³ m³, nλ³ = 2.446×10²⁵ × 4.05×10⁻³³ = 9.91×10⁻⁸

The degeneracy parameter nλ³ ≪ 1 confirms the classical (Maxwell–Boltzmann) treatment is valid. A

μ = kT ln(9.91×10⁻⁸) = 0.02585 eV × (−16.13) = −0.417 eV

Units check. nλ³ is dimensionless (m⁻³ × m³), so ln(nλ³) is a pure number and μ inherits the energy units of kT. The result is negative, as required for a dilute classical gas — the entropic gain of spreading a new particle over the huge phase-space volume outweighs any energy cost.

Problems

(Easy) A single adsorption site binds one molecule with binding energy ε = 0.15 eV. It is in contact with a gas at T = 250 K and chemical potential μ = −0.25 eV. Find the mean occupancy.
Solution
kT = 8.617×10⁻⁵ × 250 = 0.02154 eV. ε+μ = 0.15 − 0.25 = −0.10 eV. β(ε+μ) = −0.10/0.02154 = −4.643. N̅ = 1/(1 + e^4.643) = 1/(1 + 103.9) = 1/104.9 = 9.5×10⁻³, about 0.95% occupancy. Lower temperature sharpens the Fermi factor, so the site is emptier than the 300 K case despite similar ε+μ.
(Easy–medium) Starting from Ξ = Σ_s e^{−β(E_s−μN_s)}, show that N̅ = (1/β)(∂ ln Ξ/∂μ)_T,V.
Solution
Differentiate ln Ξ with respect to μ at fixed β,V: (∂ ln Ξ/∂μ) = (1/Ξ)Σ_s(βN_s)e^{−β(E_s−μN_s)}. The factor e^{−β(E_s−μN_s)}/Ξ = P_s, so the sum is βΣ_sN_sP_s = βN̅. Hence (∂ ln Ξ/∂μ) = βN̅, i.e. N̅ = (1/β)(∂ ln Ξ/∂μ) = kT(∂ ln Ξ/∂μ). Equivalently N̅ = −(∂Φ/∂μ)_T,V.
(Medium) Show that the variance of particle number satisfies σ_N² = ⟨N²⟩ − N̅² = kT(∂N̅/∂μ)_T,V, and comment on what happens near a critical point.
Solution
From problem 2, N̅ = (1/β)(∂ ln Ξ/∂μ). Differentiate again: (∂N̅/∂μ) = (1/β)(∂² ln Ξ/∂μ²). Now (∂² ln Ξ/∂μ²) = ∂/∂μ[(1/Ξ)ΣβN_se^…] = β²(⟨N²⟩−N̅²) by the standard cumulant identity (derivative of a normalised mean gives the variance). Therefore (∂N̅/∂μ) = (1/β)β²σ_N² = βσ_N², giving σ_N² = (1/β)(∂N̅/∂μ) = kT(∂N̅/∂μ)_T,V. Near a critical point the compressibility κ_T ∝ (∂N̅/∂μ) diverges, so σ_N diverges — this is critical opalescence and the breakdown of ensemble equivalence.
(Medium–hard) Compute the chemical potential of gaseous helium (M = 4.00 g/mol) at T = 300 K and P = 1 atm, and confirm the classical limit applies.
Solution
Use μ = kT ln(nλ³). Density n = P/kT = 101325/(1.381×10⁻²³×300) = 2.446×10²⁵ m⁻³ (same as argon). Mass m = 4.00×1.6605×10⁻²⁷ = 6.64×10⁻²⁷ kg. 2πmkT = 6.283×6.64×10⁻²⁷×4.143×10⁻²¹ = 1.728×10⁻⁴⁶; √{} = 1.315×10⁻²³. λ = 6.626×10⁻³⁴/1.315×10⁻²³ = 5.04×10⁻¹¹ m (larger than argon since He is lighter). λ³ = 1.28×10⁻³¹ m³. nλ³ = 2.446×10²⁵×1.28×10⁻³¹ = 3.13×10⁻⁶ ≪ 1 (classical). μ = 0.02585 eV × ln(3.13×10⁻⁶) = 0.02585 × (−12.67) = −0.328 eV. Helium's μ is less negative than argon's because its larger λ raises nλ³.
(Hard) Prove the identity Φ = −PV for a simple bulk system, and hence that ln Ξ = PV/kT. State where extensivity is used.
Solution
Energy is a first-order homogeneous (extensive) function of its extensive arguments: U(λS,λV,λN) = λU(S,V,N). Euler's theorem then gives U = S(∂U/∂S) + V(∂U/∂V) + N(∂U/∂N) = TS − PV + μN. The grand potential is defined as Φ ≡ F − μN = (U − TS) − μN. Substitute the Euler relation: Φ = (TS − PV + μN) − TS − μN = −PV. Since Φ = −kT ln Ξ, we get ln Ξ = PV/kT. Extensivity enters precisely at the Euler step; it fails for systems with dominant surface energy or long-range forces, where Φ = −PV no longer holds and the grand potential must be used directly.