The ideal gas from partition functions: translational, rotational and vibrational parts

The dilute gas is the system on which statistical mechanics was first tested and where its predictions are most precise. Because the molecules hardly interact, the partition function of the gas is built from the partition function of a single molecule, and the single-molecule partition function factorizes into translational, rotational, vibrational, electronic and nuclear parts. Each part can be computed from the energy levels of quantum mechanics, and each makes its own contribution to the heat capacity, entropy and chemical potential. This lesson carries out the program for diatomic molecules. It explains Maxwell's puzzle of the missing degrees of freedom, and it ends with one of the most striking stories in the history of the subject, the heat capacity of hydrogen, which could not be explained until the spin of the proton was taken into account.

Translation

A molecule of mass \(m\) in a box of volume \(V\) has translational energy levels so closely spaced that the sum over them is effectively an integral over phase space. With one quantum state per \(h^3\) of phase space,

where \(\lambda = h/\sqrt{2\pi mkT}\) is the thermal wavelength from lesson 3. The result has a simple meaning: \(Z_{\text{tr}}\) is the number of boxes of size \(\lambda^3\) that fit into the container, roughly the number of translational states thermally accessible. For nitrogen at room temperature \(\lambda = 19\) pm, so a liter holds about \(1.4 \times 10^{29}\) such boxes.

For \(N\) molecules, \(Z = (Z_{\text{tr}}Z_{\text{int}})^N/N!\), where \(Z_{\text{int}}\) gathers all internal motions, and the free energy is

Differentiating with respect to \(V\) gives \(P = NkT/V\) whatever the internal structure: the ideal gas law depends only on translation. Differentiating with respect to \(T\) gives the entropy, which for \(Z_{\text{int}} = 1\) is the Sackur–Tetrode result of lesson 3, now obtained without any hypersphere. And differentiating with respect to \(N\) gives the chemical potential \(\mu = -kT\ln(VZ_{\text{int}}/N\lambda^3)\), which lesson 7 uses. Because \(\ln Z\) is a sum of logarithms, the internal motions contribute additively to \(U\), \(C_V\) and \(S\), and the translational part contributes \(\frac32 Nk\) to the heat capacity at every temperature at which the gas exists.

Rotation

A diatomic molecule rotating as a rigid body with moment of inertia \(I\) has energy levels \(E_J = J(J + 1)\hbar^2/2I\), each with degeneracy \(2J + 1\), as PHYS-372 derived. It is convenient to define a rotational temperature \(\theta_{\text{rot}} = \hbar^2/2Ik\), so that

Values of \(\theta_{\text{rot}}\) are small for heavy molecules: 2.88 K for \(\mathrm{N_2}\), 2.08 K for \(\mathrm{O_2}\), 15.2 K for HCl and 85.3 K for \(\mathrm{H_2}\), the lightest. When \(T \gg \theta_{\text{rot}}\), the sum can be replaced by an integral, giving \(Z_{\text{rot}} \approx T/\theta_{\text{rot}}\); the Euler–Maclaurin formula adds a correction of \(\frac13\). Then \(\bar E_{\text{rot}} = NkT\), the equipartition value for two rotational axes, and the rotational heat capacity is \(Nk\). The first part of the program compares the exact sum with \(T/\theta + \frac13\) for HCl; at 300 K they agree to two parts in ten thousand.

For a molecule made of two identical atoms, such as \(\mathrm{N_2}\) or \(\mathrm{H_2}\), rotating by 180 degrees gives an indistinguishable configuration, and quantum mechanics requires the total wave function to have a definite symmetry under exchanging the two nuclei. At high temperature the effect is simply to halve the number of rotational states: \(Z_{\text{rot}} \approx T/\sigma\theta_{\text{rot}}\) with symmetry number \(\sigma = 2\). This halving reduces the entropy by \(R\ln 2\) per mole but does not change the heat capacity. At low temperature, the symmetry requirement has much more dramatic consequences, as hydrogen shows.

Vibration

The vibration of the bond is close to a harmonic oscillator with frequency \(\omega\). Measuring energies from the ground vibrational level and defining \(\theta_{\text{vib}} = \hbar\omega/k\),

Vibrational temperatures are large: 6332 K for \(\mathrm{H_2}\), 4152 K for HCl (from its 2,886 per cm band, the same value used in PHYS-372), 3374 K for \(\mathrm{N_2}\), and only 309 K for the weakly bound, heavy \(\mathrm{I_2}\). The last part of the program shows the consequence. At 300 K, fewer than one nitrogen molecule in 70,000 is vibrationally excited and the vibrational heat capacity is negligible; at 1000 K it has reached 0.42 of its classical value. Iodine, by contrast, is almost fully classical at room temperature. Maxwell's puzzle of 1875 is solved: air at room temperature has \(C_V = \frac52 Nk\) because translation and rotation are fully excited and vibration is frozen out.

Electronic excitations usually lie several electronvolts up, tens of thousands of kelvin, so \(Z_{\text{el}}\) is simply the degeneracy of the electronic ground state. For most molecules this is 1; for \(\mathrm{O_2}\), whose ground state is a spin triplet, it is 3, which adds \(R\ln 3 = 9.1\) J/(mol K) to its entropy. Nitric oxide has a low-lying excited electronic level only about 170 K above the ground state, and its heat capacity shows a Schottky bump on top of the rotational contribution, a direct application of lesson 5.

import math

def zrot(T, theta, parity=None):
    """Rotational partition function and energy (units of k) for a rigid rotor.
    parity: None = all J, 0 = even J only, 1 = odd J only."""
    Z = U = U2 = 0.0
    for J in range(0, 200):
        if parity is not None and J % 2 != parity: continue
        e = theta * J * (J + 1)          # energy / k, in kelvin
        w = (2 * J + 1) * math.exp(-e / T)
        Z += w; U += e * w; U2 += e * e * w
    return Z, U / Z, U2 / Z - (U / Z) ** 2

def crot(T, theta, parity=None):
    Z, U, var = zrot(T, theta, parity)
    return var / T ** 2                  # C/(N k) from the fluctuation formula

# HCl: sum versus the high-temperature formula T/theta + 1/3
th = 15.2
for T in (30.0, 100.0, 300.0):
    Z = zrot(T, th)[0]
    print(f"HCl at {T:5.0f} K: Z_rot = {Z:8.4f}, T/theta + 1/3 = {T / th + 1 / 3:8.4f}")

# Hydrogen: theta_rot = 85.3 K; para (even J) and ortho (odd J, 3 nuclear spin states)
th = 85.3
print("\nRotational heat capacity of hydrogen, C_rot/(N k)")
print("    T (K)   para    ortho   normal 1:3 mix   equilibrium   rigid-rotor (no spin rule)")
for T in (25, 50, 75, 100, 150, 200, 300, 500):
    cp, co = crot(T, th, 0), crot(T, th, 1)
    normal = 0.25 * cp + 0.75 * co
    # equilibrium: one ensemble with even J weight 1 and odd J weight 3
    Z = U = U2 = 0.0
    for J in range(200):
        g = (1 if J % 2 == 0 else 3) * (2 * J + 1)
        e = th * J * (J + 1); w = g * math.exp(-e / T)
        Z += w; U += e * w; U2 += e * e * w
    ceq = (U2 / Z - (U / Z) ** 2) / T ** 2
    frac_ortho = sum(3 * (2 * J + 1) * math.exp(-th * J * (J + 1) / T) for J in range(1, 200, 2)) / Z
    print(f"  {T:6}   {cp:6.3f}   {co:6.3f}   {normal:8.3f}        {ceq:8.3f}      {crot(T, th):8.3f}"
          f"   (ortho fraction at equilibrium {frac_ortho:.3f})")

# Vibrations: fraction of molecules out of the ground state and C_vib at 300 K and 1000 K
print("\nVibrations: theta_vib in K, excited fraction exp(-theta/T), C_vib/(N k)")
for name, tv in (('H2', 6332.0), ('HCl', 4152.0), ('N2', 3374.0), ('I2', 309.0)):
    out = []
    for T in (300.0, 1000.0):
        x = tv / T
        c = x * x * math.exp(x) / (math.exp(x) - 1) ** 2
        out.append(f"T = {T:.0f} K: excited {math.exp(-x):.2e}, C = {c:.3f}")
    print(f"  {name:3} ({tv:6.0f} K): " + "; ".join(out))

HCl at    30 K: Z_rot =   2.3446, T/theta + 1/3 =   2.3070
HCl at   100 K: Z_rot =   6.9227, T/theta + 1/3 =   6.9123
HCl at   300 K: Z_rot =  20.0736, T/theta + 1/3 =  20.0702

Rotational heat capacity of hydrogen, C_rot/(N k)
    T (K)   para    ortho   normal 1:3 mix   equilibrium   rigid-rotor (no spin rule)
      25    0.000    0.000      0.000           0.447         0.151   (ortho fraction at equilibrium 0.010)
      50    0.019    0.000      0.005           2.067         0.970   (ortho fraction at equilibrium 0.229)
      75    0.250    0.003      0.065           1.382         1.092   (ortho fraction at equilibrium 0.479)
     100    0.739    0.033      0.210           0.884         1.047   (ortho fraction at equilibrium 0.614)
     150    1.425    0.252      0.545           0.687         1.013   (ortho fraction at equilibrium 0.714)
     200    1.394    0.562      0.770           0.792         1.006   (ortho fraction at equilibrium 0.740)
     300    1.088    0.916      0.959           0.959         1.002   (ortho fraction at equilibrium 0.749)
     500    1.002    1.000      1.000           1.000         1.001   (ortho fraction at equilibrium 0.750)

Vibrations: theta_vib in K, excited fraction exp(-theta/T), C_vib/(N k)
  H2  (  6332 K): T = 300 K: excited 6.82e-10, C = 0.000; T = 1000 K: excited 1.78e-03, C = 0.072
  HCl (  4152 K): T = 300 K: excited 9.76e-07, C = 0.000; T = 1000 K: excited 1.57e-02, C = 0.280
  N2  (  3374 K): T = 300 K: excited 1.31e-05, C = 0.002; T = 1000 K: excited 3.43e-02, C = 0.418
  I2  (   309 K): T = 300 K: excited 3.57e-01, C = 0.916; T = 1000 K: excited 7.34e-01, C = 0.992

The puzzle of hydrogen

In 1912 Arnold Eucken measured the heat capacity of hydrogen gas at low temperatures. As quantum theory predicted, rotation froze out at low temperature and \(C_V\) fell from \(\frac52 Nk\) toward \(\frac32 Nk\). But the shape of the fall did not match any calculation. The straightforward rotational partition function, summing over all \(J\), predicts the curve labeled rigid rotor in the program's output, with a small overshoot above 1 near 75 K. The measured curve has no overshoot and rises much more slowly. Through the early 1920s, the best physicists of the day tried different energy levels and different statistical weights without success.

The answer, found by David Dennison in 1927, lay in the nuclear spins. Each hydrogen nucleus is a proton with spin one half, so the two nuclei together have four spin states: three symmetric (total spin 1, the triplet) and one antisymmetric (total spin 0, the singlet). Protons are fermions, so the total wave function must change sign when the two nuclei are exchanged. Exchanging the nuclei multiplies the rotational wave function by \((-1)^J\). Therefore even-\(J\) rotational states must pair with the antisymmetric nuclear state and odd-\(J\) states with the three symmetric ones. Molecules with even \(J\) are called parahydrogen, and those with odd \(J\) orthohydrogen. At high temperature, three quarters of the molecules are ortho.

The crucial fact is that converting ortho to para requires flipping a nuclear spin, which in the gas, without a catalyst, takes far longer than any heat-capacity measurement. When hydrogen is cooled in a calorimeter, its ortho-to-para ratio stays frozen at the room-temperature value of 3:1. Normal hydrogen behaves as a mixture of two separate gases, and its heat capacity is a quarter of the para curve plus three quarters of the ortho curve. This is the normal 1:3 mixture column of the program, and it matches Eucken's data. The equilibrium column, where conversion is assumed to be fast, has a large peak near 50 K because converting ortho to para releases energy; it can be observed by adding a catalyst such as activated charcoal that speeds up the conversion. Pure parahydrogen, prepared with a catalyst at low temperature, shows the para curve with its large overshoot near 150 K.

The hydrogen story had consequences far beyond heat capacities. It established that the proton has spin one half and obeys Fermi statistics, and it showed that the fundamental postulate must be applied to the microstates that are actually accessible on the time scale of the experiment. Parahydrogen conversion matters in practice: liquid hydrogen used as rocket fuel is converted to para form with catalysts during liquefaction, because otherwise the slow ortho-to-para conversion in storage would release enough heat to boil off a large fraction of the liquid.

Putting it together

For a diatomic gas, the heat capacity per molecule climbs in steps as temperature rises: \(\frac32 k\) from translation at all temperatures, an extra \(k\) once \(T\) exceeds a few times \(\theta_{\text{rot}}\), and a further \(k\) once \(T\) approaches \(\theta_{\text{vib}}\). For nitrogen, the rotational step is complete above about 10 K (in practice nitrogen condenses first), and the vibrational step occurs between roughly 1000 and 5000 K. For polyatomic molecules there are three rotational axes and several vibrational modes, each with its own \(\theta_{\text{vib}}\), and the same method applies mode by mode. These partition functions are the foundation of the thermochemical tables used in combustion and atmospheric chemistry, and lesson 12 uses them to calculate an equilibrium constant from first principles.