\nextTopic{gd_06.tex}
\previousTopic{gd_04.latex}
\listoffigures
\chapter {The heat of the matter}
{\em The Earth is a heat engine}, and as with all such devices it is doomed eventually
to die! Indeed, it is this very heat loss that is responsible for lithospheric growth
and ultimately provides the energy source driving the motions of the lithosphere.
Heat loss from the Earth's interior is accompanied both by conduction and advection,
which can be simply codified in terms of the requirement for the conservation of energy
- and which therefore can serve to illustrate one of the great conservation principles of
physics. \footnote{ Conservation principles express the simple reality of
{\em what goes in must come out!}}
\section{Conduction and advection}
Thermal energy may be transported in two fundamentally different modes; namely by
conduction and by advection. In the former heat diffuses through matter while in the latter
it is carried by the physical transport of matter.
\subsection*{Conduction}
Everyday experience witnesses the fact that thermal energy is transferred from hot bodies
to cool bodies, with the rate of thermal energy transfer, or heat flow, dependent on the
nature of the material. Fourier postulated that the heat flow, $q$, is proportional to
the negative of the temperature gradient, because heat flows down temperature gradients,
with the constant of proportionality given by a material dependent parameter known as the
conductivity, {$k$}, a postulate consistent with most subsequent experimental
investigation: \begin{equation} \label{eq:heatflow} q = -k\, {\rm grad} \, T
\end{equation} \footnote{\scriptsize{In a three-dimensional orthogonal co-ordinate system,
$x$, $y$, $z$, grad $T$ is: \[ {\rm grad} \,T = \nabla T \] \[ \, \, \, \, =
\left(\frac{\partial{T}}{\partial{z}}, \frac{\partial{T}}{\partial{y}},
\frac{\partial{T}}{\partial{z}}\right) \] That is, the gradient of T is perpendicular to
the planes of constant temperature or isotherms. }} Written in this form grad $T$
signifies the gradient of $T$ (often written in the equivalent form $ \nabla T$) and as
such is a vector operator. That is, grad defines the magnitude and direction of the
gradient of a scalar quantity such as temperature.
\subsection*{Advection}
The physical motion of material results in advective transport of heat.
Using a {\em Eulerian description} it is easy to see that the rate of heat advected
through a fixed point is dependent on the velocity of material
transport and the temperature gradient.
\footnote{\ In many problems in fluid mechanics we can distinguish between a
Eulerian and Lagrangian description. In the former, the description is in terms of a fixed
spatial reference frame through which the fluid flows while in the latter the
description is in terms of a moving reference frame attached to individual fluid particles. }
Thus, as either the tempertaure gradient or the velocity tends to zero, the rate of heat
advection must also tend to zero.
\section{The thermal energy balance}
In order to understand the thermal evolution of material subject to conductive heat
transfer it is necessary to derive the thermal energy balance that expresses the
fundamental physical constraint provided by the requirement of {\em conservation of
energy}. Ignoring the advection of heat through the movement of the medium then the thermal
energy balance in a control volume can be expressed simply as: \begin{center} {\it rate
of gain of thermal energy \\ = \\ rate of thermal energy flowing into the control
volume \\ - \\ rate\, of thermal energy flowing out of control volume \\ + \\ rate
of heat production in control volume} \end{center} \begin{equation} \label{tbalance}
\end{equation} Since the thermal energy of a volume of unit dimensions is given by: \[ C_p
\, \, \rho \, \, T \] where $C_p$ is the heat capacity, the first term in Eqn
\ref{tbalance}, the {\it rate of gain of thermal energy}, can be expressed simply as \[
\frac{\partial{\left( C_p\, \rho\, T \right)}} {\partial{t}}\] which for temperature
independent $\rho $ and $C_p$ is equivalent to \[ C_p \,\rho \,\frac {\partial{
T}}{\partial{t}} \] The difference between the second and third terms in Eqn
\ref{tbalance} gives the negative of the divergence of the heat flow across the boundary of
the control volume: \[- {\rm div}\, q \] \footnote{\scriptsize{The divergence of a
vector field such as $q$ is a scalar quantity and can be expresses equivalently as \[ {\rm
div} \,q = \nabla \cdot \, q = \frac{\partial{q_x}}{\partial{x}} +
\frac{\partial{q_y}}{\partial{y}} + \frac{\partial{q_z}}{\partial{z}} \]}} Note that the
negative sign arises because divergence actually measures the rate of loss whereas in Eqn
\ref{tbalance} we are really interested in the rate of gain (or convergence). Assuming
that $k$ is independent of $T$ then the difference between the second and third terms in
Eqn \ref{tbalance} can be expressed in terms of T in the following way by substituting the
expression for $q$ given in Eqn \ref{eq:heatflow}: \begin{eqnarray*} & k\, {\rm div}
\,\left({\rm grad} \, T \right)& = k \,\nabla \cdot \left(\nabla T \right)\\ &&= k \left(
\frac{\partial{^2T}}{\partial{x^2}} + \frac{\partial{^2T}}{\partial{y^2}} +
\frac{\partial{^2T}}{\partial{z^2}}\right) \end{eqnarray*} The final term in Eqn
\ref{tbalance} is given by: \[ \rho \, H \] where {$H $} is the heat production per unit
mass.
Thus the thermal energy balance given by Eqn \ref{tbalance} can be given in the following
ways: \begin{eqnarray} \label{diffusion} &\frac{\partial{T}}{\partial{t}} &= \kappa\left(
{\rm div} \,\left( {\rm grad}\, T\right)\right) +\frac{H}{C_p} \nonumber \\ && =
\kappa\left(\nabla \cdot \nabla T\right) + \frac{H}{C_p} \nonumber \\ && =
\kappa\left(\frac{\partial{^2T}}{\partial{x^2}} + \frac{\partial{^2T}}{\partial{y^2}} +
\frac{\partial{^2T}}{\partial{z^2}}\right) + \frac{H}{C_p} \end{eqnarray} where $\kappa$ is
the thermal diffusivity: \[\kappa = \frac{k}{\rho\, C_p} \] The various forms of the
thermal energy balance in Eqn \ref{diffusion} form the {\em diffusion equation} or,
alternatively, the heat equation.
The {\em diffusion-advection equation} appropriate to whole rock advection, is obtained by
adding an additional term to Eqn \ref{diffusion} which expresses the rate at which heat is
advected through a point. This advective contribution is given by the product of the
velocity and the temperature gradient, giving: \begin{equation} \label{difad}
\frac{\partial{T}}{\partial{t}} = \kappa \nabla^{2}T - v .\nabla \,T + \frac{H}{ C_p}
\end{equation} where the second term on the right is the {\em advective} term.
For many purposes, it is useful to recast Eqn \ref{difad} in a dimensionless form:
\begin{equation} \label{ddifad} \frac{\partial{T'}}{\partial{t'}} = \nabla^{2}T' - Pe_T
\,\nabla \,T' \end{equation} where the dimensionless variable $T'$ and $t'$ are given by
\[ T' = T/T_0\]
\[ z' = z/l
\]
\[ t' = t \kappa/ l^2 \] where $l$ is the appropriate length-scale
(for example, the thickness of the lithosphere),
and $T_0$ is the characteristic temperature difference at the appropriate length-scale.
%(see Appendix A.2 for the non-dimensionalisation procedure).
In Eqn
\ref{ddifad} $Pe_T$ is the thermal Peclet number which is a dimensionless variable \[ Pe_T =
\frac{v \, l }{\kappa} \]. The Peclet number
expresses the ratio of the advective to diffusive terms. For $Pe_T > 10$ the advective
term dominates the thermal evolution and diffusion can be largely ignored, while for $Pe_T <
1$ diffusion dominates and the advective term can be largely ignored. For $1 < Pe_T < 10 $
both diffusion and advection contribute to the thermal evolution.
\section{Thermal time constants}
It is useful to define the characteristic diffusive response time to thermal
perturbations when heat is transferred by conduction. For example consider a bunsen burner
applied to the base of a conducting plate (geologically, this is analagous to the
impingement of a mantle plume on the base of the lithosphere). We are interested in how
long does it takes for the top surface to feel the effect of the perturbed lower boundary
condition. Clearly the response time depends on the length-scale (i.e. the thickness
of the plate) and (inversely) on the thermal conductivity (or, equivalently, diffusivity).
In this problem we can ignore the effects if internal heat production (which does not
change the response time to a change in boundary conditions). From Eqn \ref{diffusion} the
term on the right, $\frac{\partial{T}}{\partial{t}}$, gives the order the temperature
divided by the timescale for thermal response, $\tau$. Similarly the frist term on the
right of Eqn \ref{diffusion}, $\kappa\left( {\rm div} \,\left( {\rm grad}\,
T\right)\right)$, is of order $\kappa T/ l^2$, where $l$ is the appropriate length scale.
Equating the two quantities gives : \[ \tau \propto \frac{l^2}{\kappa} \] The constant of
proportionality normally used is $1/\pi^2$ giving: \[ \tau = \frac{l^2}{\pi^2 \,\kappa} \]
For the lithosphere typical parameter values are $l = 125$ km and $\kappa = 10^{-6}$ m $^2$
s$^{-1}$ giving $\tau$ = 50 Ma.
An important feature of the thermal time constant is its dependence on $l^2$. For the
lithosphere this implies that short wavelength fluctuations are damped on timescales
which are very short compared to lithospheric-scale perturbations (and, incidently,
explains why a 3 m deep cellar is insensitive to seasonal atmospheric temperature
fluctuations).
\section{Continental geotherms}
In the steady-state, that is when $\frac{\partial{T}}{\partial{t}} = 0$, and when heat flow
is only in one direction, z, Eqn \ref{diffusion} reduces to the ordinary differential
equation: \begin{equation} \label{geot} 0 = k\,\frac{d^2 T}{d \, z^2} + \rho \, H
\end{equation} which can be readily integrated with appropriate boundary conditions to
yield a geotherm. For example, assuming a constant distribution of heat producing
elements with depth, Eqn \ref{geot} can be integrated using the boundary conditions
appropriate to the continental lithosphere: {$ q = - q_{s} $} at {$ z = 0 $}
(i.e., the heat flow at the surface of the lithosphere in the direction of increasing
depth is the negative of the surface heat flow), and {$ T = T_{s} $} at {$ z
= 0$} (i.e., the temperature at the top of the lithosphere is the surface temperature). The
resulting integration yields the expression \begin{equation} T_{z} = T_{s} +
\frac{q_{s}\, z}{k} - \frac{\rho\, H\, z^{2}}{2 k} \end{equation} relating the
temperature at any depth, {$z$}, to the conductivity, surface temperature, {$T_{s}$},
surface heat flow, {$q_{s}$}, and heat production, {$H$}. The base of the lithosphere is
defined by the temperature sensitive rheological properties of peridotite and can therefore
be treated as an isotherm, {$T_{l} = 1280^{o} C$}. The heat flow at the surface of the
lithosphere, {$q_{s}$}, therefore reflects the thickness of the lithosphere as well as the
distribution of heat sources within the lithosphere.
\section{Natural convection}
Natural convection \footnote{ Natural convection is distinguished from forced
convection where the convective motion is induced by "stirring".} is the motion induced
by thermally created density gradients in viscous materials, otherwise known as fluids.
\footnote{ In this sense the term {\it fluid} incorporates viscous solids, as
well as liquids and gases.} In the presence of gravitational field a compositionally
homogeneous fluid heated from below will develop an unstable density stratification. Once
convective motion is initiated and hotter, less dense material are juxtaposed against
cooler more dense material by slight upward displacement (and vice versa in a downwelling
zone) the motion will be maintained by buoyancy forces operating on the individual fluid
packets (Figure 3.1). Counteracting the tendency to rise will be the viscous drag of the
material and the tendency for the rising material to loose heat to its surroundings by
thermal diffusion (Figure 3.2 \& 3.3). Whether convective motion is initiated depends the
material properties (thermal diffusivity, coefficient of thermal expansion, viscosity).
The likelihood of convective motion will be enhanced by the fluid having a high
coefficient of thermal expansion, subject to a large temperature gradient over a large
vertical extent and diminished by a high thermal diffusivity or by high viscosity.
Initiation of the convective motion requires a finite amplitude perturbation of the
thermal, and hence density, structure in order to localize upwelling and or downwelling.
\hspace{2 cm}
\includegraphics{Figures/Chapter_05/Figure_01.jpg}
\begin{figure}[h,t,b]
%\vspace{2.5in} \hspace{0.7in} \special{epsf=:gdfigs:chap4:Fig4_1.epsscale=0.700}
\caption{ For the initiation of convection an imbalance
of forces is required. A fluid layer subject to an inverted density gradient due, for
example, to basal heating and/or surface cooling is capable of undergoing convection.
Slight upward displacement of part of the basal layer will produce a horizontal density
gradient which will lead to a net upward directed buoyancy force on the displaced fluid
parcel, further enhancing the prospect for upward motion and eventually fully developed
convective motion. Initiation of the convective motion requires a finite amplitude
perturbation of the thermal and hence density structure in order to localize upwelling and
or downwelling. } \end{figure}
\hspace{2 cm}
\includegraphics{Figures/Chapter_05/Figure_02.jpg}
\begin{figure}[h,t,b]
%\vspace{2.5in} \hspace{0.7in} \special{epsf=:gdfigs:chap4:Fig4_2.eps scale=0.700}
\caption{ The
displacement of a parcel of hot light fluid is resisted by viscous dissipation or drag in
the surrounding fluid (strictly the rate at which momentum diffuses into the surrounding
medium), a measure of which is given by the kinematic viscosity of the fluid. }
\end{figure}
The Rayleigh number, $Ra_{T}$, is a measure of the relative magnitude of the destabilizing
buoyancy forces (Figure 7.1) and the stabilizing effects of momentum and heat diffusion
(Figure 7.2 \& 7.3) and is given by: \begin{equation} Ra_{T} = \frac{gÊ\,\alpha\,\Delta
T\,ÊL^{3}}{\kappa \,Ê\nu} \end{equation} Where $g$ is the acceleration due to gravity (m
s$^{-2}$), $\alpha$ is the coefficient of thermal expansion ($^{o}K^{-1}$), $\delta T$
is the temperature difference across the fluid layer ($^{o}$K), $L$ is the vertical length
of the system (m), $\kappa$ is the thermal diffusivity (m$^2$ s$^{-1}$) and $\nu$ is the
kinematic viscosity of the fluid (m$^{2}$ s$^{-1}$). The kinematic viscosity is the ratio
of the viscosity $\mu$ to the density $\rho$. \begin{displaymath} \nu =
\frac{\mu}{\rho} \end{displaymath} A measure of the relative thickness of the viscous and
thermal boundary layers and is given by the Prandtl number, $Pr$ : \begin{equation} Pr =
\frac{\nu}{\kappa} \end{equation} Note that both the Rayleigh and Prandtl numbers are
dimensionless. When the Prandtl number is large the fluid diffuses momentum relatively
faster than it does heat. In such a fluid the slow diffusion of heat from a vertical
boundary into the fluid produces a thin layer of hot, low viscosity thermal boundary layer
adjacent to the walls. Through the rapid dissipation of momentum, the buoyancy of this
layer will drag more fluid into motion with the result that the thickness of the viscous
boundary layer is much greater than the thermal boundary layer. In a system in which the
temperature gradient is between horizontal surfaces, the critical Rayleigh number for the
onset of convection is Ra$_T$ = 10$^3$.
\hspace{2 cm}
\includegraphics{Figures/Chapter_05/Figure_03.jpg}
\begin{figure}
%[h,t,b] \vspace{2.5in} \hspace{0.7in} \special{epsf=:gdfigs:chap4:Fig4_3.eps scale=0.700}
\caption{ Displaced packets of fluids are no longer in
thermal equilibrium with the surrounds and thus gain or lose heat at a rate proportional
to the thermal diffusivity of the fluid. For large diffusivities the thermal equilibration
takes place quickly thereby annihilating the density contrasts that drive the convective
motion. } \end{figure}
\newpage
The geometrical form of convection in 3-dimensions may vary considerably with the physical
properties of the fluid.
\footnote{See Velarde and Normand, 1980, {\em Scientific American} who provide
an exccellent introduction to natural convection.} However,
in all cases continuity requirements dictate a general self-similar organization. At
initiation buoyancy forces induce upward motion of the heated fluid and due to the
continuity requirement \footnote{\scriptsize{Continuity refers to the basic requirement of
\em{conservation of mass.}}} simultaneously cooling parts of the fluid
body must begin to sink. These opposite motions obviously cannot occur at the same place
and the convective motion must organize itself into convection cells. Upwelling (and
downwelling) zones may either be planar or columnar in form.
Two types of fluid flow in convection can be defined; laminar or turbulent. Conditions for
these two modes are described by the Reynolds number, $Re$ (again a dimensionless
parameter): \begin{equation} Re = \frac{wÊ\, d}{\nu } \end{equation} where $w$ is the mean
velocity of the fluid (m s$^{-1}$), $d$ is the width or diameter of the source ($m$), and
$\nu$ the kinematic viscosity (m$^{2}$ s$^{-1}$). If not confined to a pipe then flow
becomes turbulent at values of $Re > 100$. Confined flow requires much higher values ($Re
= 2000$) to become turbulent. The nature of this flow is clearly important during any
sort of process of mixing since turbulence leads to much more efficient exchange of heat
and material between two fluids. Convective mixing in laminar flow regimes, as appropriate
to the mantle, is ultimately a process of stretching and folding. During this mixing
process initial compositional heterogeneities must thin with time, eventually to
the extent that diffusional exchange is capable of complete homogenisation and
equilibration with the host. At temperatures close to
those of the liquidi of magmas, convective flow in magma chambers will be vigorous and
turbulent and under these circumstances both efficient mixing and wall-rock exchange is to
be expected. As a magma body cools towards its solidus, it will change to laminar flow
and eventually cease to convect.
Plate tectonics is in part a consequence of thermal convection in the mantle largely driven
by radiogenic heat sources as well as residual heat generated due to the release of
gravitational energy during core-mantle separation early in earth history. In contrast,
convection does not occur in the interior of the geologically inactive planetary bodies
like the Moon and, most probably, Mars. The rate and geometry of mantle convection control
the rate of heat flow from the earth's interior to its surface, the rate at which crustal
extraction takes place and the rate at which material is recycled via subduction. Any
fully developed geodynamic hypothesis must therefore address the nature of mantle
convective geometry in the present earth. Convection requires viscous flow or
deformation. Clearly liquids or gases flow at rates observable on the scale of human
patience however the asthenosphere is also capable of viscous flow albeit with very slow
strain rates (of the order of $10^{-13} - 10^{-15}$s$^{-1}$). The very high viscosity of
the asthenosphere ($10^{20} - 10^{22}$ poise) leads to the possibility of Rayleigh numbers
appropriate for convection only when coupled with the very large length scales of
appropriate to the mantle thickness (note that the $Ra_{T}$ increases as a function of the
cube of the length scale, but decreases only to the first power of viscosity). Convection
of at least the upper part of the asthenosphere is probable on the grounds of the relative
homogeneity of MORB suggesting a quite well mixed source. This conclusion is also
suggested by the timescale and efficiency of extraction of the lithophile elements to the
crust from the depleted upper mantle.
Systems such as the Earth with substantial density stratification organise
themselves into discrete layered
systems. These layers may show internal convection with little material,
but significant thermal energy, exchange across
layer boundaries. In this case, heat transfer between a lower denser system and
the base of an overlying layer will be mediated by conduction through a thermal boundary
layer. This situation occurs in the Earth at the core-mantle boundary
at the so-called D'' (or double prime) layer. Here, the very dense
metallic core heats the base of the silicate mantle resulting
in periodic, thermally-induced, gravitational instabilities (Rayleigh-Taylor
instabilities).
When these instabilities achieve critical positive bouyancy they initiate the ascent of a
plume of hot mantle material.