BASIC MECHANISMS
4.5. HEAT TRANSFER
To close the set of equations specifying ice deformation, the temperature distribution needs to be known simultaneously in order to adjust the ice
stiffness Parameter. This is necessary, since A(T*) changes by three orders of magnitude for the ternperature range encountered in polar ice sheets (-50°
-
O°C Paterson, 1981). In the model, ice temperature is calculated in the grounded ice sheet, but prescribed in the ice shelf. Some of the model runs, in particular the time-dependent ones, also consider heat conduction in the bedrock below.
4.5.1. Ice temperature
Taking the fixed coordinate System of fig. 4.2, the temperature distribution within the ice sheet is calculated frorn the general thermodynamic equation governing the transfer of heat in a continuum. It reads:
Here T is absolute ice temperature [K], t time [Y], ki thermal conductivity [Jrn-1K- l y - I ] , p ice density [910 kgm-31, cp specific heat capacity [Jkg-1K-11, V three- dimensional ice velocity [my-11 and <I> internal frictional heating [Jm^y-11 caused by deformation. In this equation, heat transfer is considered to result frorn vertical diffusion (first term), three-dimensional advection (second term) and deformational heating. Simplifications made in (4.39) include the use of constant density and the omission of melting and refreezing processes in the variable density firn layer (which is of minor importance on the scales considered). Furthermore, horizontal conduction in an ice sheet can be disregarded as temperature gradients in the horizontal directions are usually small compared to the vertical gradient, thus V2T can safely be replaced by 32T13z2.
The calculations also include the temperature dependence of the thermal Parameters of ice. The effect is not negligible since cp and ki may change by up to 30% for temperatures ranging between O° and -50°C The following relations were adopted:
C = 2115.3
+
7.79293(T-273.15)P (Pounder, 1965)
ki = 3.101 X 10 exp 8 ( -0.0057 T ) (Ritz, 1987) (4.40)
The internal heating rate per unit volume can be expressed as:
Assuming that the deformational heating from longitudinal strain-rates is small cornpared to that from horizontal shear strain-rates (which is certainly true near the base where heating is largest), leads to:
according to (4.10), (4.1 1) and the assumptions made in (4.12), (4.13).
4.5.2. Boundary conditions
Boundary conditions are chosen as follows. At the surface, temperature is Set equal to the mean annual air-temperature at that altitude and location. This is justified because over the Antarctic ice sheet screen temperatures are generally within 1° of the firn temperature at 10 m depth, where the annual cycle fades out (Loewe, 1970). At the base, the ice sheet gains heat from both sliding friction and the geothermal heat flux. These contributions can most easily be incorporated in the basal temperature gradient:
where yg [Km-'] is the geothermal heat entering the ice expressed as a ->
ternperature gradient, two-dimensional basal shear Stress and v(h) basal sliding velocity. The geothermal heating can take two forms, depending On whether or not heat conduction is considered in the rock below:
Yg
=-
G no bedrock heat conductionki
kr
Yg
=- {
with bedrock heat conductionki
A value of -54.6 mWlrn2 or -1.72 X 106 Jm-*Y-1 is taken for the geotherrnal heat flux G (Sclater et al., 1980). This corresponds to 1.30 HFU (Heat Flow Units), which is a typical value for old Precarnbrian shields. As for the thermal conductivity of rock kr, a value of 1.041 X 108 Jm-f K-ly-1 is adopted, which is a rnean value for a number of rocks listed in Turcotte and Schubert (1982).
Phase changes at the base are incorporated in the model by keeping the basal ternperature at the pressure rnelting point whenever it is reached and using the surplus energy for melting. This may eventually lead to the formation of a temperate ice layer between z = h and z = zmpn. The basal melt rate S [rny-I] can be calculated as follows:
where ( 3 T l 3 ~ ) ~ is the basal temperature gradient after correction for pressure melting, ( 3 T 1 3 ~ ) ~ the basal temperature gradient as given in (4.42) and L the specific latent heat of fusion [3.35 X 105 ~ k g - 1 1 . The pressure rnelting point is given (e.g. Paterson, 1981, p.193) by:
with T. = 273.15 K, the triple-point of water, H ice thickness and à = 8.7 X 10-4 K m - T of ice. When the melting point is reached, the rnelt rate beneath the Antarctic ice sheet is typically of the order of a few rnmly, but may reach values of several crnly in the fast flowing outlet glaciers. The present model does not rnake use of the melt rate, but S might be a useful quantity when a relation is added that governs melt water flow. Such a process can be
mportant, for instance in the development of an advanced theory for basal sliding.
In the ice shelf, the effects of basal melting andlor basal accretion and of spatially varying density and ice stiffness were ignored. Instead, a steady state linear vertical temperature distribution was assumed, with a surface temperature of -18 ' C (representing present day ice shelf conditions and allowed to vary with climatic change) and a fixed basal temperature of -2 'C.
The flow Parameter is then taken to represent mean ice shelf temperature conditions. This irnplies that the thermal inertia of the ice shelves is ignored.
However, the assumption of stationarity seems safe to make because the associated response time scale is small compared to the one for land ice (of the order of 102 years as compared to 103-1 04 years).
4.5.3. Rock temperature
Since the thermal conductivity of rock is of a comparable magnitude to that of ice, any time-dependent model experiment should take into account the thermal inertia of the bed. As pointed out by Ritz (1987), this effect may damp the basal temperature response by up to 50% for climatic oscillations operating on the longer time scales. In the bed, only heat transfer by vertical diffusion needs to be taken into account :
where Tr is rock temperature [K], pr rock density [3300 kgm-31, C,. specific heat capacity [ I 0 0 0 Jkg-1K-11 and kr was specified above. The lower boundary condition is given by the geothermal heat flux entering through the base of the rock column. At the upper rock surface, boundary conditions follow from the basal ice temperature and by considering flux continuity at the ice-rock interface:
T (upper) = T, (basal) and
at the upper surface
at the lower surface (4.48)
where the indices r and i refer to rock and ice respectively. In the model, a rock slab of 2000 m thickness is considered, which appears to be sufficient to describe the essential features of the physical process.
Including the bedrock thermal inertia has the effect that the basal temperature gradient in the ice also depends on the temperature evolution in the rock. In effect, adding such a 'thermal buffer' gives rise to a negative feedback :
0 1 2 3 4 5 6 7 8 9 1 0
Temperature rise [K]
fig.4.4: Demonstration of the effect of heat conduction in the rock on the englacial ternperature distribution. A sudden temperature increase of 10 K is applied at time Zero at the upper ice surface. Vertical heat diffusion and advection are considered in an ice slab of 3000 rn thick that overlies a rock slab of 2000 m thickness. The upper velocity boundary condition is -5 crnly, linearly decreasing downwards. The graph shows the temperature rise throughout the cornplete ice-rock colurnn after the indicated time has elapsed. The right curve of each pair (Open squares) corresponds to an experiment where heat conduction in the bed is not considered.
whenever a warm wave is conducted into the bed, temperature gradients will decrease in the rock, in turn decreasing the amount of geotherrnal heat flowing through the ice-rock interface, thereby counteracting the original warm wave. The effect on the englacial temperature distribution is illustrated in more detail in fig. 4.4.