Chapter-04-Applications-Thermodynamics-Earth 4744 0 PDF

W. M.
White Geochemistry
Chapter 4: Applications of Thermodynamics
Chapter 4: Applications of Thermodynamics to the

Earth
4.1 Introduction
I
n the previous 2 chapters, we developed the fundamental thermodynamic relationships and
saw how they are applied to geochemical problems. The tools now in our thermodynamic
toolbox are sufficient to deal with most of the phenomena we will encounter in the second half
of this book. They are not sufficient, however, to deal with all geochemical problems. In this
chapter, will add a final few thermodynamic tools. These allow us to deal with non-ideal behavior
and exsolution phenomena in solids and silicate liquids. With that, we can use thermodynamics to
determine the pressure and temperature at which rock assemblages formed, certainly one of the
most useful applications of thermodynamics to geology. Along the way, we will see how thermo-
dynamics is related to one of the most useful tools in petrology: phase diagrams. Finally, we return
to the question of non-ideal behavior in electrolyte solutions and examine in more depth the prob-
lems of ion association and solvation and how this affects ion activities. Deviations from ideal be-
havior tend to be greater in solutions of high ionic strength, which includes such geologically im-
portant solutions as hydrothermal and ore-forming fluids, saline lake waters, metamorphic fluids,
and formation and oil field brines. We briefly examine methods of computing activity coefficients
at ionic strengths relevant to such fluids.
4.2 Activities in Non-Ideal Solid Solutions
4.2.1 Mathematical Models of Real Solutions: Margules Equations
Ideal solution models often fail to describe the behavior of real solutions; a good example is wa-
ter and alcohol, as we saw in Chapter 3. Ideal solutions fail spectacularly when exsolution occurs,
such as between oil and vinegar, or between orthoclase and albite, a phenomenon we will discuss
in more detail shortly. In non-ideal solutions, even when exsolution does not occur, more complex
models are necessary.
Power, or Maclaurin, series are often a convenient means of expressing complex mathematical
functions, particularly if the true form of the function is not known. This approach is the basis of
Margules† equations, a common method of calculating excess free energy. For example, we could
express the excess volume as a power series:
Vex = A + BX 2 + CX 22 + DX 23 + K 4.1
where X2 is the mole fraction of component 2.
Following the work of Thompson (1967), Margules equations are used extensively in geochemis-
try and mineralogy as models for the behavior of non-ideal solid solutions. We will consider two
variants of them: the symmetric and asymmetric solution models.
4.2.1.1 The Symmetric Solution Model
In some solutions, a sufficient approximation of thermodynamic functions can often be obtained
by using only a second order power series, i.e., in equ. 4.1, D = E = ... = 0. Now in a binary solu-
tion, the excess of any thermodynamic function should be entirely a function of mole fraction X 2 (or
X 1, however we wish to express it). Put another way, where X 2 = 0, we expect Vex = 0. From this
we can see that the first term in Equ. 4.1, A, must also be 0. Thus equation 4.1 simplifies to:
† Named for M. Margules, who first used this approach in 1895.
118 October 17, 2001

W. M. White Geochemistry
Vex = BX 2 + CX 22 4.2
The simplest solution of this type would be one that is symmetric about the midpoint, X 2 = 0.5;
this is called a Symmetric Solution. In essence, symmetry requires that:
BX 2 + CX 22 = BX1 + CX12 4.3
Substituting (1 – X2) for X1 and expanding the right hand side of 4.3, we have:
BX 2 + CX 22 = B - BX 2 + C - 2CX 2 + CX 22 4.4
Collecting terms and rearranging:
B(2X2 – 1) = C(1 – 2X2 ) 4.5
which reduces to B = – C. Letting WV = B in equation 4.2, we have:
V ex = WV X2 - WV X22 = WV X2 (1 - X2 ) = X1 X2 WV 4.6
W is known as an interaction parameter (recall that non-ideal behavior arises from interactions be-
tween molecules or atoms), and depends on temperature, pressure, and the nature of the solution,
but not on X. Expressions similar to 4.2–4.6 may be written for enthalpy, entropy, and free energy;
for example:
G ex = X1 X2 WG 4.7
The WG term may be expressed as: WG = WU + PWV – TWS = WH – TWS 4.8
The WH term can be written as: WH = WU + PWV
so that 4.8 may also be written: WG = WH – TWS 4.8a
The temperature and pressure dependence of WG are then
Ê ∂WG ˆ Ê ∂WG ˆ
Á ˜ = -Ws 4.9 Á ˜ = WV 4.10
Ë ∂T ¯ P Ë ∂P ¯ T
Regular solutions‡ are a special case of symmetric solutions where:
Ws = 0 and therefore W G = WH
From equation 4.9, we see that WG is independent of temperature for regular solutions. Regular so-
lutions correspond to the case where DSex = 0, i.e., where DSmixing = DS ideal. Examples of such solu-
tions include electrolytes with a single, uncoupled, anionic or cationic substitution, e.g.,
CaCl2—CaBr2, or solid solutions where there is a single substitution in just one site (e.g.,
Mg2SiO4—Fe2SiO4).
Setting equation 4.7 equal to equation 3.57, we have for binary solutions:
G ex = X1 X2 WG = RT [ X1 ln l1 + X2 ln l2 ] 4.11
For a symmetric solution we have the additional constraint that at X2 = X1, l1= l2. From this rela-
tionship it follows that:
RT ln li = X j2 WG 4.12
This leads to the relationships:
m1 = m1o + RT ln X1 + X22 WG 4.13
m 2 = m 2o + RT ln X2 + X12 WG 4.13a
The symmetric solution model should reduce to Raoult’s and Henry’s Laws in the pure substance
and infinitely dilute solution respectively. We see that as X1Æ 1 equations 4.13 and 4.13a reduce re-
spectively to:
m1 = m1o + RT ln X1 4.14
‡ The term regular solution is often used to refer to symmetric solutions. In that case, what we termed a regular
solution is called a strictly regular solution.
119 October 17, 2001

4 DG excess
m 2 = m 2o + RT ln X2 + WG 4.15
Equation 4.14 is Raoult’s Law; letting: 2
µ* = µ° + WG
DG (kJ/mol)
or WG = RT ln h
then 4.15 is Henry’s Law. Thus the 0
interaction parameter can be related DG real
to the parameters of Henry’s Law,
and activity coefficient. In the Mar- -2
gules representation, a solution that is
ideal throughout is simply the special
case where A= B = C = D = ... = 0. -4
4.2.1.2 The Asymmetric Solution
DGi ideal
Model -6
Many real solutions, for example 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
mineral solutions with asymmetric XAb
solvi, are not symmetric. This corre-
sponds to the case where D in equa- Figure 4.1. Alkali feldspar solid solution computed at 600°
tion 4.01 is nonzero; i.e., we must C and 200 MPa (2 kb) using the data of Thompson and
carry the expansion to the third order. Waldbaum (1969). DGreal = DGideal + DGexcess .
It can be shown that in this case the 0
excess free energy in binary solutions
is given by:
G ex = (WG1 X2 + WG2 X1 ) X1 X2 4.16
500∞
C
(You can satisfy yourself that this may
550∞C
-1
be written as a power-series to the
600∞C
third order of either X 1 or X 2.) The
two coefficients are related to the
DG (kJ/mol)
Henry's Law constants:

650∞
-2 C
WGi = m *i - m io = RT ln hi 4.17
700∞
and also depend on pressure. Activ- C
ity coefficients are given by:
RT ln li = (2WG j - WGi ) X j2
4.18 -3
+2(WGi - WG j ) X 3j
where j=2 when i =1 and visa versa.
800∞C
As for the symmetric solution model,
the interaction parameters of the -4
asymmetric model can be expressed 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
as the sum of the W U , WV , and WS in- XAb
teraction parameters to account for Figure 4.2. DG
real of alkali feldspar solution computed for a
temperature and pressure dependen- series of temperatures and 200 MPa.
cies.
The alkali feldspars (NaAlSi3O8
–KAlSi3O8) are an example of a solid solution exhibiting asymmetric exsolution. Figure 4.1 shows
the DG real , DG ideal, and DG excess for the alkali feldspar solid solution computed for 600° C and 200
MPa using the asymmetric solution model of Thompson and Waldbaum (1969). DG excess is com-
120 October 17, 2001

puted using equation 4.16, 2000
DG ideal is computed using equa-
tion 3.30. Figure 4.2 shows DG (kJ)
DG real computed for a series of
temperatures. Perhaps a clearer 0
picture of how DG will vary as
a function of both composition
and temperature can be ob-
tained by plotting all 3 variables
simultaneously, as in Figure 4.3.
4.3 Exsolution
Phenomena
Now consider a binary sys-
tem, such as NaAlSi3O8—
KAlSi3O8 in the example above, X Ab 0.4 0.2
of components 1 and 2, each of 0.6
which can form a pure phase, 0.8
but also together form a solu-
tion phase, which we will call 600
c. The condition for spontane-
ous exsolution of components 1
800
and 2 to form two phases a and T (K) 1000
b is simply that G a + Gb < G c.
Figure 4.5a illustrates the varia- Figure 4.3. DG surface for the alkali-feldspar solid solution as a
function of the mole fraction albite and temperature.
Example 4.1: Computing Activities Using Margules Parameters
Compute the activity of albite in an
albite (Ab) and orthoclase (Or) solid 5
solution (alkali feldspar) as a function of
the mole fraction of albite from XAb = 0
to 1 at 600° C and 200 MPa. Use the 4
asymmetric solution model and the data g
of Thompson and Waldbaum (1969) Na
given below. 3
Alkali Feldspar Margules Parameters
Ab Or
2
WV (J/MPa-mol) 3.89 4.688
WS (J/mol) 19.38 16.157 1
WH (kJ/mol) 26.485 32.105 a Na
Answer: Our first step is to calculate
W G for each end member where WG = 0
W H + Wv P – WS T. Doing so, we find
0 0.2 0.4 0.6 0.8 1
XAB
WGAb = 10.344 kJ and WGOr = 18.938 kJ.
Figure 4.4. Activity and activity coefficient of albite
We can then calculate the activity coef-
in alkali feldspar solid solution computed at 600° C
ficient as a function of XAb and XO r from and 200 MPa using the asymmetric solution model of
equ. 4.19. The activity is then com- Thompson and Waldbaum (1969).
puted from a = lX A b . The results are
plotted in Figure 4.4.
121 October 17, 2001

tion of G excess, G mixture, and G ideal in this hypothetical
system. G mixture is simply the free energy of a me- (a) G excess
chanical mixture of phases a and b (e.g., orthoclase
and albite), which are the phases of pure components
1 and 2, and Greal is equal to G excess + G ideal. So long Gmixture
as Greal is less than Gmixture, a solution is stable relative
to pure phases a and b. You can see that G ideal is al- G
ways less than Gmixture, so as long as the G excess term is Gideal
not too great, and there will be complete solution.
We can write out the equation for the (molar) free
energy in a binary solution as:
G = X1m1o + X 2 m o2 + RT ( X1 ln X1 + X 2 ln X 2 ) 0 X2 1
4.19
+ RT ( X1 ln l1 + X 2 ln l 2 ) (b)
Substituting Gexcess for the last term and differentiating
with respect to X2, we obtain: Greal
Ê ∂G ˆ X1 Ê ∂Gex ˆ
ma2 mb2
Á ˜ = m1 + m 2 + RT
o o
Á ˜ 4.20 G
Ë ∂X 2 ¯ X 2 Ë ∂X 2 ¯
This is the equation for the slope of the curve of G vs. ma1 mb1
a+b
X 2. The second derivative is:
Ê ∂2G ˆ RT Ê ∂ 2 Gex ˆ X2a X2b
Á 2˜ = X X Á 2 ˜
4.21
Ë ∂X 2 ¯ 1 2 Ë ∂X 2 ¯
This tells us how the slope of the curve changes with
0 X2 1
composition. For an ideal solution, G excess is 0, the
omogeneous soluti
second derivative is always positive, and the free en-
ergy curve is concave upward. But for real solutions
ne h o
o
n
G excess can be positive or negative. G excess depends on
the lambda terms, or interaction parameters, which in a
turn depend on composition (for non-Henry's Law T b
solutions). If for some combination of T and X (and two coexisting
P), the second derivative of G excess becomes negative solutions
and its absolute value is greater than the RT/X1X2
term, inflection points appear in the G-X curve. Exso-
lution is thermodynamically favored (Figure 4.5b) in
this situation. 0 X2 1
Suppose now that phases a and b are not pure sub- Figure 4.5. (a) Schematic isothermal, iso-
stances, but are each limited solutions of both compo- baric G-X plot for a real solution showing
nents 1 and 2, as is illustrated in Figure 4.5b. Two in- DG of mechanical mixing, ideal mixing and
flection points appear, one at X a2 and one at X 2b . We excess mixing. (b) Sum of ideal and excess
can draw a straight line that is tangent to the free en- mixing free energies shown in (a). Tangents
ergy curve at both inflection points. This line is the to the minima give the chemical potentials
free energy of a mechanical mixture of the two lim- in immiscible phases a and b. (c). T-X plot
ited solutions a and b. Phase a is mostly component for same system as in (b). Solid line is the
solvus, dashed line is the spinodal. Ex-
1, but contains X a2 of component 2. Similarly, phase
solution may not occur between the
b is mostly component 2 but contains 1 – X 2b of com- spinodal and solvus because the free energy
ponent 1. The mechanical mixture of a and b has less can locally increase during exsolution.
free energy than a single solution phase everywhere From Nordstrom and Munoz (1986).
between X a2 and X 2b . It is therefore thermodynami-
122 October 17, 2001

cally more stable, so exsolution can occur in this region.
Since the free energy of ideal mixing becomes less negative with decreasing temperature (Figure
3.6), many solids show complete solution at higher temperature and decreasing miscibility at lower
temperature. In Figure 4.2, we can see inflection points developing at about 650° C the alkali feld-
spar solution. The inflection points become more marked and occur at increasingly different values
of XAb as temperature decreases. Figure 4.5c shows a schematic drawing of a tempera-
ture–composition plot in which there is complete solution at higher temperature with a widening
two-phase region at lower temperatures. The boundary between the two-phase and one-phase re-
gions is shown as a solid line and is known as the solvus.
The analysis of exsolution above is relevant to immiscible liquids (e.g., oil and vinegar, silicate
and sulfide melts) as well as solids. There is a difference, however. In solids, exsolution must occur
through diffusion of atoms through crystal lattices, while in liquids both diffusion and advection
serve to redistribute components in the exsolving phases. As exsolution begins, the exsolving
phases begin with the composition of the single solution and must rid themselves of unwanted
components. In a solid, this only occurs through diffusion, which is very slow. This leads to a ki-
netic barrier that often prevents exsolution even though 2 exsolved phases are more stable than 1 so-
lution. This is illustrated in Figure 4.6. For example, consider a solution of composition C. It be-
gins to exsolve protophases of A and B, which initially have compositions A´ and B´. Even though
a mechanical mixture of A and B will have lower free energy than solution phase C, A´ and B´, the
initial products of exsolution, have higher free energy than C. Furthermore, as exsolution proceeds
and these phases move toward compositions A and B, this free energy excess becomes larger. Thus
exsolution causes a local increase in free energy and therefore cannot occur. This problem is not
encountered at composition C’ though, be-
cause a mixture of the exsolving proto-
phases A” and B” has lower free energy
than original solution at C´. Thus the ac- C' C
tual limit for exsolution is not inflection
S
se
points such as B but at inflection points
DGunmix 1
pha
(where ∂2G/∂X2 = 0) such as S. The locus A"
of such points is plotted in Figure 4.5c as B" A' + B'
the gray line and is known as the spinodal.
A' G
B
4.4 Thermodynamics and Phase B'
Diagrams
B
A phase diagram is a representation of A+
the regions of stability of one or more
phases as a function of two or more ther-
modynamic variables such as temperature,
pressure, or composition. In other words, X2
if we plot 2 thermodynamic variables such Figure 4.6. A small portion of a G-X plot illustrating
as temperature and pressure or tempera- the origin of the spinodal. The process of exsolution of
ture and composition, we can define an two phases from a single solid solution must overcome
area on this plot where a phase of interest an energy barrier. As exsolution from a solution of
is thermodynamically stable. Figure 4.7 is composition C begins, the two exsolving phases have
an example of a T-P phase diagram for a 1 compositions that move away from C, e.g., A’ and B’.
component system: SiO2. The diagram But the free energy of a mechanical mixture of A’ and
shows the SiO2 phase stable for a given B´ has greater free energy, by DGunmix than the original
combination of pressure and temperature. single solution phase. Exsolution will therefore be in-
Figure 4.8 is an example of a simple T-X hibited in this region. This problem does not occur if
the original solution has composition C’.
diagram for the two-component system
123 October 17, 2001

diopside-anorthite (CaMgSi2O6 or clinopy- 30
roxene and Ca-plagioclase, CaAl2Si2O8; two Coesite
of the more common igneous rock forming
minerals). In multicomponent systems we
must always be concerned with at least 3
thermodynamic variables: P, T, and X. Thus
any T-X phase diagram will be valid for P 20
b-Quartz
5 kb a-Quartz
only one pressure, 1 atm (ª 1 bar = 10 Pa) in
this case. Of course, with a three dimen- L
sional drawing it is possible to represent Cristobalite
both temperature and pressure as well as Tridymite
composition in a binary system.
It should not surprise you at this point to 10
hear that the phase relationships in a chemi- 1000 2000
cal system are a function of the thermody- T ∞C
namic properties of that system. Thus phase
Figure 4.7. P-T phase diagram for SiO2 . This system
diagrams, such as Figures 4.7 and 4.8, can
be constructed from thermodynamic data. has 1 component but 7 phases. L designates liquid,
dashed lines indicate where phase boundaries are un-
Conversely, thermodynamic information
certain. The a—b quartz transition is thought to be
can be deduced from phase diagrams.
Let’s now see how we can construct partially second order, that is, it involves only
phase diagrams, specifically T-X diagrams, stretching and rotation of bonds rather than a com-
plete reformation of bonds as occurs in first order
from thermodynamic data. Our most im- phase transitions.
portant tool in doing so will be the G – X
diagrams that we have already encountered. –
The guiding rule in constructing– phase diagrams
– from G – X diagrams is that the stable phases are those
that combine to give the lowest G . –Since a G – X diagram is valid for only one particular tempera-
ture, we will need a number of G – X dia-
grams at different temperatures to construct
a single T – X diagram (we could also con- 1600
–
struct P – X diagrams from a number of G 1553 P = 0.1 MPa
– X diagrams at different pressures). Before
we begin, we will briefly consider the ther-
modynamics of melting in simple systems.
T, ∞C
4.4.1 The Thermodynamics of Melting 1400 L

1391
One of the more common uses of phase An + L
diagrams is the illustration of melting rela-
tionships in igneous petrology. Let’s con- Di + L
sider how our thermodynamic tools can be 1274
applied to understanding melting relation- An + Di
ships. We begin with melting in a simple 1200
one component system, for example quartz. An XDi Di
At the melting point, this system will con-
sist of two phases: a solid and a melt. At Figure 4.8. Phase diagram (T-X) for the two-
the melting point, the liquid and solid are component system diopside-anorthite at 1 atm. Four
in chemical equilibrium. Therefore, ac- combinations of phases are possible as equilibrium as-
cording to equation 3.17: µl = µs . semblages: liquid (L), liquid plus diopside (L + Di),
liquid plus anorthite (L + An), and diopside plus
anorthite.
124 October 17, 2001

The Gibbs Free Energy of melting, DGm , must be 0 at the melting point (and only at the melting
point). Since
D G m = D H m - Tm D S m 4.22
and DGm = 0 at Tm , then:
DHm = Tm DSm
where DHm is the heat (enthalpy) of melting or fusion*, Tm is the melting temperature, and DSm is
the entropy change of melting. Thus the melting temperature of a pure substance is simply:
DH m
Tm = 4.23
DS m
This is a very simple, but very important, relationship. This equation tells us that temperature of
melting of a substance is the ratio of the enthaply change to entropy change of melting. Also, if we
can measure temperature and enthalpy change of the melting reaction, we can calculate the entropy
change.
The pressure dependence of the melting point is given by the Clapeyron Equation:
dT DVm
= 4.24
dP DS m
Precisely similar relationships hold for vaporization (boiling). Indeed, the temperature and pressure
boundaries between any two phases, such as quartz and tridymite, calcite and aragonite, etc., de-
pend on thermodynamic properties in an exactly analogous manner.
In equation 3.66 we found that addition
of a second component to a pure substance 1600
depresses the melting point. Assuming DSm
and DH m are independent of temperature,
we can express this effect as:
1500 1 Liquid
Ti , m R ln ail, m
Temperature, ∞C
= 1- 4.25
T DSi , m
Since enthalpies of fusion, rather than
1400 An + L
entropies, are the quantities measured,
equation 4.25 may be more conveniently
expressed as:
Di + L
Ti , m Ti , m R ln ail, m 1300
= 1- 4.26
T DHi , m An + Di
Example 4.2 shows how the approximate
phase diagram for the diopside-anorthite 1200
system (Figure 4.9) may be constructed 0 0.2 0.4 0.6 0.8 1
using this equation. X Di Di
It must be emphasized that in deriving
An
equation 3.66, and hence the equations 4.25 Figure 4.9. Computed phase diagram for the system
and 4.26, we made the assumption that the Anorthite-Diopside (CaAl2 Si2 O8 –CaMgSi2 O6 ). The eu-
solid was a pure phase. This assumption is tectic occurs at XDi = 0.7 and 1334°C. The dashed lines
a reasonably good one for ice, and for beyond the eutectic give the apparent melting points of
anorthite-diopside binary system, but it is the components in the mixture.
* The heat of fusion is often designated by DH . I have chosen to use the subscript m to avoid confusion with heat
f
of formation, for which we have already been using the subscript f.
125 October 17, 2001

Example 4.2: Calculating Melting Curves

Using the data given below and assuming (1) Tm DH m
that the melt is an ideal solution and (2) diopside °C joules/mole
and anorthite solids are pure phases, calculate a T- Diopside 1391 138100
X phase diagram for melting of an anorthite-diop- Anorthite 1553 136400
side mixture.
(Data from Stebbins et al., 1983)
Answer: Solving equation 4.25 for T, and replac-
ing activity with mole fraction (since we may
assume ideality), we have:
DHi , m
T= 4.27
DHi , m
Ti , m - R ln X
l
i
We then calculate T for every value of X An and X Di . This producestwo curves on a T-X plot, as
shown in Figure 4.09. The curvesintersect at the eutectic , or lowest point at which melt may exist
in the system.
Comparing our result with the actual phase relationships determined experimentally
(Figure 4.9), we see that while the computed phase diagram is similar to the actual one, our
computedeutectic occursat X Di = 0.70 and 1335° C and the actual eutectic occursat X Di ª 0.56 and
1274° C. The difference reflects the failure of the several assumptions we made. First, and most
importantly, silicate liquids are not ideal solutions. Second, the entropies and enthalpies of
fusion tend to decrease somewhat with decreasing temperature, violating the assumption we
made in deriving equation 4.25. Third, the diopside crystallizing from anorthite-diopside
mixtures is not pure, but contains some Al and an excess of Mg.
not generally valid. Should the solid or solids involved exhibit significant solid solution, this
assumption breaks down and these equations are invalid. In that case, melting phase diagrams can
still be constructed from thermodynamic equations, but we need to model solid solution as well as
the liquid one. Section 4.4.2.1 below illustrates an example (Anorthite-Albite) we the two solutions
can be modeled as ideal.
4.4.2 Thermodynamics of Phase Diagrams for Binary Systems
In a one component system, a phase boundary, such as the melting point, is univariant since at
that point two phases coexist and ƒ = c- p + 2 = 1 - 2 + 2 = 1. Thus specifying either temperature or
pressure fixes the other. A three-phase point, e.g., the triple point of water, is invariant. Hence
simply from knowing that three phases of water coexist (i.e., knowing we are at the triple point),
we know the temperature and pressure.
In binary systems, the following phase assemblages are possible according to the Gibbs Phase
Rule (ignoring for the moment gas phases):
Phases Free compositional variables
Univariant 2 solids + liquid, 2 liquids + solid, 3 solids or liquids 0
Divariant 1 solid + 1 liquid, 2 solids, 2 liquids 0
Trivariant 1 solid or 1 liquid 1
– –
When a G -X diagram is drawn, it is drawn for a specific temperature and pressure, i.e., G -X are
isobaric and isothermal. Thus we have already fixed two variables, and the compositions of all
phases in univariant and divariant systems are fixed by virtue of our having
– fixed T and P. Only in
trivariant systems are we free to choose the composition of a phase on a G -X diagram.
126 October 17, 2001

Figure 4.10 is schematic diagram of a two
component, one phase (trivariant) system, in which
there is complete solution between component 1 and
component 2. This phase might be either a liquid, or a
G L solid such as plagioclase. The composition of the phase
may fall anywhere on the curve. Of course, since this
diagram applies only to one temperature, we cannot
say from this diagram alone that there will be complete
solution at all temperatures.
X1 X2 Figure 4.11 illustrates four possible divariant systems.
0 X2 1 The first case (Figure 4.11a) is that of a liquid solution of
Figure 4.10. –Molar free energy vs. composition L' in equilibrium with a solid of fixed
composition (G – X2) for a one phase composition S2 (pure component 2). Because the system
system which exhibits complete solu-
tion of either a liquid or solid.
L
is divariant, there can be only one pos- GS2 S2
sible liquid composition since we have
implicitly specified P and T. As usual, the GL' L
l GL' S
equilibrium condition is described by µi
= m i (equation 3.17). For i = 2, this means
s
L' + S2 GS'
L'+S'
the tangent to the free energy curve for
the melt must intersect the X 2 = 1 line at
µ s2 as is shown. In other words, the L' L'
chemical potential of component 2 in the 0 X2 1 0 X2 1
melt must be equal to the chemical poten- (a) (b)
tial of component 2 in the solid. Again,
this diagram is valid for only one tem-
perature; at any other temperature, the GS2 S2
free energy curve for the liquid would be
S
different, but the composition of this new S1
liquid in equilibrium with solid S2 would GS1
still be found by drawing a tangent from S1 + S 2 S1 + S 2
S2 to the free energy curve of the liquid.
At sufficiently high temperature, the
tangent would always intersect below S2. S1 S2
The temperature at which this first occurs 0 X2 1 0 X2 1
is the melting temperature of S2 (because (c) (d)
it is the point at which the free energy of Figure 4.11. Plot of molar free energy vs. composition
a liquid of pure 2 is less than the solid). –
(G – X2) for two phase (divariant) systems. (a) shows a
The shaded region shows the composi- liquid solution (L) in equilibrium with a solid (S2) of
tions of systems that will have a combi- pure X2. The shaded area shows the range of composi-
nation of solid S2 and liquid L’ as their tion of systems for which L’ and S2 will be the equilib-
equilibrium phases as this temperature. rium phases as this temperature. (b) is the case of where
We can also think of the tangent line as both solid and liquid have variable composition. Equi-
defining the free energy of a mechanical librium compositions at this temperature are determined
mixture of S2 and L’. In the range of by finding a tangent to both free energy curves. L’ and
compositions denoted by the shaded re- S’ will be the equilibrium phases for systems having
gion, this mixture has a lower free energy compositions in the shaded area. (c) is the case of 2 im-
than the liquid solution, hence at equilib- miscible solids. These solids will coexist at equilibrium
over the entire compositional range at this temperature.
127 October 17, 2001

rium we expect to find this mixture rather
than the liquid solution.
Figure 4.11b illustrates a system with a
liquid plus a solid solution, each of which has S2
its own G-X curve. Again, the equilibrium
L G L
l s
condition is µ i = µ i so the compositions of
S2
the coexisting liquid and solid are given by a
tangent to both curves. Since the system is S1
divariant and we have fixed P and T the com-
positions of the solutions are fixed. All S1 L' S2 L'
0 X 2 1 0 X 2 1
system compositions in the shaded region can
be accommodated by a mixture of liquid and (a) (b)
solid. Compositions lying to the left of the re- Figure 4.12. Two univariant systems. The first is
gion would have only a liquid; compositions that of a liquid plus two solid solutions, the second
to the right of the shaded region would be ac- is that of two pure solids and a liquid. Since these
commodated by a solid solution. systems are univariant, they occur only at one fixed
Figure 4.11c illustrates the case of two T if P is fixed.
immiscible solids (pure components 1 and 2).
The molar free energy of the system is simply
that of a mechanical mixture of S1 and S2: a straight line drawn between the free energy points of the
two phases.
Figure 4.11d illustrates the case of a limited solution. We have chosen to illustrate a solid solu-
tion, but the diagram would apply equally well to the case of two liquids of limited solubility.
Figure 4.12a shows the case of two solid solutions plus one liquid. The chemical potential of
each component in each phase must be equal to the chemical potential of that component in every
other phase, so chemical potentials are given by tangents to all three phases. This is an univariant
system, specifying either temperature, pressure, or the composition of a phase fixes other variables
in the system. Because of this, if we move to a slightly higher or low temperature at fixed pressure
one of the phases must be eliminated in a phase elimination reaction. If the liquid is eliminated, the
diagram represents the eutectic, the lowest temperature at which the liquid can exist. Moving to a
higher temperature would result in elimination of one of the solids. If, alternatively, one of the sol-
ids is eliminated by moving to lower temperature, and a liquid is stable both above and below this
point, the point would be known as a peritectic. Figure 4.12b is a eutectic in a system where the
two solids are the phases of pure components 1 and 2. A line drawn between the free energies of
the pure components is also tangent to the liquid curve.
4.4.2.1 An Example of a Simple Binary System with Complete Solution: Albite–Anorthite
Phase diagrams in T-X space can be constructed by analyzing G-X diagrams at a series of temper-
atures. Let's examine how this can be done in the case of a relatively simple system of two compo-
nents albite (NaAlSi3O8) and anorthite (CaAl2Si2O8) whose solid (plagioclase) and liquid exhibit com-
plete solid solution. Figure 4.13 shows G-X diagrams for various temperatures as well as a T-X
phase diagram for this system. Since both the solid and liquid exhibit complete solution, we need
to consider G-X curves for both.
We start at the highest point at which liquid and solid coexist, Tm (T1) for anorthite. Here the
solid and liquid curves both have the same value at XAn = 1; i.e., they are at equilibrium. A G-X
plot above this temperature would show the curve for the liquid to be everywhere below that of
the solid, indicating the liquid to be the stable phase for all compositions.
At a somewhat lower temperature (T2), we see that the curves for the solid and liquid intersect at
some intermediate composition. To the right, the curve for the solid is lower than that of the liquid,
and tangents to the solid curve extrapolated to both XAb =1 and XAn =1 are always below the curve
128 October 17, 2001

for the liquid, indicating the solid is the stable phase. As we move toward Ab (left) in composition,
tangents to the solid curve eventually touch the curve for the liquid. The point where the tangent
touches each curve gives composition of the liquid and the solid stable at this temperature. In the
compositional range between the points where the tangent touches the two curves, the tangent is
below both curves, thus a mechanical mixture of solid and liquid is stable over this compositional
range at this temperature. For compositions to the left of the point where the tangent touches the
liquid curve, the liquid curve is lower than both the solid curve and a tangent to both, so it is stable
relative to both the solid and any mixture
of solid and liquid.
Plag Going to progressively lower tempera-
tures (e.g., T3), the points where a tangent
G L 1 T
intersects the two curves move toward
Ab (to the left). Eventually, at a suffi-
ciently low temperature (T4), the curve
for the solid is everywhere below that of
Plag the liquid and only solid solution is sta-
2 T
ble. By extracting information from G-X
G L curves at a number of temperatures, it is
possible to reconstruct the phase diagram
shown at the bottom of Figure 4.13.
L Since both the solid and liquid show
3 T
complete miscibility in this system, we
G Plag
can assume that both solutions are ideal
and do an approximate mathematical
treatment. We recall that the condition
for equilibrium was:
L
T4 µ ai = µ bi
We can express the chemical potential of
G Plag
each component in each phase as:
mai = m oi a + RT ln Xia 4.28
Combining these relationships, we have:
1600 Ê Xl ˆ
T1 maAbs = m oAbl + RT ln Á Abs ˜
4.29
Ë X Ab ¯
1400 T2 Ê X An
l
ˆ
as
T T3 m =m ol
+ RT ln Á s ˜ 4.30
Ë X An ¯
An An
∞C Plag + L
1200 Plag Here our standard states are the pure end
members of the melt and solid. The left
T4 side of each of these equations corre-
1000 sponds to the free energy change of melt-
ing, thus:
NaAlSi3O8 0.2 0.4 0.6 0.8 CaAl2Si2O8 Ê Xl ˆ
X An DGmAb = RT ln Á Abs ˜
4.31
Ë X Ab ¯
Figure 4.13. G-X diagrams and a T-X phase diagram
for the plagioclase-liquid system. After Richardson and Ê X An
l
ˆ
McSween (1987). DG An
= RT ln Á s ˜ 4.32
Ë X An ¯
m
129 October 17, 2001

l
Both sides of these equations reduce to 0 if and only if X i = X i = 1 and T = Tm . Rearranging:
s
DG mAb
X sAb = X Ab
l
e RT
4.33
DG mAn RT
X s
An =X e l
An 4.34
Thus the fraction of each component in the melt can be predicted from the composition of the solid
l l s s
and thermodynamic properties of the end members. Since X An = 1 – X Ab and X An = 1 – X Ab , we can
combine equations 4.33 and 4.34 to obtain:
(1 - X )e DG mAn RT
e DG m
Ab
l
Ab = 1 - X Ab
l RT
4.35
DG mAn
1- e RT
and rearranging yields: X l
Ab = 4.36
e DG m - e DG m
Ab An
RT RT
The point is that the mole fraction of any component of any phase in this system can be predicted from the
thermodynamic properties of the end-members. We must bear in mind that we have treated this as an
ideal system; i.e., we have ignored any G excess term. Nevertheless, the ideal treatment is relatively
successful for the plagioclase system.
4.5 Geothermometry and Geobarometry
An important task in geochemistry is estimating the temperature and pressure at which mineral
assemblages equilibrate. The importance extends beyond petrology to tectonics and all of geology.
Here we take a brief look at the thermodynamics underlying geothermometry and geobarometry.
Geothermometry and geobarometry involve two nearly contradictory assumptions. The first is
that the mineral assemblage of interest is an equilibrium one, the second is that the system did not
re-equilibrate during the passage through lower P and T conditions that brought the rock to the sur-
face where it could be collected. As we will see in the next chapter, reaction rates depend exponen-
tially on temperature, hence these assumptions are not quite as contradictory as they might seem. In
this section, we will focus only on “chemical” thermobarometers. In Chapter 9, we will see that
temperatures can also be deduced from the distribution of isotopes between phases.
4.5.1 Theoretical considerations
In general, geobarometers and geothermometers make use of the pressure and temperature de-
pendence of the equilibrium constant, K. In Section 3.9 we found that DG° = -RT ln K. Assuming
that DCp and DV of the reaction are independent of temperature and pressure, we can write:
DG o = DH To ,Pref - TDS To ,Pref + DVTo,Pref ( P - Pref ) = - RT ln K 4.37
where the standard state of all components is taken as the pure phase at the temperature and pres-
sure of interest, and the enthalpy, entropy and volume changes are for the temperature of interest
and a reference pressure (generally 0.1 MPa).
Solving 4.37 for ln K and differentiating the resulting equation with respect to temperature and
pressure leads to the following relations:
Ê ∂ ln K ˆ DH To ,Pref + DVTo,Pref ( P - Pref )
Á ˜ = 4.38
Ë ∂T ¯ P RT 2
Ê ∂ ln K ˆ DVTo,Pref
and Á ˜ = 4.39
Ë ∂P ¯ T RT
These equations provide us with the criteria for reactions that will make good geothermometers
and geobarometers. For a good geothermometer, we want the equilibrium constant to depend
heavily on T, but be approximately independent of P. Looking at equation 4.38, we see this means
the DH term should be as large as possible and the DV term as small as possible. A fair amount of
130 October 17, 2001

effort was devoted to development of a geothermometer based on the exchange of Fe and Mg be-
tween olivine and pyroxenes in the late 1960’s. The effort was abandoned when it was shown that
the DH for this reaction was very small. As a rule, a reaction should have a DH˚ of at least 1 kJ to
be a useful geothermometer. For a good geobarometer, we want the DV term to be as large as pos-
sible. Even though the rhodonite ([Mn,Fe,Ca]SiO3) and pyroxmangite ([Mn,Fe]SiO3) pairs com-
monly occur in metamorphic rocks, the reaction rhodonite Æ pyroxmangite does not make a useful
geobarometer because the DV of reaction is only 0.2 cc/mol. In general, a reaction should have a
DV of greater than 2 cc/mol if it is to be used for geobarometry.
The following discussion presents a few examples of useful chemical geothermometers and geo-
barometers (since most reactions are both temperature and pressure dependent, it is perhaps more
accurate to use the term “thermobarometer”). It is not an exhaustive treatment, nor should it be in-
ferred that those examples discussed are in any way superior to other geothermometers and geo-
barometers. Reviews by Essene (1982, 1989) and Bohlen and Lindsley (1987) summarize a wide
range of igneous and metamorphic thermobarometers.
4.5.2 Practical Thermobarometers
4.5.2.1 Univariant Reactions and Displaced Equilibria
We can broadly distinguish 3 main types of thermobarometers. The first is the univariant reac-
tion, in which the phases have fixed compositions. They are by far the simplest, and often make
good geobarometers as the DV of such reactions is often large. Examples include the graphite-dia-
mond transition, any of the SiO2 transitions (Figure 4.7), and the transformations of Al2SiO5, shown
in Figure 4.14. While such thermobarometers are simple, their utility for estimating temperature
and pressure is limited. This is because exact temperatures and pressures can be obtained only if
two or more phases coexist, for example, kyanite and andalusite in Figure 4.14. If kyanite and an-
dalusite are both found in a rock, we can determine either temperature or pressure if we can inde-
pendently determine the other. Where 3
phases, kyanite, sillimanite, and andalusite 800
coexist the system is invariant and P and T
are fixed. If only one phase occurs, for
example sillimanite, we can only set a
range of values for temperature and pres- 600
sure. Unfortunately, the latter case, where
only 1 phase is present, is the most likely
situation. It is extremely rare that kyanite, Ky
P, MPa
sillimanite, and andalusite occur together. 400 Sil

The term displaced equilibria refers to
variations in the temperature and pressure
of a reaction that results from appreciable
solid solution in one or more phases.
Thermobarometers based on this phe-
200
nomenon are more useful than univariant
reactions because the assemblage can co- And
exist over a wide range of P and T condi-
tions. In the example shown in Figure 200 300 400 500 600 700 800
4.15, the boundaries between garnet-bear-
ing, spinel-bearing, and plagioclase-bear- T, ∞C
ing assemblages are curved, or “dis- Figure 4.14. Phase diagram for Al SiO (kyanite-silli-
2 5
placed” as a result of the solubility of Al in manite-andalusite) as determined by Holdaway (1971).
enstatite. In addition to the experimental Due to sluggish reaction kinetics, the exact position of
calibration, determination of P and T from these phase boundaries remains somewhat uncertain.
131 October 17, 2001

displaced equilibria requires (1) careful 3.0
determination of phase composition
and (2) an accurate solution model.
Geobarometers based on the 2.5 Fo + Ga
solubility of Al in pyroxenes have been
the subject of extensive experimental 2.0 18
P, GPa
investigations for the past 25 years. The
general principal is illustrated in Figure 1.5
4.15, which shows the concentration of
Al in orthopyroxene (opx) coexisting 1.0 1 2 4 6 8 10 12 14 16 Crystals
with olivine (forsterite) and an Fo + Sp +
aluminous phase, anorthite, spinel, or 0.5 Liquid
garnet. The Al content of opx depends Fo + An
almost exclusively on pressure in the 0
presence of anorthite, is essentially in- 400 600 800 1000 1200 1400 1600
dependent of pressure in the presence T, ∞C
of spinel, and depends on both Figure 4.15. Isopleths of Al in orthopyroxene (thin red
temperature and pressure in the lines; weight percent) coexisting with forsterite plus
presence of garnet. Orthopyroxene- and aluminous phase in the CMAS (Ca-Mg-Al-Si) sys-
garnet equilibria has proved to be a tem. After Gasparik (1984).
particularly useful geobarometer.
Garnet is an extremely dense phase. So we might guess that the DV of reactions that form it will
be comparatively large, and therefore that it is potentially a good geobarometer. The concentration
of Al in opx in equilibrium with garnet may be used as a geobarometer if temperature can be
independently determined. Although there has been a good deal of subsequent work and re-
finement of this geobarometer, the underlying thermodynamic principles are perhaps best illus-
trated by considering the original work of Wood and Banno (1973).
Wood and Banno (1973) considered the following reaction:
Mg2 Si2 O6 + MgAl2 SiO6 ® Mg3 Al2 Si3 O12 4.40
Opx Solid Solution ® Pyrope Garnet
In developing a geobarometer based on this reaction, they had to overcome a number of problems.
First, the substitution of Al in orthopyroxene is a coupled substitution. For each atom of Al substi-
tuting in the M1 octahedral site, there must be another Al atom substituting for SiO2 in the tetrahe-
dral site. Second, there was a total lack of thermodynamic data on the MgAl2SiO6 phase compo-
nent. Data was lacking for a good reason: the phase does not exist and cannot be synthesized as a
pure phase. Another problem was the apparent non-ideal behavior of the system, which was indi-
cated by orthopyroxenes in Fe- and Ca-bearing systems containing less alumina than in pure MgO
systems at the same pressure.
The equilibrium constant for reaction 4.40 is:
aMg 3 Al 2 Si 3O12
K= 4.41
aMg 2 Si 2O 6 aMgAl 2 SiO 6
where the activities in the denominator represent the activities of the enstatite and the hypothetical
aluminous enstatite phase components in the enstatite solid solution. In the pure MgO system (i.e.,
no CaO, FeO, MnO, etc.), the numerator, the activity of pyrope, is 1, of course, and we may write:
DG o = RT ln( aMg 2 Si 2O 6 aMgAl 2 SiO 6 ) = DH o - TDS o + ( P - Pref ) DV o 4.42
(compare equation 4.37). For an ideal case, this may be rewritten as:
RT ln( X Mg 2 Si 2O 6 X MgAl 2 SiO 6 ) = DH o - TDS o + ( P - Pref ) DV o 4.43
132 October 17, 2001

Wood and Banno first estimated thermodynamic parameters (DH, DS, and DV for aluminous pyrox-
ene) from experimental data. They dealt with the non-ideality in two ways. First, they assumed
ideal solution behavior at 1 bar and assumed all non-ideality associated with substitution of Al in
orthopyroxene at higher pressure could be accounted for in the volume term in 4.42, which they
rewrote as:
DV o = VMg
o
3 Al 2 Si 3O12
- VMg
o opx
2 Si 2O 6
- VMgAl
o opx
2 SiO 6
4.44
As for non-ideality related to substitution of Ca and Fe in the system, they noted that non-idealities
of most silicate systems were of similar size and magnitude and hence the activity coefficients for
garnet tend to cancel those for orthopyroxene. Furthermore, the DV and DH terms are both large
and tend to reduce the errors due to non-ideal behavior.
Since equation 4.42 contains temperature as well as pressure terms, it is obvious that the tempera-
ture must be known to calculate pressure of equilibration. In the same paper, Wood and Banno
(1973) provided the theoretical basis for estimating temperature from the orthopyrox-
ene–clinopyroxene miscibility gap. Thus in a system containing garnet, orthopyroxene and clino-
pyroxene, both temperature and pressure of equilibration may be estimated from the composition
of these phases.
This geobarometer-geothermometer is commonly used to estimate the temperature and pressure
(depth) of equilibration of mantle-derived garnet lherzolite xenoliths. One of the first applications
was by Boyd (1973), who calculated P and T for a number of xenoliths in South African kimberlites,
and hence reconstructed the geotherm in the mantle under South Africa.
4.5.2.2 Solvus Equilibria
Solvus Equilibria provides a second kind of thermobarometer. Generally, these make better
geothermometers than geobarometers. A good example is the ortho- and clinopyroxene system,
illustrated in Figure 4.16. The two-pyroxene solvus has be the subject of particularly intensive ex-
perimental and theoretical work because ortho- and clinopyroxene coexist over a wide range of
conditions in Mg, Fe-rich rocks of the crust and upper mantle.
One of the inherent thermodynamic difficulties with this type of geothermometer is that since it
involves exsolution, ideal solution models will
clearly be very poor approximations. Thus
considerable effort has been made to develop En Pig 0.1 MPa
solution models for the pyroxenes. Several 1400 ss Pig + Diss
factors further complicate efforts to use the
pyroxene solvus as a thermobarometer. The
Diss
first is the existence of a third phase, pigeonite 1200
(a low-Ca clinopyroxene), at high temperatures Enss + Diss
and low pressures; the second is that the system T, ∞C
is not strictly binary: natural pyroxenes in 1000
igneous rocks are solutions of Mg, Ca, and Fe
components. The presence of iron is problem-
atic because of the experimental difficulties 800
encountered with Fe-containing systems. These
difficulties include the tendency both for iron to
dissolve in the walls of commonly used plat- 600
inum containers and for Fe2+ either to oxidize to 0 0.2 0.4 0.6 0.8 10
3+
Fe or to reduce to metallic iron, depending on X Di
the oxygen fugacity. In addition, other Figure 4.16. Phase relationships in the system
components, particularly Na and Al are often Mg2Si2O6 (enstatite) — CaMgSi2O6 (diopside) system
present in the pyroxenes, as we have just seen. (after Lindsley, 1983).
133 October 17, 2001

Di Hd
1000 900 70 0
800
1300 120 11
0 00
En Fs
Figure 4.17. Comparison of calculated (solid lines) and experimentally observed (red dashed lines)
phase relationships between clino- and orthopyroxene shown in the ‘pyroxene quadrilateral’, a part
of the CaSiO3–MgSiO3–FeSiO3 system. Di: diopside, En: enstatite, Hd: hedenbergite, Fs: ferrosilite.
Lines show the limit of solid solution at the corresponding temperatures (° C). Experiments show
solid solution to be complete in the iron end-members at 700 ° C.
Despite its complexities, the system has been modeled with some success using a symmetric solu-
tion model developed by Wood (1987). There are two octahedral sites in both ortho- and clinopy-
roxenes, generally called M1 and M2. Ca2+ occurs only in the M2 site, while Fe and Mg can oc-
cupy either site. Ignoring pigeonite and components other than Ca, Mg and Fe, we can treat mix-
ing in the M2 and M1 sites separately. Mixing in the M2 site can be treated as a ternary Mg, Fe,
and Ca solution. In a symmetric ternary solution consisting of components A, B, and C, the activi-
ties of the components may be calculated from:
(
RT ln g A = X B2WGAB + XC2WGAC + X B XC WGAB + WGAC - WGBC ) 4.45
where W GAB is the A-B binary interaction parameter, etc. Mixing of Fe and Mg between the M1 and
M2 sites was treated as a simple exchange reaction:
FeM2 + MgM1 ® FeM1 + MgM2
with DH of 29.27 kJ/mol and DS of 12.61 j/mol. Using this approach, Wood calculated the tempera-
ture dependence of the solvus in shown in Figure 4.17. The model fits experimental observation
reasonably well for the Mg-rich pyroxenes, but significant deviations occur for the Fe-rich pyrox-
enes.
4.5.2.3 Exchange Reactions
The third type of thermobarometer we will consider is the exchange reaction. These thermo-
barometers depend on the exchange of two species between phases. We will consider two exam-
ples of these.
The Roeder and Emslie olivine-liquid geothermometer is a rather simple one based on the equi-
librium between magma and olivine crystallizing from it. Consider the exchange reaction:
MgOOl + FeOliq ® MgOliq + FeOOl
where Ol denotes olivine and liq denotes liquid. We can write the equilibrium constant for this re-
action as:
liq
X Ol
FeO X MgO
KD = 4.46
X liq Ol
FeO X MgO
134 October 17, 2001

90Fo-10Fa 80Fo-20Fa
20
1350
18
16 1300
14
Mole % MgO in Liquid
1250 70Fo-30Fa
12
1200
10
1150 60Fo-40Fa
8
1100
6 1050 50Fo-50Fa
4 1000
40Fo-60Fa
2 30Fo-70Fa
20Fo-80Fa
0 4 8 12 16 20
Mole % FeO in Liquid
Fig. 4.18. Olivine saturation surface constructed by Roeder and Emslie (1970).
Recalling our criteria for a good geothermometer, we can guess that this reaction will meet at least
several of these criteria. First, olivine exhibits complete solid solution, so we might guess we can
treat it as an ideal solution, which turns out to be a reasonably good assumption. We might also
guess that the molar volumes of forsterite and fayalite and of their melts will be similar, meaning
the DV term, and hence pressure dependence, will be small. As it turns out, the DH term, which is
related to the difference in heats of fusion of forsterite and fayalite, is also relatively small, so the ex-
change reaction itself is a poor geothermometer. However, we can consider two separate reactions
here:
MgOliq Æ MgOOl and FeOliq Æ FeOOl
and we can write two expressions for KD. This was the approach of Roeder and Emslie (1970), who
deduced the following relations from empirical (i.e., experimental) results:
Ol
X MgO 3740
log liq
= - 1.87 4.47
X MgO
T
X Ol 3911
log FeO
liq
= - 2.50 4.48
X FeO T
These KD's are much more temperature dependent than for the combined exchange reaction. Sub-
tracting equation 4.47 from 4.48 yields:
135 October 17, 2001

171
log K D = - 0.63 4.49
T
where KD is defined as in equation 4.46. Note that these equations have the form of equation 3.95.
Roeder and Emslie (1970) used these equations to construct the diagram in Figure 4.18.
Example 4.3. Calculating Magma Temperatures Using the Olivine
Geothermometer
From the electron microprobe analysis of glass of a mid-ocean ridge basalt and its coexisting
olivine microphenocryst, calculate the temperature at which the olivine and liquid equilibrated:
Answer: We will answer this assuming the glass composition rep-
SiO2 50.3 resents that of the liquid and using equations 4.47 and 4.48. To use the
Al2O3 14.3 equations, we will have to convert the analysis of the glass from
SFeO 11.1 weight percent to mole fraction.
MgO 7.2 Let’s setup a spreadsheet to do these calculations. First we must
CaO 11.5 deal with the Fe analysis. The analysis reports only iron as FeO.
Na2O 2.6 Generally, about 10% of the iron in a basaltic magma will be present
K 2O 0.23 as ferric iron (Fe2O 3), so we will have to assign 10% of the total iron to
MnO 0.20 F e 2 O 3 . To do this, we get the weight percent FeO simply by
TiO2 1.71 multiplying the total FeO by 0.9. To get weight percent Fe O , we
Total 99.02 multiply total FeO (11.1%) by 0.1, then multiply by the ratio2 of3 the
Mol % Fo in Ol 82 molecular weight of Fe2O 3 to FeO and divide by 2 (since there are 2
w t % w/10% ferric Mol. wt moles mol frac. Fe atoms per ‘molecule’).
Now we are ready to cal-
SiO2 50.3 50.3 60.09 0.8371 0.5265
culate the mole fractions.
Al2O3 14.3 14.3 102 0 . 1 4 0 2 0 . 0 8 8 2 We’ll set up a column with
total FeO 1 1 . 1 11.1 molecular weights and divide
FeO 9.99 7 1 . 8 5 0 . 1 3 9 0 0 . 0 8 7 5 each weight percent by the
Fe2O3 1.22 1 5 7 . 7 0 . 0 0 7 7 0 . 0 0 4 9 molecular weight to get the
MgO 7.8 7.8 40.6 0 . 1 9 2 1 0 . 1 2 0 8 number of moles per 100
CaO 11.5 11.5 5 6 . 0 8 0 . 2 0 5 1 0 . 1 2 9 0 grams. To convert to mole
Na2O 2.6 2.6 6 1 . 9 8 0 . 0 4 1 9 0 . 0 2 6 4 fraction, we divide the
K2O 0.23 0.23 94.2 0 . 0 0 2 4 0 . 0 0 1 5 number of moles by the sum
MnO 0.2 0.2 7 0 . 9 4 0 . 0 0 2 8 0 . 0 0 1 8 of the number of moles.
TiO2 1.71 1.71 79.9 0.0214 0.0135 Since the mole fraction of
Mg in olivine is equal to the
Total 99.74 99.85 1.590 1.000
mole fraction of forsterite, we
XMgO-Ol 0.82
need only convert percent to
XFeO-Ol 0 . 1 8 fraction (i.e., divide by 100).
The mole fraction of FeO in
TMgO 1384 kelvin 1 1 1 1 °C olivine is simply 1 – XM g O .
TFeO 1390 kelvin 1 1 1 7 °C Thus XMg O (o l) = 0.82 and X F eO
(o l) = 0.18. Now we are ready to calculate temperatures. We can calculate 2 temperatues: one
from MgO, and the other from FeO. The temperature based on the FeO exchange is:
3911 3740
TFeO = and that based on MgO is: TMgO =
Ê X Ol ˆ Ê X MgO
Ol
ˆ
log Á FeO
liq ˜
+ 2.50 log Á liq ˜ + 1.87
Ë X FeO ¯ Ë X MgO ¯
We find that the temperatures of the two methods agree within 6°, which is fairly good. This in-
dicates the analyzed olivine probably was in equilibrium with the liquid.
136 October 17, 2001

The iron-titanium oxide system evaluated by
Buddington and Lindsley (1964) was one of the first
TiO 2
means of obtaining quantitative estimates of
crystallization temperatures of igneous rocks. It is
important not only because it is useful over a wide
range of temperatures and rock types, but also be-
cause it yields oxygen fugacity as well. Figure 4.19
shows the TiO2–FeO–Fe2O3 (rutile–wüstite–hematite) FeTiO 3 ilm
ternary system. The geothermometer is based on the enit
reaction: Fe 2TiO 4 e–h
ema
t ite
u lv g n e t
1
yFe 2 TiO 4 + (1 - y )Fe 3O 4 + O 2 s.s.
ma
osp ite
ine s.s.
4
4.50
l–
Ê3 ˆ
® yFeTiO 3 + Á - y˜ Fe 2 O 3
Ë2 ¯ FeO Fe 3O 4 Fe 2O 3
which describes equilibrium between the ul- Figure 4.19. The TiO2–FeO–Fe2O3 ternary
vospinel–magnetite (titanomagnetite) and ilmen- system. Phases are: FeO: wüstite; Fe2O 3:
ite–hematite solid solution series. The equilibrium hematite; TiO2: rutile; Fe2TiO4: ulvospinel;
constant expression may be written as: Fe2O 4: magnetite; FeTiO3: ilmenite. The sys-
tem also includes the FeTi2O5—Fe2TiO5 solu-
tion, which is not shown.
–5 y 3/ 2 - y
a FeTiO a Fe
K= 3 2O 3
4.51
y
a Fe 2 TiO 4
a1Fe- y3O 4 fO1/2 4
p 40
Us Usp 60 The original Buddington and
sp 20 p 80
Ups 85
p0 Lindsley geothermometer was
–10 Us U based on empirical (experimental)
Us
1 00 -
observations of compositional
m
dependence on oxygen fugacity

He
and temperature, as shown in

Figure 4.20. Having values for
p 10
Us
the compositions of the titano-

log ƒO2
–15 magnetite and ilmenite phases,

one simply read T and ƒO 2 from
the graph. To understand the
system from a thermodynamic
m15 perspective, it is better to con-
He 10 sider the two fundamental
–20 m reactions occurring separately in
He this system:
5
m m
1
He He
2Fe3 O4 + O2 ® 3Fe2 O3 4.52

2
3
Hem1 FeTiO3 +FeO®Fe2 TiO4 4.53

The first reaction is the oxida-
–25 tion of ferrous iron to ferric iron.
500 600 700 800 900 1000 1100 The second reaction is the parti-
T, ∞C tioning of Fe between the tita-
Figure 4.20. Relationship of composition of coexisting titano-
magnetite and ilmenite to temperature and oxygen fugacity.
137 October 17, 2001

nomagnetite phase and the ilmenite phase.
Several investigators have studied the iron-titanium oxides attempting to improve upon the work
of Buddington and Lindsley (1964). The approach of Spencer and Lindsley (1981) was to consider
two reactions:
Fe3O4 + FeTiO3 ® Fe2TiO4 + Fe2O3
Magnetite + Ilmenite ® Ulvospinel + Hematite
and: 4Fe3 O4 + O2 ® 6Fe2 O3
The first reaction represents a temperature dependent exchange between the titanomagnetite and ul-
vospinel solutions; the second reaction is the oxidation of magnetite to hematite. They modeled the
ilmenite as a binary asymmetric Margules solution and titanomagnetite as a binary asymmetric
Margules solution below 800° C and as an ideal binary solution above 800° C. They modeled
configurational entropy based ordering of Fe 2+ , Fe 3+ , and Ti4+ in the ilmenite lattice structure (they
assumed Fe3+ mixed randomly with Ti4+ in ‘A’ sites and Fe3+ and Fe2+ randomly in ‘B’ sites). The D G
of reactions above were written as:
DG È XUsp
a
(1 - X Ilm )a ù Èg Usp
a
g aHem ù
- = ln Í a a ú + ln Í a ú 4.54
ÍÎ (1 - XUsp ) X Ilm úû ÍÎ g Mtg Ilm úû
a
RT
DG È X 6a ù Èg Hem
6a
ù
and: - = ln Í Hem4a ú
+ ln Í 6a ú - ln fO 2 4.55
RT Î X Mt û Î g Mt û
The a parameter is related to the number of sites involved in the exchange; Spencer and Lindsley as-
sumed a was 2 for ilmenite and 1 for titanomagnetite. The excess free energy was expressed in the
usual way for an asymmetric solution (equation 4.16):
( )
Gex = WG1 X 2 + WG 2 X1 X1 X 2
for each solution series. When pressure is neglected, the free energy interaction parameter expres-
sion (equation 4.08) simplifies to:
WG = WH – TWS 4.56
Values for WH and WS were obtained from least-squares fits to experimental data. The parameters
obtained are listed in Table 4.1.
Substituting equations 4.56 and 4.16 into the free energy of solution expression (DG excess =
DGideal – DGreal), the following equation can be obtained:
AWHUsp - BWHMt - CWHIl + DWHHem + DH o

T (K ) = 4.57
AW SUsp - BW SMt - CW SIlm + DW SHem + DS o - R ln Kexch
Oxygen fugacity is determined as:
Ê12 ln(1 - Xilm ) - 4 ln(1 - XUsp ) + ˆ
Á ˜
log fO 2 = log MH + Á 1 È8 XUsp
2
(aUsp - 1)WGUsp + 4 XUsp 2
(1 - 2 XUsp )WGMt + ù˜ 2.303 4.58
Á RT Í ú˜
Ë ÍÎ12 X Ilm
2
(1 - X Ilm )WGIlm - 6 X Ilm
2
(1 - 2 X Ilm )WGHem úû¯
where:
A = 3 XUsp
2
- 4 XUsp + 1, B = 3 XUsp
2
- 2 XUsp , C = 3 X Ilm
2
- 4 X Ilm + 1, D = 3 X Ilm
2
- 2 X Ilm
(
Kexch = XUsp X Hem
2
) 2
X Mt X Ilm , DHo = 27.799 kJ/mol, DSo = 4.1920 J/K-mol
138 October 17, 2001

and MH is the magnetite-hematite buffer: log MH = 13.966 – 24634/T.
Table 4.1. Margules Parameters for Ilmenite and Titanomagnetite Solid Solutions
Usp Mag Ilm Hem
(<800 ° C) (<800° C)
WH (joules) 64835 20798 102374 36818
WS (joules) 60.296 19.652 71.095 7.7714
WG (>800° C) (joules) 0 0
o
DS Usp (joules) 4.192
o
D H Usp (joules) 27799
We have reviewed just a few of the available thermobarometers in use. These were selected to il-
lustrate the underlying principals. There are, however, many thermobarometers in use by geo-
chemists and petrologists. Some of these are listed in Table 4.2.
Example 4.4: Using the Iron-Titanium Oxide Geothermometer

An electron microprobe analysis of oxide phases in an andesite reveals that there is 68 mole per-
cent of ulvospinel in an ulvospinel–magnetite phase and 93.3% of ilmenite in an illmenite-hematite
phase. Calculate the temperature and ƒO 2 at which these phases equilibrated.
Answer: We can use equations 4.57 XUsp XIlm
and 4.58 to answer this question. The
data in Table 4.2 are relevant to the 0.68 0.933
binary asymmetric solution model for DH 27799
the system below 800° C. Above 800° DS 4.192
C, an ideal solution is assumed for the
ulvospinel-magnetite phase, so the in- R 8.314
teraction parameters for this phase go Interaction Parameters
to 0. But if we don’t know the WHU 64835 WSU 60.296
temperature, how do we know which WHM 20798 WSM 1 9 . 6 5 2
equation to use? We begin by WHI 102374 WSI 71.095
computing temperature using the
parameters for less than 800° C. If the WHH 3 6 8 1 8 WSH 7.7714
temperature computed in this way is A -0.3328
greater than 800° C (1073 K), we set the B 0.0272
W H and W S for ulvospinel and
magnetite to 0 and recompute. C - 0 .12053
Once we have temperature, we can D 0 . 7 45467
compute the WG terms using the rela- K 0.010958
tionship W G = WH –TW S , bearing in
mind that WGus p = WG Mt = 0 if the tem- T= (A*WHU-B*WHM-C*WHI+D*WHH+ D H)
perature is greater than 800° C. With (A*WSU-B*WSM–C*WSI+D*WSH+ D S-R*ln(K)
these values in hand, we can use equa-
T (<800) 1281 K 1 0 0 8 °C
tion 4.58 to calculate the ƒ O 2 . Our
T (>800) 1205 K 9 3 2 °C
spreadsheet is shown on the right.
These data we taken from one of WG=WH-T*WS
Spencer and Lindsley’s (1981) WGU - 7 8 2 9 . 5 2 WGI 16695.29
experiments, performed at 938° C and WGM - 2 8 8 5 . 2 1 WGH 2 7 4 5 2 . 4 5
log ƒO 2 = -12.76. Our calculations are MH -6.47
in good agreement with the LogƒO2 (<800) -12.58
experimental observation.
LogƒO2 (>800) -12.69
139 October 17, 2001

Table 4.2. Commonly Used Thermobarometers
Reaction Type Reference
Garnet=Biotite Fe-Mg exchange (Ferry and Spear, 1978)
(Fe,Mg)3Al2Si3O12 ® K(Mg,Fe)AlSi3O10(OH)2
Plagioclase = Garnet + Kyanite + Quartz displaced equilibria (Ghent, 1976; Koziol
3(Ca,Na)Al2Si2O8 ® (Fe,Ca)3Al2Si3O12 + 2Al2SiO5 + SiO2 and Newton, 1988)
Garnet + Quartz = Plagioclase + Wollastonite displaced equilibria (Gasparik, 1984b)
(Fe,Ca)3Al2Si3O12 + SiO2 ® (Ca,Na)Al2Si2O8 + 2CaSiO3
Calcite = Dolomite solvus equilibria Goldsmith and Newton (1978)
CaCO3 ® (Ca,Mg)CO3
Calcite = Aragonite univariant (Johannes and Puhan, 1971)
CaCO3 ® CaCO3
Ilmenite + Al2SiO5 = Garnet + Rutile + Quartz displaced equilibria (Bohlen et al., 1983)
3FeTiO3 + Al2SiO5 ® 3TiO2 + (Fe,Ca)3Al2Si3O12 + SiO2
Hercynite + Quartz = Garnet + Sillimanite displaced equilibria (Bohlen et al., 1986)
FeAl2O4 + 5SiO2 ® Fe3Al2Si3O12 + Al2SiO5
4.6 Thermodynamic Models of Magmas

Silicate liquids have played an extremely important role in the development of the Earth, as well
as other bodies in the solar system. As we shall see, the Earth’s crust formed as melts from the
mantle rose to the surface and cooled. Thus an understanding of igneous processes is an essential
part of earth science. Until the last decade or two, the primary approaches to igneous petrology
were observational and experimental. Results of melting experiments in the laboratory were used
to interpret observations on igneous rocks. This approach has proved highly successful and is
responsible for most of our understanding igneous processes. However, such an approach has
inherent limitations: virtually every magma is unique in its composition and crystallization history.
Yet the experimental database is limited: it is not practical to subject each igneous rock to melting
experiments in the laboratory. Realizing this, igneous petrologists and geochemists turned to
thermodynamic models of silicate melts as a tool to interpret their evolution. With a proper ‘model’
of the interaction of various components in silicate melts and adequate thermodynamic data, it
should be possible to predict the equilibrium state of any magma * under any given set of
conditions. The obstacles in developing proper thermodynamic models of silicate liquids,
however, have been formidable. Because they are stable only at high temperatures, obtaining basic
thermodynamic data on silicate liquids is difficult. Furthermore, silicate liquids are very complex
solutions, with 8 or more elements present in high enough concentrations to affect the properties of
the solution. Nevertheless, sufficient progress has been made on these problems that
thermodynamics is now an important tool of igneous petrology.
4.6.1 Structure of Silicate Melts
As was the case for silicate solids and electrolyte solutions, application of thermodynamics to sili-
cate liquids requires some understanding of the interactions that occur on the atomic level. Thus
we will once again have to consider the microscopic viewpoint before developing a useful
thermodynamic approach. In this section, we briefly consider the nature of silicate melts on the
atomic level.
* A magma consists not only of a liquid, but any suspended gas or crystalline phases as well.
140 October 17, 2001

Most, though not all, of our knowledge of the a
structure has come from studies of glasses rather than
melts. While the thermodynamic properties of silicate
liquids and their respective glasses differ, other
studies have confirmed the general structural simi- bridging
larities of glasses and liquids. Spectral studies of
non-bridging
glasses, which in some respects can be viewed as oxygens oxygens
supercooled liquids, have revealed that silicate liquids
have structures rather similar to those of silicate
solids. In fact, the principal difference between
silicate liquids and solids is the absence of long-range
ordering in the former; short range ordering is b
similar. As in silicate minerals, the primary structural
element of silicate liquids is the silicon tetrahedron
(see Fig. 1.11), consisting of a silicon atom
surrounded by four oxygens. As in silicate minerals,
tetrahedra may be linked by a shared oxygen, called
a bridging oxygen; not surprisingly, unshared
oxygens are termed non-bridging (Figure 4.21a).
Unlinked silica tetrahedra, that is, those with no Monomer Dimer
bridging oxygens, are termed monomers, SiO 44 -
(Figure 4.21b). Two tetrahedra linked by a single Figure 4.21. Silicate structures. a: Short
oxygen are termed dimers and have the formula range silicate structures in melts resemble
Si2O 6- those in solids. Individual tetrahedra may
7 . Tetrahedra may also be linked by two
oxygens to form infinite chains; these have a chemical be linked by bridging oxygens and linked
formula of SiO 32 - . In silicates such as quartz and to 2 silicon atoms. b. Unit in silicate melts
feldspar, the tetrahedra are all linked into a may include monomers, with no bridging
framework, and all oxygens are shared. All these oxygens, and dimers, where only 1 of 4
oxygens in each tetrahedra are ‘bridging’.
Network-forming ion
Network-modifying
Bridging oxygen ion
Nonbridging oxygen ion
(a) (b)
Figure 4.22. (a) structure of pure silica glass and (b) a silica-rich glass
with additional component ions.
141 October 17, 2001

structural can be present in silicate glasses.
The degree to which the silica tetrahedra are linked, or polymerized, in silicate liquids affects
chemical and physical properties. The degree of polymerization in turn depends on other cations
present. These may be divided into two groups, network formers and network modifiers. Relatively
small, highly charged cations such Al3+ and Fe3+ (more rarely, Ti4+, P 3+, B 3+ as well) often substitute
for silicon in tetrahedral sites and, along with Si, are termed network formers. The other common
cations of natural silicate liquids, Ca 2+, Mg 2+, K + , Na + , and H + , are network modifiers. These ions
cannot substitute for silicon in tetrahedra and their positive charges can only be balanced by non-
bridging oxygens. Addition of these ions disrupts the linkages between silica tetrahedra. Thus as
silicate melts become richer in these network modifiers they become progressively depolymerized.
This is illustrated in Figure 4.22, which compares the structure of pure silica glass (liquid) and a
silica-rich glass (liquid). Melt structure in turn affects the physiochemical properties of the melt.
For example, SiO2–rich melts tend to have low densities and high viscosities. As ions such as MgO
or CaO are added to the melt, viscosity decreases and density increases as the polymer structure is
disrupted.
4.6.2 Magma Solution Models
Advances on several fronts have moved thermodynamic modeling of magmas from an academic
curiosity to useful petrological tool. First, spectroscopic (mainly Raman and infrared spectroscopy,
both of which are sensitive to atomic and molecular vibrations) studies are revealing the structure
of silicate melts, which provides the theoretical basis for thermodynamic models. Second, more
sophisticated thermodynamic models more accurately reflect interactions in silicate melts. Third,
the thermodynamic database has become more complete and more accurate. Finally, the wide
accessibility and power of computers and appropriate programs have made the extensive matrix
calculations involved in these models possible.
Several factors complicate the task of thermodynamic modeling of magmas. First, magmas are
solutions of many components (typically 8 or more). Second, the solids crystallizing from magmas
are themselves solutions. Third, magmas crystallize over a substantial temperature range (as much
as 400-500° C, more in exceptional cases). Furthermore, crystallization may occur over a range of
pressures as a magma ascends through the Earth, and crystallization may be accompanied by
melting and assimilation of the surrounding ‘country’ rock. Despite these complications several
models that are sufficiently accurate to be useful to petrologists have been published, most notably
those of Ghiroso (Ghiorso et al., 1983; Ghiroso and Sack, 1995) and Nielsen and Dungan (1983). The
goal of these models is to describe the phases and their proportions crystallizing from a cooling
magma, and the resulting evolution of liquid composition. In the section below, we briefly
consider the model of Ghiorso.
4.6.2.1 The Regular Solution Model of Ghiorso and Others
Ghiorso (Ghiorso et al., 1983; Ghiorso, 1987; Ghiorso and Sack, 1995) noted that silicate liquids
have substantial compositional regions in which immiscibility occurs and therefore argued that sim-
plest model that might be able to describe them is the regular solution model. As we saw earlier in
the chapter, regular solution models attempt to describe excess functions with interaction, or Mar-
gules, parameters. The Margules equation for excess Gibbs Free Energy for many components is:
1
Gex = Â Â
2 i j , j πi
Xi X jWGi, j 4.59*
and the Gibbs Free Energy is:
* The 1 term arises because the sum contains both X X W ij ji ij ij

=W G terms and X jX i W G and W G
ji
=W G .
2 i j G
142 October 17, 2001

1
G = Â Xi m oi + RT Â Xi ln Xi + Â Â Xi X jWGi, j 4.60†
i i 2 i j , j πi
The chemical potentials of individual components are:
1
m i = m oi = + RT ln Xi + ÂXW j
i, j
G - Â Â X X W j, k
2 j, j π k k , k π j j k G
4.61
j , j πi
and the activity coefficients are:
1
RT ln li = ÂXW j
i, j
G - Â Â
2 j, j π k k , k π j
X j X k WGj , k 4.62
j , j πi
Having chosen a general form for the solution model, the next step is to select the components.
For practical reasons, Ghiorso et al. (1983) placed all components on an 8-oxygen basis. Ghiorso
and Sack (1995) chose oxides that are more familiar to igneous petrologist as components of the
liquid phase: SiO2, TiO2, Al2O 3, Fe 2O 3, MgCr 2O 4, Fe2SiO4, Mg 2SiO4, CaSiO3, Na2SiO3, KAlSiO4,
Ca3(PO4)2, and H 2O. For components of solid phases, they chose pure end-member phase
components (e.g., MgSiO3 in orthopyroxene).
The next task is to find values for the interaction parameters. These can be calculated from solid-
liquid equilibria experiments. The principle involved is an extension of that which we used in
constructing phase diagrams: when a solid and liquid are in equilibrium, the chemical potential of
each component in each phase must be equal. Since thermodynamic properties of the solids
involved are available (determined using standard thermodynamics techniques), the
thermodynamic properties of the co-existing liquid may be calculated.
The reaction of a solid phase, j, with the melt can be described with a set of p reactions of the
form:
j p ® Â n p ,i ci 4.63
i
where jp is the pth end member component of phase j, ci refers to the formula for the ith component
in the liquid and np,k refers to the stoichiometric coefficient of this component. Thus for reaction of
olivine with the liquid, we have two versions of 4.63:
(Mg2 SiO4 )Ol ® 2MgOl + SiO2-l 4.63a
and (Fe2 SiO4 )Ol ® 2FeOl + SiO2-l 4.63b
We can express the Gibbs Free Energy change for each of these reaction as:
DGr = DGjop + RT Â n p ,i ln ail - RT ln aj p 4.64
i
l
where a i is –the activity of the oxide component in the liquid and j p refers to phase component p in
phase j. DG r is, of course, 0 at equilibrium. For example, for reaction 4.63a above, we have:
DGr = 0 = DG Fo
o
+ RT 2 ln aMgO
l
+ ln a SiO
l
2 [
+ RT ln a Fo ]
where the subscript Fo refers to the forsterite (Mg2SiO4) component in olivine and the superscript
l refers to the liquid phase. Expanding the liquid activity term, we have:
0 = DGjop + RT Â n p ,i ln X il + RT Â n p ,i ln lli - RT ln aj p 4.65
i i
Substituting 4.62 for the activity coefficient term in 4.65 and rearranging to place the “knowns” on
the left-hand side, we have:
† For clarity, we have simplified Ghioso's equation by neglecting H O, which they treated separately.
2
143 October 17, 2001

Table 4.3. Interaction Parameters for the Ghiorso Regular Solution Model
SiO2 TiO 2 Al2 O 3 Fe 2 O 3 MgCr2 O 4 Fe 2 SiO4 Mg2 SiO4 CaSiO3 Na 2SiO3 KAlSiO4 Ca3 (PO4 )
TiO 2 26267
Al2 O 3 --39120 -29450
Fe 2 O 3 8110 -84757 -17089
MgCr2 O 4 27886 -72303 -31770 21606
Fe 2 SiO4 23661 5209 -30509 -179065 -82972
Mg2 SiO4 3421 -4178 -32880 -71519 46049 -37257
CaSiO3 -864 -35373 -57918 12077 30705 -12971 -31732
Na 2SiO3 -99039 -15416 -130785 -149662 113646 -90534 -41877 -13247
KAlSiO4 -33922 -48095 -25859 57556 75709 23649 22323 -17111 6523
Ca3 (PO4 )2 613892 25939 52221 -4214 5342 87410 -23209 37070 15572
H2 O 30967 81879 -16098 31406 28874 35634 20375 --96938 10374 43451
Values are in kJ/mol. From Ghiorso and Sack (1995).
1
-DGjop + RT ln aj p - RT Â n p ,i ln Xil = Â n p .i Â X jWGi, j - Â Â X X W k, j 4.66
i i j 2 j k,k π j k j G
The quantities on the left-hand side of the equation are terms that can be calculated from the compo-
sitions of coexisting solids and liquids and solution models of the solids. The right hand side con-
tains the unknowns. One statement of equation 4.66 can be written for each component in each
solid phase at a given temperature and pressure. With enough experiments, values for the interac-
tion parameters can be extracted from the phase relations. Ghiorso et al. (1983) and Ghiorso and
Sack (1995) used a statistical technique called least squares† to determine the interaction parameters
from a large number of published experimental data. Ghiorso and Sack (1995) also noted that the
absence of a phase in an experiment provides thermodynamic information about that phase, i.e.,
that its free energy must be higher than that of the phases that are present. Their approach made
use of this information as well, though discussion of that aspect of their method would take us too
far afield. The interaction parameters they determined are listed in Table 4.3.
With values for the interaction parameters, the model can then be used to predict the assemblage
of solids, their compositions, and the liquid composition that will be present in the system as a func-
tion of temperature and pressure. The equilibrium condition for a magma, as for any other system,
is the condition where Gibbs Free Energy is at a minimum. Thus the problem becomes finding a
composition for the liquid and coexisting solids that minimizes G at a particular temperature and
pressure. In other words, we want to find values of G l and G f1, G f2, ... G fn such that G sys is minimal
where:
Gsys = Gl + Â Gj 4.67
j
Inherent in the problem is finding which solids will be in equilibrium with the liquid for a given
bulk system composition at specified temperature and pressure. In Ghiorso’s approach, an initial
guess is made of the state of the system. This is done by taking the liquid composition as equal to
the system composition and estimating what phases are likely to be in equilibrium with this liquid.
† ‘Least squares’ is a numerical technique that attempts to minimize the square of the difference between
calculated a nd observed value of some parameter. The square is taken to give greater weight to large deviations.
Thus least squares techniques yield results where there are relatively few large deviations between the calculated
and observed value of the parameter of interest. We discuss this technique further in Chapter 8.
144 October 17, 2001

Then G is expanded as a 3-term Taylor Se-
ries‡ about that initial point, N´, where N´
FeO
is the composite vector containing the
mole fractions describing the compositions
of all phases in the system. The second
term in the expansion is the first derivative
of G with respect to ni , the moles of
JJ
component i, which is simply the chemical JJJJ JJJ
J
potential. A minimum of G occurs where J
J
the first derivative is 0. Thus the second J J J
term, the chemical potentials, is set to 0 J
and solution sought by successive it-

erations. After each iteration N´ is reset to
the composition found in the most recent
iteration. This approach clearly involves
90 Tridymite
Plagioclase
Percent Crystallized
70 Na2O MgO
Liquid + K2O
50 Augite
Orthopyroxene Figure 4.23. Calculated compositional evolution of
30 magmas from Thingmuli, Iceland using the regular
10
Pigeonite solution model of Ghiorso compared with actual
Spinel analyses plotted on a FeO–MgO–Na2O+K2O (“AFM”)
ternary diagram. Solid line is for 0.1 MPa pressure,
SiO2 ¥ 13 Shergotty dashed line is for 200 MPa pressure. Both assume ƒO 2
WSeight Percent
20 on the FMQ buffer. From Ghiorso and Carmichael

GGG G G G G G
FeO CaO (1985).
Al2O3
10 BBB J J J
JJJ B B
JJ
TiO2 B O O
B B MgO repetitive matrix calculations and would not be
Na2O O
O
OO
F F FF
Fe
FF 3
2 O practical without a computer, but they can easily be
FOH
OF
HH H HH HH performed on the current generation of computers.
1000 1100 1200 Figure 4.23 and 4.24 illustrate application of the
Temperature, ∞C Ghiorso model. Figure 4.23 compares the calculated
compositional evolution during fractional crystalliza-
Figure 4.24. (a) Cumulative mineral pro- tion of magmas from Thingmuli, eastern Iceland,
portions and liquid composition predicted with actual analyses. During fractional crystalliza-
by the “MELTS” model of Ghiorso and
tion, solid phases are isolated from the liquid and
Sack (1975) compared with experimental
therefore do not continue to equilibrate with the
determined values of Stolper and
magma after they have formed. Figure 4.24a shows
McSween (1979) for the Shergotty
the calculated cumulative proportions of spinel,
meteorite. Shergotty is one of the SNC
olivine, pigeonite, plagioclase, and clinopyroxene
class of meteorites believed to have come
that crystallized from the meteorite Shergotty as well
from Mars (see Chapter 10). From
as the composition of the remaining liquid as a
Ghiorso et al. (1994).
‡ A Taylor series expansion of a function ƒ(z) in the vicinity of some point z = a has the form:
(z - a) (z - a)2
f ( z) = f ( a ) + '(
f a ) + f ''( a ) + K
1! 2!
where ƒ’ and ƒ” are the first and second derivatives of ƒ with respect to z.
145 October 17, 2001

function of temperature. This calculation assumes that all solids remain in equilibrium with the
melt as it continues to cool, a process referred to as equilibrium crystallization.
The latest version of the Ghiorso model, embodied in a computer program called “MELTS” that
runs of a variety of UNIX workstations (and as JAVA applets on personal computers), has recently
been made available to the scientific community through the World Wide Web (WWW) (URL:
“http://www.geology.washington.edu/~ghiorso/MeltsWWW/Melts.html”). The MELTS program
also considers the effects of pressure, allowing for modeling of melting and crystallization up to 1
GPa.
4.6.2.2 The Two Lattice Model of Nielsen and Dungan
Nielsen and Dungan (1983) approached the problem of modeling the behavior of silicate melts
quite differently. Their model is based upon a structural model of silicate melts developed by Bot-
tinga and Weill (1972) to predict melt viscosities. The Bottinga and Weill model is based upon the
observation that components in silicate melts can be divided into network formers and network
modifiers. Nielsen and Dungan modeled silicate melts as a mixture of two solutions, or quasi-
lattices: one lattice, or solution, consists of the network-forming components and the other of the
network-modifying components, hence the name “two lattice model”. Within each solution,
components are assumed to mix ideally. The two lattices, or solutions are assumed to exist entirely
independently of each other, so all effects related to the mixing of the two lattices are ignored.
Having established this framework, the next task is to decide upon the components and assign
them to either the network forming solution of the network modifying one. SiO2 is, of course,
assigned exclusively to the network-forming solution, and MgO, FeO, CaO, TiO2, MnO, and CrO1.5
are assigned to the network-modifying solution. However, experimental evidences suggests that Al
and Fe3+ can be both network modifiers and network formers, and that Na and K can complex with
network-forming Al. Thus Nielsen and Dungan combined Na and K with aluminum to form
components NaAlO2 and KAlO2. Any excess Al is assigned to the network modifying solution as
AlO1.5 . Any excess Na and K are combined with Fe 3+ to form components NaFeO 2 and NaFeO 2;
Fe 3+ is otherwise assigned to the network-modifying solution as FeO1.5 . Activities are then
calculated assuming ideality, e.g.:
l
X Na
a l
=X NF
= 4.68
l
+ X Kl + X Sil
NaAlO 2 NaAlO 2
X Na
l
X MgO
a l
=X NM
= 4.69
MgO MgO
ÂX l
NM
where NF and NM refer to network-formers and network-modifiers respectively. Equilibrium con-

stants for reaction between components in the melt and in various possible minerals were extracted
from experimental data using linear regression. Their expressions for the equilibrium constant have
the form:
a
ln K = +b 4.70
T
with parameters a and b being the slope and intercept of the regression results. From Chapter 3, we
can identify parameter a with –DH r /R and parameter b with DSr /R. Nielsen and Dungan did not
consider the effects of pressure, and their model is restricted to mineral-melt equilibria at or close to
atmospheric pressure.
The Nielsen and Dungan Two Lattice Model is substantially simpler than Ghiorso and has a
number of theoretical weaknesses. The assumption of ideality is one such weakness, given that
silicate liquids can exhibit immiscibility. A more seriously weakness is its inability to predict true
activities. In the Nielsen and Dungan model, for example, adding MgO to the liquid does not
146 October 17, 2001

decrease the activity of SiO2, which is certainly incorrect. Despite these theoretical shortcomings,
the model has been widely applied with considerable success, with many, though not all magma
compositions (the Ghiroso model also does not work for all compositions). Furthermore, their
model predicts trace element behavior as well, as we shall see in Chapter 7. A computer program,
MIXFRAC, for the Nielsen and Dungan Model that runs on personal computers are available on
the Geochemical Earth Reference Model (GERM) Web Site (http://earthref.org/GERM/).
4.7 Reprise: Thermodynamics of Electrolyte Solutions
We discussed the nature of electrolyte solutions and introduced one approach to dealing with
their non-ideality, namely the Debye-Hückel activity coefficients, in Section 3.7. We also noted a
number of theoretical weaknesses in the Debye-Hückel approach and that this approach is
restricted to fairly dilute solutions (ionic strengths less than 0.1 M). In this section we will return to
the problem of electrolyte solutions and examine the causes of non-ideal behavior in high ionic
strength solutions in more detail. Before doing so, however, we need to introduce a new variation
on our now-familiar thermodynamic parameters, namely mean ionic quantities.
4.7.1 Mean Ionic Quantities
Consider an aqueous NaCl solution. In Chapter 3 we saw that the thermodynamic properties of
a salt are related to those of its component ions by:
YAB ∫ nA YA + nBYB (3.73)
So, for example, the chemical potential of NaCl in solution is:
m NaCl = m Na+ + m Cl – m NaCl = m Na + + m Cl –
which we can express as:
(
m NaCl = m oNa + + m oCl - + RT ln aNa + + ln aCl - ) 4.71
or: ( ) (
m NaCl = m oNa + + m oCl - + RT ln m Na + + ln mCl - + RT ln g Na + + ln g Cl - )
Though we can certainly determine the concentrations of Na and Cl in solution, how do we
independently determine their activity coefficients? Since we cannot create a pure Na + solution or
a pure Cl– one, we cannot say what part of the non-ideality of NaCl solution is due to Na+ and what
part is due to Cl–. The practical solution then is to assign all non-ideality equally to both ions. This
leads to the concept of the mean ion activity coefficient:
g ± = (g Na + g Cl - )1/ 2 4.72
Thus the mean activity coefficient of a salt is the multiplicative mean of the activity coefficients of its
component ions. Equation 4.71 then becomes:
(
m NaCl = m oNa + + m oCl - + RT ln m Na + + ln mCl - + ln g ±2 )
Equation 4.72 is valid for 1:1 salts (i.e., 1 cation for each anion). A general expression for the mean
activity coefficient of a salt of composition An+ Bn– is:
+ -
g ± = (g n+ g n- )1/n 4.73
where n is the sum of the component positive and negative ions:
n = n+ + n– 4.74
Mean activity coefficients have the advantage that they are readily measurable (through electro-
chemical means or solubility, for example). Given a well-behaved salt, such as KCl, where the rela-
147 October 17, 2001

tionship g– =g+ appears to hold, it is then Example 4.5: Calculating Single Ion Activity
possible to determine single ion activity Coefficients from Mean Ionic Activity
coefficients. For example, we can obtain Coefficients
gNa+ in our NaCl solution by first
The measured mean ionic activity coefficient of KCl in
determining gCl – in KCl*:
a solution of 1.0 m ionic strength is 0.604; that of CaCl2
gCl – = gK – = g± KCl in a solution of the same ionic strength is 0.449. What
then determining the mean ion activity is the activity coefficient of Ca2+? Assume gCl– = gK+ .
coefficient of NaCl experimentally in a Answer: We begin by noting that gCl– = gK+ = g±KCl
solution of the same ionic strength and and therefore that gCl– = 0.604. According to equ. 4.73,
calculating gNa+ as: the mean ion activity coefficient for CaCl2 is related to
g ±2 NaCl the single ion activity coefficients as:
g Na + =
g Cl - ( )
1/ 3
g ± CaCl 2 = g Ca + g Cl2 -
We can extend the concept of mean Solving this for gCa 2+ we have:
ionic quantities to other thermodynamic
g 3± CaCl2 0.449 3
variables as well. The mean ionic po- g Ca + = = = 0.248
tential , µ±, is defined as: g Cl2 - 0.604 2
n+ m n– m
m = + – 4.75
Thus the mean ionic potential is simply the arithmetic mean of the potential of the individual ions
weighted by their stoichiometric coefficients. We could also express the mean ionic potential as:
m ± = m o± +
(
RT ln an+ + ln an-
+ -
) 4.76
n
Rearranging once more, we obtain:
( )
+ - 1/n
m ± = m o± + RT ln an+ an- 4.77
Comparing this relationship with equation 4.71, we define a mean ionic activity such that:
( )
+ - 1/n
a± = an+ an- 4.78
We can also define mean ionic molalities such that a± = g± m±. Substituting a + = g + m + , and a– =
g–m –, we find the mean ionic molality is then:
( )
+ - 1/n
m ± = mn+ mn- 4.79
Mass balance requires that:
m+ = n+ m and m – = n –m 4.80
Substituting this into equation 4.79, we see that:
( )
+ - 1/n
m ± = m nn+ nn- 4.81
Let’s returning to our NaCl example. Dissociation is essentially complete and n+ and n– are
unity, so that:
mNa+ = mNaCl and mCl – = mNaCl
Since n = 2: m ± NaCl = m NaCl
2
= m NaCl
* The use of KCl as a reference for determining mean ion activity coefficients is based on the observation that K +
and Cl– have about the same effective radius and ion mobility and is known as the MacInnes Convention. Like
that of Debye-Hückel, h owever, this approach breaks down at high ionic strength.
148 October 17, 2001

Mean ionic molality is simply equal to molality for a 0.6
completely dissociated salt consisting of monovalent ions (b)
such as NaCl. 0.4
The mean ionic activity coefficient, or the stoichiometric
aNaCl
activity coefficient as it’s sometimes referred to, of NaCl
would be the square root of the product of the
0.2
component activity coefficients according to equation
4.78, as would the mean ionic activity. The individual 0.2 0.4 0.6 0.8 1.0
ion activities can be measured in a number of ways. 2
Therefore, the above relationships allow calculation of mNaCl
0.6
the mean ionic activity coefficient from measurable (a)
quantities.
0.4
aNaCl
For strong electrolytes, i.e., salts that completely
dissociate, it can also be shown that mean activity
coefficient and mean activity of the salt are related to its
0.2
activity coefficient and activity by:
g = g n± 4.82 0.1 0.2 0.3 0.4 0.5 0.6
and a = an± mNaCl
4.83
We can modify the Debye-Hückel equations to obtain Figure 4.25. (a) Relationship between
activity and molality of NaCl in
mean ion activity coefficients as follows:
aqueous solution. The activity is very
Debye-Hückel Extended Law: low and the “Henry’s Law Slope” is
almost 0 at low concentrations. (b)
Relationship between activity and the
square of molality of NaCl in aqueous
Example 4.6: Mean Ionic Parameters for solution.
a fully dissociated electrolyte
If the molality of a CaCl2 solution is 0.3 M - Az+ z- I
and the activity coefficients of Ca2+ and Cl– are log10 g ± = 4.84
0.5 and 0.7 respectively, calculate the activity 1 + Bå I
and mean ionic molality of CaCl2 in the solu- Limiting Law:
tion. Assume that CaCl2 fully dissociates. log10 g ± = - Az+ z- I 4.85
Answer: For CaCl2, n + = 1, n– = 2, and n = 3. where å is taken as the sum of the radii of the
So we can use equation 3.58 to calculate mean anion and cation, i.e., å = å+ + å– .
ionic molality:
1/3 4.7.1.1 Relationship between Activity and
m ± CaCl2 = mCaCl2 11 22 = 41/4mCaCl2 Molality of a Salt
Substituting 0.3 for m, we find that m± = 0.4762 Let’s consider the relationship between ac-
M. tivity and molality of a salt in an electrolyte
We then use equation 4.81 to calculate the solution such as a NaCl solution. Figure 4.27a
mean ionic activity coefficient: illustrates this relationship. What we immedi-
+ – 1/n 1/3
n n ately notice is that the slope in the Henry’s
g±= g+ g–
( ) = (0.510.72) = 0.625 Law region is essentially zero, which is not at
The mean ionic activity is then: all what we expect for Henry’s Law behavior.
a ± = g±m ± = 0.625 ¥ 0.4762 = 0.2980 It can easily be shown that the
relationship in 4.25a is a simple consequence
and the activity of CaCl2 is: of the dissociation of the NaCl into Na + and
n 3 Cl - ions. From 3.75 we have:
acaCl = a± = g±m± = 0.2980 = 0.0263 M
2
149 October 17, 2001

m NaCl = m Na + + m Cl - 4.86
aq aq
Substituting this into equation 3.45, we obtain:

m NaCl = m oNa + + m oCl - + RT ln aNa + + RT ln aCl -
In the reference state of infinitely dilute solution, mi = ai, so that:
m NaCl = m oNa + + m oCl - + RT ln m Na + + RT ln mCl - 4.87
Furthermore, charge balance requires that:
mNa+ = mCl- = mNaCl 4.88
Substituting 4.88 into 4.87 and rearranging:
m NaCl = m oNa + + m oCl - + 2 RT ln m NaCl = m oNa + + m oCl - + RT ln m NaCl
2
4.89
Comparing this equation with equation 3.45, we see that
aNaCl µ m NaCl
2
When we plot activity versus the square of molality, we obtain a linear relationship (Fig. 4.25b).
Generalizing this result for dissociation of a substance into a positive ion A and negative ion B,
such as:
+ -
An+ Bn- ® n+ A z+ + n– Bz– A + B - ® n + An + n + Bn
n n
the relationship between activity of a salt and its molality is:
a AB µ mnAB 4.90
For example, n is 3 for CaCl2, 4 for FeCl3, etc. 1.0
Now let’s see what happens if we substi-
tute the mean ion activity for activity. Since: 0.9
n
a = a AB
± 0.8
We have: a = mnAB or
n
a± µ m AB 25∞C
±
This is the relationship that we observed in 0.7
Figure 4.25, so we see that the mean ionic 100∞C
activity accounts for the effects of dissocia- 0.6
tion. g±
0.5
4.7.2 Activities in High Ionic Strength
Solutions 0.4 200∞C
As we noted in Chapter 3, there are two
phenomena that are not considered in the 0.3
Debye-Hückel approach that result in inaccu-
racies in the Debye-Hückel equation at ionic 0.2 300∞C
strengths above about 0.1 m. This is illus-
trated in Figure 4.26, which shows the exper- 0.1
imentally determined mean ion activity coef-
ficient for NaCl as a function of ionic strength 0
0.1 0.10 1.0 10.0
and temperature. At low temperatures, the Ionic Strength, mol/kg
activity begins to increase about ionic
Figure 4.26. Observed mean ion activity coeffi-
strengths of 1 m, whereas Debye-Hückel
cient, g± , of NaCl as a function of ionic strength and
predicts continual decrease. The activities of
many electrolytes eventually exceed 1 at high temperature (solid lines; data from Helgeson, 1981)
compared with value predicted by the Debye-
concentrations. The difference between the
Hückel Law (computed as (gNa + gCl– )1/2 ).
150 October 17, 2001

observed activity coefficients and those predicted by the Debye-Hückel equation are due to the
effects of ion association and solvation. Debye and Hückel explicitly assumed complete
dissociation, i.e., no ion associations, and while their treatment included in a general way the
dielectric properties of water, it neglected the effects of solvation. As we noted in Chapter 3, the
effects of both ion association and solvation become increasingly important with increasing ionic
strength. It should be no surprise then that the Debye-Hückel treatment breaks down at high ionic
strength. Here we will consider these effects in greater detail.
4.7.2.1 Correction for the Concentration of Water
At low and moderate ionic strength, we can assume that the mole fraction of water in solution is
1 (1 k/l). For example in seawater, with an ionic strength of 0.7, the mole fraction of water about
0.99. Generally, activity coefficients and equilibrium constants are not known within 1%, so the
error introduced by this assumption is still small compared to other errors. In higher ionic
strengths, however, this assumption is increasingly invalid (for example, at a molality of 3, the mole
fraction of water has decreased below 0.95), and this must be taken into account. A convenient
way to do this is to incorporate it into the activity coefficient. The corrected activity coefficient is:
g
g corr =
(1 + 0.018Â m )
4.91
i i
4.7.2.2 Effects of Solvation

Water molecules bound to ions in solvation shells have lost their independent translational
motion and move with the ion as a single entity. These water molecules are effectively unavailable
for reaction, hence solvation has the effect of reducing the activity of water, which increases the
apparent concentration, or activity, of the solutes. In addition to solvation, i.e., the direct
association of some water molecules with the ion, the charge of the ion causes collapse of the water
structure beyond the solvation shell.
For a solution consisting of a single salt, Robinson and Stokes (1959) proposed the contribution
of solvation to the mean ion activity coefficient could be expressed as:
h
log g ±solv = - log aw - log(1 - 0.018 hm ) 4.92
n
solv
where g± is the solvation contribution to the mean ion activity coefficient, h is the number of
moles of water molecules bound to each mole of salt, aw is the activity of water, m is the
concentration of the salt in solution, and n is a defined in equation 4.74 (i.e., total moles of ions
produced upon dissolution of a mole of salt). Table 4.3 listed estimated values for the solvation
number, i.e., number of water molecules in the solvation shell of each ion. From these, the value
of h for equ. 4.91 can be calculated. The activity of water can be adequately estimated as:
aw = 1 – 0.04m
Table 4.3. Ion Solvation Numbers
Figure 4.27 illustrates the effect of solvation on the
Species h Species h activity coefficient. As may be seen, solvation sub-
Li+ 2.3 OH – 7.6 stantially affects the activity coefficient at ionic
Na + 3.3 F– 6.7 strengths above about 0.5 m.
K+ 2.3 Cl– 2.7
Rb + 2.3 Br_ 1.7 4.7.2.2 Effects of Ion Association
Mg 2+ 8.9 CO 32 - 14.4 An ion pair can be considered to have formed
Ca2+ 8.9 SO 24 - 10.4 when ions approach closer than some critical dis-
Cd2+ 6.3 tance, rc, where the electrostatic energy, which tends
Ba2+ 9.2 to bind them, exceeds twice the thermal energy,
151 October 17, 2001

6 which tends to move them apart. When this
happens, the ions are electrostatically bound
5 and their motions are linked. They are said
to form an ion pair. The thermal energy of
4 an ion is kT and electrostatic interaction
energy is:
q1q2
3 Uelectro. = 4.93
g± 4per
The ratio of these two energies when the dis-
2
t ion only tance is less than the critical one is then:
solva c
+ e lec t rostati Uelectro
=
z1z2 e 2
>2 4.94
1 solvation Utherm 4pe0er rkT
electrostatic only We can use this equation to solve for the
0 critical distance rc:
0 1 2 3 4 5 6 z1z2 e 2
rc = 4.95
I, m 8pe0er rkT
Figure 4.27. Comparison of the electrostatic For two singly charged ions, the critical dis-
contribution to the mean ion activity coefficient of tance is 3.57 Å. In a 1 molar solution, the
NaCl (calculated by the Debye-Hückel Extended average separation between ions is about 12
Law), the solvation contribution (calculated from Å, so even in such a relatively concentrated
equation 4.92 assuming h = 4) and the sum of the solution, ion pairs will not form between
two. singly charged ions. Indeed, the critical
distance is smaller than the combined
Debye-Hückel radii of all pairs of singly charged ions. Thus we do not expect ion associations to
form from pairs of singly charged ions under most circumstances. In contrast, the critical distance
for ion association between a singly and a doubly charged ion is 70 Å, considerably greater than
the sum of their Debye-Hückel radii. It also exceeds the average separation of ions in a 0.01 m
solution (about 55 Å), so that even in dilution solutions, we would expect significant ion pair
formation for multiply charged ions.
As we saw earlier, all ions in solution are surrounded by a solvation shell of water molecules.
This solvation shell may or may not be disrupted when ion pair formation occurs (Fig. 4.28). If it is
not disrupted, and the two solvation shells remain intact, an outer sphere ion pair (also called an
outer sphere complex) is said to have formed. If water molecules are excluded from the space
between the ions, an inner sphere ion pair (or complex) is said to have formed.
For some purposes, ion pairs can be treated as distinct species having charge equal to the
algebraic sum of the charge of the ions involved. These can be included, for example, in
calculation of ionic strength to obtain a somewhat more accurate estimate of activities. On the
other hand, ion pairs, including neutral ones, can be highly dipolar and may behave as charge-
separated ions.
Ion associations affect activities in two ways. First, associated ions are less likely to participate in
reactions, thus reducing the activity of the ions involved. Second, ion association reduces the ionic
strength of the solution, and hence reduces the extent of electrostatic interactions among ions. This
has the effect of increasing activity. To understand the first effect, consider the case where if a cer-
tain fraction of the free ions re-associates to form ion pairs, e.g.:
n + A z + + n - B z - ®( An + Bv - )0aq
where the superscript 0 indicates neutrality and the subscript aq a dissolved aqueous species. A salt
that only partially dissociates in solution is called a weak electrolyte. Let a be the fraction of the
152 October 17, 2001

Outer Sphere Ion Pair Inner Sphere Ion Pair
+ – – +–
+
Solvation shells intact Partial disruption of solvation shells Disruption of solvation shells
Figure 4.28. In formation of ion pairs, the solvation shells may remain intact or be partially or
totally disrupted. The former results in an outer sphere ion pair, the latter results in an inner
sphere ion pair.
ions that associate to form ion pairs or complexes. The associate of this fraction of ions as ion pairs
will be thermodynamically equivalent to that fraction of the substance not dissociating to begin
with. The fraction of free ions is then 1 – a. Equation 4.80 becomes:
m+ = (1 – a)n+ m and m– = (1 – a)n–m 4.96
where m is the molality of the solute. We can rewrite equation 4.78 as:
[ ]
1/n
a± = [(g+m+ )n (g–m– )n ]1/n a± = (g +m + )n (g -m - )n
+ – + -
4.97
Substituting 4.95 into 4.96 and rearranging, we obtain:
1/n
( ) [ Ï
] [(1 - a )n m ] ¸
1/n
Ì (1 - a )n m
n+ n-
n+ n- + -
a± = g g + - ý
Ó þ
A little more rearranging and we have:
1/ n
( ) [ Ï
] (n ) (n ) ¸
1/n
Ì (1 - a )m
(n + +n - ) n+ n-
n+ n- + -
a± = g g + - ý
Ó þ
Finally, since n = n+ + n– , we obtain:
1/n
( ) (1 - a )m ÏÌÓ n + ( ) (n ) ¸
1/n n+ n-
n+ n- -
a± = g g + - ý 4.98
þ
We can recognize the last term as m ± . Since a ± = g±m ±, we see that the mean ionic activity
coefficient will be
( )
1/n
g ± = (1 - a ) g n+ g n-
+ -
4.99
for an incompletely dissociated electrolyte. Thus the mean ion activity coefficients are reduced by a
factor of 1 – a. Provided we have appropriate stability constants for the ion pairs or complexes, a
can be calculated and an appropriate correction applied.
Now consider a CaSO4 solution of which some fraction of the Ca2+ and SO 24 - ions, a, associate to
form CaSO o4 . The ionic strength of the this solution would be
I=
(1 - a )
2
(
4mCa 2+ + 4m SO 2-
4
)
153 October 17, 2001
Thus the ionic strength is reduced by a factor
of 1 - a as well.
0.8
Ion pairs and complexes need not be 0.7
neutral species (AlCl2+, for example). When
they are not, they will contribute to ionic 0.6
strength. A general expression for ionic
strength taking account of ion associations 0.5
must include charged ion pairs and
complexes: g ± 0.4
a max = 0
1È ù
I = ÍÂ (1 - a )mi zi2 + Â c n zn2 ú 4.100 0.3
2Î i n û
where ai is the fraction of each ion involved
0.2
in ion associations, and cn is the concentra- a max = 0.4
tion of each ion pair or complex and zn is its
0.1 a max= 0.2
charge. We could use this result directly in
0
the Debye-Hückel equation to make an im- 0 1 2 3 4 5
proved estimate of ionic strength, and hence
of the single ion activity coefficient.
I (apparent)
Figure 4.29. Effects of ion association on the activity
Figure 4.29 illustrates the effect of ion pair
coefficient. Mean ion activity coefficient of CaCl2 for
formation for a hypothetical CaCl2 solution in 2+
which some fraction of the ions combine to varying extents of– ion association. Fraction of Ca
ions forming CaCl was assumed to increase linearly
form ion pairs. The fraction of Ca2+ ions
with ionic strength up to a maximum value (amax) at
forming CaCl– was assumed to increase
linearly with ionic strength up to the maxi- = 5 m. Solid line shows electrostatic term (Debye-
I
Hückel) after correction for ion association, dashed
mum value shown.
line shows the combined electrostatic and solvation
If the formation of ion pairs depends on term.
the ratio of thermal to electrostatic energy,
we might expect that ion pair formation will
decrease with temperature. However, the relative permittivity of water decreases with
temperature, allowing increased electrostatic interaction between ions, and this effect dominates
over the increased thermal energy of ions. As a result, the extent of ion association increases with
temperature. Increasing pressure, on the other hand, favors dissociation of ions.
4.7.2.2 Alternative Expressions for Activity Coefficients
There have been a number of attempts to develop working equations that account for all the
effects on activity coefficients at high ionic strength. Many of these are ultimately based on the
specific ion interaction theory of Brønsted (1922). Brønsted proposed an equation of the form:
log g i = am1/ 2 + bi m 4.101
where a is a constant that is independent of the solute ions and b is the “specific ion interaction pa-
rameter” and is different for each ionic species. Guggenheim (1935) replaced the first term on the
right with a simplified form of the Debye-Hückel equation and the second term with the
summation of ion-ion interaction parameters:
- zi2 A I
log g i = + 2 Â bi , k m k 4.102
1+ I k
where bik is parameter describing the interactions between ions i and k. For natural waters with
many species, the Guggenheim equation becomes complex. Also starting from Debye-Hückel,
Truesdell and Jones (1974) proposed the following simpler equation:
154 October 17, 2001

Species Conc
Example 4.7. Activity Coefficients in a Brine
g/kg
The following concentrations were measured in a shield brine from Sudbury,
Canada at 22° C. Calculate the activity coefficients of these species using the Na 18.9
Truesdell-Jones equation. K 0.43
Answer: Our first task is to convert g/kg to molal concentrations. We do this by Ca 63.8
dividing by molecular weight. Next, we need to calculate ionic strength Mg 0.078
(equation 3.73) which we find to be 5.9 m. Calculation of activity coefficients is SO 4 0.223
then straightforward using the parameters in Tables 3.2 and 4.6. Finally, we HCO3 0.042
apply a correctionfor the decreased concentrationof water (equation 4.91). Our Cl 162.7
final spreadsheet is shown below.
‰ m z å_TJ b_TJ log (gamma) ganma gamma cor

Na 18.9 0.822 1 5 0.165 0.728 5.341 4.741
K 0.43 0.017 1 3.5 0.015 -0.238 0.579 0.514
Ca 63.8 1.595 2 5 0.165 -0.017 0.963 0.855
Mg 0.078 0.003 2 5.5 0.2 0.264 1.836 1.630
SO4 0.223 0.002 2 5 -0.04 -1.229 0.059 0.052
HCO3 0.058 0.001 1 5.4 0 -0.233 0.585 0.519
Cl 162.7 4.590 1 3.5 0.015 -0.238 0.579 0.514
m 7.030 A 0.5092
I 5.913 B 0.3283
- zi2 A I
1.25 log g i = + bi I 4.103
1.05 1 + Båi I
MgCl2 The first term on the right is identical in
0.85
ies
form to Debye-Hückel; the second term is

es
Dav
0.65 n similar to the Brønsted specific ion

0.45 ell-Jo interaction term. Truesdell and Jones
uesd
g± 0.25 Tr determined parameters å and b empiri-
cally. Table 4.5 lists these parameters for
0.05 some common ions. Figure 4.30
compares mean activity coefficient of
-0.15 calculated with the Debye-Hückel,
-0.35 Debeye-Hückel Davies, and Truesdell-Jones equations
-0.55 with the actual measured values. The
Truesdell-Jones equations fit these
-0.75 observations very well. This is not al-
0 5 10 15 ways the case, however. The fit for
I, m Na2CO 3, for example is little better than
Figure 4.30. Measured mean ionic activity for Debye-Hückel.
coefficients in MgCl 2 solution as a function of ionic Other equations include those
strength compared with values calculated from the developed by Pitzer (1979) and the
Debye-Hückel, Davies and Truesdell-Jones National Bureau of Standards. While
equations. these equations are generally more
accurate than the above, their complexity
places them beyond the scope of this book. The interested reader is referred to any of several texts
155 October 17, 2001

Table 4.5. Truesdell-Jones Parameters on geochemical thermodynamics that treat them

(Nordstrom and Munoz, 1986; Flectcher, 1993;
Ion å b Anderson and Crerar, 1993) as well as the original
Na + 4.0 0.075 literature.
K+ 3.5 0.015
Mg 2+ 5.5 0.20 References and Suggestions For Further
Ca2+ 5.0 0.165 Reading
Cl– 3.5 0.015 Anderson, G. M. and D. A. Crerar. 1993.
SO 24 - 5.0 –0.04 Thermodynamics in Geochemistry. New York:
CO 32 - 5.4 0 Oxford Univ. Press.
HCO 3 5.4 0 Bohlen, S. R., W. A. Dollase and V. J. Wall. 1986.
Calibration and application of spinel equilibria in
the system FeO-Al2O3-SiO2. J. Petrol. 27: 1143-1156.
Bohlen, S. R. and D. H. Lindsley, 1987. Thermometry and barometry of igneous and metamorphic
rocks. Ann. Rev. Earth Planet. Sci. 15:397-420.
Bohlen, S. R., V. J. Wall and A. L. Boettcher. 1983. Experimental investigations and geologica appli-
cations of equilibria in the system FeO-TiO2-Al2O3-SiO2-H2O. Am. Mineral. 68: 1049-1058.
Bottinga, Y. and D. Weill. 1972. The viscosity of magmatic silicate liquids: a model for calculation.
Am. J. Sci. 272: 438-475.
Boyd, F.R. 1973. A pyroxene geotherm. Geochim. Cosmochim. Acta, 37, 2533-2546.
Bradley, R. S. 1962. Thermodynamic calcuations on phase equilibria involving fused salts. Part II.
Solid solutions and applications to the olivines. Am. J. Sci., 260, 550-554.
Brønsted, J. N. 1922. Studies on solubility, IV. The principle of specific interaction of ions. J. Am.
Chem. Soc. 44: 877-898.
Buddington, A. F. and D. H. Lindsley. 1964. Iron-titanium oxide minerals and synthetic equivalents,
J. Petrol. 5: 310-357.
Carmichael, I. S. E. and H. P. Eugster (eds). 1987. Thermodynamic Modeling of Geological Materials:
Minerals, Fluids and Melts (Reviews in Mineralogy, vol 17). Washington: Mineral. Soc. Am.
Essene, E. 1982. Geologic thermometry and barometry. In Characterization of metamorphism through
mineral equilibria, Reviews in Mineralogy Vol. 10, J. M. Ferry. ed., pp. 153-206. Washington: Min-
eral. Soc. Am.
Essene, E. 1989. The current status of thermobarometry in metamorphic rocks. In The Evolution of
Metamorphic Belts, Vol. J. S. Daly, R. A. Cliff and B. W. D. Yardley. ed., pp. 1-44. London: Black-
well Sci. Pub.
Ferry, J. M. and F. S. Spear. 1978. Experimental calibration of the partitioning of Fe and Mg
between biotite and garnet. Contrib. Mineral. Petrol. 25: 871-893.
Fletcher, P. 1993. Chemical Thermodynamics for Earth Scientists, Essex: Logman Scientific and Tech-
nical.
Garrels, R. M., and Christ, C. L. 1965. Solutions, Minerals, and Equilibria. New York: Harper and
Row.
Gasparik, T. 1984a. Two-pyroxene thermobarometry with new experimental data in the system
CaO-MgO-Al2O 3-SiO2. Contrib. Mineral. Petrol. 87: 87-97.
Gasparik, T. 1984b. Experimental study of subsolidus phase relations and mixing properties of py-
roxene and plagioclase in the system CaO-Al2O 3-SiO2. Geochim. Cosmochim. Acta. 48: 2537-2446.
Ghent, E. D. 1976. Plagioclase-garnet-Al2O5-quartz: a potential geobarometer-geothermometer. Am.
Mineral. 61: 710-714.
Ghiorso, M. S. 1987. Modelling magmatic systems: Thermodynamic relations. in Thermodynamic
Modelling of Geological Materials: Minerals, Fluids, and Melts. eds. I. S. E. Carmichael and H. P.
Eugster. 443-465. Washington: Mineral. Soc. Am.
156 October 17, 2001

Ghiorso, M. S. and I. S. E. Carmichael. 1985. Chemical mass transfer in magmatic processes: II.
Applications in equilibrium crystallization, fractionation crystallization, and assimilation, Contrib.
Mineral. Petrol. 90: 121-141.
Ghiorso, M. S., I. S. E. Carmichael, M. L. Rivers and R. O. Sack. 1983. The Gibbs free energy of
mixing of natural silicate liquids; an expanded regular solution approximation for the calculation
of magmatic intensive variables. Contrib. Mineral. Petrol. 84: 107-145.
Ghiorso, M. S., M. M. Hirshmann and R. O. Sack. 1994. New software models thermodynamics of
magmatic systems. EOS. 75: 571-576.
Ghiorso, M. S. and R. O. Sack. 1995. Chemical mass transfer in magmatic processes IV. A revised
and internally consistent thermodynamic model for the interpolation and extrapolation of liquid-
solid equilibra in magmatic systems at elevated temperatures and pressures. Contrib. Mineral.
Petrol. 119: 197-212.
Goldsmith, J. R. and R. C. Newton. 1969. P-T-X relations in the system CaCO3-MgCO3 at high tem-
peratures and pressures. Am. J. Sci. 267A: 160-190.
Grover, J., 1977. Chemical mixing in multicomponent solutions: An introduction to the use of
margules and other thermodynamic excess functions to represent non-ideal behavior, in
Thermodynamics in Geology. ed. D. G. Fraser. 67-97. Dordrecht: D. Reidel.
Guggenheim. 1935. The specific thermodynamic properties of aqueous solutions of strong
electrolytes. Phil. Mag. 19: 588.
Helgeson, H. C. 1969. Thermodynamics of hydrothermal systems at elevated temperatures and pres-
sures. Am. J. Sci. 267: 729-804.
Helgeson, H. C. and D. H. Kirkham. 1974. Theoretical prediction of the thermodynamic behavior
of aqueous electrolytes at high pressures and temperatures, Debye-Hückel parameters for activity
coefficients and relative partial molar quantities. Am J. Sci. 274: 1199-1261.
Helgeson, H. C., Kirkham, D. H., and G. C. Flowers. 1981. Theoretical predicitons of the
thermodynmaic behavior of electrolytes at high pressures and temperatures: IV. Calculation of
activity coefficients, osmotic coefficients, and apparent molar and standard and relative partial
molar properties to 600°C and 5 kbar. Am. J. Sci. 281: 1249-1316.
Holdaway, M. J. 1971. Stability of andalusite and the aluminum silicate phase diagram. Am. J. Sci.
271: 91-131.
Holland, T. J. B., A. Navrotsky, and R. C. Newton, 1979. Thermodynamic parameters of
CaMgSi2 O6 – Mg2 Si2 O6 pyroxenes based on regular solution and cooperative disordering
models. Contrib. Mineral. Petrol. 69: 337-344.
Johannes, W. and D. Puhan. 1971. The calcite-aragonite equilibrium reinvestigated. Contrib. Min-
eral. Petrol. 31: 28-38.
Koziol, A. M. and R. C. Newton. 1988. Redetermineation of the garnet breakdown reaction and
improvement of the plagioclase-garnet-Al2O 5-quartz geobarometer. Am. Mineral. 73: 216-223.
Navrotsky, A. 1994. Physics and Chemistry of Earth Materials. Cambridge: Cambridge University
Press.
Nicholls, J. and J. K. Russell. 1990. Modern Methods in Igneous Petrology, Reviews in Mineralology,
vol. 24, Washington: Min. Soc. Am.
Nielsen, R. L. and M. A. Dungan. 1983. Low pressure mineral-melt equilibria in natural anhydrous
mafic systems, Contrib. Mineral. Petrol., 84: 310-326.
Nordstrom, , D. K. and J. L. Munoz, 1986. Geochemical Thermodynamics, Palo Alto: Blackwell Scien-
tific Publ.
Pitzer, K. S. 1979. Theory: ion interaction approach. in Activity Coefficients in Electolytes, ed. R. M .
Pytkowicz. 157-208. Boca Raton: CRC Press.
Robinson, R. A. and R. H. Stokes. 1959. Electrolyte Solutions. London: Butterworths.
Roeder, P.L. and R.F. Emslie. 1970. Olivine-liquid equilibrium. Contrib. Mineral. Petrol. 29: 275-289.
157 October 17, 2001

Spencer, K. and D. H. Lindsley. 1981. A solution model for coexisting iron-titanium oxides. Am.
Mineral. 66: 1189-1201.
Stolper, E. and H. Y. McSween. 1979. Petrology and origin of the shergottite meteorites. Geochem.
Cosmochim. Acta. 43: 1475-1498.
Thompson, J. B. 1967. Thermodynamic –
V S
properties of simple solutions, in Researches in f (cm3) (J/K-mol)
Geochemistry, vol. 2, ed. P. H. Abelson. New
kyanite 44.09 242.30
York: Wiley and Sons.
andelusite 51.53 251.37
Thompson, J. B. and D. R. Waldbaum, 1969.
sillmanite 49.90 253.05
Mixing properties of sanidine crystalline
solutions: III. Calculations based on two-phase data. Am. Mineral. 54: 811–838.
Truesdell, A. H. and B. F. Jones. 1974. WATEQ, a computer program for calculating chemical
equilibria in natural waters. J. Res. U. S. Geol. Surv. 2: 233-248.
Wood, B. J., 1987. Thermodynamics of multicomponent systems containing several solutions, in
Thermodynamic Modeling of Geological Materials: Minerals, Fluids and Melts , eds. I. S. E.
Carmichael and H. P. Eugster. 71-96. Washington: Mineral. Soc. Am.
Wood, B.J. and S. Banno. 1973. Garnet-orthopyroxene and orthopyroxene-clinopyroxene relation-
ships in simple and complex systems. Contrib. Mineral. Petrol. 42:109-124.
Problems
1. Kyanite, andelusite, and silimanite (all polymorphs of Al2SiO5) are all in equilibrium at 500°C and
376 MPa. Use this information and the adjacent table to construct an approximate temperature-pres-
sure phase diagram for the system kyanite-sillamanite-andelusite. Assume DV and DS are indepen-
dent of temperature and pressure. Label each field with the phase present.
2. Show that: Gexcess = WG1X2 + WG2X1 X1X2 may be written as a 4 term power expansion, i.e.:
G ex = A + BX 2 + CX 22 + DX 32
3. Construct G-bar–X diagrams for a a regular solution with W= 12 kJ (W is the interaction
perameter in a non-ideal solution) at 100° temperature intervals from 200 to 700° C. Sketch the
corresponding phase diagram.
4. Interaction parameters for the enstatite-diopside solid solution have been determined as follows:
WH-En = 34.0 kJ/mol, WH-Di = 24.74 kJ/mol (assume WV and WS are 0).
a.) Use the asymmetric solution model to calculate DGreal as a function of X2 (let diopside be com-
ponent 2) curves for this system at 100 K temperature from 1000 K to 1500 K. Label your curves.
b.) What is the maximum mole fraction of diopside that can dissolve in enstatite in this tempera-
ture range:?
c.) Sketch the corresponding T-X phase diagram.
5. Sketch G-bar–X diagrams for 1600° C, 1500° C, 1300° C, and 1250° C for the system Diopside-
Anorthite (Figure 4.8). Draw tangents connecting the equilibrium liquids and solids.
6. Suppose you conduct a 1 atm melting experiment on a plagioclase crystal. Predict the mole
fractions of anorthite in the liquid and solid phases at a temperature of 1425° C. Assume both the
liquid and solid behave as ideal solutions. Albite melts at 1118° C, anorthite at 1553°C. DH m for
albite is 54.84 kJ/mol; DHm for anorthite is 123.1 kJ/mol.
158 October 17, 2001

8. Given the following 2 analyses of basaltic glass and coexisting olivine phenocrysts, determine
the K D for the MgO ® FeO exchange reaction, and calculate the temperatures at which the olivine
crystallized using both MgO and FeO. Assume Fe2O 3 to be 10 mole% of total iron (the analysis
below includes only the total iron, calculated as FeO; you need to calculate from this the amount of
FeO by substracting an appropriate amount to be assigned as Fe2O 3). Note that the mole % Fo in
olivine is equivalent to the mole % Mg or MgO. (HINT: you will need to calculate the mole fraction of
MgO and FeO in the liquid).
Glass (liquid) composition:

Sample TR3D-1 DS-D8A
(wt % oxide) (wt % oxide)
SiO2 50.32 49.83
Al3 O2 14.05 14.09
SFe as FeO 11.49 11.42
MgO 7.27 7.74
CaO 11.49 10.96
Na 2 O 2.3 2.38
K2 O 0.10 0.13
MnO 0.17 0.20
TiO2 1.46 1.55
olivine
Mole % Fo (=mole % Mg) 79 81
8. Determine the temperature and oxygen fugacity of equilibration for the following set of
coexisting iron-titanium oxides in lavas from the Azores:
titanomagnetite s.s. phase ilmenite s.s. phase
mole percent magnetite mole percent hematite
G-4 groundmass 29.0 10.3

SJ-8 phenocrysts 41.9 13.0
SM-28 microphenocrysts 54.5 7.0
T-8 groundmass 33.7 8.1
F-29 microphenocrysts 36.2 6.0
Make a plot of ƒO2 vs. temperature using your results and compare with Fig. 3.21. What buffer
do the data fall near?
7. Starting from equations 4.23 and 4.62, use the fundamental relationships between free energy, en-
tropy, enthalpy, and the equilibrium constant to derive the temperature
dependence of the titan–omagnetite–ilmentite exchange (equation 4.63). Species Conc
g/kg
10. In a melt having a composition, in wt %, of: Na + 63.00
SiO2 58.12% TiO2 0.92% K+ 6.15
Al2O3 16.47% Fe2O3 1.82% Mg 2+ 2.77
MgO 5.62% FeO 9.94% CaO 7.11% Ca2+ 44.6
Use the Ghiorso regular solution model and the interaction parameters in Cl– 200.4
Table 4.3 to: 2-
– – SO4 0.13
a.) calculate the G ex and DG mixing for this composition at 1300°C. HCO 3
–
0.03
159 October 17, 2001

b.) calculate the activity of Si4O8 at this temperature.
11. An analysis of an oil field brine from Mississippi with a temperature of 150° C is shown to the
right. Calculate the activities of these species at that temperature using the Truesdell-Jones equation.
12. Show that for a strong electrolyte, i.e., one in which dissociation is complete and:
m– = n–m and m + = n+m
where m is the molality of the solute component An+Bn-, that:
1/n
m ± = m nn+ n–n
+ –
13. Mean ionic activity coefficients were measured for the following solutions at an ionic strength
of 3: gKCl = 0.569, gNaCl = 0.734, gNa2 CO 3 = 0.229. Assuming gCl – = gK – = g±KCl,
what is the activity coefficient of CO 32 - ? g±
I, m observed
14. Calculate the electrostatic and solvation contributions to the mean ionic 0.001
activity coefficient of MgCl2 at concentrations of 0.0033, 0.01, 0.033, 0.05, 0.1, 0.003 0.887
0.33, 0.5, and 1 using the Debye-Hückel equation and Robinson and Stokes 0.006 0.847
(equ. 4.92) equ. respectively. Plot results, as well as gelect+solv = gelect *gsolv as 0.01
a function of ionic strength. 0.015 0.78
0.03 0.716
15. Calculate the mean ionic activity coefficient for NaCO3 using the Debye-
0.06 0.644
Hückel and Truesdell-Jones equations and compare your results with the ob-
0.1
served values to the right. Overall, which fits the data better?
0.15 0.541
0.3 0.462
0.6 0.385
1
1.5 0.292
3 0.229
6 0.182
160 October 17, 2001

Chapter-04-Applications-Thermodynamics-Earth 4744 0 PDF

Uploaded by

Copyright:

Available Formats

Chapter-04-Applications-Thermodynamics-Earth 4744 0 PDF

Uploaded by

Document Information

Original Title

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

Chapter-04-Applications-Thermodynamics-Earth 4744 0 PDF

Uploaded by

Copyright:

Available Formats

W. M.

Chapter 4: Applications of Thermodynamics to the

† Named for M. Margules, who first used this approach in 1895.

118 October 17, 2001

119 October 17, 2001

Henry's Law constants:

120 October 17, 2001

121 October 17, 2001

122 October 17, 2001

123 October 17, 2001

4.4.1 The Thermodynamics of Melting 1400 L

124 October 17, 2001

125 October 17, 2001

Example 4.2: Calculating Melting Curves

126 October 17, 2001

127 October 17, 2001

128 October 17, 2001

129 October 17, 2001

130 October 17, 2001

sillimanite, and andalusite occur together. 400 Sil

131 October 17, 2001

132 October 17, 2001

133 October 17, 2001

134 October 17, 2001

135 October 17, 2001

136 October 17, 2001

dependence on oxygen fugacity

and temperature, as shown in

the compositions of the titano-

–15 magnetite and ilmenite phases,

2Fe3 O4 + O2 ® 3Fe2 O3 4.52

Hem1 FeTiO3 +FeO®Fe2 TiO4 4.53

137 October 17, 2001

AWHUsp - BWHMt - CWHIl + DWHHem + DH o

138 October 17, 2001

Example 4.4: Using the Iron-Titanium Oxide Geothermometer

139 October 17, 2001

4.6 Thermodynamic Models of Magmas

140 October 17, 2001

141 October 17, 2001

and the Gibbs Free Energy is:

* The 1 term arises because the sum contains both X X W ij ji ij ij

142 October 17, 2001

143 October 17, 2001

144 October 17, 2001

and solution sought by successive it-

20 on the FMQ buffer. From Ghiorso and Carmichael

145 October 17, 2001

where NF and NM refer to network-formers and network-modifiers respectively. Equilibrium con-

146 October 17, 2001

147 October 17, 2001

148 October 17, 2001

149 October 17, 2001

Substituting this into equation 3.45, we obtain:

the relationship between activity of a salt and its molality is:

150 October 17, 2001

4.7.2.2 Effects of Solvation

151 October 17, 2001

152 October 17, 2001