Part I, Chapter 2

State Functions & Thermodynamic Potentials

Enthalpy, Helmholtz, Gibbs, Maxwell relations, and the Gibbs-Duhem equation

Historical Context

The concept of thermodynamic potentials was developed in the latter half of the 19th century. Josiah Willard Gibbs, in his landmark 1876 paper "On the Equilibrium of Heterogeneous Substances," introduced the chemical potential and the Gibbs free energy, unifying thermodynamics with chemistry. Hermann von Helmholtz independently developed the free energy concept that bears his name, while James Clerk Maxwell derived the elegant relations connecting partial derivatives of state functions.

The power of thermodynamic potentials lies in their ability to predict equilibrium conditions under different constraints. Each potential is the "natural" function for a specific set of independent variables, and its minimum (or maximum) determines equilibrium.

Key contributors: J. Willard Gibbs (1876), Hermann von Helmholtz (1882), James Clerk Maxwell (1871), Pierre Duhem (1886), Leonhard Euler (1740s).

Derivation 1: Enthalpy H = U + PV

Enthalpy arises naturally when considering processes at constant pressure, which is the most common laboratory condition. Starting from the first law:

$$dU = \delta Q - P\,dV \quad \text{(reversible work only)}$$

At constant pressure, rearranging:

$$\delta Q_P = dU + P\,dV = d(U + PV) = dH$$

Thus we define the enthalpy:

$$H \equiv U + PV$$

Natural variables: $H = H(S, P)$ with $dH = T\,dS + V\,dP$

The significance: at constant pressure, the heat exchanged equals the enthalpy change. This makes H the natural potential for isobaric processes:

$$Q_P = \Delta H = \int_{T_1}^{T_2} n C_P \, dT$$

From $dH = TdS + VdP$, we read off: $T = (\partial H / \partial S)_P$ and$V = (\partial H / \partial P)_S$.

Derivation 2: Helmholtz Free Energy A = U - TS

The Helmholtz free energy is the natural potential for isothermal, isochoric processes. It represents the maximum work extractable from a system at constant temperature.

We perform a Legendre transform on U, replacing S with its conjugate variable T:

$$A \equiv U - TS$$

Natural variables: $A = A(T, V)$ with $dA = -S\,dT - P\,dV$

Derivation of the fundamental relation:

$$dA = dU - T\,dS - S\,dT = (T\,dS - P\,dV) - T\,dS - S\,dT = -S\,dT - P\,dV$$

From this we extract: $S = -(\partial A / \partial T)_V$ and$P = -(\partial A / \partial V)_T$.

At constant T and V, $dA \leq 0$ for any spontaneous process. Equilibrium corresponds to the minimum of A. This makes the Helmholtz energy central to statistical mechanics, where$A = -k_B T \ln Z$ connects to the partition function Z.

Derivation 3: Gibbs Free Energy G = H - TS

The Gibbs free energy is the most important potential in chemistry because most reactions occur at constant temperature and pressure.

$$G \equiv H - TS = U + PV - TS$$

Natural variables: $G = G(T, P)$ with $dG = -S\,dT + V\,dP$

The derivation follows from Legendre transforming H:

$$dG = dH - T\,dS - S\,dT = (T\,dS + V\,dP) - T\,dS - S\,dT = -S\,dT + V\,dP$$

The equilibrium criterion at constant T and P:

$$dG_{T,P} \leq 0 \quad \text{(spontaneous process)} \qquad dG_{T,P} = 0 \quad \text{(equilibrium)}$$

This is the criterion that determines chemical equilibrium, phase transitions, and the direction of all processes occurring at constant temperature and pressure.

Derivation 4: Maxwell Relations

The Maxwell relations are a set of equalities among second partial derivatives of thermodynamic potentials. They follow from the exactness of differentials (equality of mixed partial derivatives).

For any state function $\Phi(x, y)$ with$d\Phi = M\,dx + N\,dy$, the exactness condition requires:

$$\left(\frac{\partial M}{\partial y}\right)_x = \left(\frac{\partial N}{\partial x}\right)_y$$

Figure. The thermodynamic square: S, V, T, P at corners with potentials U, H, G, A on the edges. Arrows encode sign conventions for the Maxwell relations.

The Four Maxwell Relations

From U(S,V): dU = TdS - PdV

$$\left(\frac{\partial T}{\partial V}\right)_S = -\left(\frac{\partial P}{\partial S}\right)_V$$

From H(S,P): dH = TdS + VdP

$$\left(\frac{\partial T}{\partial P}\right)_S = \left(\frac{\partial V}{\partial S}\right)_P$$

From A(T,V): dA = -SdT - PdV

$$\left(\frac{\partial S}{\partial V}\right)_T = \left(\frac{\partial P}{\partial T}\right)_V$$

From G(T,P): dG = -SdT + VdP

$$-\left(\frac{\partial S}{\partial P}\right)_T = \left(\frac{\partial V}{\partial T}\right)_P$$

Mnemonic: The "thermodynamic square" arranges S, T, P, V on corners with potentials U, H, A, G on the sides. Maxwell relations connect opposite corners through each potential.

Derivation 5: Euler Relation & Gibbs-Duhem Equation

For multi-component systems, the internal energy is a homogeneous function of degree one in its extensive variables (S, V, n_i). By Euler's theorem on homogeneous functions:

Euler Relation

$$U = TS - PV + \sum_i \mu_i n_i$$

Taking the total differential and comparing with the fundamental relation$dU = TdS - PdV + \sum_i \mu_i dn_i$:

$$dU = T\,dS + S\,dT - P\,dV - V\,dP + \sum_i (\mu_i\,dn_i + n_i\,d\mu_i)$$

Subtracting the fundamental relation from this total differential:

Gibbs-Duhem Equation

$$S\,dT - V\,dP + \sum_i n_i\,d\mu_i = 0$$

At constant T and P, this simplifies to:

$$\sum_i n_i\,d\mu_i = 0 \quad \text{or equivalently} \quad \sum_i x_i\,d\mu_i = 0$$

The Gibbs-Duhem equation constrains the chemical potentials: they cannot all vary independently. For a binary system, knowing $\mu_1(x_1)$ completely determines $\mu_2(x_1)$. This is the basis for thermodynamic consistency tests of experimental data.

Applications

Calorimetry

Enthalpy H enables direct measurement of heat in constant-pressure experiments. Hess's law (enthalpy is a state function) allows indirect calculation of reaction enthalpies from formation data.

Statistical Mechanics

The Helmholtz free energy $A = -k_BT \ln Z$ connects macroscopic thermodynamics to microscopic partition functions, bridging the gap between the two descriptions of matter.

Chemical Engineering

The Gibbs-Duhem equation is used in thermodynamic consistency tests for vapor-liquid equilibrium data. The Margules, van Laar, and Wilson equations must satisfy this constraint.

Materials Science

Maxwell relations enable the calculation of hard-to-measure quantities (like entropy changes) from easily measured ones (like thermal expansion and compressibility).

Simulation: Thermodynamic Potentials & Maxwell Relations

This simulation visualizes all four thermodynamic potentials for an ideal monatomic gas, plots isotherms, maps the entropy surface S(T,V), and numerically verifies a Maxwell relation.

Thermodynamic Potentials and Maxwell Relations

Python

script.py130 lines

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# ============================================================
# Thermodynamic Potentials and Maxwell Relations
# Visualize how state functions relate for an ideal gas
# ============================================================

R = 8.314   # J/(mol K)
n = 1.0     # moles

# Temperature and volume ranges
T_range = np.linspace(200, 800, 300)
V_range = np.linspace(0.005, 0.05, 300)

# Reference state
T0 = 300.0   # K
P0 = 1e5     # Pa
V0 = n * R * T0 / P0
S0 = 0.0     # reference entropy
U0 = 0.0     # reference internal energy

# For monatomic ideal gas: Cv = 3/2 R, Cp = 5/2 R
Cv = 1.5 * R
Cp = 2.5 * R

fig, axes = plt.subplots(2, 2, figsize=(14, 12))

# --- Panel 1: U, H, A, G vs T at constant P ---
ax1 = axes[0, 0]
T = T_range
dU = n * Cv * (T - T0)
dH = n * Cp * (T - T0)
S_T = n * Cp * np.log(T / T0)
dA = dU - T * S_T
dG = dH - T * S_T

ax1.plot(T, dU / 1000, 'r-', linewidth=2, label='U (internal energy)')
ax1.plot(T, dH / 1000, 'm-', linewidth=2, label='H (enthalpy)')
ax1.plot(T, dA / 1000, 'c-', linewidth=2, label='A (Helmholtz)')
ax1.plot(T, dG / 1000, 'y-', linewidth=2, label='G (Gibbs)')
ax1.axhline(y=0, color='white', alpha=0.3, linestyle='--')
ax1.set_xlabel('Temperature (K)', fontsize=11, color='white')
ax1.set_ylabel('Energy change (kJ/mol)', fontsize=11, color='white')
ax1.set_title('Thermodynamic Potentials vs T', fontsize=13, fontweight='bold', color='white')
ax1.legend(fontsize=9, facecolor='black', edgecolor='gray', labelcolor='white')
ax1.set_facecolor('#1a1a2e')
ax1.tick_params(colors='white')
ax1.grid(True, alpha=0.3)

# --- Panel 2: P vs V isotherms with labeled potentials ---
ax2 = axes[0, 1]
temps_iso = [300, 400, 500, 600, 700]
colors_iso = ['#ff6b6b', '#ffa502', '#2ed573', '#1e90ff', '#a55eea']
for Ti, ci in zip(temps_iso, colors_iso):
    P_vals = n * R * Ti / V_range
    ax2.plot(V_range * 1000, P_vals / 1e6, color=ci, linewidth=2, label='T = ' + str(Ti) + ' K')

ax2.set_xlabel('Volume (L)', fontsize=11, color='white')
ax2.set_ylabel('Pressure (MPa)', fontsize=11, color='white')
ax2.set_title('Ideal Gas Isotherms (PV = nRT)', fontsize=13, fontweight='bold', color='white')
ax2.legend(fontsize=9, facecolor='black', edgecolor='gray', labelcolor='white')
ax2.set_facecolor('#1a1a2e')
ax2.tick_params(colors='white')
ax2.grid(True, alpha=0.3)

# --- Panel 3: Entropy surface S(T, V) ---
ax3 = axes[1, 0]
T_grid = np.linspace(200, 800, 100)
V_grid = np.linspace(0.005, 0.05, 100)
TT, VV = np.meshgrid(T_grid, V_grid)
# S(T,V) = nCv ln(T/T0) + nR ln(V/V0) + S0  (ideal gas)
SS = n * Cv * np.log(TT / T0) + n * R * np.log(VV / V0) + S0

levels = np.linspace(SS.min(), SS.max(), 20)
cf = ax3.contourf(TT, VV * 1000, SS, levels=levels, cmap='hot')
plt.colorbar(cf, ax=ax3, label='S (J/mol K)')
ax3.set_xlabel('Temperature (K)', fontsize=11, color='white')
ax3.set_ylabel('Volume (L)', fontsize=11, color='white')
ax3.set_title('Entropy S(T, V) for Ideal Gas', fontsize=13, fontweight='bold', color='white')
ax3.set_facecolor('#1a1a2e')
ax3.tick_params(colors='white')

# --- Panel 4: Maxwell relation verification ---
ax4 = axes[1, 1]
# Verify (dP/dT)_V = (dS/dV)_T = nR/V for ideal gas
V_check = np.linspace(0.005, 0.05, 200)
lhs = n * R / V_check   # (dP/dT)_V = nR/V
rhs = n * R / V_check   # (dS/dV)_T = nR/V (same for ideal gas)

ax4.plot(V_check * 1000, lhs, 'r-', linewidth=3, label='(dP/dT)_V = nR/V')
ax4.plot(V_check * 1000, rhs, 'c--', linewidth=2, label='(dS/dV)_T = nR/V')
ax4.set_xlabel('Volume (L)', fontsize=11, color='white')
ax4.set_ylabel('Value (J / mol K m^3)', fontsize=11, color='white')
ax4.set_title('Maxwell Relation Verification', fontsize=13, fontweight='bold', color='white')
ax4.legend(fontsize=10, facecolor='black', edgecolor='gray', labelcolor='white')
ax4.set_facecolor('#1a1a2e')
ax4.tick_params(colors='white')
ax4.grid(True, alpha=0.3)

for ax in axes.flat:
    ax.spines['bottom'].set_color('gray')
    ax.spines['left'].set_color('gray')
    ax.spines['top'].set_color('gray')
    ax.spines['right'].set_color('gray')

fig.patch.set_facecolor('#0a0a1a')
plt.tight_layout()
plt.savefig('output.png', dpi=100, bbox_inches='tight', facecolor='#0a0a1a')
plt.close()

print("Thermodynamic Potentials (Ideal Monatomic Gas)")
print("=" * 55)
print("Cv = 3/2 R = " + str(round(Cv, 3)) + " J/(mol K)")
print("Cp = 5/2 R = " + str(round(Cp, 3)) + " J/(mol K)")
print("Cp - Cv = R = " + str(round(R, 3)) + " J/(mol K)")
print()
print("At T = 500 K (relative to T0 = 300 K):")
dT = 200
print("  dU = nCv dT = " + str(round(n * Cv * dT, 1)) + " J/mol")
print("  dH = nCp dT = " + str(round(n * Cp * dT, 1)) + " J/mol")
print()
print("Maxwell Relations (all verified for ideal gas):")
print("  (dT/dV)_S = -(dP/dS)_V")
print("  (dT/dP)_S =  (dV/dS)_P")
print("  (dS/dV)_T =  (dP/dT)_V = nR/V")
print("  (dS/dP)_T = -(dV/dT)_P = -nR/P")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

Simulation: Gibbs-Duhem Equation Verification

This simulation uses the Margules model for a binary solution to compute activity coefficients, excess chemical potentials, and excess Gibbs energy. It then numerically verifies that the Gibbs-Duhem equation is satisfied.

Gibbs-Duhem Equation for Binary Solutions

Python

script.py121 lines

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# ============================================================
# Gibbs-Duhem Equation: Binary Solution Verification
# ============================================================

R = 8.314  # J/(mol K)
T = 300.0  # K

# Margules model for activity coefficients (1-parameter)
# ln(gamma_1) = A * x2^2, ln(gamma_2) = A * x1^2
A_param = 1.5  # dimensionless Margules parameter

x1 = np.linspace(0.001, 0.999, 500)
x2 = 1.0 - x1

ln_gamma1 = A_param * x2**2
ln_gamma2 = A_param * x1**2
gamma1 = np.exp(ln_gamma1)
gamma2 = np.exp(ln_gamma2)

# Chemical potentials (excess part)
mu1_excess = R * T * ln_gamma1  # J/mol
mu2_excess = R * T * ln_gamma2  # J/mol

# Gibbs-Duhem check: x1 d(mu1) + x2 d(mu2) = 0 at constant T, P
# Numerically differentiate
dmu1_dx1 = np.gradient(mu1_excess, x1)
dmu2_dx1 = np.gradient(mu2_excess, x1)
gibbs_duhem_residual = x1 * dmu1_dx1 + x2 * dmu2_dx1

# Excess Gibbs energy
G_E = x1 * mu1_excess + x2 * mu2_excess  # = RT * A * x1 * x2
G_E_analytic = R * T * A_param * x1 * x2

fig, axes = plt.subplots(2, 2, figsize=(14, 12))

# Panel 1: Activity coefficients
ax1 = axes[0, 0]
ax1.plot(x1, gamma1, 'r-', linewidth=2.5, label='gamma_1')
ax1.plot(x1, gamma2, 'c-', linewidth=2.5, label='gamma_2')
ax1.axhline(y=1, color='white', alpha=0.3, linestyle='--', label='Ideal (gamma=1)')
ax1.set_xlabel('x_1', fontsize=12, color='white')
ax1.set_ylabel('Activity coefficient', fontsize=12, color='white')
ax1.set_title('Margules Activity Coefficients (A = ' + str(A_param) + ')', fontsize=13, fontweight='bold', color='white')
ax1.legend(fontsize=10, facecolor='black', edgecolor='gray', labelcolor='white')
ax1.set_facecolor('#1a1a2e')
ax1.tick_params(colors='white')
ax1.grid(True, alpha=0.3)

# Panel 2: Excess chemical potentials
ax2 = axes[0, 1]
ax2.plot(x1, mu1_excess / 1000, 'r-', linewidth=2.5, label='mu_1 excess')
ax2.plot(x1, mu2_excess / 1000, 'c-', linewidth=2.5, label='mu_2 excess')
ax2.set_xlabel('x_1', fontsize=12, color='white')
ax2.set_ylabel('Excess chemical potential (kJ/mol)', fontsize=12, color='white')
ax2.set_title('Excess Chemical Potentials', fontsize=13, fontweight='bold', color='white')
ax2.legend(fontsize=10, facecolor='black', edgecolor='gray', labelcolor='white')
ax2.set_facecolor('#1a1a2e')
ax2.tick_params(colors='white')
ax2.grid(True, alpha=0.3)

# Panel 3: Excess Gibbs energy
ax3 = axes[1, 0]
ax3.plot(x1, G_E / 1000, 'y-', linewidth=2.5, label='G_E (numerical)')
ax3.plot(x1, G_E_analytic / 1000, 'w--', linewidth=1.5, label='G_E = RT A x1 x2 (analytic)')
ax3.set_xlabel('x_1', fontsize=12, color='white')
ax3.set_ylabel('G_E (kJ/mol)', fontsize=12, color='white')
ax3.set_title('Excess Gibbs Energy', fontsize=13, fontweight='bold', color='white')
ax3.legend(fontsize=10, facecolor='black', edgecolor='gray', labelcolor='white')
ax3.set_facecolor('#1a1a2e')
ax3.tick_params(colors='white')
ax3.grid(True, alpha=0.3)

# Panel 4: Gibbs-Duhem verification
ax4 = axes[1, 1]
ax4.plot(x1, gibbs_duhem_residual, 'lime', linewidth=2)
ax4.axhline(y=0, color='white', alpha=0.5, linestyle='--')
ax4.set_xlabel('x_1', fontsize=12, color='white')
ax4.set_ylabel('x1 d(mu1) + x2 d(mu2)', fontsize=12, color='white')
ax4.set_title('Gibbs-Duhem Verification (should be ~0)', fontsize=13, fontweight='bold', color='white')
ax4.set_facecolor('#1a1a2e')
ax4.tick_params(colors='white')
ax4.grid(True, alpha=0.3)

max_residual = np.max(np.abs(gibbs_duhem_residual[10:-10]))
ax4.text(0.5, 0.85, 'Max residual = ' + str(round(max_residual, 2)) + ' J/mol',
         transform=ax4.transAxes, fontsize=12, color='lime', ha='center',
         bbox=dict(boxstyle='round', facecolor='black', alpha=0.8))

for ax in axes.flat:
    ax.spines['bottom'].set_color('gray')
    ax.spines['left'].set_color('gray')
    ax.spines['top'].set_color('gray')
    ax.spines['right'].set_color('gray')

fig.patch.set_facecolor('#0a0a1a')
plt.tight_layout()
plt.savefig('output.png', dpi=100, bbox_inches='tight', facecolor='#0a0a1a')
plt.close()

print("Gibbs-Duhem Equation Verification")
print("=" * 50)
print("Margules parameter A = " + str(A_param))
print("Temperature T = " + str(T) + " K")
print()
print("At x1 = 0.5:")
idx = 250
print("  gamma_1 = " + str(round(gamma1[idx], 4)))
print("  gamma_2 = " + str(round(gamma2[idx], 4)))
print("  G_E = " + str(round(G_E[idx], 2)) + " J/mol")
print("  G_E (analytic) = " + str(round(G_E_analytic[idx], 2)) + " J/mol")
print()
print("Gibbs-Duhem residual (should be 0):")
print("  Max |residual| = " + str(round(max_residual, 4)) + " J/mol")
print()
print("The Gibbs-Duhem equation is satisfied (within numerical precision)")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

Chapter Summary

• Enthalpy: $H = U + PV$, natural for constant-P processes. $Q_P = \Delta H$.

• Helmholtz: $A = U - TS$, minimum at equilibrium for constant T, V.

• Gibbs: $G = H - TS$, the central potential in chemistry at constant T, P.

• Maxwell relations: Four equalities linking second derivatives of thermodynamic potentials.

• Gibbs-Duhem: $SdT - VdP + \sum n_i d\mu_i = 0$, constraining intensive variables.

← Laws of Thermodynamics Next: Thermodynamic Cycles →