D.57 The par­ti­cle en­ergy dis­tri­b­u­tions

This note de­rives the Maxwell-Boltz­mann, Fermi-Dirac, and Bose-Ein­stein en­ergy dis­tri­b­u­tions of weakly in­ter­act­ing par­ti­cles for a sys­tem for which the net en­ergy is pre­cisely known.

The ob­jec­tive is to find the shelf num­bers $\vec{I}$ $\vphantom0\raisebox{1.5pt}{$=$}$ $(I_1,I_2,I_3,\ldots)$ for which the num­ber of eigen­func­tions $Q_{\vec{I}}$ is max­i­mal. Ac­tu­ally, it is math­e­mat­i­cally eas­ier to find the max­i­mum of $\ln(Q_{\vec{I}})$, and that is the same thing: if $Q_{\vec{I}}$ is as big as it can be, then so is $\ln(Q_{\vec{I}})$. The ad­van­tage of work­ing with $\ln(Q_{\vec{I}})$ is that it sim­pli­fies all the prod­ucts in the ex­pres­sions for the $Q_{\vec{I}}$ de­rived in de­riva­tion {D.56} into sums: math­e­mat­ics says that $\ln(ab)$ equals $\ln(a)$ plus $\ln(b)$ for any (pos­i­tive) $a$ and $b$.

It will be as­sumed, fol­low­ing de­riva­tion {N.24}, that if the max­i­mum value is found among all shelf oc­cu­pa­tion num­bers, whole num­bers or not, it suf­fices. More dar­ingly, er­rors less than a par­ti­cle are not go­ing to be taken se­ri­ously.

In find­ing the max­i­mum of $\ln(Q_{\vec{I}})$, the shelf num­bers can­not be com­pletely ar­bi­trary; they are con­strained by the con­di­tions that the sum of the shelf num­bers must equal the to­tal num­ber of par­ti­cles $I$, and that the par­ti­cle en­er­gies must sum to­gether to the given to­tal en­ergy $E$:

$$\sum_s I_s = I \qquad \sum_s I_s {\vphantom' E}^{\rm p}_s = E.$$

Math­e­mati­cians call this a con­strained max­i­miza­tion prob­lem.

Ac­cord­ing to cal­cu­lus, with­out the con­straints, you can just put the de­riv­a­tives of $\ln(Q_{\vec{I}})$ with re­spect to all the shelf num­bers $I_s$ to zero to find the max­i­mum. With the con­straints, you have to add penalty terms that cor­rect for any go­ing out of bounds, {D.48}, and the cor­rect func­tion whose de­riv­a­tives must be zero is

$$F = \ln(Q_{\vec I}) - \epsilon_1\bigg(\sum_s I_s - I\bigg)... ...\epsilon_2\bigg(\sum_s I_s {\vphantom' E}^{\rm p}_s - E\bigg)$$

where the con­stants $\epsilon_1$ and $\epsilon_2$ are un­known penalty fac­tors called the La­grangian mul­ti­pli­ers.

At the shelf num­bers for which the num­ber of eigen­func­tions is largest, the de­riv­a­tives $\partial{F}$$\raisebox{.5pt}{$/$}$​$\partial{I}_s$ must be zero. How­ever, that con­di­tion is dif­fi­cult to ap­ply ex­actly, be­cause the ex­pres­sions for $Q_{\vec{I}}$ as given in the text in­volve the fac­to­r­ial func­tion, or rather, the gamma func­tion. The gamma func­tion does not have a sim­ple de­riv­a­tive. Here typ­i­cal text­books will flip out the Stir­ling ap­prox­i­ma­tion of the fac­to­r­ial, but this ap­prox­i­ma­tion is sim­ply in­cor­rect in parts of the range of in­ter­est, and where it ap­plies, the er­ror is un­known.

It is a much bet­ter idea to ap­prox­i­mate the dif­fer­en­tial quo­tient by a dif­fer­ence quo­tient, as in

$$\begin{eqnarray*} \lefteqn{0 =\frac{\partial F}{\partial I_s} \approx \frac{\D... ...- F(I_1,I_2,\ldots,I_{s-1},I_s,I_{s+1},\ldots)}{I_s + 1 - I_s}. \end{eqnarray*}$$

This ap­prox­i­ma­tion is very mi­nor, since ac­cord­ing to the so-called mean value the­o­rem of math­e­mat­ics, the lo­ca­tion where $\Delta{F}$$\raisebox{.5pt}{$/$}$​$\Delta{I}_s$ is zero is at most one par­ti­cle away from the de­sired lo­ca­tion where $\partial{F}$$\raisebox{.5pt}{$/$}$​$\partial{I}_s$ is zero. Bet­ter still, $I_s+\frac12$ $\vphantom0\raisebox{1.5pt}{$\equiv$}$ $I_{s,{\rm {best}}}$ will be no more that half a par­ti­cle off, and the analy­sis al­ready had to com­mit it­self to ig­nor­ing frac­tional parts of par­ti­cles any­way. The dif­fer­ence quo­tient leads to sim­ple for­mu­lae be­cause the gamma func­tion sat­is­fies the con­di­tion $(n+1)!$ $\vphantom0\raisebox{1.5pt}{$=$}$ $(n+1) n!$ for any value of $n$, com­pare the no­ta­tions sec­tion un­der !.

Now con­sider first dis­tin­guish­able par­ti­cles. The func­tion $F$ to dif­fer­en­ti­ate is de­fined above, and plug­ging in the ex­pres­sion for $Q^{\rm {d}}_{\vec{I}}$ as found in de­riva­tion {D.56} pro­duces

$$F = \ln(I!) + \sum_s\left[I_s\ln(N_s)-\ln(I_s!)\right] - \... ...epsilon_2 \bigg(\sum_s I_s {\vphantom' E}^{\rm p}_s - E\bigg)$$

For any value of the shelf num­ber $s$, in the limit $I_s\downarrow-1$, $F$ tends to neg­a­tive in­fin­ity be­cause $I_s!$ tends to pos­i­tive in­fin­ity in that limit and its log­a­rithm ap­pears with a mi­nus sign. In the limit $I_s\uparrow+\infty$, $F$ tends once more to neg­a­tive in­fin­ity, since $\ln(I_s!)$ for large val­ues of $I_s$ is ac­cord­ing to the so-called Stir­ling for­mula ap­prox­i­mately equal to $I_s\ln(I_s)-I_s$, so the $-\ln(I_s!)$ term in $F$ goes to mi­nus in­fin­ity more strongly than the terms pro­por­tional to $I_s$ might go to plus in­fin­ity. If $F$ tends to mi­nus in­fin­ity at both ends of the range $\vphantom{0}\raisebox{1.5pt}{$-$}$1 $\raisebox{.3pt}{$<$}$ $I_s$ $\raisebox{.3pt}{$<$}$ $\infty$, there must be a max­i­mum value of $F$ some­where within that range where the de­riv­a­tive with re­spect to $I_s$ is zero. More specif­i­cally, work­ing out the dif­fer­ence quo­tient:

$$\frac{\Delta F}{\Delta I_s} = \ln(N_s) - \ln(I_s+1) - \epsilon_1 -\epsilon_2 {\vphantom' E}^{\rm p}_s = 0$$

and $-\ln(I_s+1)$ is in­fin­ity at $I_s$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vphantom{0}\raisebox{1.5pt}{$-$}$1 and mi­nus in­fin­ity at $I_s$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\infty$. Some­where in be­tween, $\Delta{F}$$\raisebox{.5pt}{$/$}$​$\Delta{I}_s$ will cross zero. In par­tic­u­lar, com­bin­ing the log­a­rithms and then tak­ing an ex­po­nen­tial, the best es­ti­mate for the shelf oc­cu­pa­tion num­ber is

$$I_{s,{\rm best}}= I_s + {\textstyle\frac{1}{2}} = \frac{N_... ...\vphantom' E}^{\rm p}_s+\epsilon_1}} - {\textstyle\frac{1}{2}}$$

The cor­rect­ness of the fi­nal half par­ti­cle is clearly doubt­ful within the made ap­prox­i­ma­tions. In fact, it is best ig­nored since it only makes a dif­fer­ence at high en­er­gies where the num­ber of par­ti­cles per shelf be­comes small, and surely, the cor­rect prob­a­bil­ity of find­ing a par­ti­cle must go to zero at in­fi­nite en­er­gies, not to mi­nus half a par­ti­cle! There­fore, the best es­ti­mate $\iota^{\rm {d}}$ $\vphantom0\raisebox{1.5pt}{$\equiv$}$ $I_{s,{\rm {best}}}$$\raisebox{.5pt}{$/$}$​$N_s$ for the num­ber of par­ti­cles per sin­gle-par­ti­cle en­ergy state be­comes the Maxwell-Boltz­mann dis­tri­b­u­tion. Note that the de­riva­tion might be off by a par­ti­cle for the lower en­ergy shelves. But there are a lot of par­ti­cles in a macro­scopic sys­tem, so it is no big deal.

The case of iden­ti­cal fermi­ons is next. The func­tion to dif­fer­en­ti­ate is now

$$\begin{eqnarray*} F & = & \sum_s\left[\ln(N_s!)-\ln(I_s!)-\ln((N_s-I_s)!)\right... ...- \epsilon_2\bigg(\sum_s I_s {\vphantom' E}^{\rm p}_s - E\bigg) \end{eqnarray*}$$

This time $F$ is mi­nus in­fin­ity when a shelf num­ber reaches $I_s$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vphantom{0}\raisebox{1.5pt}{$-$}$1 or $I_s$ $\vphantom0\raisebox{1.5pt}{$=$}$ $N_s+1$. So there must be a max­i­mum to $F$ when $I_s$ varies be­tween those lim­its. The dif­fer­ence quo­tient ap­prox­i­ma­tion pro­duces

$$\frac{\Delta F}{\Delta I_s} = -\ln(I_s+1) + \ln(N_s-I_s) - \epsilon_1 -\epsilon_2 {\vphantom' E}^{\rm p}_s = 0$$

which can be solved to give

$$I_{s,{\rm best}}= I_s + {\textstyle\frac{1}{2}} = \frac{N_... ...silon_1}}{1+e^{\epsilon_2{\vphantom' E}^{\rm p}_s+\epsilon_1}}$$

The fi­nal term, less than half a par­ti­cle, is again best left away, to en­sure that 0 $\raisebox{-.3pt}{$\leqslant$}$ $I_{s,{\rm {best}}}$ $\raisebox{-.3pt}{$\leqslant$}$ $N_s$ as it should. That gives the Fermi-Dirac dis­tri­b­u­tion.

Fi­nally, the case of iden­ti­cal bosons, is, once more, the tricky one. The func­tion to dif­fer­en­ti­ate is now

$$\begin{eqnarray*} F & = & \sum_s\left[\ln((I_s+N_s-1)!)-\ln(I_s!)-\ln((N_s-1)!)... ...- \epsilon_2\bigg(\sum_s I_s {\vphantom' E}^{\rm p}_s - E\bigg) \end{eqnarray*}$$

For now, as­sume that $N_s$ $\raisebox{.3pt}{$>$}$ 1 for all shelves. Then $F$ is again mi­nus in­fin­ity for $I_s$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vphantom{0}\raisebox{1.5pt}{$-$}$1. For $I_s\uparrow\infty$, how­ever, $F$ will be­have like $-(\epsilon_1+\epsilon_2{\vphantom' E}^{\rm p}_s)I_s$. This tends to mi­nus in­fin­ity if $\epsilon_1+\epsilon_2{\vphantom' E}^{\rm p}_s$ is pos­i­tive, so for now as­sume it is. Then the dif­fer­ence quo­tient ap­prox­i­ma­tion pro­duces

$$\frac{\Delta F}{\Delta I_s} = \ln(I_s+N_s) - \ln(I_s+1) - \epsilon_1 -\epsilon_2 {\vphantom' E}^{\rm p}_s = 0$$

which can be solved to give

$$I_{s,{\rm best}}= I_s + {\textstyle\frac{1}{2}} = \frac{N_... ...hantom' E}^{\rm p}_s+\epsilon_1}-1} - {\textstyle\frac{1}{2}}.$$

The fi­nal half par­ti­cle is again best ig­nored to get the num­ber of par­ti­cles to be­come zero at large en­er­gies. Then, if it is as­sumed that the num­ber $N_s$ of sin­gle-par­ti­cle states on the shelves is large, the Bose-Ein­stein dis­tri­b­u­tion is ob­tained. If $N_s$ is not large, the num­ber of par­ti­cles could be less than the pre­dicted one by up to a fac­tor 2, and if $N_s$ is one, the en­tire story comes part. And so it does if $\epsilon_1+\epsilon_2{\vphantom' E}^{\rm p}_s$ is not pos­i­tive.

Be­fore ad­dress­ing these nasty prob­lems, first the phys­i­cal mean­ing of the La­grangian mul­ti­plier $\epsilon_2$ needs to be es­tab­lished. It can be in­ferred from ex­am­in­ing the case that two dif­fer­ent sys­tems, call them $A$ and $B$, are in ther­mal con­tact. Since the in­ter­ac­tions are as­sumed weak, the eigen­func­tions of the com­bined sys­tem are the prod­ucts of those of the sep­a­rate sys­tems. That means that the num­ber of eigen­func­tions of the com­bined sys­tem $Q_{\vec{I}_A\vec{I}_B}$ is the prod­uct of those of the in­di­vid­ual sys­tems. There­fore the func­tion to dif­fer­en­ti­ate be­comes

$$\begin{eqnarray*} F & = &\ln(Q_{\vec I_A}Q_{\vec I_B}) \\ && - \epsilon_{1,A}... ..._A} + \sum_{s_B} I_{s_B} {\vphantom' E}^{\rm p}_{s_B}- E \bigg) \end{eqnarray*}$$

Note the con­straints: the num­ber of par­ti­cles in sys­tem $A$ must be the cor­rect num­ber $I_A$ of par­ti­cles in that sys­tem, and sim­i­lar for sys­tem $B$. How­ever, since the sys­tems are in ther­mal con­tact, they can ex­change en­ergy through the weak in­ter­ac­tions and there is no longer a con­straint on the en­ergy of the in­di­vid­ual sys­tems. Only the com­bined en­ergy must equal the given to­tal. That means the two sys­tems share the same La­grangian vari­able $\epsilon_2$. For the rest, the equa­tions for the two sys­tems are just like if they were not in ther­mal con­tact, be­cause the log­a­rithm in $F$ sep­a­rates, and then the dif­fer­en­ti­a­tions with re­spect to the shelf num­bers $I_{s_A}$ and $I_{s_B}$ give the same re­sults as be­fore.

It fol­lows that two sys­tems that have the same value of $\epsilon_2$ can be brought into ther­mal con­tact and noth­ing hap­pens, macro­scop­i­cally. How­ever, if two sys­tems with dif­fer­ent val­ues of $\epsilon_2$ are brought into con­tact, the sys­tems will ad­just, and en­ergy will trans­fer be­tween them, un­til the two $\epsilon_2$ val­ues have be­come equal. That means that $\epsilon_2$ is a tem­per­a­ture vari­able. From here on, the tem­per­a­ture will be de­fined as $T$ = 1/${\epsilon_2}k_{\rm B}$, so that $\epsilon_2$ $\vphantom0\raisebox{1.5pt}{$=$}$ 1/${k_{\rm B}}T$, with $k_{\rm B}$ the Boltz­mann con­stant. The same way, for now the chem­i­cal po­ten­tial $\mu$ will sim­ply be de­fined to be the con­stant $-\epsilon_1$$\raisebox{.5pt}{$/$}$​$\epsilon_2$. Chap­ter 11.14.4 will even­tu­ally es­tab­lish that the tem­per­a­ture de­fined here is the ideal gas tem­per­a­ture, while de­riva­tion {D.61} will es­tab­lish that $\mu$ is the Gibbs free en­ergy per atom that is nor­mally de­fined as the chem­i­cal po­ten­tial.

Re­turn­ing now to the nasty prob­lems of the dis­tri­b­u­tion for bosons, first as­sume that every shelf has at least two states, and that $({\vphantom' E}^{\rm p}_s-\mu)$$\raisebox{.5pt}{$/$}$​${k_{\rm B}}T$ is pos­i­tive even for the ground state. In that case there is no prob­lem with the de­rived so­lu­tion. How­ever, Bose-Ein­stein con­den­sa­tion will oc­cur when ei­ther the num­ber den­sity is in­creased by putting more par­ti­cles in the sys­tem, or the tem­per­a­ture is de­creased. In­creas­ing par­ti­cle den­sity is as­so­ci­ated with in­creas­ing chem­i­cal po­ten­tial $\mu$ be­cause

$$I_s = \frac{N_s-1}{e^{({\vphantom' E}^{\rm p}_s-\mu)/{k_{\rm B}}T}-1}$$

im­plies that every shelf par­ti­cle num­ber in­creases when $\mu$ in­creases. De­creas­ing tem­per­a­ture by it­self de­creases the num­ber of par­ti­cles, and to com­pen­sate and keep the num­ber of par­ti­cles the same, $\mu$ must then once again in­crease. When $\mu$ gets very close to the ground state en­ergy, the ex­po­nen­tial in the ex­pres­sion for the num­ber of par­ti­cles on the ground state shelf $s$ $\vphantom0\raisebox{1.5pt}{$=$}$ 1 be­comes very close to one, mak­ing the to­tal de­nom­i­na­tor very close to zero, so the num­ber of par­ti­cles $I_1$ in the ground state blows up. When it be­comes a fi­nite frac­tion of the to­tal num­ber of par­ti­cles $I$ even when $I$ is macro­scop­i­cally large, Bose-Ein­stein con­den­sa­tion is said to have oc­curred.

Note that un­der rea­son­able as­sump­tions, it will only be the ground state shelf that ever ac­quires a fi­nite frac­tion of the par­ti­cles. For, as­sume the con­trary, that shelf 2 also holds a fi­nite frac­tion of the par­ti­cles. Us­ing Tay­lor se­ries ex­pan­sion of the ex­po­nen­tial for small val­ues of its ar­gu­ment, the shelf oc­cu­pa­tion num­bers are

$$\begin{eqnarray*} && I_1 = \frac{(N_1-1){k_{\rm B}}T}{{\vphantom' E}^{\rm p}_1-... ...vphantom' E}^{\rm p}_3-{\vphantom' E}^{\rm p}_2)} \\ && \vdots \end{eqnarray*}$$

For $I_2$ to also be a fi­nite frac­tion of the to­tal num­ber of par­ti­cles, ${\vphantom' E}^{\rm p}_2-{\vphantom' E}^{\rm p}_1$ must be sim­i­larly small as ${\vphantom' E}^{\rm p}_1-\mu$. But then, rea­son­ably as­sum­ing that the en­ergy lev­els are at least roughly equally spaced, and that the num­ber of states will not de­crease with en­ergy, so must $I_3$ be a fi­nite frac­tion of the to­tal, and so on. You can­not have a large num­ber of shelves each hav­ing a fi­nite frac­tion of the par­ti­cles, be­cause there are not so many par­ti­cles. More pre­cisely, a sum roughly like $\sum_{s=2}^\infty{\rm {const}}/s\Delta{E}$, (or worse), sums to an amount that is much larger than the term for $s$ $\vphantom0\raisebox{1.5pt}{$=$}$ 2 alone. So if $I_2$ would be a fi­nite frac­tion of $I$, then the sum would be much larger than $I$.

What hap­pens dur­ing con­den­sa­tion is that $\mu$ be­comes much closer to ${\vphantom' E}^{\rm p}_1$ than ${\vphantom' E}^{\rm p}_1$ is to the next en­ergy level ${\vphantom' E}^{\rm p}_2$, and only the ground state shelf ends up with a fi­nite frac­tion of the par­ti­cles. The re­main­der is spread out so much that the shelf num­bers im­me­di­ately above the ground state only con­tain a neg­li­gi­ble frac­tion of the par­ti­cles. It also fol­lows that for all shelves ex­cept the ground state one, $\mu$ may be ap­prox­i­mated as be­ing ${\vphantom' E}^{\rm p}_1$. (Spe­cific data for par­ti­cles in a box is given in chap­ter 11.14.1. The en­tire story may of course need to be mod­i­fied in the pres­ence of con­fine­ment, com­pare chap­ter 6.12.)

The other prob­lem with the analy­sis of the oc­cu­pa­tion num­bers for bosons is that the num­ber of sin­gle-par­ti­cle states on the shelves had to be at least two. There is no rea­son why a sys­tem of weakly-in­ter­act­ing spin­less bosons could not have a unique sin­gle-par­ti­cle ground state. And com­bin­ing the ground state with the next one on a sin­gle shelf is surely not an ac­cept­able ap­prox­i­ma­tion in the pres­ence of po­ten­tial Bose-Ein­stein con­den­sa­tion. For­tu­nately, the math­e­mat­ics still partly works:

$$\frac{\Delta F}{\Delta I_1} = \ln(I_1+1) - \ln(I_1+1) - \epsilon_1 -\epsilon_2 {\vphantom' E}^{\rm p}_1 = 0$$

im­plies that $\epsilon_1-\epsilon_2{\vphantom' E}^{\rm p}_1$ $\vphantom0\raisebox{1.5pt}{$=$}$ 0. In other words, $\mu$ is equal to the ground state en­ergy ${\vphantom' E}^{\rm p}_1$ ex­actly, rather than just ex­tremely closely as above.

That then is the con­densed state. With­out a chem­i­cal po­ten­tial that can be ad­justed, for any given tem­per­a­ture the states above the ground state con­tain a num­ber of par­ti­cles that is com­pletely un­re­lated to the ac­tual num­ber of par­ti­cles that is present. What­ever is left can be dumped into the ground state, since there is no con­straint on $I_1$.

Con­den­sa­tion stops when the num­ber of par­ti­cles in the states above the ground state wants to be­come larger than the ac­tual num­ber of par­ti­cles present. Now the math­e­mat­ics changes, be­cause na­ture says “Wait a minute, there is no such thing as a neg­a­tive num­ber of par­ti­cles in the ground state!” Na­ture now adds the con­straint that $I_1$ $\vphantom0\raisebox{1.5pt}{$=$}$ 0 rather than neg­a­tive. That adds an­other penalty term, $\epsilon_3I_1$ to $F$ and $\epsilon_3$ takes care of sat­is­fy­ing the equa­tion for the ground state shelf num­ber. It is a sad story, re­ally: be­low the con­den­sa­tion tem­per­a­ture, the ground state was awash in par­ti­cles, above it, it has zero. None.

A sys­tem of weakly in­ter­act­ing he­lium atoms, spin­less bosons, would have a unique sin­gle-par­ti­cle ground state like this. Since be­low the con­den­sa­tion tem­per­a­ture, the el­e­vated en­ergy states have no clue about an im­pend­ing lack of par­ti­cles ac­tu­ally present, phys­i­cal prop­er­ties such as the spe­cific heat stay an­a­lyt­i­cal un­til con­den­sa­tion ends.

It may be noted that above the con­den­sa­tion tem­per­a­ture it is only the most prob­a­ble set of the oc­cu­pa­tion num­bers that have ex­actly zero par­ti­cles in the unique ground state. The ex­pec­ta­tion value of the num­ber in the ground state will in­clude neigh­bor­ing sets of oc­cu­pa­tion num­bers to the most prob­a­ble one, and the num­ber has nowhere to go but up, com­pare {D.61}.