Method of the Most Probable Distribution

How does a gas "decide" how to spread its energy among molecules? It does not decide at all — it simply lands, overwhelmingly often, in the arrangement that can be built the greatest number of ways. That single idea, made quantitative, is Boltzmann's combinatorial route to equilibrium statistical mechanics: the method of the most probable distribution. Partition the available single-particle states into groups labelled by energy, assign occupation numbers $n_i$ , and count the microstates $W(\{n_i\})$ that realize each choice. By the fundamental postulate of equal a priori probabilities every microstate is equally likely, so the most probable macrostate is the one with the largest $W$ . Equilibrium is then, by definition, the distribution that maximizes $W$ — or equivalently $\ln W = S/k_{\mathrm B}$ , which is why maximizing multiplicity is maximizing entropy. The physical constraints enter through Lagrange multipliers. Conservation of particle number $\sum_i n_i=N$ brings a multiplier $\alpha$ ; conservation of energy $\sum_i n_i\varepsilon_i=E$ brings a multiplier $\beta$ . The stationarity condition $\delta(\ln W-\alpha\sum n_i-\beta\sum n_i\varepsilon_i)=0$ delivers the equilibrium occupation numbers in closed form, and comparison with thermodynamics fixes $\beta=1/k_{\mathrm B}T$ and $\alpha=-\mu/k_{\mathrm B}T$ . Feed in three different counting rules and the same maximization produces the three great distributions of statistical physics: Maxwell–Boltzmann for distinguishable (or dilute) particles, Bose–Einstein for indistinguishable particles with unrestricted occupancy, and Fermi–Dirac for particles obeying the Pauli exclusion principle. The legitimacy of replacing an entire ensemble by a single distribution rests on the astonishing sharpness of the maximum. Because $\ln W$ is extensive ( $\sim N$ ) while the width of its peak grows only as $\sqrt N$ , the relative fluctuation about the most probable distribution vanishes as $N^{-1/2}$ . For $N\sim10^{23}$ that is a part in $10^{11}$ : the most probable distribution, the mean distribution, and essentially every distribution carrying non-negligible weight all coincide. This is precisely why Boltzmann's combinatorics reproduces, term for term, the canonical and grand-canonical ensemble results that Gibbs would later derive more formally — the method is both the historical foundation and a transparent, physically motivated shortcut.

Key takeaways

The method identifies equilibrium as the set of occupation numbers $\{n_i\}$ that maximizes the multiplicity $W(\{n_i\})$ — the number of microstates realizing a macrostate — subject to fixed total particle number and total energy.
Maximizing $\ln W$ with Lagrange multipliers $\alpha$ (for $N$ ) and $\beta$ (for $E$ ) yields the equilibrium distribution; the multipliers acquire thermodynamic meaning $\beta=1/k_{\mathrm B}T$ , $\alpha=-\mu/k_{\mathrm B}T$ .
The same maximization, with three different counting rules for $W$ , produces Maxwell–Boltzmann ( $n_i=g_ie^{-\alpha-\beta\varepsilon_i}$ ), Bose–Einstein ( $-1$ in the denominator), and Fermi–Dirac ( $+1$ ).
The maximum is overwhelmingly sharp: the curvature $\partial^2\ln W/\partial n_i^2=-1/n_i$ makes fluctuations Gaussian with relative width $1/\sqrt{n_i}\sim N^{-1/2}\to0$ , so one distribution stands in for the whole ensemble.
The result coincides with the canonical/grand-canonical ensemble averages, providing the historical bridge from Boltzmann's combinatorics to Gibbs's ensemble theory.
Key idea: Equilibrium is the macrostate of maximum multiplicity. Maximize $\ln W$ over occupation numbers under $\sum_i n_i=N$ and $\sum_i n_i\varepsilon_i=E$ using Lagrange multipliers $\alpha,\beta$ ; the stationarity condition delivers the equilibrium distribution, and the sharpness of the peak (relative width $\sim N^{-1/2}$ ) makes that single distribution effectively the entire ensemble.

The core insight

Equilibrium is the macrostate of maximum multiplicity. Maximize $\ln W$ over occupation numbers under the constraints $\sum_i n_i=N$ and $\sum_i n_i\varepsilon_i=E$ using Lagrange multipliers $\alpha,\beta$ ; the stationarity condition delivers the equilibrium distribution, and the sharpness of the peak (relative width $\sim N^{-1/2}$ ) makes that single distribution effectively the entire ensemble.

Where this lives today

Every time a computational chemist evaluates a partition function, or a semiconductor engineer computes the carrier density in a doped crystal, they are using the equilibrium distributions this method produces. The Fermi–Dirac occupancy is the workhorse of solid-state electronics — it sets the position of the Fermi level, the carrier statistics in transistors, and the thermionic emission from cathodes. The Bose–Einstein occupancy, applied to photons ( $\mu=0$ ), is Planck's blackbody law, and applied to trapped ultracold atoms it predicts the Bose–Einstein condensate first realized in 1995 (Nobel Prize 2001). The method also illuminates exotic states. For a system with a bounded energy spectrum — such as $N$ nuclear spins in a magnetic field — maximizing $\ln W$ at high energy yields population inversion and a negative absolute temperature ( $\beta<0$ ), hotter than any positive temperature. Negative-temperature spin states, and more recently negative-temperature motional states of cold atoms in optical lattices, are direct experimental realizations of what happens when $W(E)$ decreases with energy — a scenario the most probable distribution method makes completely transparent.

The physical picture

Imagine distributing $N$ identical coins among many boxes labelled by energy. Configurations that pile every coin into one box are rare; configurations that spread the coins broadly are overwhelmingly numerous. Equilibrium is simply the spreading pattern realizable in the largest number of ways — it is multiplicity, not energy, that is being maximized, with energy entering only as a constraint that tilts the balance toward lower-energy boxes. The Lagrange multipliers are just the "prices" of the conserved quantities. The multiplier $\beta$ is the marginal entropy cost of adding energy — it has to equal $1/k_{\mathrm B}T$ — and $\alpha$ is the cost of adding a particle, which fixes the chemical potential. When two systems touch and exchange energy, they equalize $\beta$ : that is the zeroth law. When they also exchange particles, they equalize $\alpha$ : that is chemical equilibrium. The abstract multipliers are the temperature and chemical potential in disguise.

Figure

Same maximization, three counting rules: MB, BE, FD occupation

The only difference between Maxwell–Boltzmann, Bose–Einstein, and Fermi–Dirac statistics is the counting rule fed into $W$ . Bose counting lets any number of indistinguishable particles share a state; Fermi counting caps occupancy at one, enforcing Pauli exclusion; Maxwell–Boltzmann treats particles as distinguishable (or the states as dilutely occupied). Same maximization, three combinatorial inputs, three distributions — a remarkable economy of principle.

Figure

The maximum is overwhelmingly sharp: ln W peaks with Gaussian width √N

The core equations

Start by counting. For $N$ distinguishable particles distributed with $n_i$ in each group $i$ (each group holding $g_i$ single-particle states), the number of ways is a multinomial times the ways of placing particles among the states within each group.

W(\{n_i\}) = \frac{N!}{\prod_i n_i!}\prod_i g_i^{\,n_i} \qquad (1)

where: $W$ — number of microstates realizing the macrostate $\{n_i\}$ ; $N$ — total particles; $n_i$ — occupation of group $i$ ; $g_i$ — degeneracy (number of single-particle states) in group $i$ . Reading it physically: the multinomial counts which particles land in which group; the $g_i^{n_i}$ counts which state within the group each takes. The physics enters as two conservation laws that the distribution must respect.

\sum_i n_i = N,\qquad \sum_i n_i\,\varepsilon_i = E \qquad (2)

where: $\varepsilon_i$ — single-particle energy of group $i$ ; $E$ — fixed total energy. Reading it physically: these constraints are what prevent the trivial answer (pile everything into the lowest state) — energy conservation forces some spread. Because the $n_i$ are not independent, we cannot just set $\partial\ln W/\partial n_i=0$ . We attach the constraints with undetermined multipliers and extremize the combination freely.

\delta\!\left(\ln W - \alpha\sum_i n_i - \beta\sum_i n_i\varepsilon_i\right)=0 \qquad (3)

where: $\alpha$ — multiplier for particle number, later $\alpha=-\mu/k_{\mathrm B}T$ ; $\beta$ — multiplier for energy, later $\beta=1/k_{\mathrm B}T$ . Reading it physically: $\alpha$ and $\beta$ are the marginal "prices" of a particle and a unit of energy, respectively. To differentiate we need $\ln W$ as a smooth function, which Stirling's approximation provides by turning factorials into logarithms.

\ln W \approx N\ln N - \sum_i n_i\ln n_i + \sum_i n_i\ln g_i \qquad (4)

where: the linear $-N$ terms from Stirling on $N!$ and each $n_i!$ cancel because $\sum_i n_i=N$ . Reading it physically: this is $S/k_{\mathrm B}$ as a function of the occupation numbers — the object whose maximum defines equilibrium. Extremizing (4) under the constraints gives the Maxwell–Boltzmann distribution — an exponential in energy.

n_i = g_i\, e^{-\alpha-\beta\varepsilon_i} \qquad (5)

where: $e^{-\alpha}$ — normalization fixing $N$ (related to fugacity); $e^{-\beta\varepsilon_i}$ — the Boltzmann factor. Reading it physically: occupancy falls off exponentially with energy, at a rate set by $\beta=1/k_{\mathrm B}T$ . The multipliers are not free parameters; thermodynamics pins them down through $S=k_{\mathrm B}\ln W$ .

\beta=\frac{1}{k_{\mathrm B}T},\qquad \alpha=-\frac{\mu}{k_{\mathrm B}T} \qquad (6)

where: $T$ — absolute temperature; $\mu$ — chemical potential. Reading it physically: $\beta=(\partial S/\partial E)/k_{\mathrm B}$ is inverse temperature, and $\alpha=-(\partial S/\partial N)/k_{\mathrm B}$ carries the chemical potential — equality of these across subsystems is thermal and diffusive equilibrium. Change the counting rule and the same method yields quantum statistics. Bose–Einstein allows unlimited occupancy; Fermi–Dirac allows at most one per state.

W_{\mathrm{BE}}=\prod_i \frac{(n_i+g_i-1)!}{n_i!\,(g_i-1)!},\quad W_{\mathrm{FD}}=\prod_i \frac{g_i!}{n_i!\,(g_i-n_i)!} \qquad (7)

where: BE counts indistinguishable particles in "bins" with repetition allowed; FD counts choices of $n_i$ occupied states out of $g_i$ (each at most once). Reading it physically: the two combinatorial forms encode the two ways nature lets identical particles share states. Maximizing (7) under the same constraints collapses both into a single master formula.

\langle n_i\rangle = \frac{g_i}{e^{\alpha+\beta\varepsilon_i}\pm 1} \qquad (8)

where: $+1$ — Fermi–Dirac; $-1$ — Bose–Einstein. Reading it physically: the lone sign in the denominator is the entire mathematical difference between fermions and bosons, and it produces degeneracy pressure on one side and Bose condensation on the other.

Derivations

Derivation AMaxwell–Boltzmann distribution from maximizing ln W

We maximize the Boltzmann-corrected multiplicity under the number and energy constraints to obtain the exponential distribution $n_i\propto e^{-\beta\varepsilon_i}$ . The strategy is to write $\ln W$ via Stirling as a smooth function of the occupation numbers, attach the two constraints with Lagrange multipliers so the $n_i$ can be varied independently, differentiate level by level, and solve; the payoff is that the same three-line calculation, with a different $W$ , will later give the quantum distributions too.

Step 1. From equation (4), $\ln W = N\ln N + \sum_i n_i\ln g_i - \sum_i n_i\ln n_i$ . Form the constrained functional $\Phi=\ln W-\alpha(\sum_i n_i-N)-\beta(\sum_i n_i\varepsilon_i-E)$ .

Step 2. Differentiate with respect to a single $n_j$ : $\partial\Phi/\partial n_j=\ln g_j-\ln n_j-1+1-\alpha-\beta\varepsilon_j=0$ . Why the $\pm1$ cancels: $\partial(n_j\ln n_j)/\partial n_j=\ln n_j+1$ supplies a $-1$ , while $\partial(N\ln N)/\partial n_j=1$ (since $\sum n_i=N$ ) supplies a $+1$ ; they annihilate, leaving $\ln(g_j/n_j)=\alpha+\beta\varepsilon_j$ .

Step 3. Solve: $n_j=g_j\,e^{-\alpha-\beta\varepsilon_j}$ .

Step 4. Fix the multipliers. Imposing $\sum_i g_i e^{-\alpha-\beta\varepsilon_i}=N$ defines the single-particle partition function $Z_1=\sum_i g_i e^{-\beta\varepsilon_i}$ , so $e^{\alpha}=Z_1/N$ ; the thermodynamic identification gives $\beta=1/k_{\mathrm B}T$ , hence $n_j=(N/Z_1)g_j e^{-\varepsilon_j/k_{\mathrm B}T}$ . Verified: symbolic solve of the stationarity condition returns exactly $n_i=g_i e^{-\alpha-\beta\varepsilon_i}$ ; the result equals the canonical single-particle Boltzmann weight $p_i\propto g_i e^{-\beta\varepsilon_i}$ (Pathria & Beale Ch. 3).

Go deeper · Derivation B · 1 minSharpness of the maximum: Gaussian fluctuations and the 1/√N law

Derivation BSharpness of the maximum: Gaussian fluctuations and the 1/√N law

We expand $\ln W$ to second order about its maximum to prove that the relative width of the peak scales as $N^{-1/2}$ — the structural result that justifies replacing the entire ensemble by the single most probable distribution. The strategy is to write $n_i=n_i^*+\delta n_i$ , use the fact that the first-order term vanishes on the constraint surface, compute the curvature explicitly, and read off the variance of the fluctuations; the conclusion is that for macroscopic $N$ the peak is a needle so thin it may as well be a delta function.

Step 1. Taylor expand: $\ln W=\ln W^*+\sum_i(\partial\ln W/\partial n_i)|_*\,\delta n_i+\tfrac12\sum_i(\partial^2\ln W/\partial n_i^2)|_*(\delta n_i)^2$ . The first-order coefficient is $\alpha+\beta\varepsilon_i$ , and the constrained variations satisfy $\sum_i\delta n_i=0$ and $\sum_i\varepsilon_i\delta n_i=0$ , so the linear term vanishes.

Step 2. Compute the curvature. From the $-n_i\ln n_i$ piece, $\partial^2\ln W/\partial n_i^2=-1/n_i$ . Why this is the whole story: only the diagonal curvature survives at leading order, and it is uniformly negative, confirming a genuine maximum.

Step 3. Hence $W(\{n_i^*+\delta n_i\})/W(\{n_i^*\})=\exp[-\tfrac12\sum_i(\delta n_i)^2/n_i^*]$ — a Gaussian with variance $\langle(\delta n_i)^2\rangle=n_i^*$ , so $\Delta n_i/n_i^*=1/\sqrt{n_i^*}$ .

Step 4. Since typical $n_i^*\sim N$ , relative fluctuations $\sim N^{-1/2}$ . For $N=10^{23}$ this is $\sim10^{-11}$ : the most probable distribution dominates so completely that it is equilibrium. Verified: symbolic curvature $\partial^2(-n\ln n)/\partial n^2=-1/n$ ; numeric binomial widths confirm relative width $=1/\sqrt N$ for $N=10^2$ – $10^8$ ; consistent with $\langle(\delta n)^2\rangle=\langle n\rangle$ Poisson statistics in the classical limit (Reif §6).

Go deeper · Derivation C · 1 minBose–Einstein and Fermi–Dirac distributions from the same method

Derivation CBose–Einstein and Fermi–Dirac distributions from the same method

We apply the identical variational procedure to the quantum multiplicities $W_{\mathrm{BE}}$ and $W_{\mathrm{FD}}$ to obtain the $\mp1$ distributions, showing that the only thing that changes is the combinatorial input. The strategy mirrors the MB derivation — Stirling, Lagrange, differentiate, solve — but with the Bose "bins-with-repetition" count and the Fermi "at most one per state" count, and the punchline is that a single sign flip in the final denominator separates all of quantum statistics.

Step 1 (Bose). With $g_i\gg1$ , $\ln W_{\mathrm{BE}}=\sum_i[(n_i+g_i)\ln(n_i+g_i)-n_i\ln n_i-g_i\ln g_i]$ . Extremizing $\ln W_{\mathrm{BE}}-\alpha\sum n_i-\beta\sum n_i\varepsilon_i$ gives $\ln\frac{n_i+g_i}{n_i}-\alpha-\beta\varepsilon_i=0$ .

Step 2 (Bose result). Solve: $n_i/g_i=1/(e^{\alpha+\beta\varepsilon_i}-1)$ — the Bose–Einstein distribution, with the telltale $-1$ .

Step 3 (Fermi). With occupancy capped at $g_i$ , $\ln W_{\mathrm{FD}}=\sum_i[g_i\ln g_i-n_i\ln n_i-(g_i-n_i)\ln(g_i-n_i)]$ . Extremizing gives $\ln\frac{g_i-n_i}{n_i}-\alpha-\beta\varepsilon_i=0$ .

Step 4 (Fermi result). Solve: $n_i/g_i=1/(e^{\alpha+\beta\varepsilon_i}+1)$ — the Fermi–Dirac distribution, with the $+1$ that enforces $\langle n_i\rangle\le g_i$ . Why one sign matters so much: the $+1$ caps occupancy (Pauli), giving degeneracy pressure and the linear electronic heat capacity; the $-1$ lets occupancy diverge as $\varepsilon\to\mu$ , giving Bose–Einstein condensation.

Step 5. When $e^{\alpha+\beta\varepsilon_i}\gg1$ (dilute/non-degenerate limit) both reduce to $n_i=g_ie^{-\alpha-\beta\varepsilon_i}$ , recovering Maxwell–Boltzmann. Verified: symbolic solve returns $n/g=1/(e^{a}\mp1)$ for BE/FD; the classical limit $a\ge5$ agrees with MB to under 1%; the FD occupancy entropy $-[f\ln f+(1-f)\ln(1-f)]$ peaks at $\ln2$ at $f=\tfrac12$ (Landau–Lifshitz §53–55).

Limiting cases

The master result $\langle n_i\rangle=g_i/(e^{\alpha+\beta\varepsilon_i}\pm1)$ specializes cleanly:

Case	Result	Meaning
Maxwell–Boltzmann	$\langle n_i\rangle=g_i e^{(\mu-\varepsilon_i)/k_{\mathrm B}T}$	Classical / dilute limit, valid when $e^{\alpha+\beta\varepsilon_i}\gg1$ , i.e. occupancy $\ll1$ per state.
Bose–Einstein	$\langle n_i\rangle=\dfrac{g_i}{e^{(\varepsilon_i-\mu)/k_{\mathrm B}T}-1}$	Indistinguishable, unrestricted occupancy; requires $\mu<\varepsilon_0$ . Bose condensation appears as $\mu\to\varepsilon_0$ .
Fermi–Dirac	$\langle n_i\rangle=\dfrac{g_i}{e^{(\varepsilon_i-\mu)/k_{\mathrm B}T}+1}$	Pauli exclusion, occupancy $\le1$ ; at $T\to0$ a step at the Fermi energy $\varepsilon_F=\mu(0)$ .
Thermodynamic limit	$\dfrac{\Delta n_i}{\langle n_i\rangle}\sim\dfrac{1}{\sqrt{\langle n_i\rangle}}\to0$	Most probable distribution equals the mean; fluctuations negligible, validating the whole method.

Common trapThe most probable distribution is just one of many — why trust it?

Because the maximum is not a gentle hill but a needle. It is tempting to worry that singling out one distribution ignores all the others, but the curvature $\partial^2\ln W/\partial n_i^2=-1/n_i$ makes the peak Gaussian with relative width $1/\sqrt{n_i}\sim N^{-1/2}$ . For $N\sim10^{23}$ that is $10^{-11}$ — every distribution with non-negligible weight is indistinguishable from the most probable one. The separate caution: Stirling's approximation assumes $n_i\gg1$ , so the smooth maximization can quietly fail for sparsely populated high-energy levels; there one must keep more terms or treat the levels discretely.

Worked examples

Example 1 — The equilibrium of a two-level system by direct counting. Take $N=40$ distinguishable particles, each either in a ground state ( $\varepsilon=0$ ) or an excited state ( $\varepsilon>0$ ). Find the most probable number in the excited state (i) with no energy constraint and (ii) with a Boltzmann weight $e^{-\beta\varepsilon}=e^{-0.8}$ , then read off the single-particle occupancy.

Step 1. Multiplicity of " $k$ excited" is $W(k)=\binom{N}{k}$ . With no energy constraint, $\ln W$ is maximized where $\binom{40}{k}$ peaks, i.e. at $k=N/2=\mathbf{20}$ — equal populations, the infinite-temperature ( $\beta=0$ ) state.

Step 2. Now weight each configuration by energy: the equilibrium mean is $\langle k\rangle=\sum_k k\binom{N}{k}e^{-\beta\varepsilon k}\big/\sum_k\binom{N}{k}e^{-\beta\varepsilon k}$ . Computing with $\beta\varepsilon=0.8$ gives $\langle k\rangle=\mathbf{12.40}$ .

Step 3. For $N$ independent two-level particles the single-particle excited occupancy is $f=1/(e^{\beta\varepsilon}+1)$ , so $\langle k\rangle=N/(e^{0.8}+1)=40/(2.2255+1)=40/3.2255=\mathbf{12.40}$ — the direct count and the closed-form agree exactly. The most probable excited population drops from 20 (no energy cost) to 12.40 (with $\beta\varepsilon=0.8$ ), and the two independent calculations match to all figures — the method delivers a real number, not just a formula.

Example 2 — Sharpness of the peak for a macroscopic gas. Estimate the relative fluctuation of the occupation of a typical energy level for $N=10^{20}$ particles, and interpret it.

Step 1. From the Gaussian result, $\Delta n_i/n_i^*=1/\sqrt{n_i^*}$ . For a typical level with $n_i^*\sim N/(\text{number of populated levels})$ ; take $n_i^*\sim10^{18}$ as a representative macroscopic occupancy.

Step 2. Then $\Delta n_i/n_i^*=1/\sqrt{10^{18}}=\mathbf{10^{-9}}$ .

Step 3. Interpret: a level holding $10^{18}$ particles fluctuates by only about $10^{9}$ particles — one part in a billion. No measurement of a macroscopic gas could ever resolve a deviation from the most probable distribution. The relative fluctuation is $\sim10^{-9}$ , confirming that for macroscopic $N$ the most probable distribution is equilibrium to a precision far beyond any experiment.

History band

Origins · 1 minLudwig Boltzmann and how this was discovered

Ludwig Boltzmann

1844–1906 · Austrian physicist who founded statistical mechanics. His combinatorial counting of microstates and the identification of entropy with

\log

of multiplicity —

S=k\log W

, engraved on his tombstone in Vienna — is exactly the method of the most probable distribution. Boltzmann derived the equilibrium velocity distribution and the H-theorem, defending the atomic hypothesis against fierce opposition; the Lagrange-multiplier maximization of

\ln W

is his enduring route from microscopic counting to thermodynamics, later reformulated by Gibbs as ensemble theory.

History

1860

— Maxwell derives the equilibrium velocity distribution of gas molecules

1868

— Boltzmann generalizes Maxwell's distribution to include potential energy

1872

— Boltzmann's transport equation and H-theorem describe the approach to equilibrium

1877

— Boltzmann relates entropy to microstate multiplicity, S = k log W

1902

— Gibbs recasts the combinatorial method as ensemble theory

1924

— Bose and Einstein introduce the counting rule for indistinguishable bosons

1926

— Fermi and Dirac introduce the exclusion-respecting counting for fermions

Going further

References. Pathria & Beale, Statistical Mechanics, 3rd ed., Ch. 3 & 6 · Kardar, Statistical Physics of Particles, Ch. 4 & 7 · Reif, Fundamentals of Statistical and Thermal Physics, Ch. 6 & 9 · Sethna, Statistical Mechanics: Entropy, Order Parameters, and Complexity, Ch. 3 & 7 · Landau & Lifshitz, Statistical Physics Part 1, §7, §40, §53–55.

Provenance. All distributions, multiplicities, and the Lagrange-multiplier derivation are standard and recomputable: verify the MB derivation against Pathria & Beale Ch. 3, the BE/FD multiplicities and distributions against Pathria Ch. 6 and Landau–Lifshitz §53–55, the sharpness-of-maximum argument against Reif Ch. 6, and Stirling's approximation $\ln N!\approx N\ln N-N$ . The identifications $\beta=1/k_{\mathrm B}T$ , $\alpha=-\mu/k_{\mathrm B}T$ follow from $S=k_{\mathrm B}\ln W$ with $1/T=(\partial S/\partial E)_N$ . The two-level worked example ( $N=40$ , $\beta\varepsilon=0.8$ ) gives $\langle k\rangle=12.40$ by both direct summation and the closed form $N/(e^{\beta\varepsilon}+1)$ ; this and all other numbers are independently recomputed. No non-standard numbers are asserted.

Builds on: Microstates and Macrostates · Fundamental Postulate of Equal A Priori Probabilities · Accessible Phase Space and Energy Shells · Boltzmann Entropy Leads to: Statistical Ensembles · Gibbs Entropy

Why maximize ln W instead of W itself?

They have the same maximizer — $\ln$ is monotonic — but $\ln W$ is vastly easier to work with. $W$ is an astronomically large product of factorials ( $\sim10^{10^{23}}$ ), impossible to differentiate directly. Taking the log turns the product into a sum, and Stirling's approximation then turns each $\ln n_i!$ into the smooth $n_i\ln n_i-n_i$ , giving a differentiable function you can extremize with calculus. Physically, $\ln W$ is also the more natural object: it equals $S/k_{\mathrm B}$ , so maximizing $\ln W$ is literally maximizing entropy.

Where does the Boltzmann factor

e^{-\beta\varepsilon}

actually come from here?

It falls straight out of the maximization. Differentiating $\ln W$ under the energy constraint gives $\ln(g_i/n_i)=\alpha+\beta\varepsilon_i$ , and exponentiating yields $n_i=g_ie^{-\alpha-\beta\varepsilon_i}$ . The $e^{-\beta\varepsilon_i}$ is not assumed — it is forced by the requirement that the distribution maximize multiplicity while conserving energy. The multiplier $\beta$ multiplies $\varepsilon_i$ precisely because energy is the constrained quantity, and thermodynamics then identifies $\beta=1/k_{\mathrm B}T$ . So the exponential decay of occupancy with energy is a theorem about counting, not a postulate.

Do I use the same

\alpha

and

\beta

for bosons, fermions, and classical particles?

The roles are identical — $\beta=1/k_{\mathrm B}T$ always, and $\alpha=-\mu/k_{\mathrm B}T$ always — because they are fixed by the same thermodynamic relations $1/T=(\partial S/\partial E)_N$ and $-\mu/T=(\partial S/\partial N)_E$ . What differs is only the multiplicity $W$ fed into the maximization, and therefore the functional form of the resulting occupancy: $g_ie^{-\alpha-\beta\varepsilon_i}$ (MB), $g_i/(e^{\alpha+\beta\varepsilon_i}-1)$ (BE), $g_i/(e^{\alpha+\beta\varepsilon_i}+1)$ (FD). Same multipliers, same constraints, three counting rules.

Is the "most probable distribution" the same as the canonical ensemble average?

In the thermodynamic limit, yes — and that equivalence is the whole point. The canonical ensemble computes averages over all distributions weighted by their probability; the most probable distribution method uses only the single peak. Because the peak's relative width vanishes as $N^{-1/2}$ , the peak and the weighted average coincide to a part in $\sqrt N$ . So Boltzmann's combinatorial method and Gibbs's ensemble method give identical results for macroscopic systems; they differ only in bookkeeping, and the ensemble method is cleaner for computing fluctuations, which the peak method must handle by expanding around the maximum.

Why does energy conservation need a multiplier but nothing else does?

Only conserved (constrained) quantities get multipliers. In a closed system with fixed $N$ and $E$ , exactly two things are held fixed, so exactly two multipliers appear — $\alpha$ for $N$ and $\beta$ for $E$ . If the system could also exchange volume or another conserved charge, each would bring its own multiplier (pressure, another chemical potential). Quantities that are not conserved — like the number of photons, which can be created and destroyed — carry no multiplier, which is why photon gases have $\alpha=0$ , i.e. $\mu=0$ . The count of multipliers is simply the count of independent conservation laws.

How can a temperature be negative in this framework?

Through $\beta=(\partial\ln W/\partial E)$ . For most systems $W(E)$ increases with energy, so $\beta>0$ and $T>0$ . But if the energy spectrum is bounded above — as for $N$ two-level spins, whose maximum energy is all-excited — then $W(E)$ rises, peaks, and falls as $E$ approaches its maximum. On the falling side $\partial\ln W/\partial E<0$ , so $\beta<0$ and the absolute temperature is negative. Such states are genuinely hotter than $T=+\infty$ : put them in contact with any positive-temperature system and heat flows out of them. The method makes this transparent because $\beta$ is manifestly the slope of $\ln W(E)$ .