Featured post

Welcome to Physics Pages

This blog consists of my notes and solutions to problems in various areas of mainstream physics. An index to the topics covered is contained in the links in the sidebar on the right, or in the menu at the top of the page.

This isn’t a “popular science” site, in that most posts use a fair bit of mathematics to explain their concepts. Thus this blog aims mainly to help those who are learning or reviewing physics in depth. More details on what the site contains and how to use it are on the welcome page.

Despite Stephen Hawking’s caution that every equation included in a book (or, I suppose in a blog) would halve the readership, this blog has proved very popular since its inception in December 2010. Details of the number of visits and distinct visitors are given on the hit statistics page.

Many thanks to my loyal followers and best wishes to everyone who visits. I hope you find it useful. Constructive criticism (or even praise) is always welcome, so feel free to leave a comment in response to any of the posts.

Before leaving a comment, you may find it useful to read the “Instructions for commenters“.

Mixed systems: effect of interaction energy and the solubility gap

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.58.

So far when we’ve examined the Gibbs free energy and entropy of a mixture of two pure substances, we’ve assumed that the interaction energy is the same between all molecules in the mixture, that is, that a species {A} molecule has the same interaction energy with another {A} molecule as it does with a {B} molecule (and likewise for {B} molecules).

To get an idea of how a difference in the interaction energies affects the phase transition of a mixture, we’ll consider a simple model as follows. We have our usual two species {A} and {B}, and the potential interaction energy between two nearest-neighbour {A}s is the same as between two {B}s, and we’ll call this energy {u_{0}}. Between a nearest neighbour {A} and {B}, however, the energy is different, say {u_{AB}}. Further, suppose each molecule has {n} nearest neighbours, and any given molecule interacts only with its nearest neighbours. For a system composed of unmixed samples of {A} and {B}, with total number {N=N_{A}+N_{B}}, we can work out the total potential energy as follows.

Any particular molecule has {n} other molecules with which it interacts, so its total interaction energy is {nu_{0}}. There are a total of {N} molecules, but if we simply multiply {nu_{0}} by {N}, we are counting each interaction twice, so the total interaction energy is

\displaystyle  U_{unmixed}=\frac{1}{2}Nnu_{0} \ \ \ \ \ (1)

Now suppose we mix the two populations thoroughly. Let {x} be the proportion of molecules that are of species {B}, so that {1-x} is the proportion for {A}. On average, an {A} molecule now has {xn} nearest neighbours of type {B} and {\left(1-x\right)n} of type {A}, so its interaction energy is now

\displaystyle  U_{A}=xnu_{AB}+\left(1-x\right)nu_{0} \ \ \ \ \ (2)

Similarly for type {B}:

\displaystyle  U_{B}=\left(1-x\right)nu_{AB}+xnu_{0} \ \ \ \ \ (3)

To get the total energy, we multiply 2 by {N_{A}=\left(1-x\right)N} and 3 by {N_{B}=xN}, add and divide by 2 to avoid double-counting:

\displaystyle   U_{mixed} \displaystyle  = \displaystyle  \frac{1}{2}\left[xnu_{AB}+\left(1-x\right)nu_{0}\right]\left(1-x\right)N+\frac{1}{2}\left[\left(1-x\right)nu_{AB}+xnu_{0}\right]xN\ \ \ \ \ (4)
\displaystyle  \displaystyle  = \displaystyle  x\left(1-x\right)nNu_{AB}+\frac{1}{2}\left(\left(1-x\right)^{2}+x^{2}\right)nNu_{0}\ \ \ \ \ (5)
\displaystyle  \displaystyle  = \displaystyle  x\left(1-x\right)nNu_{AB}+\frac{1}{2}\left[1-2x\left(1-x\right)\right]nNu_{0} \ \ \ \ \ (6)

We can now subtract 1 from this to get the energy change resulting from mixing:

\displaystyle   \Delta U \displaystyle  = \displaystyle  U_{mixed}-U_{unmixed}\ \ \ \ \ (7)
\displaystyle  \displaystyle  = \displaystyle  x\left(1-x\right)nN\left(u_{AB}-u_{0}\right) \ \ \ \ \ (8)

Depending on the sign of {u_{AB}-u_{0}}, a plot of {\Delta U} looks like this:

The red line is for {u_{AB}-u_{0}>0} and the blue line for {u_{AB}-u_{0}<0}.

Note that the slopes of {\Delta U} are finite at the endpoints:

\displaystyle  \frac{d\Delta U}{dx}=\left(1-2x\right)nN\left(u_{AB}-u_{0}\right) \ \ \ \ \ (9)

which is finite when {x=0} or {x=1}.

As we saw earlier, the Gibbs free energy of an unmixed system is

\displaystyle  G=\left(1-x\right)G_{A}^{\circ}+xG_{B}^{\circ} \ \ \ \ \ (10)

and the change in entropy due to mixing is

\displaystyle  \Delta S_{mixing}=-Nk\left[x\ln x+\left(1-x\right)\ln\left(1-x\right)\right] \ \ \ \ \ (11)

Combining this with 8 we get the Gibbs energy for our system, taking account of both the difference in interaction energy and the entropy of mixing:

\displaystyle  G=\left(1-x\right)G_{A}^{\circ}+xG_{B}^{\circ}+x\left(1-x\right)nN\left(u_{AB}-u_{0}\right)+NkT\left[x\ln x+\left(1-x\right)\ln\left(1-x\right)\right] \ \ \ \ \ (12)

In order to plot the Gibbs energy, we need to plug in a few numbers. I’ll take {G_{A}^{\circ}=100\mbox{ J}}, {G_{B}^{\circ}=50\mbox{ J}}, {N=6.02\times10^{23}}, {n=4} and {u_{AB}-u_{0}=10^{-21}\mbox{ J}}. With these values, we can plot {G} for several temperatures:

From the top down, the red curve is for {T=50\mbox{ K},} dark green for {T=100\mbox{ K}}, blue for {T=144.928\mbox{ K}} (the reason for this peculiar value will be explained below), light green for {T=175\mbox{ K}} and yellow for {T=200\mbox{ K}}.

Note that the shape of the curve changes from concave down for low temperatures to concave up for high temperatures, with a combination of concavities for the intermediate temperature of {T=100\mbox{ K}} (dark green curve). It’s worth looking at this latter curve in more detail.

Suppose we expand the plot of this curve and draw the straight line (in grey) that is the lowest line that is tangent to the curve at two points:

The grey line represents a system consisting of two unmixed substances, one of which has {x=0.081} (the value of {x} at the left-most tangent point) and the other of which has {x=0.919} (the right-most tangent). For values of {x} between these points, the grey line represents a system which consists of two separate subsystems, one with mostly {A} molecules and the other with mostly {B} molecules. Since the Gibbs energy of the separated system is lower than that of a fully mixed system (represented by the green curve), the system will be found in this relatively unmixed state for these values of {x}. This is known as a solubility gap. By the way, note that these two unmixed subsystems are not pure systems, since they both contain molecules of both species. It’s just that one subsystem consists of about 92% {A} molecules while the other consists of about 92% {B} molecules. As we vary {x}, the proportion of each of these two subsystems varies, but the composition within each subsystem remains the same: there is no mixing between them.

Returning now to the compound plot of {G} for various temperatures above, we note that the central maximum that is present in the dark green curve gets shallower until at some critical temperature it disappears and then as the temperature is increased further, it turns into a minimum. We can find this critical temperature as follows.

The two tangent points (where the grey line intersects the green curve above) can be found by taking the derivative of 12:

\displaystyle  \frac{dG}{dx}=G_{B}^{\circ}-G_{A}^{\circ}+\left(1-2x\right)nN\left(u_{AB}-u_{0}\right)+NkT\ln\frac{x}{1-x} \ \ \ \ \ (13)

The tangent points are defined by the condition that their slope is equal to that of an unmixed system, which is {G_{B}^{\circ}-G_{A}^{\circ}}, so the tangent points must satisfy

\displaystyle  \left(1-2x\right)nN\left(u_{AB}-u_{0}\right)+NkT\ln\frac{x}{1-x}=0 \ \ \ \ \ (14)

One solution of this equation is always {x=0.5}, since this makes both terms on the LHS zero. This corresponds to the maximum point on the green curve above. To find the two tangent points we’d need to solve the transcendental equation for {x} (which I did using Maple to get the plot above), but we actually don’t need to do that to find the critical temperature. The key is to look at the second derivative at the maximum point. The second derivative is negative for a maximum and positive for a minimum, so the critical temperature (the highest temperature at which a solubility gap exists) is that value of {T} where the second derivative is zero. That is

\displaystyle   \frac{d^{2}G}{dx^{2}} \displaystyle  = \displaystyle  -2nN\left(u_{AB}-u_{0}\right)+NkT\left.\left(\frac{1}{x}+\frac{1}{1-x}\right)\right|_{x=0.5}=0\ \ \ \ \ (15)
\displaystyle  T_{c} \displaystyle  = \displaystyle  \frac{n\left(u_{AB}-u_{0}\right)}{2k} \ \ \ \ \ (16)

For the values we used above, {T_{c}=144.928\mbox{ K}}, which explains why I chose this value for the blue curve above.

Finally, we can generate a phase diagram which shows the temperature at which a solubility gap appears. From 14, if a tangent occurs at a given value of {x}, the corresponding temperature must be

\displaystyle  T=\frac{\left(1-2x\right)n\left(u_{AB}-u_{0}\right)}{k\ln\frac{1-x}{x}} \ \ \ \ \ (17)

For the parameters we used above, a plot of this is as follows:

For temperatures above the curve, full mixing occurs, while for temperatures below the curve, the system separates into two unmixed subsystems. The maximum on this curve occurs at {x=0.5} and {T=144.928\mbox{ K}}, as above.

Entropy of mixing in a small system

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.57.

As an example of the entropy changes when two pure substances are mixed, consider a system of 100 molecules, which may vary in composition from 100% of species {A} through a mixture of {A} and {B} to 100% pure {B}. The entropy of mixing is given by

\displaystyle  \Delta S_{mixing}=-Nk\left[x\ln x+\left(1-x\right)\ln\left(1-x\right)\right]

where {N=N_{A}+N_{B}} is the fixed total number of molecules (100 here) and {x=N_{A}/N}.

For a small system such as this, we can generate an array of {\Delta S_{mixing}/k} values for each value of {N_{A}} from 0 to 100. Plotting this as a bar chart, we get

Starting from {N_{A}=0} where {\Delta S/k=0} (since there is only one species at this point, there is no mixing), we see that the entropy increase per molecule as we convert successive molecules from {A} to {B} decreases. The changes in {\Delta S/k} for the first few steps are:

add molecule number change in {\frac{\Delta S}{k}}
{1} {5.60}
{2} {4.20}
{3} {3.67}
{4} {3.32}
{5} {3.06}

The rate at which the entropy increases declines as we convert more molecules from {A} to {B}. If we add a slight impurity into an initially pure mixture of 100% {B}, this generates a larger increase in entropy than if we add a bit more impurity to an already mixed system.

Gibbs free energy of a mixture of two ideal gases

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.56.

To study phase changes of mixtures of substances, rather than pure substances on their own, it’s best to start by looking at the Gibbs free energy of the mixture. If we start with a collection of molecules of two types, {A} and {B}, that are initially separated but whose total number {N_{A}+N_{B}=N} is fixed. Since the two populations are separated, the total Gibbs energy is just the sum of the energies of the two individual populations. If population {B} makes up a fraction {x} and {A} a fraction {1-x} of the total, then

\displaystyle  G=\left(1-x\right)G_{A}^{\circ}+xG_{B}^{\circ} \ \ \ \ \ (1)

What happens if we now mix the two populations, but keep the pressure and temperature constant? The Gibbs energy is defined as

\displaystyle  G\equiv U+PV-TS \ \ \ \ \ (2)

It’s possible that the energy {U} will change due to interactions between the two species being different than interactions between molecules of the same species. It’s also possible that the volume will change, if the two species either attract or repel each other differently than molecules of the same species. The biggest change, however, is likely to come from a change in entropy, because with the two populations mixed there are now a great many more ways that the molecules can be arranged.

As we saw earlier, if the two substances are ideal gases at the same pressure and temperature and we allow them to mix so that the total volume is unchanged, the change in entropy is

\displaystyle   \Delta S_{mixing} \displaystyle  = \displaystyle  -Nk\left[x\ln x+\left(1-x\right)\ln\left(1-x\right)\right]\ \ \ \ \ (3)
\displaystyle  \displaystyle  = \displaystyle  -nR\left[x\ln x+\left(1-x\right)\ln\left(1-x\right)\right] \ \ \ \ \ (4)

where {n} is the total number of moles of both species and {R} is the gas constant.

If we ignore changes in {U} and {V}, then the Gibbs energy after mixing is

\displaystyle  G=\left(1-x\right)G_{A}^{\circ}+xG_{B}^{\circ}+nRT\left[x\ln x+\left(1-x\right)\ln\left(1-x\right)\right] \ \ \ \ \ (5)

The change in entropy looks like this:

The slopes at the two ends are actually infinite, as we can see by taking the derivative of 4:

\displaystyle  \frac{1}{nR}\frac{d\Delta S_{mixing}}{dx}=-\ln x+\ln\left(1-x\right)=\ln\frac{1-x}{x} \ \ \ \ \ (6)

As {x\rightarrow0}, {\frac{1-x}{x}\rightarrow+\infty} so the log also tends to {+\infty}. As {x\rightarrow1}, {\frac{1-x}{x}\rightarrow0} and the log tends to {-\infty}.

A comparison of 1 and 5 looks like this:

Here we’ve taken {G_{A}^{o}=100}, {G_{B}^{o}=75}, {n=1} and {T=100}. The yellow line shows the Gibbs energy from 1, where the two populations are unmixed. The green curve shows 5, where the two populations are mixed. Since the energy for the mixture is lower (because of the increase in entropy), it is the favoured state.

The Gibbs energy also has infinite slopes at {x=0} and 1:

\displaystyle  \frac{dG}{dx}=-G_{A}^{o}+G_{B}^{o}+nRT\ln\frac{x}{1-x} \ \ \ \ \ (7)

Because the numerator and denominator of the log have swapped places from 6, the slope at {x=0} is now {-\infty} and the slope at {x=1} is now {+\infty}.

The infinite slopes show that even if the proportion of one of the species is very small, there is a very large increase in entropy when the two species are allowed to mix, so there is a strong tendency for this to occur.

Van der Waals fluid at the critical point

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.55.

We can use the van der Waals equation of state to investigate the behaviour of a fluid near the critical point (the point where the distinction between the liquid and gas phases disappears). In the problem statement in Schroeder’s book, he recommends that we use reduced variables, but all parts of the problem use ordinary variables, so we’ll use reduced variables and then translate the results back into ordinary variables at the end.

We start with the van der Waals equation in reduced variables

\displaystyle p=\frac{8tv^{2}-9v+3}{\left(3v-1\right)v^{2}} \ \ \ \ \ (1)


where {v\equiv V/V_{c}}, {t\equiv T/T_{c}} and {p\equiv P/P_{c}}, with:

\displaystyle V_{c} \displaystyle = \displaystyle 3Nb\ \ \ \ \ (2)
\displaystyle T_{c} \displaystyle = \displaystyle \frac{8}{27}\frac{a}{bk}\ \ \ \ \ (3)
\displaystyle P_{c} \displaystyle = \displaystyle \frac{a}{27b^{2}} \ \ \ \ \ (4)

Near the critical point all of {p}, {v} and {t} are close to 1, so we can rewrite 1 in terms of the deviation of these quantities from their values at the critical point. Using

\displaystyle A\equiv v-1=\frac{1}{V_{c}}\left(V-V_{c}\right) \ \ \ \ \ (5)

we have

\displaystyle p=\frac{8t\left(A+1\right)^{2}-9\left(A+1\right)+3}{\left(3\left(A+1\right)-1\right)\left(A+1\right)^{2}} \ \ \ \ \ (6)

Near the critical point {A\approx0} so we can expand this equation in a Taylor series about {A=0} to get (I used Maple, but you can grind through the derivatives by hand if you like):

\displaystyle p=4t-3+6\left(1-t\right)A+9\left(t-1\right)A^{2}+\left(12-\frac{27}{2}t\right)A^{3}+\mathcal{O}\left(A^{4}\right) \ \ \ \ \ (7)


Near the critical point {t\approx1} and {A\approx0} so the term in {A^{2}} will be much smaller than the other terms. [OK, the term in {A} is also the product of two near-zero terms, but for small {A}, {A^{2}\ll A} so we’ll discard the {A^{2}} term and keep the {A} term.] Our approximation for the pressure near the critical point is then

\displaystyle p\left(v\right) \displaystyle = \displaystyle 4t-3+6\left(1-t\right)\left(v-1\right)+\left(12-\frac{27}{2}t\right)\left(v-1\right)^{3}\ \ \ \ \ (8)
\displaystyle P\left(V\right) \displaystyle = \displaystyle P_{c}p\ \ \ \ \ (9)
\displaystyle \displaystyle = \displaystyle P_{c}\left[4t-3+6\left(1-t\right)\left(v-1\right)+\left(12-\frac{27}{2}t\right)\left(v-1\right)^{3}\right] \ \ \ \ \ (10)

For {t=0.95}, a plot of {p} looks like this:

Because 8 is antisymmetric about the the vertical axis, if we draw a horizontal line passing through the point {p\left(v-1=0\right)} (the yellow line above), the areas between this line and the curve on either side of the axis are equal, so we have effectively performed a Maxwell construction. This shows that the pressure at which the blue curve crosses the vertical axis is the vapour pressure for the given temperature. From 8 with {v=1}, we have

\displaystyle p_{v}=4t-3 \ \ \ \ \ (11)

Comparing this to our earlier results, for {t=0.95} our approximation gives {p_{v}=0.8} (the earlier result was {p_{v}=0.8119}), while for {t=0.8}, we get {p_{v}=0.2} compared to the earlier result of {p_{v}=0.3834}. As we’d expect, as {t} gets farther from 1, the approximation deteriorates, but for {t=0.95}, it’s not too bad.

Using this approximation, the slope of the phase boundary is approximately

\displaystyle \frac{dp_{v}}{dt} \displaystyle = \displaystyle 4\ \ \ \ \ (12)
\displaystyle \frac{dP}{dT} \displaystyle = \displaystyle \frac{P_{c}}{T_{c}}\frac{dp_{v}}{dt}=\frac{4k}{8b}=\frac{k}{2b} \ \ \ \ \ (13)

We can also use our approximation to get an estimate of the difference in volumes between the gas and liquid phases. As we can see from the above plot, there are three values of {v-1} that give the same pressure. The largest volume is {v_{g}}, the volume of the gas phase, while the smallest volume is {v_{l}}, the volume of the liquid. When {p=p_{v}}, 8 reduces to

\displaystyle 6\left(1-t\right)\left(v-1\right)+\left(12-\frac{27}{2}t\right)\left(v-1\right)^{3}=0 \ \ \ \ \ (14)

which has 3 roots: {v=1} and

\displaystyle v-1=\pm\sqrt{\frac{12\left(1-t\right)}{27t-24}} \ \ \ \ \ (15)

The operand of the square root is positive for {\frac{24}{27}\approx0.889<t<1}, so the approximation should be reasonably valid in this region. The difference {v_{g}-v_{l}} is therefore

\displaystyle v_{g}-v_{l}=2\sqrt{\frac{12\left(1-t\right)}{27t-24}} \ \ \ \ \ (16)

We can expand this in a series (I used Maple, but again you’re welcome to do the binomial expansion by hand):

\displaystyle v_{g}-v_{l}\approx4\sqrt{1-t}+\mathcal{O}\left(\left(1-t\right)^{3/2}\right) \ \ \ \ \ (17)

Converting back to ordinary variables, we have

\displaystyle V_{g}-V_{l}\approx\frac{4V_{c}}{\sqrt{T_{c}}}\sqrt{T_{c}-T} \ \ \ \ \ (18)

That is, {V_{g}-V_{l}\propto\left(T_{c}-T\right)^{1/2}}. As Schroeder notes in the problem, the experimental value of the critical exponent {\beta} in {V_{g}-V_{l}\propto\left(T_{c}-T\right)^{\beta}} is closer to {\frac{1}{3}}.

From this and 13, we can use the Clausius-Clapeyron equation to get an estimate of the latent heat of vapourization:

\displaystyle L \displaystyle = \displaystyle T\Delta V\frac{dP}{dT}\ \ \ \ \ (19)
\displaystyle \displaystyle = \displaystyle T\frac{4V_{c}}{\sqrt{T_{c}}}\sqrt{T_{c}-T}\frac{k}{2b}\ \ \ \ \ (20)
\displaystyle \displaystyle = \displaystyle \frac{2kV_{c}}{b\sqrt{T_{c}}}T\sqrt{T_{c}-T} \ \ \ \ \ (21)

The shape of the curve as a function of temperature is determined by the {T\sqrt{T_{c}-T}} factor, and looks like this:

The latent heat decreases to zero at the critical temperature, where there is no difference between the liquid and gas phases, so no transformation takes place.

Another critical exponent called {\delta} is defined by the relation between pressure and volume near the critical point:

\displaystyle \left(P-P_{c}\right)\propto\left(V-V_{c}\right)^{\delta} \ \ \ \ \ (22)

The value of {\delta} in our approximation is {\delta=3} since we’re saving terms up to the cubic term in the Taylor expansion 8. If we saved higher order terms, we’d get a higher value for {\delta}, and as Schroeder points out, experimental values are around 4 or 5.

There is one final critical exponent known as {\gamma} which is defined by the behaviour of the isothermal compressibility

\displaystyle \kappa\equiv-\frac{1}{V}\left(\frac{\partial V}{\partial P}\right)_{T} \ \ \ \ \ (23)

The derivative is the reciprocal of the slope of the {pv} curve in the plot above. As the temperature approaches the critical temperature, {t\rightarrow1} and the {pv} plot looks like this:

There is now only one phase possible and it occurs when {v-1=0}. At this point, the curve has an inflection point and its slope is zero, so {\partial V/\partial P\rightarrow\infty}. For temperatures higher than {T_{c}}, the curve looks like this:

The slope is now non-zero at {v-1=0} so {\partial V/\partial P} is again finite.

We’d like an estimate of how rapidly the derivative {\partial V/\partial P} goes to infinity, that is, we’d like to find {\gamma} in the proportion

\displaystyle \left(\frac{\partial V}{\partial P}\right)_{T}\propto\left(T_{c}-T\right)^{-\gamma} \ \ \ \ \ (24)

As the shape of the {pv} curve is qualitatively different on each side of {t=1} (for {t<1}, there are two distinct phases and three possible volumes at the vapour pressure; for {t\ge1}, there is only one phase and corresponding volume), we need to look at the two cases separately.

Using our approximation 8, we have

\displaystyle \frac{\partial p}{\partial v}=6\left(1-t\right)+\left(v-1\right)^{2}\left(36-\frac{81}{2}t\right) \ \ \ \ \ (25)


For {t>1}, there is only one value of {v} that corresponds to a physical state of the fluid, and that is {v=1}, so in this case

\displaystyle \left(\frac{\partial V}{\partial P}\right)_{T}\propto\left(T_{c}-T\right)^{-1} \ \ \ \ \ (26)

For {t<1}, however, two physical phases exist, and they correspond to the smallest and largest values of {v}. It can be seen from the top plot above that the slopes at these two points are equal (since the curve is antisymmetric). We can rewrite 25 as (this is just a Taylor expansion about {t=1}):

\displaystyle \frac{\partial p}{\partial v}=-\frac{9}{2}\left(v-1\right)^{2}+\left(6+\frac{81}{2}\left(v-1\right)^{2}\right)\left(1-t\right) \ \ \ \ \ (27)

As we approach the critical temperature from below, the volumes of the two phases both approach {v=1}, so again we see that

\displaystyle \frac{\partial p}{\partial v}\rightarrow6\left(1-t\right) \ \ \ \ \ (28)

Thus we have, for {t\rightarrow1} from below

\displaystyle \left(\frac{\partial V}{\partial P}\right)_{T}\propto\left(T_{c}-T\right)^{-1} \ \ \ \ \ (29)

The value of {\gamma} is thus {\gamma=1} on both sides of the critical temperature.

Helmholtz free energy of a van der Waals fluid

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.54.

For the van der Waals equation of state, we’ve looked at how the Gibbs free energy predicts a phase transition. We can also analyze a van der Waals fluid using the Helmholtz free energy {F}. We begin with the thermodynamic identity

\displaystyle dF=-S\;dT-P\;dV+\mu\;dN \ \ \ \ \ (1)

At constant temperature and number, this reduces to

\displaystyle dF=-P\;dV \ \ \ \ \ (2)

so we can write this as a derivative:

\displaystyle \left(\frac{\partial F}{\partial V}\right)_{T,N}=-P \ \ \ \ \ (3)

We can use the van der Waals equation of state to eliminate {P}:

\displaystyle \left(\frac{\partial F}{\partial V}\right)_{T,N}=-P=-\frac{NkT}{V-Nb}+a\frac{N^{2}}{V^{2}} \ \ \ \ \ (4)


Integrating, we get

\displaystyle F\left(V\right)=-NkT\ln\left(V-Nb\right)-\frac{aN^{2}}{V}+c\left(T\right) \ \ \ \ \ (5)


where {c\left(T\right)} is a function of temperature only, and can be ignored in what follows since we’ll be specifying a particular temperature in each example, so {c\left(T\right)} is effectively a constant in each case.

This equation is easier to work with if we rewrite it in terms of the reduced variables, defined earlier as {v\equiv V/V_{c}}, {t\equiv T/T_{c}} and {p\equiv P/P_{c}}, with:

\displaystyle V_{c} \displaystyle = \displaystyle 3Nb\ \ \ \ \ (6)
\displaystyle T_{c} \displaystyle = \displaystyle \frac{8}{27}\frac{a}{bk}\ \ \ \ \ (7)
\displaystyle P_{c} \displaystyle = \displaystyle \frac{a}{27b^{2}} \ \ \ \ \ (8)

If we divide both sides of 5 by {NkT_{c}} and simplify using the above equations, we get

\displaystyle F_{c}\equiv\frac{F}{NkT_{c}}=-t\ln\left(3v-1\right)-\frac{9}{8v}+c^{\prime}\left(t\right) \ \ \ \ \ (9)


where {c^{\prime}\left(t\right)=c\left(T\right)/NkT_{c}-t\ln\left(Nb\right)} is another function that depends only on temperature, so can be ignored in what follows.

A plot looks like this, with the red curve being {F_{c}}:

We saw in the last post that for {t=0.8}, the vapour pressure (the pressure at which the phase transition occurs) is {p_{0}=0.3834}. In terms of reduced variables, the van der Waals equation of state is

\displaystyle p=\frac{8tv^{2}-9v+3}{\left(3v-1\right)v^{2}} \ \ \ \ \ (10)


For given values of {p} and {t}, we can solve this equation to get three values for {v} (because this is a cubic equation in {v}). Using Maple for this, we find

\displaystyle v_{1} \displaystyle = \displaystyle 0.5174\ \ \ \ \ (11)
\displaystyle v_{2} \displaystyle = \displaystyle 1.2083\ \ \ \ \ (12)
\displaystyle v_{3} \displaystyle = \displaystyle 4.1725 \ \ \ \ \ (13)

The smallest volume {v_{1}} corresponds to the liquid state; the largest {v_{3}} to the gas (the other state is unstable). We can locate these volumes on the plot above and join them with a straight line (the green line). The first thing to notice is that this line appears to be the lowest line we can draw that is tangent to the curve in two points, with these two points being {v_{1}} and {v_{3}}. (You can check this is true by calculating {dF_{c}/dv} at these two points; you’ll find that the slope of the line is indeed equal to the tangents to the curve at both these points.) Thus it would appear that we can use a plot of {F_{c}} to find the two volumes at which the phase transition occurs by drawing the lowest line that is tangent to the curve in two places. The vapour pressure can then be found by plugging either of these volumes into 10 (and as a check, both volumes should give the same pressure).

[As an aside, the value of {dF_{c}/dv} at {v_{2}} is also the same as at {v_{1}} and {v_{3}}, so the tangent at that point is parallel to the green line but lies on the top side of the red curve. I don’t think this is physically significant, though.]

As for a proof that for volumes between {v_{1}} and {v_{3}}, the stable form of the fluid is a mixture of gas and liquid phases (rather than some homogenous state that is neither liquid nor gas), we can observe that if the total volume of the system is constant, {F} satisfies the condition

\displaystyle \left(\frac{\partial F_{\ell}}{\partial V}\right)_{T,N}=\left(\frac{\partial F_{g}}{\partial V}\right)_{T,N} \ \ \ \ \ (14)

That is, the tangent lines to the curve are the same for the liquid and gas phases. This condition must hold for any mixture of the two phases, so the value of {F} for such a mixture must be the straight line between {v_{1}} and {v_{3}}. Since this line is below the red curve, the Helmholtz energy is lower, so the mixture is more stable than the homogeneous state. [I’m not sure this constitutes a rigorous proof, but it seems to make sense.]

Van der Waals equation of state: Maxwell construction

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problems 5.52 – 5.53.

To get a better feel for the van der Waals equation of state, we can look at how the Gibbs free energy varies with the volume and pressure. Schroeder derives an expression for {G} in his equation 5.56, which I will just quote:

\displaystyle G=-NkT\ln\left(V-Nb\right)+\frac{\left(NkT\right)\left(Nb\right)}{V-Nb}-\frac{2aN^{2}}{V}+c\left(T\right) \ \ \ \ \ (1)


where {c\left(T\right)} is a function of temperature only, and can be ignored in what follows since we’ll be specifying a particular temperature in each example, so {c\left(T\right)} is effectively a constant in each case.

This equation is easier to work with if we rewrite it in terms of the reduced variables, defined earlier as {v\equiv V/V_{c}}, {t\equiv T/T_{c}} and {p\equiv P/P_{c}}, with:

\displaystyle V_{c} \displaystyle = \displaystyle 3Nb\ \ \ \ \ (2)
\displaystyle T_{c} \displaystyle = \displaystyle \frac{8}{27}\frac{a}{bk}\ \ \ \ \ (3)
\displaystyle P_{c} \displaystyle = \displaystyle \frac{a}{27b^{2}} \ \ \ \ \ (4)

If we divide both sides of 1 by {NkT_{c}} and simplify using the above equations, we get

\displaystyle G_{c}\equiv\frac{G}{NkT_{c}}=-t\ln\left(3v-1\right)+\frac{t}{3v-1}-\frac{9}{4v}+c^{\prime}\left(t\right) \ \ \ \ \ (5)


where {c^{\prime}\left(t\right)=c\left(T\right)/NkT_{c}-t\ln\left(Nb\right)} is another function that depends only on temperature, so can be ignored in what follows.

To plot {G_{c}} as a function of reduced pressure, we can use the equation derived in the last post:

\displaystyle p=\frac{8tv^{2}-9v+3}{\left(3v-1\right)v^{2}} \ \ \ \ \ (6)


We can specify a temperature {t} and then vary {v}, using it to calculate {p} and the corresponding {G_{c}} from the last two equations, and thus generate a plot of {G_{c}} versus {p}. For {t=0.95} the isotherm looks like this:

If we plot the Gibbs energy, we get:

As we reduce the volume, we move from right to left on the {pv} plot, and the corresponding path on the Gibbs plot is the line coming up from the bottom. At point {C}, we can continue to decrease {v}, which causes {G} to carry on along the same path above point {C}, until it reaches a sharp corner, at which point it decreases until it reaches point {B}, which is directly above (that is, it has the same pressure) as point {C}. However, note that {B} does not have the minimum Gibbs energy for that pressure, since point {C} is below it. Thus {B} is unstable and does not correspond to a physical state of the substance. Back on the {pv} plot, as we decrease {v} below point {B}, the Gibbs path decreases until it reaches the lower sharp corner. This corner corresponds to the minimum between {A} and {B} on the {pv} plot. Again, note that this isn’t the minimum Gibbs energy for that pressure, so this is also unstable. In fact the whole section of the {pv} plot where the pressure decreases with decreasing volume (between the minimum between {A} and {B} and the maximum between {B} and {C}) lies on a Gibbs path that is not the minimum energy, so this whole range of volumes is unstable and is not found in the real world.

As we continue to decrease {v} from the minimum point down to point {A}, the Gibbs energy returns from the lower sharp corner until it reaches point {A}, where the Gibbs energy is equal to the value it had at point {C}. Thus points {A} and {C} have the same Gibbs energy at the same pressure, but point {A} has a much lower volume than point {C}, so it corresponds to a condensed state (a liquid) while {C} corresponds to a gas.

In fact, this phase transformation occurs for only one pressure, as we can see from the following argument. At constant temperature, the thermodynamic identity states that

\displaystyle \left(\frac{\partial G}{\partial P}\right)_{T}=V \ \ \ \ \ (7)

As we go around the triangular loop in the Gibbs plot above, the overall change in {G} is zero, so if we integrate it around the loop, we must get zero:

\displaystyle 0=\oint dG=\oint\left(\frac{\partial G}{\partial P}\right)_{T}dP=\oint V\;dP \ \ \ \ \ (8)

If we divide this equation by {V_{c}P_{c}} it changes nothing, since both these quantities depend only on {T}, which is held constant. Thus

\displaystyle \oint v\;dp=0 \ \ \ \ \ (9)

To do this integral directly requires inverting 6 to get {v} as a function of {p}. However, this requires solving a cubic equation in {v}, which gives 3 different roots. Just for completeness, I did actually do this using Maple, and if we plot the result we get, for {t=0.95}:

The red, blue and green curves represent contributions from each of the three roots. As we go round the triangular loop above, we go from point {A} to {B} to {C}. Requiring this integral to be zero is equivalent (as Schroeder shows) to the area bounded by the vertical line and curve between {B} and {C} being equal to the area bounded by the vertical line and curve between {A} and {B}. This is known as the Maxwell construction (yes, it’s the same Maxwell that devised Maxwell’s equations in electrodynamics).

The problem is then to find the value {p_{0}} of {p} that makes these two areas equal. Doing this requires moving the yellow line in the plot left or right until the areas are equal. Solving the problem analytically is a bit trickier, of course, since we’re dealing with a cubic equation.

It’s a lot easier to solve if we revert to the traditional {pv} plot above. The area bounded by the vertical line and curve between {B} and {C} is then given by

\displaystyle A_{BC}=\int_{v_{B}}^{v_{C}}p\;dv-p_{0}\left(v_{C}-v_{B}\right) \ \ \ \ \ (10)


That is, it’s the area under the {p} curve between {v_{B}} and {v_{C}} minus the area of the rectangle under the horizontal line between the same two volumes. By the same reasoning, the area in the {vp} plot bounded by the vertical line and curve between {A} and {B} is

\displaystyle A_{AB}=p_{0}\left(v_{B}-v_{A}\right)-\int_{v_{A}}^{v_{B}}p\;dv \ \ \ \ \ (11)


Surprisingly, the integral turns out to be fairly simple, using Maple:

\displaystyle \int p\;dv \displaystyle = \displaystyle \int\frac{8tv^{2}-9v+3}{\left(3v-1\right)v^{2}}dv\ \ \ \ \ (12)
\displaystyle \displaystyle = \displaystyle \frac{3}{v}+\frac{8}{3}t\ln\left(3v-1\right) \ \ \ \ \ (13)

The problem now becomes:

  1. Specify the temperature.
  2. Choose a starting value for the pressure {p_{0}}. This can be estimated from the {pv} plot by judging by eye what value of {p_{0}} makes the two areas roughly equal.
  3. Solve 6 to find the 3 volumes {v_{A}}, {v_{B}} and {v_{C}} for this value of {p_{0}}. This requires solving a cubic equation in {v}.
  4. Calculate the two areas using 10 and 11.
  5. Adjust {p_{0}} and repeat from step 3 until we get convergence.

I did this using a Maple procedure, with the best result (for {t=0.95}):

\displaystyle p_{0}=0.8119 \ \ \ \ \ (14)

As a check, we can compare this with the Gibbs energy plot above and note that the point {A,C} does indeed occur at this pressure.

We can repeat the procedure for {t=0.8} and get the result:

\displaystyle p_{0}=0.3834 \ \ \ \ \ (15)

The Gibbs plot for this case is:

Again, we can see that the crossover between the two phases occurs at around {p_{0}=0.3834} so it looks like all is well.

Van der Waals equation of state

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problems 5.48 – 5.51.

The van derWaals model of a substance is able to predict (qualitatively) the existence of the liquid-gas phase transition and the critical point (where there is no clear distinction between the liquid and gas phases). Schroeder gives a good explanation of the van der Waals model in section 5.3 of his book so I won’t go over the whole thing here, but it’s worth stating the main points.

The model is a refinement of the ideal gas equation of state, {PV=NkT}, which looks like this:

\displaystyle \left(P+\frac{aN^{2}}{V^{2}}\right)\left(V-Nb\right)=NkT \ \ \ \ \ (1)


where {a} and {b} are constants whose values depend on the particular substance we’re describing. The correction to the volume term is due to the fact that in a real substance, it is not possible to reduce the volume to zero since the molecules have a size below which they cannot be compressed further. Thus the minimum volume of a system of {N} molecules is {Nb}, where {b} depends on the nature of the substance.

The correction to the pressure is obtained from the following argument. Due to electric interactions, all molecules exhibit a long term attraction to each other. The potential energy due to the interaction between one molecule and its neighbours depends on the density of molecules {N/V}; a higher density means that a given molecule has more neighbours, and thus a higher potential energy due to electrical interactions. The total potential energy of the whole system due to these interactions is then proportional to the number {N} of molecules times the potential energy per molecule, which is {-aN^{2}/V}, where {a} is the constant of proportionality, which depends on the substance. [The potential energy is negative since we’re dealing with an attractive interaction. Thus two molecules separated by some finite distance require positive work done on them to pull them apart to an infinite distance, at which point the potential energy is zero. That is, the work is required to pull them out of a potential well, so their potential energy is negative.]

The effect of this on the pressure can be obtained from the thermodynamic identity {dU=T\;dS-P\;dV}. If we consider a system in which all the molecules are frozen in place then the entropy will not change. Increasing the volume of this system by {dV} will reduce the potential energy by an amount {dU} (since the neighbours will be further apart), so that

\displaystyle dU \displaystyle = \displaystyle -P\;dV\ \ \ \ \ (2)
\displaystyle P \displaystyle = \displaystyle -\left(\frac{\partial U}{\partial V}\right)_{S}\ \ \ \ \ (3)
\displaystyle \displaystyle = \displaystyle -\frac{\partial}{\partial V}\left[-\frac{aN^{2}}{V}\right]\ \ \ \ \ (4)
\displaystyle \displaystyle = \displaystyle -a\frac{N^{2}}{V^{2}} \ \ \ \ \ (5)

Thus the attractive interaction produces a negative correction to the overall pressure, as we would expect (if the molecules are mutually attracted to each other, this pulls them together and reduces the pressure they exert on their container. The net pressure therefore becomes

\displaystyle P=\frac{NkT}{V-Nb}-a\frac{N^{2}}{V^{2}} \ \ \ \ \ (6)


which gives the van der Walls equation 1.

To get a feel for the model, we can plot some isotherms on a {PV} diagram:

[These are actually plots of {p\equiv P/P_{c}} versus {v\equiv V/V_{c}} where the critical pressure {P_{c}} and volume {V_{c}} will be defined in a minute (see 17 below). At this point, we’re interested only in the shapes of the curves.] The blue curve is for a temperature {0.8T_{c}} (where the critical temperature {T_{c}} will also be defined in a minute; have patience!), the red curve is for {T_{c}} and the green curve is for {1.2T_{c}}. For {T<T_{c}}, the curve shows a minimum pressure as we reduce the volume, which seems to indicate that as we compress the substance, its pressure actually decreases. The red curve shows that at the critical temperature {T=T_{c}}, the minimum in the pressure curve becomes an inflection point, where both {\frac{dP}{dV}} and {\frac{d^{2}P}{dV^{2}}} are zero (this is a point like that in the graph of {y=x^{3}} at {x=0}). As we increase {T} above {T_{c}}, the minimum in the pressure disappears and we are left with a curve that gets closer to that for an ideal gas, where {P=NkT/V} (which is a hyperbola).

So does the van der Waals model actually predict a decrease in pressure as we compress the substance? The answer is ‘no’ as we can see if we examine the Gibbs free energy. Schroeder gives a complete description of how this works, which I won’t repeat here. However, it’s worth looking at one more graph to see what actually happens. We’ll examine the blue curve (for {T=0.8T_{c}}) above in a bit more detail:

Suppose we follow the yellow horizontal line from right to left. Look first at the right-most intersection between the line and the blue curve, at around {v=4.8}. At this point, the substance is a gas. According to the van der Waals equation, this same pressure is also obtained when the volume is around {v=1.1}, and again when the volume has decreased to about {v=0.5}. However, if we look at the Gibbs energies for these three points, we find that the Gibbs energy is the same for {v=0.5} as it is for {v=4.8}, but that the Gibbs energy for the point {v=1.1} is higher than the other two points. In other words, when the pressure is that given by the yellow line, the system will be found either at {v=4.8} or at {v=0.5} but not at {v=1.1}, since the latter state is unstable. That is, there are two stable states at that pressure and temperature: a high-volume state (a gas) and a low-volume state (a liquid).

The critical point occurs at the lowest temperature for which there is only a single intersection between a horizontal line and the {PV} curve, which occurs when the curve has an inflection point, as is shown by the red curve in the first graph above. We can find this by taking the first two derivatives of {P} in 6 and setting them to zero:

\displaystyle \frac{\partial P}{\partial V} \displaystyle = \displaystyle -\frac{NkT}{\left(V-Nb\right)^{2}}+\frac{2aN^{2}}{V^{3}}=0\ \ \ \ \ (7)
\displaystyle \frac{\partial^{2}P}{\partial V^{2}} \displaystyle = \displaystyle \frac{2NkT}{\left(V-Nb\right)^{3}}-\frac{6aN^{2}}{V^{4}}=0 \ \ \ \ \ (8)

We can solve these 2 equations for {V} and {T}:

\displaystyle T \displaystyle = \displaystyle 2aN\frac{\left(V-Nb\right)^{2}}{kV^{3}}\ \ \ \ \ (9)
\displaystyle 0 \displaystyle = \displaystyle \frac{4aN^{2}}{\left(V-Nb\right)V^{3}}-\frac{6aN^{2}}{V^{4}}\ \ \ \ \ (10)
\displaystyle V_{c} \displaystyle = \displaystyle 3Nb\ \ \ \ \ (11)
\displaystyle T_{c} \displaystyle = \displaystyle \frac{8}{27}\frac{a}{bk} \ \ \ \ \ (12)

Substituting these back into 6 gives {P_{c}}:

\displaystyle P_{c}=\frac{a}{27b^{2}} \ \ \ \ \ (13)

We can rewrite the van der Waals equation 6 in terms of dimensionless reduced variables {v\equiv V/V_{c}}, {t\equiv T/T_{c}} and {p\equiv P/P_{c}}. We get

\displaystyle V \displaystyle = \displaystyle 3Nbv\ \ \ \ \ (14)
\displaystyle T \displaystyle = \displaystyle \frac{8}{27}\frac{a}{bk}t\ \ \ \ \ (15)
\displaystyle P \displaystyle = \displaystyle \frac{a}{27b^{2}}p \ \ \ \ \ (16)

Making these substitutions and simplifying, we get

\displaystyle p=\frac{8tv^{2}-9v+3}{\left(3v-1\right)v^{2}} \ \ \ \ \ (17)


which is conveniently independent of {a} and {b}. It is this equation that was used to generate the plots above.

Using the critical parameters above, we can calculate the compression factor of the fluid at its critical point. The compression factor is defined as the ratio {PV/NkT} and since this is 1 for an ideal gas, its actual value can be used as a measure of how much the fluid differs from an ideal gas. Plugging in the critical values above, we find

\displaystyle \frac{P_{c}V_{c}}{NkT_{c}}=\frac{3}{8}=0.375 \ \ \ \ \ (18)

As Schroeder points out, the compression factors of real fluids at their critical points are usually lower than this. Water, for example, has a compression factor of 0.227.

Finally, we can work out the critical parameters for a few substances, assuming that the van der Waals model is accurate for them (which it probably isn’t). Rather than use Schroeder’s values for {a} and {b}, I’ll use these from Wikipedia. The values for {a} are given in strange units ({\mbox{L}^{2}\cdot\mbox{bar}\cdot\mbox{mol}^{-2}}) where the L stands for ‘litre’. To convert these to {\mbox{J}\cdot\mbox{m}^{3}} for a single molecule, we divide by {10N_{A}^{2}}. The values for {b} are given in {\mbox{L\ensuremath{\cdot}\ensuremath{\ensuremath{\mbox{mol}^{-1}}}}}; to convert this to {\mbox{m}^{3}}, divide by {1000N_{A}}.

{a\;\left(\mbox{J}\cdot\mbox{m}^{3}\right)} {b\;\left(\mbox{m}^{3}\right)} {T_{c}\;\left(\mbox{K}\right)} {P_{c}\;\left(\mbox{bar}\right)} {\frac{V_{c}}{N}\;\left(\mbox{m}^{3}\right)}
{\mbox{N}_{2}} {3.78\times10^{-49}} {6.43\times10^{-29}} {126.2} {33.9} {1.93\times10^{-28}}
{\mbox{H}_{2}\mbox{O}} {1.53\times10^{-48}} {5.06\times10^{-29}} {649.2} {221} {1.58\times10^{-28}}
{\mbox{He}} {9.55\times10^{-51}} {3.95\times10^{-29}} {5.19} {2.27} {1.19\times10^{-28}}

Surface tension and cloud formation

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.46.

At the boundary between two phases, there are some molecules that are in a sort of “transition phase”, that is, in neither one phase nor the other. To describe the phase boundary more fully, then, we need to add a term to the Gibbs free energy for these boundary molecules. It turns out that the thickness of the boundary layer is more or less constant, independent of the area of the boundary. The extra Gibbs free energy is therefore proportional to the area of the boundary:

\displaystyle G_{boundary}=\sigma A \ \ \ \ \ (1)

where the proportionality coefficient {\sigma} is called the surface tension. It represents the minimum work (per unit area) required to reshape a portion of a substance into a form with a larger surface area. For liquid water at {20^{\circ}\mbox{ C}}, {\sigma=0.073\mbox{ J m}^{-2}}.

As an example, suppose we have a spherical droplet of water of radius {r} containing {N_{l}} molecules surrounded by {N-N_{l}} molecules of water vapour. If we write the Gibbs free energy in terms of the chemical potential, the Gibbs energy of this system, initially neglecting surface tension, is

\displaystyle G=N_{l}\mu_{l}+\left(N-N_{l}\right)\mu_{g} \ \ \ \ \ (2)

where the subscripts {l} and {g} stand for ‘liquid’ and ‘gas’.

Given the volume per molecule {v_{l}} in the liquid, we can write

\displaystyle v_{l} \displaystyle = \displaystyle \frac{4\pi r^{3}}{3N_{l}}\ \ \ \ \ (3)
\displaystyle G \displaystyle = \displaystyle -\frac{4\pi r^{3}}{3v_{l}}\left(\mu_{g}-\mu_{l}\right)+N\mu_{g} \ \ \ \ \ (4)

Adding in the energy due to surface tension gives

\displaystyle G=-\frac{4\pi r^{3}}{3v_{l}}\left(\mu_{g}-\mu_{l}\right)+N\mu_{g}+4\pi r^{2}\sigma \ \ \ \ \ (5)

Equilibrium values of the sphere’s radius {r} are obtained when {\partial G/\partial r=0}, so that there is no change in Gibbs energy when {r} is changed by an amount {dr}. This gives

\displaystyle \frac{\partial G}{\partial r} \displaystyle = \displaystyle -\frac{4\pi r^{2}}{v_{l}}\left(\mu_{g}-\mu_{l}\right)+8\pi r\sigma=0\ \ \ \ \ (6)
\displaystyle r_{c} \displaystyle = \displaystyle \frac{2\sigma v_{l}}{\mu_{g}-\mu_{l}} \ \ \ \ \ (7)

In order for {r_{c}>0}, we need {\mu_{g}>\mu_{l}}. The Gibbs energy for this case looks qualitatively like this [units are abitrary; it’s just the shape of the graph we’re interested in]:

Because {G} is a maximum at the equilibrium point, this is an unstable equilibrium, meaning that if {r} is nudged slightly lower, the tendency is for the droplet to continue to decrease in size due to evaporation, while if {r} is nudged slightly higher, the tendency is for the droplet to continue to increase in size.

If {\mu_{g}<\mu_{l}} we get:

There is a minimum at {r=0}, but this isn’t a physical solution.

To convert 7 into a form involving the relative humidity, we can start with the thermodynamic identity for the Gibbs energy:

\displaystyle dG=-S\;dT+V\;dP+\mu\;dN \ \ \ \ \ (8)

At constant temperature and number we have

\displaystyle dG=N\;d\mu=V\;dP \ \ \ \ \ (9)

So for the liquid and gas phases, a slight change in the vapour pressure {P} results in

\displaystyle d\mu_{g} \displaystyle = \displaystyle \frac{V_{g}}{N_{g}}dP\ \ \ \ \ (10)
\displaystyle \displaystyle = \displaystyle v_{g}dP\ \ \ \ \ (11)
\displaystyle d\mu_{l} \displaystyle = \displaystyle v_{l}dP \ \ \ \ \ (12)

Therefore the change in the difference between the chemical potentials is given by

\displaystyle d\left(\mu_{g}-\mu_{l}\right) \displaystyle = \displaystyle \left(v_{g}-v_{l}\right)dP\ \ \ \ \ (13)
\displaystyle \displaystyle \approx \displaystyle v_{g}dP\ \ \ \ \ (14)
\displaystyle \displaystyle = \displaystyle \frac{kT}{P}\;dP \ \ \ \ \ (15)

where in the second line we used the fact that {v_{g}\gg v_{l}}, since the volume per gas molecule is much greater than the volume per liquid molecule. The last line assumes that the vapour is an ideal gas, so {v_{g}=V_{g}/N_{g}=kT/P}.

If the liquid and gas are in diffusive equilibrium, {\mu_{g}=\mu_{l}} and the vapour is saturated so that {P=P_{v}}, so we can integrate to get

\displaystyle \mu_{g}-\mu_{l} \displaystyle = \displaystyle \int_{0}^{\mu_{g}-\mu_{l}}d\left(\mu_{g}-\mu_{l}\right)^{\prime}\ \ \ \ \ (16)
\displaystyle \displaystyle = \displaystyle kT\int_{P_{v}}^{P}\frac{dP^{\prime}}{P^{\prime}}\ \ \ \ \ (17)
\displaystyle \displaystyle = \displaystyle kT\ln\frac{P}{P_{v}}\ \ \ \ \ (18)
\displaystyle \displaystyle = \displaystyle kT\ln h \ \ \ \ \ (19)

where {h} is the relative humidity. Plugging this into 7 we get

\displaystyle r_{c}=\frac{2\sigma v_{l}}{kT\ln h} \ \ \ \ \ (20)

The snag with this equation is that the relative humidity is usually between 0 and 1, so {\ln h<0} and {r_{c}<0}, which is clearly nonsense. In order for this equation to give a physical radius, we need {h>1}, so the air is ‘supersaturated’. In this case, a plot of {r_{c}} versus {h} looks like this:

In order for droplets to form, the initial droplet size has to be very large (relatively speaking; a radius of a few nanometres still contains a very large number of water molecules) so it would seem unlikely that such droplets would form spontaneously.

Wet adiabatic lapse rate

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.45.

Earlier, we worked out the dry adiabatic lapse rate, which is the temperature gradient in dry air at which convection begins, that is, it is the point at which air near the surface achieves a density low enough that it starts to rise adiabatically. However, in the real atmosphere there will always be some water vapour, and as the air rises and cools, some of this water vapour will condense out, forming clouds. When water vapour condenses to liquid, it releases an amount of heat given by the latent heat of vapourization {L}, causing a decrease in the rate at which cooling occurs with height. We can get an estimate of this revised rate of cooling, called the wet adiabatic lapse rate (sometimes the moist adiabatic lapse rate) by combining some of our earlier results.

First, we can use the first law of thermodynamics (conservation of energy) to include the effect of condensation. Suppose we have a mass of air containing {n} moles of air molecules and {n_{w}} moles of water vapour (we’ll assume that {n_{w}\ll n} in what follows). If an amount {dn_{w}} of water vapour condenses, it releases an amount of heat {L\;dn_{w}}. As the air rises, it also expands by a volume {dV}, doing work {P\;dV}, so the net energy change of the air mass is

\displaystyle dU=-P\;dV-L\;dn_{w} \ \ \ \ \ (1)


The first negative sign is because the gas does work (and hence loses energy) as it expands, while the second negative sign is because a decrease in {n_{w}} (that is, {dn_{w}<0}) releases heat, thus increasing {dU}.

We want to convert this to a formula for the change in temperature {dT} in terms of {dP} and {dn_{w}}. From the equipartition theorem, the energy of the air mass is

\displaystyle U=\frac{f}{2}nRT \ \ \ \ \ (2)

where {f} is the number of degrees of freedom of an air molecule, which is usually taken to be {f=5}. We therefore have

\displaystyle dU=\frac{5}{2}nR\;dT \ \ \ \ \ (3)

We now need to convert the {P\;dV} term. My first thought was that, since the gas is supposed to be expanding adiabatically, we can use the relation

\displaystyle PV^{\gamma}=K \ \ \ \ \ (4)


where {\gamma=\left(f+2\right)/f=7/5} and {K} is a constant. Taking differentials gives

\displaystyle dV=-\frac{V}{\gamma P}dP=-\frac{5V}{7P}dP \ \ \ \ \ (5)

and plugging all this back into 1 gives

\displaystyle \frac{5}{2}nR\;dT \displaystyle = \displaystyle \frac{5V}{7}dP-L\;dn_{w}\ \ \ \ \ (6)
\displaystyle dT \displaystyle = \displaystyle \frac{2V}{7nR}dP-\frac{2}{5}\frac{L}{nR}dn_{w} \ \ \ \ \ (7)

Applying the ideal gas law {PV=nRT} to the first term, we get

\displaystyle dT=\frac{2}{7}\frac{T}{P}dP-\frac{2}{5}\frac{L}{nR}dn_{w} \ \ \ \ \ (8)

This does not agree with Schroeder’s result, which has {\frac{2}{7}} instead of {\frac{2}{5}} in the last term.

We can get Schroeder’s answer if instead of using the adiabatic relation 4, we use the ideal gas law and take differentials:

\displaystyle P\;dV+V\;dP \displaystyle = \displaystyle nR\;dT\ \ \ \ \ (9)
\displaystyle P\;dV \displaystyle = \displaystyle nR\;dT-V\;dP \ \ \ \ \ (10)

We now get

\displaystyle \frac{5}{2}nR\;dT \displaystyle = \displaystyle -nR\;dT+V\;dP-L\;dn_{w}\ \ \ \ \ (11)
\displaystyle \frac{7}{2}nR\;dT \displaystyle = \displaystyle V\;dP-L\;dn_{w}\ \ \ \ \ (12)
\displaystyle \displaystyle = \displaystyle \frac{nRT}{P}dP-L\;dn_{w}\ \ \ \ \ (13)
\displaystyle dT \displaystyle = \displaystyle \frac{2}{7}\frac{T}{P}dP-\frac{2}{7}\frac{L}{nR}dn_{w} \ \ \ \ \ (14)

I’m not entirely sure why the adiabatic relation gives the wrong answer, although it could have something to do with the fact that, strictly speaking, the process isn’t adiabatic, since the condensing vapour injects some heat into the air mass.

In any case, we can proceed with the derivation of the wet adiabatic lapse rate. The next step is to assume that the air remains saturated (100% humidity) at all stages of the expansion. In this case, the ratio {n_{w}/n} satisfies

\displaystyle \frac{n_{w}}{n}=\frac{P_{v}\left(T\right)}{P} \ \ \ \ \ (15)

where {P_{v}} is the vapour pressure at which condensation occurs at temperature {T}. Thus, {n_{w}=n_{w}\left(P,T\right)} is a function of total air pressure {P} and temperature {T} at a given height {z}. We can now take the total derivative to get

\displaystyle \frac{dn_{w}}{dz}=\frac{\partial n_{w}}{\partial P}\frac{dP}{dz}+\frac{\partial n_{w}}{\partial T}\frac{dT}{dz} \ \ \ \ \ (16)

The partial derivatives are

\displaystyle \frac{\partial n_{w}}{\partial P} \displaystyle = \displaystyle -\frac{nP_{v}}{P^{2}}\ \ \ \ \ (17)
\displaystyle \frac{\partial n_{w}}{\partial T} \displaystyle = \displaystyle \frac{n}{P}\frac{dP_{v}}{dT} \ \ \ \ \ (18)

So we get

\displaystyle \frac{dn_{w}}{dz}=-\frac{nP_{v}}{P^{2}}\frac{dP}{dz}+\frac{n}{P}\frac{dP_{v}}{dT}\frac{dT}{dz} \ \ \ \ \ (19)


From the Clausius-Clapeyron equation, we have

\displaystyle \frac{dP_{v}}{dT}=\frac{L}{T\Delta V} \ \ \ \ \ (20)


where {L} is the latent heat of vapourization per mole, and {\Delta V} is the difference in volume between the gas and liquid phases. Since the volume of the liquid is negligible relative to the volume of the gas, we can approximate {\Delta V} using the ideal gas law:

\displaystyle \Delta V=\frac{n_{w}RT}{P_{v}}=\frac{RT}{P_{v}} \ \ \ \ \ (21)

[Note that {n_{w}=1} here since {L} in 20 is the latent heat per mole.] This gives

\displaystyle \frac{dP_{v}}{dT}=\frac{LP_{v}}{RT^{2}} \ \ \ \ \ (22)

Therefore, 19 becomes

\displaystyle \frac{dn_{w}}{dz}=-\frac{nP_{v}}{P^{2}}\frac{dP}{dz}+\frac{nLP_{v}}{PRT^{2}}\frac{dT}{dz} \ \ \ \ \ (23)


Returning to 14, we can divide both sides by {dz} and substitute in 23 to get

\displaystyle \frac{dT}{dz} \displaystyle = \displaystyle \frac{2}{7}\frac{T}{P}\frac{dP}{dz}-\frac{2}{7}\frac{L}{nR}\left[-\frac{nP_{v}}{P^{2}}\frac{dP}{dz}+\frac{nLP_{v}}{PRT^{2}}\frac{dT}{dz}\right]\ \ \ \ \ (24)
\displaystyle \left[1+\frac{2}{7}\frac{P_{v}}{P}\left(\frac{L}{RT}\right)^{2}\right]\frac{dT}{dz} \displaystyle = \displaystyle \frac{2}{7}\left[\frac{T}{P}+\frac{LP_{v}}{RP^{2}}\right]\frac{dP}{dz} \ \ \ \ \ (25)

We can now use the barometric equation on the RHS:

\displaystyle \frac{dP}{dz} \displaystyle = \displaystyle -\frac{mg}{kT}P\ \ \ \ \ (26)
\displaystyle \displaystyle = \displaystyle -\frac{Mg}{RT}P \ \ \ \ \ (27)

where {m} is the mass of one air molecule and {M} is the mass of a mole of air molecules. This gives

\displaystyle \left[1+\frac{2}{7}\frac{P_{v}}{P}\left(\frac{L}{RT}\right)^{2}\right]\frac{dT}{dz} \displaystyle = \displaystyle -\frac{2}{7}\frac{Mg}{R}\left[1+\frac{P_{v}}{P}\frac{L}{RT}\right]\ \ \ \ \ (28)
\displaystyle \frac{dT}{dz} \displaystyle = \displaystyle -\frac{2}{7}\frac{Mg}{R}\left[1+\frac{P_{v}}{P}\frac{L}{RT}\right]\left[1+\frac{2}{7}\frac{P_{v}}{P}\left(\frac{L}{RT}\right)^{2}\right]^{-1} \ \ \ \ \ (29)

If we plug in some numbers we can see the effect of condensation on the lapse rate. We’ll use the earlier expression for {P_{v}} and other values from Schroeder’s book:

\displaystyle P_{v} \displaystyle = \displaystyle Ke^{-L/RT}\ \ \ \ \ (30)
\displaystyle K \displaystyle = \displaystyle 1.63\times10^{11}\mbox{ Pa}\ \ \ \ \ (31)
\displaystyle m \displaystyle = \displaystyle 4.81\times10^{-26}\mbox{ kg}\ \ \ \ \ (32)
\displaystyle M \displaystyle = \displaystyle N_{A}m=0.029\mbox{ kg}\ \ \ \ \ (33)
\displaystyle L \displaystyle = \displaystyle 4.399\times10^{4}\ \ \ \ \ (34)
\displaystyle R \displaystyle = \displaystyle 8.314\ \ \ \ \ (35)
\displaystyle g \displaystyle = \displaystyle 9.8 \ \ \ \ \ (36)

[The last 3 quantities are in SI units.]

At a temperature of {T=25^{\circ}\mbox{ C}=298\mbox{ K}} and pressure of 1 bar ({10^{5}\mbox{ Pa})}, we have

\displaystyle \frac{dT}{dz}=-0.00395\mbox{ K m}^{-1} \ \ \ \ \ (37)

or about 4 K per kilometre. The dry lapse rate is around 9.75 K per kilometre, so condensation has a considerable effect.

At the lower temperature of {0^{\circ}\mbox{ C}=273\mbox{ K}} and pressure of 1 bar we have

\displaystyle \frac{dT}{dz}=6.5\mbox{ K km}^{-1} \ \ \ \ \ (38)

Presumably this is because at lower temperatures the air can hold less water vapour, so the air is dryer, which means less vapour is available for condensation.

Cloud formation

Reference: Daniel V. Schroeder, An Introduction to Thermal Physics, (Addison-Wesley, 2000) – Problem 5.44.

As warm, moist air at the Earth’s surface rises, it expands and cools until the vapour pressure reaches the dew point at a certain height. When the air reaches that height, the vapour starts to condense and form clouds. We can get a rough estimate of the height at which clouds should form by combining a few of our earlier results.

First, we’ll assume that the rising air expands adiabatically and that the temperature gradient in the atmosphere is at the adiabatic lapse rate. This rate was calculated earlier to be

\displaystyle  \frac{dT}{dz}=-0.00976\mbox{ deg m}^{-1} \ \ \ \ \ (1)

This is the rate at which the temperature decreases as we go higher in the atmosphere. Since this rate is a constant, the temperature as a function of height {z} above the Earth’s surface is

\displaystyle  T\left(z\right)=T_{0}-0.00976z \ \ \ \ \ (2)

where {T_{0}} is the temperature at the surface.

Secondly, if we assume that the atmosphere is just at the critical point where convection begins, then its pressure still obeys the barometric equation:

\displaystyle  P\left(z\right)=P_{0}e^{-mgz/kT} \ \ \ \ \ (3)

where {P_{0}} is the pressure at the surface. The values we used earlier are {m=4.81\times10^{-26}\mbox{ kg}} for the mass of a typical air molecule, {g=9.8\mbox{ m s}^{-2}} and {k=1.38\times10^{-23}} in SI units.

Finally, we can use the vapour pressure equation which gives the vapour pressure at which an ideal gas condenses into liquid:

\displaystyle  P_{v}\left(T\right)=Ke^{-L/RT} \ \ \ \ \ (4)

where {K=1.63\times10^{11}\mbox{ Pa}} is a constant, {L=4.399\times10^{4}} is the latent heat of vapourization and the gas constant is {R=8.314}, both in SI units. Using 2 we can write this as a function of {z}:

\displaystyle  P_{v}\left(z\right)=Ke^{-L/R\left(T_{0}-0.00976z\right)} \ \ \ \ \ (5)

To apply all this to cloud formation, we’ll throw in some specific numbers. Let’s assume that the surface temperature is {T_{0}=25^{\circ}\mbox{ C}=298\mbox{ K}} and the surface humidity is 50%. This gives a vapour pressure at {z=0} of

\displaystyle  P_{v}\left(0\right)=0.5Ke^{-L/R\times298}=1585\mbox{ Pa} \ \ \ \ \ (6)

If we use this as {P_{0}} in 3 (that is, we’re assuming that the vapour pressure follows the same exponential relation as the overall pressure), and take the temperature as given by 2, then we have an equation for the partial pressure of water vapour as a function of height {z}.

\displaystyle  P\left(z\right)=P_{0}e^{-mgz/k\left(T_{0}-0.00976z\right)} \ \ \ \ \ (7)

The height at which clouds will form is the value of {z} such that {P\left(z\right)=P_{v}\left(z\right)}, since then the actual vapour pressure is equal to the vapour pressure where condensation occurs. That is, we want {z} such that

\displaystyle  P_{0}e^{-mgz/k\left(T_{0}-0.00976z\right)}=Ke^{-L/R\left(T_{0}-0.00976z\right)} \ \ \ \ \ (8)

Although this equation may look scary, it’s actually quite easy to solve for {z} by taking the logs:

\displaystyle  \ln\frac{P_{0}}{K}-\frac{mgz}{k\left(T_{0}-0.00976z\right)}=-\frac{L}{R\left(T_{0}-0.00976z\right)} \ \ \ \ \ (9)

Rearranging, we get

\displaystyle   \left(T_{0}-0.00976z\right)\ln\frac{P_{0}}{K} \displaystyle  = \displaystyle  \frac{mgz}{k}-\frac{L}{R}\ \ \ \ \ (10)
\displaystyle  z \displaystyle  = \displaystyle  \left[T_{0}\ln\frac{P_{0}}{K}+\frac{L}{R}\right]\left[0.00976\ln\frac{P_{0}}{K}+\frac{mg}{k}\right]^{-1} \ \ \ \ \ (11)

Plugging in the numbers, we get

\displaystyle  z=1416\mbox{ m} \ \ \ \ \ (12)

Thus we’d expect clouds to form at around 1.4 km above the surface.