7.1.1 Equilibrium conditions

We start by asking whether there is a limit to the number of phases that can coexist. An elegant answer is provided by the Gibbs phase rule, Eq. (7.13). The chemical potential of substances in coexisting phases has the same value in each of the phases in which coexistence occurs.⁵ Consider two phases of a substance, I and II. Because matter and energy can be exchanged between phases in physical contact, equilibrium is achieved when T and P are the same in both phases, and when the chemical potentials are equal, μI=μII (see Section 1.12). We know from the Gibbs-Duhem equation,⁶ (P1.1), that μ=μ(T,P), and thus chemical potential can be visualized as a surface μ=μ(T,P) (see Fig. 7.2). Two phases of the same substance coexist when

μI(T,P)=μII(T,P) .

(7.1)

The intersection of the two surfaces defines the locus of points P=P(T) for which Eq. (7.1) is satisfied—the coexistence curve (see Fig. 7.2). Three coexisting phases (I,II,III) would require the equality of three chemical potential functions,

μI(T,P)=μII(T,P)=μIII(T,P) .

(7.2)

Equation (7.2) implies two equations in two unknowns and thus three phases can coexist at a unique combination of T and P, the triple point. By this reasoning, it would not be possible for four phases of a single substance to coexist (which would require three equations in two unknowns). Coexistence of four phases of the same substance is not known to occur.

Multicomponent phases have more than one chemical species. Let μjγ denote the chemical potential of species j in the γ phase (we use Roman letters to label species and Greek letters to label phases). Assume k chemical species, 1≤j≤k,, and π phases, 1≤γ≤π,, where π is an integer. The first law for multiphase, multicomponent systems is the generalization of Eq. (1.21):

dU=TdS−PdV+∑γ=1π∑j=1kμjγdNjγ .

(7.3)

Extensivity implies the scaling property U(λS,λV,λNjγ)=λU(S,V,Njγ). By Euler’s theorem,

U=S∂U∂SV,Njγ+V∂U∂VS,Njγ+∑γ=1π∑j=1kNjγ∂U∂NjγS,V=TS−PV+∑γ=1π∑j=1kμjγNjγ ,

(7.4)

where the derivatives follow from Eq. (7.3) and Njγ¯ indicates to hold fixed particle numbers except Njγ. Equation (7.4) implies

G=∑γ=1π∑j=1kμjγNjγ ,

(7.5)

the generalization of Eq. (1.54). Taking the differential of Eq. (7.4) and using Eq. (7.3), we have the multicomponent, multiphase generalization of the Gibbs-Duhem equation, (P1.1),

∑γ=1π∑j=1kNjγdμjγ=−SdT+VdP .

(7.6)

By taking the differential of G in Eq. (7.5), and making use of Eq. (7.6),

dG=−SdT+VdP+∑γ=1π∑j=1kμjγdNjγ ,

(7.7)

and thus we have the alternate definition of chemical potential,⁷ μjγ=∂G/∂NjγT,P,Njγ¯, the energy to add a particle of type j in phase γ holding fixed T, P, and the other particle numbers.

The Gibbs energy is a minimum in equilibrium.⁸ From Eq. (7.7),

dGT,P=∑γ=1π∑j=1kμjγdNjγT,P=0 .

(7.8)

If the particle numbers Njγ could be independently varied, one would conclude from Eq. (7.8) that μjγ=0. But particle numbers are not independent. The number of particles of each species spread among the phases is a constant,⁹ ∑γNjγ= constant, and thus there are k equations of constraint

∑γ=1πdNjγ=0 .  j=1,⋯,k

(7.9)

Constraints are handled through the method of Lagrange multipliers¹⁰ λ_j, which when multiplied by Eq. (7.9) and added to Eq. (7.8) leads to

∑j=1k∑γ=1πμjγ+λjdNjγ=0 .

(7.10)

We can now treat the particle numbers as unconstrained, so that Eq. (7.10) implies

μjγ=−λj .

(7.11)

The chemical potential of each species is independent of phase. Equation (7.11) is equivalent to

μj1=μj2=⋯=μjπ .    j=1,⋯,k

(7.12)

There are k(π−1) equations of equilibrium for k chemical components in π phases.

7.1.2 The Gibbs phase rule

How many independent state variables can exist in a multicomponent, multiphase system? In each phase there are Nγ≡∑j=1kNjγ particles, and thus there are k−1 independent concentrations cjγ≡Njγ/Nγ, where ∑j=1kcjγ=1. Among π phases there are π(k−1) independent concentrations. Including P and T, there are 2+π(k−1) independent intensive variables.

There are k(π−1) equations of equilibrium, Eq. (7.12). The variance of the system is the difference between the number of independent variables and the number of equations of equilibrium,

f≡2+π(k−1)−k(π−1)=2+k−π .

(7.13)

Equation (7.13) is the Gibbs phase rule.[11, p96] It specifies the number of intensive variables that can be independently varied without disturbing the number of coexisting phases ( f≥0).

• k=1,π=1⇒f=2: a single substance in one phase. Two intensive variables can be independently varied; T and P in a gas.

• k=2,π=1⇒f=3: two substances in a single phase, as in a mixture of gases. We can independently vary T, P, and one mole fraction.

• k=1,π=2⇒f=1: a single substance in two phases; a single intensive variable such as the density can be varied without disrupting phase coexistence.

• k=1,π=3⇒f=0: a single substance in three phases; we cannot vary the conditions under which three phases coexist in equilibrium. Unique values of T and P define a triple point.

One should appreciate the generality of the phase rule, which doesn’t depend on the type of chemical components, only that the Gibbs energy is a minimum in equilibrium.

7.1.3 The Clausius-Clapeyron equation

The latent heat, L, is the heat released or absorbed during a phase change, and is measured either as a molar quantity (per mole) or as a specific quantity (per mass). One has the latent heat of vaporization (boiling), fusion (melting), and sublimation.¹¹ Latent heats are also called the enthalpy of vaporization (or fusion or sublimation).¹² At a given T,

L(T)=hv−hl=T sv(T,P(T))−sl(T,P(T)) ,

(7.14)

where v and l refer to vapor and liquid, lower-case quantities such as s≡S/n indicate molar values, and P(T) describes the coexistence curve.¹³ The difference in molar entropy between phases is denoted Δs≡sv−sl; likewise with Δh. Equation (7.14) is written compactly as L=Δh=TΔs.

As we move along a coexistence curve, T and P vary in such a way as to maintain the equality μI(T,P)=μII(T,P). Variations in T and P induce changes δμ in the chemical potential, and thus along the coexistence curve, δμI=δμII. For a single substance, dμ=−sdT+vdP (Gibbs-Duhem equation). At a coexistence curve, −sIdT+vIdP=−sIIdT+vIIdP, implying

dPdTcoexist=sI−sIIvI−vII=ΔsΔv=ΔhTΔv=LTΔv .

(7.15)

Equation (7.15) is the Clausius-Clapeyron equation; it tells us the local slope of the coexistence curve. If one had enough data for L as a function of T, P and the volume change Δv, Eq. (7.15) could be integrated to obtain the coexistence curve in a P–T diagram.

7.3.1 The van der Waals critical point

Thus, we pass from the low temperature, high density system (low specific volume), in which Eq. (7.18) has one real root, to the high temperature, low density system (large specific volume) which has one nontrivial root (ideal gas law). We shouldn’t be surprised if in between there is a regime for which Eq. (7.18) has multiple roots (which hopefully can be interpreted physically). Detailed studies show there is a unique set of parameters P_c, T_c, and (V/N)c (which as the notation suggests are the pressure, temperature, and specific volume of a critical point¹⁸), at which Eq. (7.18) has a triple root. For (V/N)c≡vc a triple root of Eq. (7.18), the cubic equation

v−vc3=v3−3vcv2+3vc2v−vc3=0

(7.20)

must be equivalent to Eq. (7.18) evaluated at P_c and T_c. Comparing coefficients of identical powers of v in Eqs. (7.18) and (7.20), we infer the correspondences

3vc=b+kTcPc    3vc2=aPc    vc3=abPc .

These equations imply the values of the critical parameters (show this):

vc=3b    Pc=a27b2    kTc=8a27b .

(7.21)

The relations in Eq. (7.21) predict P_c and T_c reasonably well, but overestimate v_c (Exercise 7.5). Equation (7.21) also predicts Pcvc/(kTc)=38 (show this), a universal number, the same for all gases. In actuality, Pcvc/(kTc)≈0.2 –0.3 for real gases.¹⁹ Should we be discouraged by this disagreement? Not at all! The van der Waals model is the simplest model of interacting gases—if we want better agreement with experiment, we have to include more relevant physics in the model.

7.3.2 The law of corresponding states and the Maxwell equal area rule

Express P, T, and v in units of the critical parameters. Let

P¯≡P/Pc    T¯≡T/Tc    v¯≡v/vc .

(7.22)

The dimensionless quantities P¯,T¯,v¯ are known as the reduced state variables of the system. With P=PcP¯, T=TcT¯, and v=vcv¯, the van der Waals equation of state can be written (show this)

(P¯+3v¯2) (3v¯−1)=8T¯.

(7.23)

Equation (7.23) is a remarkable development: It’s independent of material-specific parameters, indicating that when P, T, v are scaled in units of P_c, T_c, v_c (which are material specific), the equation of state is the same. Systems having their own values of (a,b), yet the same values of P¯, T¯, v¯, are said to be in corresponding states. Equation (7.23) is the law of corresponding states: Systems in corresponding states behave the same. Thus, argon at T=300 K, for which Tc=151 K, will behave the same as ethane (C₂H₆) at T=600 K, for which Tc=305 K (both have T¯≈2).²⁰ Figure 7.3 is a compelling illustration of the law of corresponding states. For ρ_l (ρ_g) the density of liquid (vapor) in equilibrium at T<Tc, by the law of corresponding states we expect ρl/ρc and ρg/ρc to be universal functions of T/Tc, where ρ_c denotes the critical density. Figure 7.3 shows data for eight substances (of varying complexity²¹ from noble-gas atoms Ne, Ar, Kr, and Xe, to diatomic molecules N₂, O₂, and CO, to methane, CH₄), which, when plotted in terms of reduced variables, ostensibly fall on the same curve![83] The solid line is a fit to the data points of Ar (solid circles) assuming the critical exponent β=13 (see Section 7.4). Figure 7.3 shows the coexistence region—values of temperature and density for which phase coexistence occurs.²²

Because the equation of state is predicted to be the same for all systems (when expressed in reduced variables), any thermodynamic properties derived from it are also predicted to be the same. Universality is the term used to indicate that the behavior of systems near critical points is largely independent of the details of the system.²³ Universality has emerged as a key aspect of modern theories of critical phenomena, as we’ll see in Section 8.4.²⁴

Figure 7.4 shows the isotherms calculated from Eq. (7.23). For T≥Tc, there is a unique specific volume for every pressure—the equation of state is single valued.²⁵ For T<Tc, however, there are three possible volumes associated with a given pressure, such as the values of v/vc shown at A,B,C for T=0.9Tc. How do we interpret these multiple solutions?

Figure 7.5 shows an expanded view of an isotherm for T<Tc. We note that (∂P/∂v)T is positive along the segment bcd, which is unphysical. An increase in pressure upon an isothermal expansion would imply a violation of the second law of thermodynamics. The root c in Fig. 7.5 (or B in Fig. 7.4) is unphysical and can be discarded.

What to make of the other two roots? Figure 7.6 indicates what we expect phenomenologically and hence what we want from a model of phase transitions. As a gas is compressed, pressure rises.²⁶ When (for T<Tc) the gas has been compressed to point B, condensation first occurs (coalescence of gas-phase atoms into liquid-phase clusters). At B, the specific volume v_g represents the maximum density a gas can have at T<Tc. Upon further isothermal compression, the pressure remains constant along the horizontal segment AB (at the saturation pressure, Psat(T)) as more gas-phase atoms are absorbed into the liquid phase.²⁷ At point A, when all gas-phase atoms have condensed into the liquid, the specific volume v=vl is the minimum density a liquid can have at T<Tc. Further compression at this point results in a rapid rise in pressure.

At a point in the coexistence region, such as point D in Fig. 7.6, there are N_l atoms in the liquid phase and N_g in the vapor phase, where N=Nl+Ng is the total number of atoms. Let xl≡Nl/N and xg≡Ng/N denote the mole fractions of the amount of material in the liquid and gas phases. Clearly xl+xg=1. The volume occupied by the gas (liquid) phase at this temperature is Vg=Ngvg ( Vl=Nlvl), implying V=Ngvg+Nlvl. The specific volume v_d associated with point D is:²⁸

vd=VN=Ngvg+NlvlN=xgvg+xlvl=vd (xg+xl),

(7.24)

implying

xlvd−vl=xgvg−vd .

(7.25)

Equation (7.25) is the lever rule.²⁹ As vd→vg, xl→0, and vice versa. See Exercise 7.7.

Note that points A and B in Fig. 7.6, which occur at separate volumes in a P–v phase diagram, are, in a P–T phase diagram (such as Fig. 7.1), located at the same point on a coexistence curve. Along AB in Fig. 7.6, μl(T,P,vl)=μg(T,P,vg), and therefore from Ndμ=−SdT+VdP (Gibbs-Duhem), dμ=0 (phase coexistence) and dT=0 (isotherm) implies dP=0. A horizontal, subcritical isotherm in a P–v diagram indicates liquid-gas coexistence. From the Gibbs phase rule, for a single substance in two phases, one intensive variable (the specific volume) can be varied and not disrupt phase coexistence.

The van der Waals model clearly does not have isotherms featuring horizontal segments. In 1875, Maxwell suggested a way to reconcile undulating isotherms with the requirement of phase coexistence [85]. From the Gibbs-Duhem equation applied to an isotherm ( dT=0), Ndμ=VdP, and thus for an isotherm to represent phase coexistence, we require

0=μl−μg=∫lgdμ=∫lgvdP=∫lg[d (Pv)−Pdv]=(Pv)|lg−∫lgPdv ,

(7.26)

implying

Psat(T) (vg−vl)=∫lgP(v)dv .

(7.27)

Equation (7.27) is the Maxwell equal area rule. It indicates that v_l, v_g, and Psat(T) should be chosen so that the area under the curve of an undulating isotherm is the same as the area under a horizontal isotherm.³⁰ Maxwell in essence modified the van der Waals model so that it describes phase coexistence:

P(v)=kTv−b−av2T≥Tc,all v;T<Tc,v>vg,v<vlPsat(T)T<Tc,vl<v<vg .

(7.28)

The Maxwell rule is an ad hoc fix that allows one to locate the coexistence region in a P–v diagram. One could object that it uses an unphysical isotherm to draw physical conclusions. It’s not a theory of phase transitions in the sense that condensation is shown to occur from the partition function.³¹ In an important study, Kac, Uhlenbeck, and Hemmer showed rigorously that a one-dimensional model with a weak, long-range interaction (in addition to a hard-core, short-range potential), has the van der Waals equation of state together with the Maxwell construction.[88] This work was generalized by Lebowitz and Penrose to an arbitrary number of dimensions.[89] Models featuring weak, long-range interactions are known as mean field theories (Section 7.9); the van der Waals-Maxwell model is within a class of theories, mean field theories. The energy per particle of the van der Waals gas is modified from its ideal-gas value ( 32kT) to include a contribution that’s proportional to the density, −an, the hallmark of a mean-field theory; see Exercise 7.8.

7.3.3 Heat capacity of coexisting liquid-gas mixtures

At a point in the coexistence region (such as D in Fig. 7.6), the total internal energy of the system,

U=Nlul+Ngug ,

(7.29)

where u_l (u_g) is the specific energy, the energy per particle u≡U/N, of the liquid (gas) phase at the same pressure and temperature as the coexisting phases.³² Dividing Eq. (7.29) by N,

u=xlul+xgug .

(7.30)

Compare Eq. (7.30) with Eq. (7.24), which has the same form.

To calculate the heat capacity it suffices to calculate the specific heat, c≡C/N. Starting from CV=∂U/∂TV,N (Eq. (1.38)), the specific heat cv≡CV/N=∂u/∂Tv, where v=V/N. One might think it would be a simple matter to find c_v by differentiating Eq. (7.30) with respect to T, presuming that the mole fractions xl,xg are temperature independent. The mole fractions, however, vary with temperature at fixed volume. As one increases the temperature along the line associated with fixed volume v_d in Fig. 7.6, xl,xg vary because the shape of the coexistence region changes with temperature; consider the lever rule, Eq. (7.25). Of course, dxl=−dxg because of the constraint xl+xg=1. We have, using Eq. (7.30),

cv=(∂u∂T)v=xl(∂ul∂T)coex+xg(∂ug∂T)coex−(ug−ul)(∂xl∂T)coex ,

(7.31)

where ∂/∂Tcoex indicates a derivative taken at the sides of the coexistence region.³³ Evaluating these derivatives is relegated to Exercise 7.10. The final result is:

cv=xlcvl−T∂P∂vlT∂vl∂Tcoex2+xgcvg−T∂P∂vgT∂vg∂Tcoex2 .

(7.32)

Equation (7.32) is a general result, not specific to the van der Waals model.

7.4.1 Shape of the critical coexistence region, the exponent β

The Maxwell rule can be used to infer an important property of the volumes v_l, v_g in the vicinity of the critical point. From its form in Eq. (7.26), ∫lgvdP=0,

0=∫lgv¯dP¯=∫lg(v¯−1+1) dP¯=∫lgϕdP¯+∫lgdP¯=0∫lgϕ(∂P¯∂ϕ)tdϕ=−∫ϕlϕgϕ (6t+92ϕ2) dϕ ,

(7.35)

where ∫PlPgdP=0 for coexisting phases, we’ve used Eqs. (7.33) and (7.34), and we’ve retained terms up to first order in |t| in ∂P¯/∂ϕt ( ϕ~|t|). The integrand of the final integral in Eq. (7.35) is an odd function of ϕ under ϕ→−ϕ, and the simplest way to ensure the vanishing of this integral for any |t|≪1 is to take ϕl=−ϕg, i.e., the limits of integration are symmetrically placed about ϕ=0. We can make this conclusion only sufficiently close to the critical point where higher-order terms in Eq. (7.34) can be neglected; close to the critical point the coexistence region is symmetric about v_c (which is not true in general—see Fig. 7.3).³⁵

The symmetry P(t,−ϕg)=P(t,ϕg) applied to Eq. (7.34) implies 3ϕg−4|t|+ϕg2=0, or

ϕg=2|t|~T−Tc1/2 .    |t|→0

(7.36)

In the van der Waals model, therefore, vg→vc with the square root of T−Tc, and thus it predicts β=12. The measured value of β for liquid-gas critical points is β≈0.32 –0.34[90, 91]. We see from Fig. 7.3 that the coexistence region is well described using β=13.

7.4.2 The heat capacity exponent, α

How does the specific heat behave in the critical region? Let’s approach the critical point along the critical isochore, a line of constant volume, v=vc. As T→Tc−, cvg,cvl→cvc, and xl,xg→12 (use the lever rule and the symmetry of the critical region, vc−vl=vg−vc). From Eq. (7.32),

cv(Tc−)=cv(Tc+)−Tc2limT→Tc− [(∂P∂vl)T(∂vl∂T)coex2+(∂P∂vg)T(∂vg∂T)coex2] ,

(7.37)

where we’ve set cv(Tc+)≡12 (cvl+cvg) as the limit T→Tc+, i.e., the terms in square brackets disappear for T>Tc. For the van der Waals gas, c_v has the ideal-gas value for T>Tc, cv0≡32k (see Exercise 7.8). Equation (7.37) is a general result, and is not specific to the van der Waals model.

For the van der Waals model, we find using Eq. (7.36),

∂v∂Tcoex=∓vcTc1|t| ,

(7.38)

where the upper (lower) sign is for gas (liquid). For the other derivative in Eq. (7.37), we can use Eq. (7.34),

∂P∂vT=Pcvc∂P¯∂ϕt=Pcvc6|t|−18|t|ϕ−92ϕ2=−12Pcvc|t|+O(|t|3/2) ,

(7.39)

where we’ve used Eq. (7.36). From Eq. (7.37), taking the limit,

cv(Tc−)−cv0=−Tc−12Pcvc|t|vcTc21|t|=12PcvcTc=92k ,

where we’ve used Eqs. (7.38) and (7.39), and that Pcvc/Tc=38k (see Eq. (7.21)). Thus,

cv−cv0=92kT→Tc−0T→Tc+ .

(7.40)

The specific heat is discontinuous at the van der Waals critical point—a second-order phase transition in the Ehrenfest classification scheme.

For many systems C_V does not show a discontinuity at T=Tc, it diverges as T→Tc (from above or below),³⁶ such as in liquid ⁴He at the superfluid transition,³⁷ the so-called “λ-transition.” To cover these cases, a heat capacity critical exponent α is introduced,³⁸

CV~T−Tc−α .    T→Tc

(7.41)

The value α=0 is assigned to systems showing a discontinuity in C_V at T_c, such as for the van der Waals model. The measured value of α for liquid-gas phase transitions is α≈0.1. Extracting an exponent from measurements can be difficult. If a thermodynamic quantity f(t) shows singular behavior as t→0 in the form f(t)~Atx, the exponent is obtained from the limit

x≡lim⁡t→0ln⁡f(t)ln⁡t .

(7.42)

If x<0, f(t) diverges at the critical point; if x>0, f(t)→0 as T→Tc.

From thermodynamics, CP>CV, Eq. (1.42). The heat capacity C_V for the van der Waals gas is the same as for the ideal gas (Exercise 7.8). C_P, however, has considerable structure for T>Tc; see Eq. (P7.4). In terms of reduced variables, Eq. (P7.4) is equivalent to (show this):

CP−CVNk=(1+ϕ)3(1+t)34ϕ2+ϕ3+t(1+ϕ)3 .

(7.43)

Because ϕ→0 as |t|, Eq. (7.36), we see that C_P diverges as T→Tc+, with

CP~T−Tc−1 .    (T→Tc+)

(7.44)

Equation (7.44) is valid for T>Tc; C_P isn’t well defined in the coexistence region.³⁹

7.4.3 The compressibility exponent, γ

We’ve introduced two exponents that characterize critical behavior: α for the heat capacity and β for vg−vl. Is there a γ? One can show that at the van der Waals critical point,

∂P∂v|T=Tc,v=vc=∂2P∂v2|T=Tc,v=vc=0 .

(7.45)

The critical isotherm therefore has an inflection point at v=vc, as we see in Figs. 7.4 or 7.6. The two conditions in Eq. (7.45) define critical points in fluid systems. One can show that

(∂P∂v)t=−2a27b3 [34ϕ2+ϕ3+t(1+ϕ)3(1+ϕ)3(1+32ϕ)2] .

(7.46)

Thus, (∂P/∂v)t~t as t→0. Its inverse, the compressibility βT=−(1/v)(∂v/∂P)T therefore diverges at the critical point,⁴⁰which is characterized by the exponent γ,

βT~T−Tc−γ .

(7.47)

The van der Waals model predicts γ=1; experimentally its value is closer to γ=1.25.

7.4.4 The critical isotherm exponent, δ

The equation of state on the critical isotherm is obtained by setting t=0 in Eq. (7.34), P¯≈1−32ϕ3, implying

P¯−1=P−PcPc=−32v−vcvc3~v−vcvcδ ,

(7.48)

where δ is the conventional symbol for the critical isotherm exponent. The van der Waals model predicts δ=3, whereas for real fluids [91], δ≈4.7 –5.

The definition of the critical exponents α,β,γ,δ and their values are summarized in Table 7.1.

7.5.1 The magnetization exponent, β

To find β is easy—it’s implied by Eq. (7.51) for T<Tc. In the limit T→Tc−,

M2=(53N2μ2) (Tc−TTc) ~ (Tc−TTc)2β ,    (T→Tc−)

(7.54)

implying β=12 for the Weiss model. Experimental values of β for magnetic systems are in the range 0.30–0.325[94, 95].

7.5.2 The critical isotherm exponent, δ

At T=Tc, M→0 as B→0; the vanishing of M in characterized by the critical exponent δ, defined as B~Mδ at T=Tc. To extract the exponent, we need to “dig out” B from the equation of state, Eq. (7.49), for small M. Combining Eq. (7.52) with Eq. (7.49), where y=M/(Nμ),

y=LβμB+3TcTy .

(7.55)

We can isolate the field term B in Eq. (7.55) by invoking the inverse Langevin function L−1(y):

βμB=L−1(y)−3TcTy .

(7.56)

A power series representation of L−1 can be derived from the Lagrange inversion theorem [50, p14], where for −1<y<1,

L−1(y)=∑k=1∞ykk!dk−1dxk−1xL(x)kx=0≡∑k=1∞ckyk=3y+95y3+297175y5+O(y7) .