Stacey F.D., Davis P.M. Physics of the Earth

Подождите немного. Документ загружается.

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

317

– [314–336] 13.3.2008 2:45PM

high temperature for kT > h, so that the energy

levels are effectively blurred into a continuum,

allowing the classical approximation (Eq. (19.3))

to be used. However, the mode theory is the

starting point for an understanding of the

Gr ¨uneisen parameter, as in the following

section.

The level to which a mode is excited ther-

mally is determined by the average thermal

energy, kT, of the other modes with which it

interacts. The probability of excitation to the

nth level is

pðnÞ¼expðnh=kTÞ=

n¼0

expðnh=kTÞ: (19:5)

This is the Boltzmann energy distribution, with a

denominator that is a normalizing factor to

make the sum of all p(n ) equal to unity. The

state n, when it occurs, contributes energy nh,

so that its contribution to the average energy of

the mode,

E,isnhp(n). Summing over all states,

n ¼0to1,

E ¼ h

n¼0

n expðnh=kTÞ=

n¼0

expðnh=kTÞ:

(19:6)

Writing x ¼exp(h/kT), the denominator of

Eq. (19.6) is recognized as a geometric progression,

n¼0

¼ 1=ð1  xÞ; (19:7)

and the numerator is

h

n¼0

¼ h

n¼0

xdðx

Þ=dx

¼ hxd

n¼0

=dx ¼ hx =ð1  xÞ

(19:8)

Using Eqs. (19.7) and (19.8) in (19.6), with

substitution for x,

E ¼ hx=ð1  xÞ¼h=½expðh=kTÞ1: (19:9)

This is the Einstein model of specific heat,

applicable to harmonic oscillators all of the same

frequency, .Athightemperatures,kT h,it

reduces to the classical situation,

E  kT,andat

low temperatures, kT h, there is negligible

excitation above the ground state,

h.Thecor-

responding contribution to specific heat is

obtained by differentiating Eq. (19.9) with respect

to T at constant V, with the assumption that  is

independent of T if V is held constant,

¼ð@E=@TÞ

¼ kðh=kTÞ

expðh=kTÞ=½expðh=kTÞ1

(19:10)

This reduces to k for h/kT 1 (the classical

limit). The assumption that  is constant at con-

stant V would be valid if a lattice mode were a

perfect harmonic oscillator, but, as in the classi-

cal theory referred to above, it must be seen as an

approximation when applied to an anharmonic

oscillator, with a potential function having a

form represented by Fig. 18.4.

Crystals have vibrational modes with a wide

range of frequencies and can be represented as

collections of Einstein oscillators with different

values of . There is a corresponding wide range

of temperatures, below which the modes are

‘frozen out’, becoming thermally inactive. Thus

the transition from low to high temperature

regimes is more spread than in the Einstein

model with a single mode frequency. With

detailed knowledge of atomic bonding, and

resort to numerical methods, realistic vibra-

tional spectra can be calculated, even for quite

complicated crystal structures. Figure 19.1 gives

an example. It is compared with a mathemati-

cally simple, widely used reference model by

P. Debye, dating from the early twentieth century.

The Debye theory treats a material as an elas-

tic continuum, with standing waves of all orien-

tations and wavelengths that can be fitted into it

down to a limit that gives the required total

number of modes, three times the number of

atoms. The number of standing waves in the

frequency range  to ( þd) that can be fitted

into a crystal is proportional to 

d and this is

the spectrum of frequencies that the Debye

model assumes, with frequency related to wave-

length by an averaged acoustic velocity. Since

the total number of modes in an N atom crystal

is 3N, the spectrum is cut off at the frequency, 

that gives this number. Thus the Debye model

19.2 SPECIFIC HEAT 317

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

318

– [314–336] 13.3.2008 2:45PM

has a characteristic temperature, the Debye

temperature, 

, which is related to 

k

¼ h

: (19:11)

At temperatures much higher than 

, all of the

modes are fully excited and C

approaches the

classical, high temperature limit, C

, given by

Eq. (19.3). At lower temperatures the high fre-

quency modes become inactive and C

has a

reduced value. Integration of Eq. (19.10),

weighted by the 

d spectrum, gives the Debye

function, which is the sum of specific heats of a

collection of Einstein oscillators (Eq. (19.10))

with a Debye spectrum of frequencies

¼ 3







ðe

 1Þ

dx: (19:12)

Although the Debye theory cannot accurately

represent the fine details of thermal properties,

because the spectrum of mode frequencies is not

very realistic, in many situations it is a conven-

ient approximation, especially its identification

of a characteristic (Debye) temperature. The

spectra in Fig. 19.1 look very different, but

when both are used to calculate specific heat

curves, as in Fig. 19.2, the Debye theory is seen

to reproduce the essential features. The mineral

represented in this figure, MgSiO

, perovskite,

has a high Debye temperature, about 950 K, and

at laboratory temperature C

(290 K)  0.65 C

For all higher temperatures the discrepancy

between the curves of Fig. 19.2 is quite small.

There is no analytical solution for the Debye

function (Eq. (19.12)) but for T >

thermal

energy and specific heat can be expanded as

polynomial functions of (

/T):

EðDebyeÞ¼C

T 



1680



þ



;

(19:13)

ðDebyeÞ¼C

1 





560





þ

(19:14)

In estimating the total stored energy in the Earth

(Chapter 21), the (3/8)

term in Eq. (19.13) is

important, especially as 

increases systemati-

cally with pressure. This term is represented by

the area between the Debye curve and the clas-

sical limit extrapolated to T ¼0 in Fig. 19.2.

It is useful to keep in mind the accuracy of the

approximation that C

is constant at high tem-

peratures. It is not valid for metals, as discussed

below, but for insulators we can differentiate

Eq. (19.14) for an adequate indication:

ð@ ln C

=@ ln TÞ

¼ð1=10Þð

=TÞ

ð3=1400Þð

=TÞ

þ:

(19:15)

Debye model

Full calculation

Number of lattice modes per

unit frequency interval

0 9 18 27 36

Frequency (THz)

FI G U R E 19.1 The spectrum of lattice modes of

MgSiO

perovskite, compared with the Debye theory.

Redrawn from Oganov et al. (2000).

Debye model

Full calculation

500

1000

1500

2000

Temperature (K)

Classical limit

(J K

–1

mol

–1

)

FI G U R E 19.2 Specific heat of MgSiO

perovskite

calculated from the lattice mode spectrum in Fig. 19.1,

compared with Debye theory. Redrawn from Oganov

et al. (2000).

318 THERMAL PROPERTIES

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

319

– [314–336] 13.3.2008 2:45PM

In the deep Earth (T/

) 1.6, giving (@lnC

@lnT)

0.04. This is no more than various uncer-

tainties, including the mean atomic weight and

the effect of anharmonicity (Section 19.8). For

the Debye model (@lnC

/@lnV)

¼0, as demonstra-

ted following Eq. (19.29), so by Eq. (E.6) in

Appendix E (@lnC

/@lnV)

0.05.

Up to this point the discussion concerns only

the lattice heat capacity. This is the total for

insulators, including the crust and mantle, but

for metals there is an addition al c ontribution by

conduction electrons. In the core it accounts for

about 30% of the total specific heat. In solids the

electron energy levels are spread by interac-

tions in to bands with wide energy ranges.

These bands are occupied by electrons up to a

level, known as the Fermi energy, set by t he

number of electrons. The characteristic of met-

als is that they have electron energy bands that

overlap and are partly filled, so that electrons

with energies close to the Fermi level may

change states in response to an electric field

(giving electrical conduction) or to thermal agi-

tation. Some are thermally excited to energies

of order kT above the Fermi level, leaving vacant

states below that level. The number of electrons

so excited and their average excitation energy

are each proportional to T, so that t he electron

thermal energy is proportional to T

and elec-

tron specific heat is proportional to T.Thisis

represented by

¼ T: (19:16)

For laboratory iron,  ¼4.98 mJ K

2

mol

1

0.0892 J K

2

1

. This means that a tem-

perature of 5000 K would make C

equal to the

lattice heat capacity by Eq. (19.3), but compres-

sion increases the spread of energy bands and

so decreases the density of energy states at

the Fermi level, reducing the value of .

Boness et al. (1986) reported a numerical cal-

culation of the band structure of iron with diffe-

rent crystal structures over the core pressure

range and summarized their results in a sim-

ple analytical approximation for  as a function

of density, , relative to the zero pressure

value, 

With substitution of the atomic weight of

pure iron, this is

 ¼ð6:113 

=  1:144ÞmJ K

2

mol

1

¼ð0:1094 

=  0:0205ÞJK

2

1

(19:17)

The band structure and the value of  hardly

depend on the crystal form and at core pressures

we assume that Eq. (19.17) is a good approxima-

tion also for liquid iron. Minor amounts of

lighter solute would have only a minor effect

on band structure. This equation is used to esti-

mate the electron contribution to core heat

capacity in the thermal model in Appendix G,

noting that Eqs. (19.16) and (19.17) give the elec-

tron contribution to C

and that Eq. (19.4) must

be used to obtain C

. It is important also to elec-

trical conductivity (Section 24.4) and conse-

quently also to the thermal conductivity of the

core (Section 19.6).

19.3 Thermal expansion and the

Gr¨uneisen parameter

A particular reason for geophysical interest in

thermal expansion is that the dilation of heated

material makes it buoyant relative to cooler,

surrounding material. In an extended medium

this causes convection, with an upward transfer

of heat by rising hot material. Thermal convec-

tion generates the mechanical energy driving

plate tectonics and the geological activity that

results from it. A uniform medium in which the

temperature gradient exceeds the adiabatic

one, that is the rate at which temperature

increases with depth due to compression, is

convectively unstable. Thus convection tends

to reduce any steeper gradient to the adiabatic

oneandthroughmuchoftheEarththetem-

perature gradient is believed to be close to

adiabatic. Calculations of this gradient

(Section 19.5) and of the mechanical energy of

convection (Chapter 22) rely on knowledge of

the volume coefficient of thermal expansion, ,

but there is no direct observation of  for the

deep Earth. It is estimated indirectly from the

Gr ¨uneisen parameter, .

The physical reason why materials expand

when heated is seen in Fig. 18.4. A bond between

neighbouring atoms that oscillates in length

19.3 THERMAL EXPANSION AND THE GRU

NEISEN PARAMETER 319

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

320

– [314–336] 13.3.2008 2:45PM

(between A and B) is extended more than it is

compressed and it spends more time in the

extended state because the restoring force (repre-

sented by the gradient of the curve) is then

weaker. The time-averaged extension is the

thermal expansion, which is a consequence of

the bond asymmetry. This has another effect. As

a material is compressed, the externally applied

pressure can be represented in Fig. 18.4 by an

added positive gradient, pushing the energy mini-

mum to smaller r and sharpening it. This means

increasing the elastic moduli. The increase in elas-

tic moduli with pressure is well observed by seis-

mology, which provides a direct measure of the

bond asymmetry. Thus thermal expansion and

the pressure dependences of elastic moduli have

a common cause and can be related by an appro-

priate theory. The Gr

uneisen parameter, , makes

this theoretical connection, which is the reason

for the effort by geophysicists to understand it.

The parameter originally established in the

solid state physics literature by E. Gr ¨uneisen, and

bearing his name, was derived from the equa-

tions for an Einstein oscillator (Eqs. (19.9) and

(19.10)). For a vibrational mode of a crystal lat-

tice, identified as an oscillator of frequency 

the Gr ¨uneisen definition is



¼ð@ ln 

=@ ln VÞ

: (19:18)



is the mode Gr ¨uneisen parameter, or mode

gamma, and is not in general identical to the 

of other modes. The interest in it arises from the

fact that, as shown below, an appropriate aver-

age of all of the 

is equivalent to the dimension-

less combination of thermodynamic parameters

in Eq. (19.1). In common usage, mention of the

Gr ¨uneisen parameter means the definition in

Eq. (19.1), but if this needs to be distinguished

from Gr ¨uneisen’s definition it is referred to as

the thermodynamic gamma. Its usefulness in

thermodynamic relationships is particularly

obvious in two frequently used identities (see

Table E.2 in Appendix E):

ð@ ln T=@ ln Þ

¼ð@ ln T=@ ln VÞ

¼ K

ð@ ln T=@PÞ

¼ ; (19:19)

ð@P=@TÞ

¼ K

¼ C

: (19:20)

Equation (19.19) gives the adiabatic temperature

variation in terms of either density or pressure

and Eq. (19.20) is the differential form of the

Mie–Gr ¨uneisen equation, relating (thermal) pres-

sure to thermal energy of a material heated at

constant volume,

Thermal

¼ 

C

dT  E

Thermal

=V; (19:21)

where the approximate equality invokes the

assumption that  is independent of T at con-

stant V. Over large depth ranges it varies rather

little and the assumption of constant  may be a

useful approximation. Then Eq. (19.19) integra-

tes to give

ðV



: (19:22)

This equation is familiar in ideal gas physics, for

which  ¼(C

1), because T ¼1, simplifying

Eq. (19.4). It illustrates the point that by using 

we introduce some of the simple insights of ideal

gas physics to solids.

To relate  by Eq. (19.1) to Gr ¨uneisen’s defini-

tion (Eq. (19.18)), we differentiate Eq. (19.9) with

respect to V at constant T and substitute for the

resulting two terms by Eqs. (19.9) and (19.10),

with introduction of subscript i to indicate that

we are considering a single mode. This gives the

volume dependence of mean mode energy,

=@VðÞ

@ ln 

@ ln V



h

exp h

=kTðÞ1





h

=kTðÞ

exp h

=kTðÞ

exp h

=kTðÞ1½

)

@ ln 

@ ln V









¼

@ ln 

@ ln V





: ð19:23Þ

Now E

is the internal energy attributable to

mode i and we make use of a general identity

for internal energy, U (see entries for @U

and @P

in Table E.2 of Appendix E),

ð@U=@VÞ

¼ð@P=@TÞ

T  P ¼ T

½@ðP=TÞ=@T

: (19:24)

320 THERMAL PROPERTIES

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

321

– [314–336] 13.3.2008 2:45PM

Identifying this expression with (@E

/@T)

Eq. (19.23) and dividing both sides by T

½@ðP

=TÞ=@T

¼ð@ ln 

=@ ln VÞ

ð1=VÞ½@ðE

=TÞ=@T

(19:25)

By making the assumption that 

and therefore



¼(@ln

/@lnV)

are independent of T at con-

stant V, as in the Einstein theory (the starting

point of the Gr ¨uneisen theory), we can integrate

Eq. (19.25) with respect to T at constant V to give

the relationship between thermal pressure due

to mode i and its thermal energy:

¼ð@ ln 

=@ ln VÞ

ðE

=VÞ: (19:26)

This has the form of the Mie–Gr ¨uneisen equation

(Eq. (19.21)), requiring Gr ¨uneisen’s definition of



by Eq. (19.18). Note that the Einstein assump-

tion, that 

is independent of T at constant V,

neglects anharmonicity (Section 19.8).

The value of  for materials having modes

with different 

is an average over all modes,

weighted according to the heat capacities of the

modes. Then, noting that  is the coefficient

relating thermal pressure to thermal energy, as

in Eq. (19.12),

 ¼



; (19:27)

where the C

are the various mode C

by Eq. (19.10).

If we restrict interest to high temperatures then

k (Boltzmann’s constant) for each mode and 

is then a simple average of the 

. A temperature-

dependence arises at low temperatures (addi-

tional to the small anharmonic effect of vibra-

tion amplitude on the 

) when modes of

different 

contribute in different proportions

because of their reduced C

by Eq. (19.10).

However, even in this situation the effect is com-

monly small. This is evident from a thermody-

namic identity (see Tables E.3 and E.1 in

Appendix E),

ð@=@ ln TÞ

¼ð@ ln C

=@ ln VÞ

; (19:28)

which vanishes identically in the Debye approx-

imation. This situation is represented by a

Debye version of Gr ¨uneisen’s definition of 

(Eq. (19.18)),



¼ð@ ln 

=@ ln VÞ

¼ð@ ln 

=@ ln VÞ

;

(19:29)

because 

and 

are related by Eq. (19.11). With

the assumption that, as for all 

in Gruneisen’s

theory, 

is independent of T at constant V, the

isothermal derivative in Eq. (19.29) is the same as

an adiabatic derivative, so that in comparing

Eqs. (19.19) and (19.29) we see that (T/

) and

therefore C

are constant on an adiabat. Of

course Debye theory is only an approximation,

but it is useful in understanding why C

approximately constant on an adiabat.

Equation (19.18) is the starting point for a

derivation of the most widely used formula for

estimating  for the Earth’s interior. This is the

‘acoustic gamma’, 

, which is based on the

assumption that there are only two kinds of lat-

tice mode, corresponding to compressional and

shear waves, with two shear modes for each

compressional mode. This is better than the

Debye method of averaging wave speeds. It

assumes an isotropic medium but makes the

best possible use of the seismologically observed

P- and S-wave speeds, being a weighted sum of

gammas for compressional (P) and shear (S)

modes:



¼ð1=3Þ

þð2=3Þ

: (19:30)

Taking mode i, controlled by elastic modulus X

to be a standing wave of speed V

¼(X

/)

1/2

, with a

wavelength l

that is a fixed number of lattice

spacings and therefore proportional to V

1/3

, the

mode frequency is



¼ V

/ X

1=2

1=6

; (19:31)

and so, by differentiation,



¼

@ ln 

@ ln V



¼









;

(19:32)

where K

¼V(@P/@V)

is the isothermal bulk

modulus. For compressional waves X

¼(K

þ(4/

3)) and for shear waves X

¼, so that by

Eq. (19.30)

19.3 THERMAL EXPANSION AND THE GRU

NEISEN PARAMETER 321

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

322

– [314–336] 13.3.2008 2:45PM







@







@





: ð19:33Þ

The temperature variation through most of

the Earth is believed to be close to adiabatic, so

that, to the extent that seismological models give

reliable gradients of the elastic moduli, they give

(@K

/@P)

and (@/@P)

. To convert these to the

isothermal derivatives required by Eq. (19.33)

we need the temperature derivatives, which are

written in dimensionless forms, as defined in

Appendix E,



¼ð1=K

Þð@K

=@TÞ

; (19:34)

" ¼ð1=Þð@=@TÞ

: (19:35)

By writing a derivative identity for any para-

meter X,

ð@X=@PÞ

¼ð@X=@P Þ

ð@X=@TÞ

ð@T=@PÞ

; (19:36)

with (@T/@P)

¼T/K

from entries in Table E.2 of

Appendix E, taking X to be either K

or , and

using Eqs. (19.34) and (19.35) to substitute for the

temperature derivatives, we obtain the required

relationships:

ð@K

=@PÞ

¼ð@K

=@PÞ

þ T

; (19:37)

ð@=@PÞ

¼ð@=@PÞ

þ T"=K

: (19:38)

These expressions may be obtained from entries

in Table E.3 of Appendix E.

The calculation of 

and " (Eqs. (19.34) and

(19.35)) is a subject of Section 19.7, which is con-

cerned with the interpretation of deep seismic

velocity anomalies. For the lower mantle, for

which there are sufficient data for detailed calcu-

lations, numerical values are given in Appendix G.

The temperature terms in Eqs. (19.37) and (19.38)

add 0.1 to 0.3 to the estimate of 

at deep Earth

temperatures and pressures, although they have

often been neglected.

A critical assessment of the assumptions

and approximations implicit in the derivation

of Eq. (19.33) is given by Stacey and Davis

(2004), who conclude that it survives the doubts

very well, whereas alternative formulae face

unresolved difficulties. However, care is

required in applying seismological estimates

of dK/dP and d/dP for the deep Earth. Models

such as PREM (Appendix F) represent seismic

velocities and density as polynomials in radius

over each of several radius ranges. This is a

mathematically convenient but unphysical rep-

resentation of material properties, so that differ-

entiation to obtain dK/dP and d/dP yields

satisfactory averages but implausible depth var-

iations. The difficulty is overcome by fitting the

PREM model to a finite strain theory (Chapter 18)

for which derivative properties have been prop-

erly constrained by thermodynamic principles.

The Gr ¨uneisen theory and derivation of

the acoustic formula (Eq. (19.33)) are specific to

solid insulators, but the thermodynamic defini-

tion of  (Eq. (19.1)) is applicable to every material.

We need to consider its application to metals and

to liquids because  is as important to the core as

it is to the mantle. There are other approaches

that we now consider, but they require calibra-

tion by the acoustic calculation.

The ‘free’ conduction electrons in metals

respond to heat as well as electric fields. Their

contribution to specific heat is discussed in

Section 19.2 and thermal conductivity in

Section 19.6. In the core the electron heat

capacity is approximately 30% of the total and

there is a corresponding electron component of

thermal pressure. This is related to thermal

energy by an electronic Gr ¨uneisen parameter,



. The total  for a metal is the sum of lattice

and electron contributions, weighted according

to the contributions to specific heat in the same

way as the lattice contributions are added in

Eq. (19.27). In principle it is possible to calculate



with a sufficiently accurate model of electron

band structure, but then we do not have an inde-

pendent estimate of the lattice  because the

conduction electrons contribute to elasticity. In

any case most of the core is liquid, precluding

use of Eq. (19.33), in which  is prominent, so a

different approach is required.

An alternative class of formulae, equally appli-

cable to solids and liquids, has the general form

 ¼½K

=2  1=6 ðf =3Þð1  P=3KÞ=½1 ð2 f =3ÞP=K;

(19:39)

322 THERMAL PROPERTIES

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

323

– [314–336] 13.3.2008 2:45PM

with different values of f according to what is

assumed about the thermal motions of atoms

(or electrons). We refer to these formulae as the

free volume type because Eq. (19.39) with f ¼2

was derived by Vashchenko and Zubarev (1963)

from free volume theory. A linear model by

Dugdale and MacDonald (1953), giving f ¼1, is

still sometimes used. A special case, Slater’s

gamma, for which f ¼0 (Eq. (18.54)), is used in

Section 18.6 for the infinite pressure extrapola-

tion. Slater’s formula is derived directly from

Gr ¨uneisen’s definition (Eq. (19.18)) with the

assumption that all 

are equal because all

mode frequencies have the same volume

dependence. As in Eq. (19.32), this means that

all moduli have the same volume dependence,

that is Poisson’s ratio,  (not to be confused with

mode frequency), is independent of pressure.

Although  is observed to increase with compres-

sion, for reasons discussed in Section 18.8, caus-

ing Slater’s formula to overestimate , the final

paragraph of Section 18.8 points out that in the

infinite pressure extrapolation Slater’s assump-

tion becomes valid, making Eq. (18.55) rigorous.

It provides an important constraint on finite

strain theories, as outlined in Section 18.6. A

difficulty with Eq. (19.39) is that we have no

satisfactory theoretical value of f. An empirical

solution to this problem, adopted by Stacey and

Davis (2004), is to match Eq. (19.39) at P ¼0 to the

zero pressure extrapolation of the acoustic for-

mula (Eq. (19.33)). Fitted to lower mantle data, it

gives f ¼1.44. A justification for this is that it

gives values of  indistinguishable from those

obtained by a third theory, based on thermo-

dynamic requirements on derivatives of ,as

discussed by Stacey and Davis (2004).

It is interesting to note that  appears to be

little affected by melting. This can be understood

by recognizing thermal motion as vibrations of

lattice modes that are dominated by very high

frequencies (10

to 10

Hz). As we discuss in the

following section, the fluidity of a liquid is due to

the mobility of the crystal dislocations with

which it is saturated, but dislocations cannot

respond to stress cycles with these frequencies.

So, in the consideration of lattice vibrations, the

elasticities of liquid and solid are not very

different.

To use Eq. (19.39) for the core we need to

justify its application to metals with strong con-

tributions to  by conduction electrons. The

point is that this equation is calculated by apply-

ing the Mie–Gr ¨uneisen equation (Eq. (19.21)),

which is a model-independent identity if  is

independent of T. It simply depends on a calcu-

lation of thermal pressure and its relationship to

the volume dependence of bulk modulus, and in

this respect electron pressure is no different

from lattice pressure. While the electron and

lattice gammas may be very different, they con-

tribute to the total in the manner of Eq. (19.27).

The essential feature of Eq. (19.39) is that it

makes no appeal to  but only to K and its pres-

sure dependence.

19.4 Melting

The standard thermodynamic identity for the

variation of melting point, T

, with pressure, P,

is the Clausius–Clapeyron equation,

=dP ¼ V=S; (19:40)

where DV and DS are the volume and entropy

increments of the melting process. The deriva-

tion of this equation relies on the fact that for a

substance in thermodynamic equilibrium at

fixed pressure the Gibbs free energy, G (defined

in Table E.1 of Appendix E), is a minimum. (At

constant V the Helmholtz free energy, F,isa

minimum and at constant S it is the enthalpy,

H, that is a minimum.) The coexistence of solid

and liquid at the melting point means that G has

the same value for both states at all points on the

melting curve. By considering an increment dG

in G due to simultaneous increments in P and T

and substituting for derivatives by entries in

Table E.2 (Appendix E), we have

dG ¼ð@G =@PÞ

dP þð@G=@TÞ

¼ VdP  SdT:

(19:41)

Identifying T with T

, so that the increments

follow the melting curve, we require that dG be

the same for solid (subscript S) and liquid (sub-

script L) and therefore

19 . 4 M E LT I N G 323

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

324

– [314–336] 13.3.2008 2:45PM

dP  S

dT ¼ V

dP  S

dT; (19:42)

which rearranges to the form of Eq. (19.40). This

is an identity, applicable to all phase transitions,

and is referred to also in Section 22.3, in connec-

tion with convection through solid–solid phase

transitions in the mantle.

In the case of melting, DS ¼S

S

is always

positive because latent heat L ¼DST

must be

applied to cause melting. In what we refer to

here as ‘normal’ or ‘simple’ melting DV is also

positive, that is the liquid is less dense than the

solid, so that T

increases with P. The best

known exception is water for which DV is neg-

ative (at low pressure), and it is useful to keep

this case in mind in considering ‘normal’ or

‘simple’ melting and why there are exceptions.

A successful general theory of melting is that it is

a free proliferation of crystal dislocations (see

Fig. 14.4) and that it occurs when the free ener-

gies of the undislocated crystal and one satu-

rated with dislocations (identified with the

liquid) are equal. Accepting this as a theory of

‘simple’ melting we can see why water does not

fit in. The theory assumes that there is no major

change in atomic coordination between the solid

and liquid states, because it is not changed much

by the introduction of dislocations, but common

ice is structurally quite different from liquid

water. As in water, strongly oriented polar

bonds occur in silicates, which also do not fit

well with the concept of ‘simple’ melting, but

most metals (bismuth excepted) do so. When

very high pressure is applied both solid and

liquid structures become more close-packed

and, regardless of low pressure behaviour,

assume greater structural similarity, allowing

simple melting theory to be applied with increas-

ing confidence. The principal geophysical appli-

cation is to the melting point of iron and

solidification of the inner core.

Equation (19.40) is basic to the theory of melt-

ing, but it can be extrapolated to high pressures

only with assumptions about both DV and DS and

this is not directly useful, so that several alter-

native theories of melting have arisen. The one

with the strongest influence on modern ideas

originated with F. A. Lindemann, who deduced

that melting point varied with volume and

Debye temperature as V

2/3



. Using Eq. (19.29)

and assuming  ¼

(Eq. (19.29)), this differenti-

ates to give (Problem 19.6)

ð1=T

ÞdT

=dP ¼ 2ð  1=3Þ=K

: (19:43)

Equation (19.43) is often attributed to Gilvarry

(1956), although comparison with his equations

(38) and (31) shows that his theory gives an addi-

tional factor [1 þ2( 1/3)T

]

1

on the right-

hand side. The basis of its derivation was an

assumption that melting occurs when the ampli-

tude of atomic vibration reaches a critical frac-

tion of the atomic spacing. There are now

derivations of the same or very similar formulae

with stronger fundamental foundations, one of

which has a thermodynamic basis and is pre-

sented here. But there is an underlying, implicit

assumption that the melting process is at least

very similar to a proliferation of dislocations,

and a simple model by Stacey and Irvine (1977)

shows why. Atoms in a dislocation are displaced

from their equilibrium positions so that some

bonds are stretched and others compressed, but

the forces are balanced. The Stacey and Irvine

calculation simulated this situation with two

lines of atoms locked together at their ends but

having unequal numbers of atoms, so that one

line was compressed and the other stretched.

Bond asymmetry, as in Fig. 18.4, causes an ave-

rage extension, which is identified with DV in

Eq. (19.40), and the energy increment, relative

to the undislocated lines, is identified with latent

heat, T

DS. This allows the equation to be writ-

ten in terms of derivatives of the atomic poten-

tial function and when K, P and K

are substituted

by Eqs. (18.18) to (18.20), the equation has the

form of Eq. (19.43) with the one-dimensional

(Dugdale–MacDonald) formula for  (Eq. (19.39)

with f ¼1). This general conclusion is fundamen-

tal and important because it shows that an equa-

tion with the form of Eq. (19.43) is independent

of assumptions about dislocation structure.

We present here a more rigorous, thermo-

dynamic derivation of the Gilvarry-type melting

law. Adopting an argument by Stacey et al.

(1989), if we consider a mass m of solid to be

heated at constant volume (without melting)

then the applied heat mC

DT is (neglecting

324 THERMAL PROPERTIES

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

325

– [314–336] 13.3.2008 2:45PM

anharmonicity) equipartitioned between increases

in atomic kinetic energy DE

and bond potential

energy DE

,sothat

E

¼ð1=2ÞmC

T: (19:44)

When the same temperature increment is

applied at constant pressure, heat mC

DT is

required, but since the temperature increment

is the same, DE

is still given by Eq. (19.44) and

therefore

E

¼ mðC

 C

=2ÞT: (19:45)

With the relationship between C

and C

(Eq. (19.4)) this becomes

E

¼ð1=2ÞmC

Tð1 þ 2TÞ: (19:46)

This is the increment in average bond potential

energy caused by heating at constant P. At the

same time there is thermal expansion,

V ¼ VT; (19:47)

and the ratio is

V=E

¼ 2=ðm=VÞC

ð1 þ 2TÞ

¼ 2=½C

ð1 þ 2TÞ; (19:48)

which, with the definition of  (Eq. (19.1)), is

V=E

¼ 2=½K

ð1 þ 2TÞ: (19:49)

Now melting (at constant pressure) is an appli-

cation of heat that causes no temperature

rise so that the latent heat is entirely bond

potential energy and Eq. (19.49) can be identi-

fied with Eq. (19.40) by writing DE

¼T

DS,so

that

ð1=T

ÞdT

=dP ¼ 2=½K

ð1 þ 2T

Þ: (19:50)

This is similar to Eq. (19.43), which, as mentioned

in the discussion of this equation, is often attrib-

uted to Gilvarry (1956) but does not correctly quote

him. His theory leads to an equation with the form

of Eq. (19.50) but with ( 1/3) replacing .

Equation (19.50) can be rewritten to repre-

sent the variation of T

with density at melting

by applying a melting curve modulus, derived by

Stacey et al. (1989),

¼ ð@P=@Þ

T¼T

¼ K

ð1 þ K

=dPÞ: (19:51)

Combined with Eq. (19.50) it yields a simpler

result:

dlnT

=dln ¼ 2: (19:52)

Both Eqs. (19.50) and (19.52) slightly overesti-

mate the pressure dependence of T

. They rely

on the assumption that bonds are stretched and

compressed but not broken. In the melting of

atomic close-packed structures roughly 2% of

atomic bonds are broken and this is a measure

of the errors in these equations.

At first sight comparison of Eqs. (19.19) and

(19.52) suggests a factor of two for the ratio of

melting point and adiabatic gradients, but it is

somewhat less than this because the variations

of density with pressure for the two equations

are represented by different moduli. This point is

elaborated in the following section.

Equation (19.50) or (19.52) can be applied to

iron under core conditions but are not satisfac-

tory for extrapolation from zero pressure. The

assumptions in the derivation break down in the

vicinity of a triple point and there are two (at

least) on the iron melting curve. The low pres-

sure  (body centred cubic) structure converts

first to  (face centred cubic) and then to " (hexa-

gonal close packed) or to something very similar.

As understood in the dislocation theory of melt-

ing, the structure of a melt in equilibrium with

solid crystals is a dislocated version of the crystal

structure. Under conditions remote from triple

points this is unambiguous, but in the vicinity of

a triple point there is a progressive transition in

liquid structure from a dislocated low pressure

form to a dislocated high pressure form. It is not

a sharp transition, as in the corresponding sol-

id–solid phase transition, but smeared over a

range in pressure and varying also with temper-

ature, as demonstrated by Sanloup et al. (2000)

for liquid iron in the vicinity of its  ––liquid

triple point.

The melting point of iron at inner core boun-

dary (ICB) pressure is experimentally inaccessi-

ble. It is beyond the range of diamond anvils

and shock waves cause melting before that pres-

sure is reached. However, observations of the

melting point of iron in its " phase can be

extrapolated to ICB pressure because no further

phase transitions are involved. There is still a

19 . 4 M E LT I N G 325

//FS2/CUP/3-PAGINATION/SDE/2-PROOFS/3B2/9780521873628C19.3D

–

326

– [314–336] 13.3.2008 2:45PM

problem that diamond anvil measurements (e.g.

Boehler, 1993) have given consistently lower

estimates of T

than shock waves, but in view

of difficulties with both temperature measure-

ment and pressure calibration in shock wave

experiments, we extrapolate the diamond

anvil results to ICB pressure by Eq. (19.52),

obtaining T

¼5750 K for iron. Depression of

the melting point by solutes is also uncertain,

but we allow 750 K, giving 5000 K as the ICB

temperature. This is a conveniently rounded

number that acknowledges the uncertainty,

which is about 500 K.

The density increment by freezing is also of

geophysical interest. The conventional approach

to ‘simple’ melting, as discussed by Poirier

(2000), gives an entropy increment for n moles:

S ¼ nR ln 2 þ K

V; (19:53)

which can be used with Eqs. (19.40) and (19.50) to

give

V ¼ð2T

ÞnR ln 2: (19:54)

Applied to the core at ICB conditions this gives

a density increment of 200 kg m

3

(1.6%). The

difference between this and the observed

820 kg m

3

density contrast (Masters and

Gubbins, 2003) is a measure of the compositional

difference and the consequent gravitational

energy release by progressive freezing of the

inner core (Section 22.6).

19.5 Adiabatic and melting point

gradients

Global geological processes at all levels are con-

trolled or caused by convection, either thermal

convection, as in the mantle, or with a strong

compositional effect, as in the core. In the case of

thermal convection, heat sources must maintain

a temperature gradient steeper than the adiabat,

even if only slightly so, and mechanical stirring,

as by compositional convection, also establishes

or maintains an adiabatic gradient. The exis-

tence of adiabatic gradients is an essential fea-

ture of an active planet.

To see why thermal convection can occur

only in a medium with a temperature gradient

exceeding the adiabatic value, consider a homo-

geneous medium with an arbitrary temperature

gradient in which a small volume of material is

displaced vertically upwards without allowing

any heat transfer to or from it. Its temperature

falls by adiabatic decompression as it rises. If it is

then at the same temperature as its new sur-

roundings the medium has an adiabatic gra-

dient. If the gradient in the medium is less than

this then the element is cooler and therefore

denser than its new surroundings and tends to

sink back to its original level. The medium is

stable, with no tendency to convect. On the

other hand, if the gradient in the medium is

steeper than the adiabat, then the displaced vol-

ume element is hotter and less dense than its

new surroundings and so tends to rise further.

The medium is then unstable and may convect

spontaneously. The steeper is the temperature

gradient, the stronger is the convection. Since

hot, rising material is displaced by cooler falling

material there is a net upward transfer of heat,

and if the heat sources are not maintained the

temperature gradient falls and convection ceases

when it reaches the adiabatic value. In the core

the temperature gradient is believed to be very

close indeed to adiabatic and the PREM model of

the lower mantle indicates that it is only slightly

steeper there (Section 17.5).

The adiabatic gradient can be written in

terms of Eq. (19.19),

ð@T@zÞ

Adiabatic

¼ð@T@PÞ

dP=dz ¼ðT=K

Þg;

(19:55)

where g is gravity. With substitution for  by

Eq. (19.1) we have an alternative form in terms of ,

ð@T@zÞ

Adiabatic

¼ Tg=C

; (19:56)

but, for the purpose of obtaining temperature

differences, integration by Eq. (19.55) or (19.19)

is more convenient. This is particularly obvious

if constant  is a sufficient approximation, so

that Eq. (19.22) can be used. As pointed out in

Chapter 22, it is the adiabatic temperature ratios

between heat sources and sinks in a convecting

medium, calculated by Eq. (19.55) or (19.56), that

326 THERMAL PROPERTIES