Francesco Dell'isola

Kolumban Hutter

What are the dominant thermomechanical processes in the basal sediment layer of large ice sheets?

Keywords: glacier-ice sheet flow, sliding on soft bed, soil-water mixture, creeping flow, thermodynamic balances, interface theory

The till-layer below the temperate base of an ice sheet is modelled as a saturated viscous solid-fluid mixture. Its solid mass grows by abrading material from the underlying rock-bed, and the water balance is affected by the draining of water into the rock-bed and melting of ice at the interface of the overlying ice sheet. This latter process is influenced by the dissipation due to the power expended by: (i) the stresses; and (ii) Darcy-type interaction forces; as well as (iii) the sliding of the ice over the till-layer and of the latter over the rock-bed. All these mechanisms are likely to influence the thermodynamics of the layer, a fact which is corroborated by a scaling analysis.

Introduction

In the article Irregular oscillations of the West Antartic Ice Sheet (WAIS) [START_REF] Macayeal | Irregular oscillations of the West Antartic Ice Sheet[END_REF] writes that concern that the WAIS might collapse in response to future climate warming is motivated by the geological evidence suggesting its collapse some time during the last million years. He explains that his model computations suggest that, during the last million years, sporadic, perhaps chaotic collapses of the WAIS may have occurred. He relates this irregular behaviour to the long relaxation times of the sediment, the so-called till, that is responsible for the amount of lubrication provided to the bottom boundary of the ice sheet. His computations also indicate that a possible prediction of the collapse of WAIS under a future Greenhouse warming crucially depends on the past history of the thermomechanical state of the till and in particular the amount of water present in the layer.

The very simple existing models for the 5-10 m thick water saturated sediment layer beneath WAIS assume that the sediment-water mixture is nonlinearly viscous with a fluidity that increases with increasing interstitial (pore) pressure. MacAyeal's analysis, that was preceded by simpler models of Alley et al. (1987a, b), [START_REF] Blankenship | Till beneath ice stream B. 1. Properties derived from seismic travel times[END_REF], [START_REF] Rooney | Till beneath ice stream B. 2. Structure and continuity[END_REF], [START_REF] Engelhardt | Physical conditions at the base of a fast moving Antarctic ice stream[END_REF] essentially agrees with that of Kamb (1991) in that it accounts for the following processes. (1) The till obeys a layer integrated mass balance in which the layer-parallel tillflux divergence is balanced by the time rate of change of the layer thickness and the sediment production rate at the rock-bed.

(2) The dilatant pore volume is filled by the melt water from the ice sheet above.

(3) The melting is supported by the geothermal heat and the dissipation (mainly due to the shearing of the till in the layer).

(4) The sediment concentration itself affects the amount of frictional dissipation between the ice and the till-layer at their common boundary.

This list of dominant physical processes and their realization in a model is impressive, but falls short in: (i) not clearly presenting a thermomechanical model for the saturated mixture of sediment and water; (ii) not formulating appropriate boundary conditions for mass, momentum and energy at the ice-till-layer interface and the till-layer-rock-bed interface; and (iii) not deducing within the framework of the presented model that the above listed dominant processes are the only relevant ones. For instance, it is a priori not clear why the Darcy-interaction force should not contribute to the dissipation, and why the sensible and latent heat flows as well as the diffusively convective constituent energies across the ice-sediment interface should not affect the melting of the basal ice.

In the present paper a careful analysis is performed in order to strongly support the supposition that the last-mentioned processes, ignored by the previous analyses, are likely to be of comparable importance as those accounted for. The emerging nonlinear model will, as a dynamical system, most likely behave in a more complex fashion, and could as such give rise to either more or else fewer disintegrations of the WAIS than MacAyeal's original model did. In short, the question, whether the WAIS may be stable or become unstable under a future Greenhouse warming, is not yet settled. Indeed, we move in this paper only a first step towards a complete answer of the above question. However, we hope to make clear that a more detailed analysis than so far performed is essential to a proper understanding; future papers will attempt consecutively more complete solutions.

Let the system under consideration consist of three disjoint subregions which possess partly common boundaries as shown in figure 1. The ice sheet is polythermal, i.e. it consists of cold regions with temperatures below the melting point and temperate regions with temperatures at the pressure corrected melting point. This ice sheet rests on a layer of sediments, usually fine-grained minerals of abraded rock, called till, of which the pore volume is filled with ice where the basal ice is cold, but filled with water where it is temperate. The two regions are separated by a phase change surface. At cold contact points, the ice adheres to the interface, at temperate ones it slides over it; the corresponding dissipation will contribute to the melting of the basal ice. The soft sediment layer itself is bounded from below by a hard rock-bed. Along the common boundary the sediment will be assumed to slide, abrade till from the rock-bed and thereby generate heat. The geothermal heat from below, the dissipation due to sliding at the two interfaces and due to the deformation rate of the sediment-water layer will all contribute to the melting rate at the base of the ice sheet. A simple model of bedrock sinking will also be incorporated.

We shall be concerned here only with that portion of the ice sheet that reaches the melting point at the interface and the layer below it †. The 'ice' will be considered as Figure 1. Sketch of a polythermal ice sheet, consisting of cold and temperate ice regions and resting on a thin layer of sediment and water or ice. This layer itself rests on a rock-bed. The existence of the sediment-water layer is restricted to those regions where the basal ice reaches the melting temperature. Below cold portions the layer consists of sediments and ice. a mixture of pure ice plus water in which two mass balances for the water and the mixture as a whole, one momentum balance and one energy balance for the entire mixture will be formulated. The 'ice' is treated as a nonlinear viscous fluid according to [START_REF] Hutter | Thermomechanically coupled ice sheet response, cold, polythermal, temperate[END_REF] or Greve (1997a, b), and properties are regarded as known. The sediment-water mixture is a two-phase saturated continuum, in which both phases are nonlinear heat conducting interpenetrating continua with constant true densities, for which a consistent thermodynamic theory was presented by [START_REF] Svendsen | On the thermodynamics of a mixture of isotropic materials with constraints[END_REF]; resulting equations are quoted without derivation. The rock-bed, finally, is a heat-conducting deformable body of approximately 5 km thickness. Its treatment is not directly related to the problem of the till-layer, but to accomodate for its thermal inertia, heat conduction is usually accounted for, while dissipation due to deformation is ignored. The response of the combined lithosphere-asthenosphere to the ice load is grossly taken into account by a relaxation-type relation for the bedrock sinking with a relaxation time in the order of a few thousand years. Its role is to vertically adjust the position of the ice sheet in order to locate it in the appropriate thermal environment of the climate.

In what follows, we shall in § 2 outline the field equations for the saturated binary mixture modelling the till-layer consisting of two density preserving constituents. The analysis will be limited to constitutive models for a linearly viscous fluid and a nonlinearly viscous isotropic solid with an interaction force consisting of a quasistatic and a Darcy-type dynamic contribution. The governing field equations turn out to be nine scalar partial differential equations for the volume fraction of the solid, the saturation pressure, the constituent velocity fields and the common temperature of the mixture.

Section 3 presents the boundary conditions at the two interfaces which separate the sediment-water layer from the ice sheet above, and the rock-bed below. The transition conditions involve both kinematic equations that describe the motions of the respective interfaces and jump conditions of mass, momentum and energy. In these jump conditions account is taken of the production of water due to melting of ice at the top interface and production of till due to abrasion of rock at the bottom interface. Sliding plays a significant role, for both the mechanical and thermal fields. The ice is sliding over its bed and the sediment-water layer over its rock-bed; sliding laws are proposed that make the corresponding dissipation quadratically dependent upon the jump of the tangential velocity component of the mixture. Adjusting the drag coefficient accordingly makes both perfect sliding and no-slip possible limiting scenarios. The rates of the abraded mass of, and the drainage of, water into the rock must also be phenomenologically prescribed. They are given here as functional relations depending on the shear and normal tractions exerted at the sedimentlayer-rock-bed interface. Finally, bedrock sinking is parametrized by postulating that the normal speed of the bottom interface into the rock-bed is proportional to the deviation of its actual position from the position it would have when overlain by the actual ice sheet for an infinitely long time. The constant of proportionality is the inverse of a relaxation time.

Section 4 provides an order of magnitude analysis of the various terms arising in the field equations and jump conditions. At first, estimates for the various process quantities arising in the theory are given. They are suggested by 'physical intuition' and 'derive' from experience with other related models, but provide information about possible values of phenomenological coefficients (such as drag). Then, using these estimates as well as the known and tabulated values for the physical constants, it is shown that in the jump conditions of energy at the two interfaces all terms, i.e. the jump in internal energy times the mass flow through the interface, the jump in energy flow and the stress power, are of the same order of magnitude. In particular, the power expended by the stresses on the diffusive motion and the diffusively convected energy cannot be ignored in comparison to the conductive heat flow. Furthermore, in the field equation for the internal energy, the mechanical contributions to the energy flux and the power of the partial stresses contribute to the energy balance as does the conductive heat. Thus, the energy equation is not simply a heat conduction equation.

Section 5 is devoted to the construction of the simplest steady shear flow problem, the intention being to gain a first qualitative physical insight into the performance of the model, but perhaps equally important to test the mathematical adequacy of the proposed field equations and boundary conditions. For a sediment-water layer of constant thickness we prescribe a steady and uniform traction, temperature and energy flow from above and the geothermal heat flow from below; we ignore abrasion but incorporate a water flow from above. These assumptions still contain some ingredients of an interaction with the ice above the sediment-water layer, yet are simple enough to make the emerging problem amenable to analysis. On the basis of a plane flow assumption we determine the cross layer profiles of the solid volume fraction, saturation pressure, sediment and water velocities in the layer-parallel and layer-normal directions, respectively, and demonstrate that the water flow from above, which crucially determines these fields, is itself determined by a second-order differential equation for the temperature, whose integration is subject to the top and bottom thermal boundary conditions. This shows that it is the thermal regime at the base that governs the sediment layer flow.

In § 6 we briefly summarize our findings and indicate further applications.

Field equations

The till-layer, situated at the bottom of the ice sheet, is modelled here as a saturated structured mixture of a granular material whose pore volume is filled with water. The following fields define this mixture. ν s , ν f volume fractions of solid, fluid ρ s , ρ f , ρ mass densities of solid, fluid and of whole mixture v s , v f , v velocities of the solid, fluid, barycentre T s , T f , T Cauchy stresses of solid, fluid and mixture f s , f f specific body forces of solid, fluid m s , m f mixture interaction forces of solid, fluid ε s , ε f , ε internal energy for solid, fluid, whole mixture q s , q f , q energy fluxes for solid, fluid, whole mixture ϑ common temperature for solid and fluid constituents r specific energy supply D = sym grad v barycentric strain rate D a = sym grad v a constituent strain rate.

To describe its behaviour, two balance laws of mass and linear momentum for each constituent and one energy balance for the mixture as a whole are needed. In the local form they read

∂ρ s ∂t + div(ρ s v s ) = 0, (2.1) ∂ρ f ∂t + div(ρ f v f ) = 0, (2.2) ρ s vs = div T s + ρ s f s + m s , (2.3) ρ f vf = div T f + ρ f f f + m f , (2.4) ρ ε = -div q + tr(T D) + ρr, (2.5)
in which production rates of mass have been assumed to vanish as no mass exchange between the constituents occurs.

As is common in mixture theory, the mixture fields can be related to the peculiar fields of the constituents. Thus,

ρ = a ρ a , ρv = a ρ a v a =⇒ D = a [ξ a D a + sym(grad ξ a ⊗ v a )], T = T I - a ρ a u a ⊗ u a , ε = a (ξ a ε a + 1 2 ξ a u 2 a )
,

q := q I + a {T T a u a -ρ a (ε a + 1 2 u 2 a )u a },                 
(2.6) where the superscript (•) T denotes transposition, a = s, f for the solid and the fluid, respectively, and quantities indexed by I, (•) I , are the inner parts of the respective field quantities, i.e. the sum of the constituent quantities. Moreover, the notation

ξ a := ρ a /ρ, a ξ a = 1 , constituent mass fractions, u a := v a -v, a ξ a u a = 0 , constituent diffusion velocities,          (2.7)
is employed, and it is required that the interaction forces conform with Newton's third law so that m s = -m f =: m, (2.8) and that the sources are f s = f f = 0, r= 0, (2.9) implying that boundary loads are the only driving force mechanisms and radiation is absent. Finally, it will also be assumed that the constituents do not exchange internal angular momenta so that the constituent Cauchy stress tensors may be assumed to be symmetric T a = T T a .

(2.10) In what follows the constituent solid and fluid components will be assumed to be true density preserving, i.e. their true mass densities ρs , ρf are supposed to remain constant under deformation: consequently (2.11) This also implies that the peculiar densities are no longer field variables as they can be replaced by the corresponding volume fractions. We also impose the saturation condition, which states that the pores are always filled by the water so that we have

ρ s = ν s ρs , ρ f = ν f ρf .
a ν a = 1. We shall henceforth use ν := ν s , ν f = 1 -ν.
Imposing the saturation condition is tantamount to reducing the independent fields by one variable. In the theory by [START_REF] Svendsen | On the thermodynamics of a mixture of isotropic materials with constraints[END_REF] this apparent incompatibility is remedied by introducing a constraint field, the saturation pressure, p, as a new variable. The fields ν, p, v s , v f and ϑ are then regarded as the (nine) independent fields characterizing the state of the considered mixture, and constitutive equations are formulated for C := {T s , T f , m, ε s , ε f , q I , ψ, η} in which ψ and η are, respectively, the Helmholtz free energy and entropy. The constituent materials are assumed to be heat conducting, viscous and isotropic so that

C = C(ν, p, v s -v f , ϑ, grad ϑ, grad ν, D s , D f ).
The following explicit expressions are in conformity with the model developed by [START_REF] Svendsen | On the thermodynamics of a mixture of isotropic materials with constraints[END_REF].

(1) Constitutive equations for partial stresses:

T s = -ν(β s + p)I + t * s , T f = -(1 -ν)pI + t f , t f = 2µ f (ν)D * f + 1 3 λ f (ν) tr(D f )I, D * s = A s (ν)F (II(t * s ))t * s ,
(2.12)

where

β s = β 0 s + β 1 s βs (ν) + β 2 s (ϑ -ϑ 0 )
(2.13) and β i s (i = 0, 1, 2) are constants. (2) Constitutive equation for interaction force:

m = [p + (1 -ξ s )β s ] grad ν + m d (v f -v s ), m d = ν(1 -ν)ρ f g/K.
(2.14)

(3) Constitutive equation for internal energy:

ε ε I = (1-ξ s )ε f +ξ s ε s = (1-ξ s )ε f (ϑ 0 )+ξ s ε s (ϑ 0 )+[(1-ξ s )c f +ξ s c s ](ϑ-ϑ 0 ). (2.15)
(4) Constitutive equation for the heat flux:

q I = -(k f + k s ) grad ϑ, (2.16) in which β s = ρ ∂ψ I ∂ν , η = - ∂ψ I ∂ϑ , ψ I = ψ I (ϑ, ν), F(II(t * s )) = k + k 1 (II(t * s )) (n-1)/2 , n > 1,
(2.17) and where for any generic tensor S we used the notation II(S) = 1 2 tr(S 2 ), (S) * = S -1 3 (tr S)I.

(2.18)

Note that by prescribing ε f (ϑ) and ε s (ϑ) as well as β s , ε I (ϑ, ν) is known, and the inner free energy can, if so desired, be deduced from

∂ψ I ∂ν = β s ρ , ∂ ∂ϑ 1 ϑ ψ I = - ε I (ϑ, ν) ϑ 2 ,
simply by integration. In the above, the spherical stress contributions are the thermodynamic equilibrium stresses, while t s and t f are the corresponding viscous non-equilibrium contributions to the partial stresses; t s has no bulk viscous part, t s = t * s , but t f may have. In particular, the effective fluidity A s (ν)F (II(t * s )) has been multiplicatively decomposed into a volume-fraction-dependent creep-rate factor A s and a stress-dependent creep response function F . This last function is a linear combination of powers of II(t * s ), (see (2.17) 2 ) with k and k 1 being the Newtonian and nonlinear rheological moduli, n being the power exponent which is approximately between three and five †. Moreover, † For order of magnitude arguments, it is sometimes advantageous to use a Newtonian constitutive relation for the viscous part of the partial stress t * s = 2µsD * s where the shear viscosity µs is interpreted as an effective viscosity. λ f and µ f are the bulk and shear viscosities of the fluid, both of which may also depend on the volume fraction of the solid. As evident, the solid viscous stress does not entail dissipation in dilatant deformation while that of the fluid does. To ignore it we shall henceforth set λ f = 0. There is experimental support for the postulate of a polynomial constitutive relation for the dissipative stress in till-layer deformation. Vulliet & Hutter (1988) have demonstrated this for soils under gliding conditions of slope movements and [START_REF] Boulton | Sediment deformation beneath glaciers: rheology and geological consequences[END_REF] have analysed the till found in the front area of a glacier in Iceland. Power or polynomial creep behaviour is also suggested by [START_REF] Clarke | Subglacial till: a physical framework for its properties and processes[END_REF].

The constitutive relation for the interaction force also consists of an equilibrium contribution-the term involving grad ν-and a non-equilibrium term, here taken in the simplest possible form of a Darcy-type force which is linear in the relative velocity with the coefficient m d being proportional to the porosity (1 -ν) and inversely proportional to the permeability K [m s -1 ]. The equilibrium term is given as deduced by [START_REF] Svendsen | On the thermodynamics of a mixture of isotropic materials with constraints[END_REF] via the entropy principle. The coefficient of the grad ν term is additively decomposed into the saturation pressure and the thermodynamic term derived from the free energy ψ I , which is given here as a function of ν and (ϑ -ϑ 0 ) with phenomenological coefficients β i s (i = 0, 1, 2). We will choose β 2 s = 0, since expected temperature variations are very small. The inner part of the internal energy is expanded about a reference temperature ϑ 0 (the melting temperature at normal pressure, i.e. 273.15 K); thus c f and c s are the heat capacities of the fluid and the solid constituents, respectively, which in principle can depend on the volume fraction ν. Finally k s and k f are the thermal conductivities of the solid and fluid, respectively.

Physically the above system has two volume dissipation mechanisms: the first is due to the viscous part of the partial stresses and is contained in the term tr(T D) in the energy equation, the second is due to the interaction force m which spends power on the relative velocity field v s -v f . In the energy equation this dissipation is accounted for by the term a div T a • u a contained in the non-thermal contribution to the divergence of energy flux q (recall also balances of linear momenta).

In realistic motions arising during creeping flows in the considered till, both the solid and fluid velocity fields are small at all times. As a consequence the accelerations in both momentum equations (2.3), (2.4) can be neglected and Stokes conditions prevail. In conformity with this, nonlinear contributions in the diffusion velocity (quadratic and cubic terms) can be dropped, implying that

T T I = a T a , ε ε I = a ε a , q q I + a {T T a u a -ρ a ε a u a }. (2.19)
The mixture stress and the mixture internal energy are in this approximation given by their inner parts, but in the energy flux the purely thermal contribution is complemented by the power of the partial stresses on the diffusive motion and the diffusively convective transport of internal energies.

With the above-mentioned reductions implemented, the balance laws of mass, momentum and energy take the forms

∂ν ∂t + div(νv s ) = 0, - ∂ν ∂t + div((1 -ν) v f ) = 0, div(-ν(β s + p)I + t s ) + [p + (1 -ξ s )β s ] grad ν + m d (v f -v s ) = 0, div(-(1 -ν)pI + t f ) -[p + (1 -ξ s )β s ] grad ν -m d (v f -v s ) = 0, ρ ε = -div q + tr(T D).              (2.20)
Complemented by the constitutive relations for t s , t f , q I (equations (2.12), (2.16)) and considering equation (2.19) these are a set of nine scalar partial differential equations for the unknown fields ν, p, v s , v f and ϑ. When complementing them by suitable boundary and initial conditions and suitable conditions determing the motion of the boundaries of the region occupied by the sediment the resulting initial boundary-value problem should be well posed.

In summary, the model equations involve the Stokes flow assumption of creeping motion for the solid and the fluid; the solid is nonlinearly viscous with polynomiallaw rheology but no bulk viscous response, and the fluid is linearly viscous with shear and bulk viscosity, (the latter later being ignored). The non-equilibrium interaction force is simply Darcy's law. The equilibrium stresses are isotropic with saturation pressure distributed among the constituents according to 'pressure equilibrium', but the solid equilibrium stress also contains a thermodynamic pressure. Similarly, the thermodynamic equilibrium part of the interaction force which is proportional to the solid-volume-fraction gradient involves both the saturation pressure and a thermodynamic pressure derivable from the free energy. Finally, owing to the Stokes flow assumption the nonlinear diffusive contributions to the mixture stress, mixture energy flux and mixture internal energy are ignored.

Boundary conditions

The above specified binary mixture models the till-layer beneath an ice sheet. Thus it is bounded from above by the ice which is assumed to reach the melting point (i.e. to be temperate) at the base. The interface between the sediment-water layer and the ice sheet is assumed to be material for the sediment but non-material for the water. As lower boundary for the till layer we assume a rock-bed, across which water may be drained, because of the likely presence-in the rock-of interconnected cracks; along this interface fine grained till may in addition be abraded via the viscous sliding of the sediment-water mixture. Both these interfaces can also be moving due to the lithosphere-asthenosphere deformation imposed by the weight of the ice. Thus, the lower bounding interface is neither material for the sediment nor for the fluid.

(a ) The ice-till-layer interface Let (x, y, z) be Cartesian coordinates; (x, y) are horizontal and z is vertical opposite to the direction of gravity. Suppose that the ice-till layer interface, henceforth called the top interface Ω t , can uniquely be described by an equation z = f t (x, y, t).

(3.1)

Then the kinematic equation for its motion is given by

∂f t /∂t 1 + |∇f t | 2 + w t • n = 0, (3.2)
in which the vector field w t is defined by

w t (x, y) := v s (x, y, f t (x, y, t)); (3.3)
n is the unit normal to Ω t pointing into the ice and ∇ denotes the two-dimensional horizontal gradient operator. Equation (3.3) identifies the velocity of the interface with that of the sediment component at the interface; this assumption ignores diffusion of sediment into the ice.

Since the interface is a singular surface, across which several fields may suffer (finite) jumps, jump conditions of constituent mass density and mixture total momentum and energy must hold as follows:

[[ρ f (v f -v s ) • n]] = P t f , [[T n]] -[[ρv(v -v s ) • n]] = 0, [[ρ(ε + 1 2 v 2 )(v -v s ) • n]] + [[q • n -vT • n]] = 0.      on Ω t . (3.4)
In these equations, (3.3) has been used for identifying w t ; this is why the jump condition for the mass of the sediment is identically satisfied. Furthermore, we have only written the jump condition of momentum for the mixture as a whole, since its formulation for the sediment alone is not meaningful. Thirdly, P t f on the right-hand side of (3.4) 1 is the production rate of water at the interface Ω t due to melting of ice from the ice sheet above the interface. Defining by

M ± t := ρ ± f (v s -v ± f )
• n, the flow of fluid mass from the ice sheet onto the interface (plus sign) and the corresponding fluid mass flow away from it into the layer (minus sign), (3.4) 1 simply reads

M - t = M + t + P t f (3.5
) and states that the surplus flow of fluid mass into the till-layer must be due to the melting processes at the interface.

A more detailed interpretation of these processes can be gained from the jump conditions of total mass

[[ρ(v -v s ) • n]] = -[[M t ]] = 0, M t := -ρ(v -v s ) • n. (3.6)
This statement includes now the ice and water from the ice sheet, the water and the sediment in the layer (which however does not contribute to the formula). On the negative side we have

ρ -(v --v s ) • n = ρ - f (v - f -v s ) • n = -M - t , (3.7)
and on the positive side we obtain

ρ + (v + -v s ) • n = ρ + i (v + i -v s ) • n + ρ + f (v + f -v s ) • n = -P t f -M + t , (3.8)
where we defined

ρ + i (v + i -v s ) • n =: -P t f , ( 3
.9) and where ρ + i and v + i denote the ice density and ice velocity immediately above the interface.

In the theory of polythermal ice of [START_REF] Hutter | Thermomechanically coupled ice sheet response, cold, polythermal, temperate[END_REF] and Greve (1997a) the following decomposition is used:

-M + t = j + • n + ρ + ξ + f (v + -v s ) • n, (3.10)
where

j + := ρ + f (v + f -v)
is the vector of moisture flux and ξ + f := ρ + f /ρ + is the moisture content of ice. The fluid mass flow M + t into the interface thus consists of two contributions, a diffusive and a convective flow of moisture. In summary, (3.7), (3.8) prove that (3.6) implies (3.5), and we can write (3.12) As the theory of polythermal ice only operates with the total stress T + , a postulate is needed to distribute this stress between the partial stress tensors T - s and T - f . This is achieved by the following postulate.

M t = M - t = M + t + P t f . ( 3
Postulate 3.1. Let T -n be the total traction exerted on the mixture of sediment and water at the interface; let, moreover, ν h be the solid volume fraction at the upper interface and α a real parameter with α ∈ [0, (1/ν h ) 2/3 ]. Then we require that

T - s n = αν 2/3 h T -n, T - f n = (1 -αν 2/3 h )T -n. (3.13)
For α = 1 the total stress is distributed among the constituents according to their volume fractions; α = (1/ν h ) 2/3 or α = 0 distributes it solely over the sediment or fluid traction, respectively. In reality α is a constitutive quantity to be prescribed. An alternative and less ad hoc approach would be to introduce constituent surface stresses or, equivalently, surface momentum productions and to postulate constitutive relations for these. This is the view taken by [START_REF] Hutter | On interfacial transition conditions in two phase gravity flow[END_REF] and [START_REF] Wu | On the role of the mechanical interaction in a gravity-driven shear flow of an ice-till mixture[END_REF]. The above postulate could be phrased in such a context. With this postulate and the condition (3.12), the following stress conditions are obtained:

αν 2/3 h T + n = -ν(β s + p)n + t s n, (1 -αν 2/3 h )T + n = -(1 -ν)pn + t f n, on Ω t .
(3.14)

Consider next the energy jump condition (3.4) 3 which we now write as (higherorder terms in the velocity are dropped)

[[ε I ]]M t = [[q • n]] -[[v • T n]].
(3.15)

The contributions on the left-hand side are †

ε + I = (1 -ξ + f )ε + i + ξ + f ε + f = ε + i + ξ + f (ε + f -ε + i ) = ε + i (ϑ 0 ) + c + i (ϑ -ϑ 0 ) + ξ + f L, (ε + f -ε + i ) =: L + (1/ρ f -1/ρ i )(n • T + n) ∼ = L, ε - I = (1 -ξ - s )ε - f + ξ - s ε - s = (1 -ξ - s )ε - f (ϑ 0 ) + ξ - s ε - s (ϑ 0 ) + [(1 -ξ - s )c - f + ξ - s c - s ](ϑ -ϑ 0 ),                 
(3.16) in which use has been made of the fact that the ice is at the melting point; ϑ 0 is the melting temperature at normal pressure and ϑ that at the local pressure. Furthermore, a distinction is made between ε - f and ε + f even though this is strictly not needed, since the fluid is the same material on both sides of the interface. Using the approximation ε - f = ε + i + L the above formulae imply:

[[ε I ]] = (ξ + f -1)L + ξ - s (ε - f (ϑ 0 ) -ε - s (ϑ 0 )) + [c + i -(1 -ξ - s )c - f -ξ - s c - s ](ϑ -ϑ 0 ), (3.17)
which together with M t defines the left-hand side of (3.15). In much the same way we may also show that

a ρ - a ε - a (v a -v) • n = ρ -ξ - s (1 -ξ - s )[ε - f (ϑ 0 ) -ε - s (ϑ 0 ) + (c - f -c - s )(ϑ -ϑ 0 )](v f -v s )
• n, (3.18) an expression which is used in the evaluation of [[q • n]]. Using the definition (2.6) 5 for q -and q + = q + therm + Lj + • n, (see Greve 1997a), we may readily deduce that

[[q • n]] = (q + therm -q - I ) • n + Lj + • n - a (T T a u a -ρ a ε a u a ) = (q + therm -q - I ) • n + Lj + • n + a ρ - a ε - a (v a -v) • n +{ν(β s + p)(v s -v) --t T s (v s -v) -} • n +{(1 -ν)p(v f -v) --t T f (v f -v) -} • n. (3.19)
The terms in the second line are the jump in the conductive heat flow, plus the corresponding jump of the diffusive latent heat flow, plus the diffusively transported peculiar energies explicitly given in (3.18). The terms in the next two lines express the corresponding peculiar stress powers on the diffusive normal velocities of the solid and the fluid.

With the definitions

τ := (T n) = (I -n ⊗ n)T n, σ := -(T n) ⊥ = -n • T n,
the jump of the stress power in (3.15) has the decomposition (3.20) into the dissipative power due to sliding and the power due to the jump in the normal component of the barycentric velocity. The former requires formulation of a sliding law σ), between the jump of the tangential component of the barycentric velocity and the shear traction with a drag coefficient C t which may itself be a function of the shear traction |τ | and the pressure σ. C t = 0 corresponds to continuity of the tangential component of the barycentric velocity, while C t = ∞ to perfect sliding (with vanishing shear traction at the interface). Following Calov & Hutter (1996a) and slightly adjusting their formula to the notation at hand, the power law representation

[[v • T n]] = [[v]] • T n = [[v ]]τ -[[v ⊥ ]]σ,
[[v ]] = C t τ , C t = Ĉ t (|τ |,
C t = C t (|τ |) pt ρg(σ) qt , p t = q t = 2, (3.21)
will be chosen; a parameter study, in which the optimal shape of the Greenland ice sheet was sought (Calov & Hutter 1996b) showed that C t was approximately 10 3 -10 4 [a -1 ] with an optimal value C t = 6 × 10 4 [a -1 ].

On the other hand, we have

[[v ⊥ ]] = [[v • n]] = [[(v -v s ) • n]] = - 1 ρ M t , so that [[v • T n]] = C t [[v ]] 2 + 1 ρ M t σ, (3.22)
is an alternative expression for (3.20). The second term on the right will be shown to be negligible in comparison to the first. This completes the formal analysis of the jump conditions at the top interface.

(b ) The till-layer-rock-bed interface Let z = f b (x, y, t), (3.23) be the equation of the till-layer-rock-bed interface, henceforth called the bottom interface Ω b , and let w b (x, y, t) be a vector field such that w b • n is its normal speed into the rock-bed, where n is the unit normal pointing into the rock-bed. Then,

∂f b /∂t 1 + |∇f b | 2 + v + r • n = (v + r -w b ) • n =: -a ⊥ =: - M b s ρ + r , (3.24)
is the evolution equation of the bottom interface; v + r denotes the material velocity of the point immediately below the interface in the rock-bed. For a motionless rock-bed v + r = 0 and w b • n = a ⊥ . Thus, it is apparent that a ⊥ is the abrasion rate, i.e. the volume per unit interface area of till torn away from the rock-bed base. Alternatively, M b s is the corresponding abrading mass and ρ + r is the mass density of abraded rock. A constitutive relation determining the abrasion rate will be given below.

As the bottom interface is a singular surface, the following jump conditions of mass:

[

[ρ f (v f -w b ) • n]] = 0, [[ρ s (v s -w b ) • n]] = 0, ( 3
.25) must hold, which, alternatively, can be written as

(1 -ν)(v - f -w b ) • n = M b f ρf , ν(v - s -w b ) • n = ρ + r ρs (v + r -w b ) • n = - ρ + r ρs a ⊥ = - M b s ρs .        (3.26)
Here M b f is the mass flow of fluid per unit area at the bottom interface, i.e., the drainage of the fluid into the rock-bed, also prescribed by a constitutive equation given below.

The jump conditions for stresses are not needed, as there will be no detailed analysis of the deformation of the rock-bed. However, the energy jump condition

[[ρ(ε + 1 2 v 2 )(v -w b ) • n]] + [[q • n]] -[[v]] • T n = 0, (3.27)
is again of significance †. Since no phase changes occur at this boundary, the first term may be ignored, so that the jump in heat flow is balanced by the stress power due to the jump in barycentric velocity. If it is assumed that the rock-bed is a diffusive mixture of a solid rock through which the water is transported by diffusion alone, † In this energy jump condition we use the fact that the second term in the jump condition for momentum,

[ [T n] ] -[ [ρv(v -w b )
• n] ] = 0, for the mixture is small so that T n can be regarded as continuous across the interface.

then-with appropriate interpretation-formula (3.19) also applies at the bottom interface, namely

[[q • n]] = (q + r -q - I ) • n + {ν(β s + p)(v s -v) --t T s (v s -v) -} • n +{(1 -ν)p(v f -v) --t T f (v f -v) -} • n + a ρ - a ε - a (v a -v) • n, (3.28)
in which q + r = -k r grad ϑ, is the conductive heat flux in the rock-bed and where there is no latent heat flux. Furthermore, the contribution due to diffusively transported peculiar energies (last term above) is given by (3.18).

The jump in stress power, [[v • T n]], is decomposed into tangential and normal contributions as in (3.20). With the sliding law

[[v ]] = C b τ , C b = Ĉ b (τ , σ),
and the representation

[[v ⊥ ]] = 1 ρ (M b f -M b s ) = - 1 ρ M b , M b := (M b f -M b s ),
we obtain, analogously to (3.22),

[[v • T n]] = C b [[v ]] 2 + 1 ρ M b σ. (3.29)
The second term on the right-hand side can again be ignored in comparison to the first. Because the adjacent materials at the interface are rock and gravel, it may not be justified to take C b as a constant. A power-law representation

C b = C b 1 ρg (|τ |) p b (σ) q b , (3.30)
with constant C b and suitably chosen exponents p b and q b may be appropriate; alternatively, a stick-slip law of the form

C b =    0, |τ | < τ 0 (σ), C b 1 ρg (|τ |) p b (σ) q b , |τ | τ 0 (σ), (3.31)
where τ 0 (σ) is a threshold shear traction, which itself is a function of normal pressure, may account for the friction effects. We shall only consider the power law (3.30) and follow Calov & Hutter (1996b) in choosing

p b = 2, q b = 2, C b = λC t , λ ∈ [5, 10 2 ].
This guarantees that, for a given stress state, sliding of the till-layer over the rockbed is λ-times slower than that of the ice over the till-layer. Summarizing these specializations, the energy jump condition reduces to

[[q • n]] = C b ρg |τ | σ 2 [[(I -n ⊗ n)v]] 2 .
This formulation of the bottom interface conditions would be complete, were it not for the abrasion rate and drainage functions which must be prescribed by constitutive relations. It is reasonable to suppose that M b f and M b s are both functions of the shear and normal tractions

M b s = M s (|τ |, σ), M b f = M f (|τ |, σ f
), and the power laws

M b s = m b s (|τ |) m s 1 (σ) m s 2 , M b f = m b f (|τ |) m f 1 (σ f ) m f 2 , (3.32)
where σ f is the normal fluid pressure, with constant coefficients m b s and m b f and constant exponents m s 1 , m s 2 , m f 1 , m f 2 , may be appropriate; if the stick-slip law (3.31) is used, then the above power law for M b s must correspondingly be replaced by

M b s = 0, |τ | < τ 0 (σ), m b s (|τ |) m s 1 (σ) m s 2 , |τ | τ 0 (σ). (3.33)
In this study we choose the exponents in (3.32) as

m s 1 = m s 2 = 2, m f 1 = 0, m f 2 = 4, (3.34)
(intuitively they must be positive, and certainly m s 1 > 1, m s 2 > 1, m f 2 > 1 are expected behaviours). With these choices we shall-in § 4-give estimates about the coefficients m b s and m b f . This completes the formulation of the model equations for the sediment-water layer.

(c ) Treatment of the rock-bed

As it is well known, the combined lithosphere-asthenosphere deforms under the weight of an ice sheet and may respond to this weight in a relaxation-type fashion; i.e. the base of a former ice sheet may still move upwards many thousand years after the disintegration of the ice sheet. We propose here the simplest possible model for its motion. To this end, let z = f 0 b (x, y) be the equation of the relaxed steady lithosphere without ice load, and let ρ a be the mass density of the asthenosphere fluid. If we envisage that the resistance to bedrock sinking is only due to the buoyancy force between the diving lithosphere into the asthenosphere, then

f 0 b - ρ i ρ a H
is the steady z-coordinate of the rock-bed interface in the presence of an ice sheet with total thickness H. The difference

f b -f 0 b - ρ i ρ a H(x, y, t)
is a measure of how much the actual position of the bottom interface deviates from its steady-state position. It will now be assumed that the velocity, v + r , of a material point immediately below the bottom interface is vertical and its magnitude is proportional to this distance with relaxation time τ v (of approximately 5000 years); thus

v + r = - 1 τ v f b -f 0 b - ρ i ρ a H(x, y, t) e z ,
or when combining this with the kinematic equation (3.24),

∂f b ∂t + 1 τ v f b -f 0 b - ρ i ρ a H(x, y, t) = - M b s ρ + r 1 + |∇f b | 2 . (3.35)
As the thermal response of the rock-bed to climate changes is about three times slower than that of ice (see Calov & Hutter 1996a, b), it is necessary to solve in a rock layer of a few km thickness the energy balance equation (3.36) subject to the boundary (transition) conditions

ρ r c r dϑ dt = ρ r c r ∂ϑ ∂t + ∂ϑ ∂z v + r = k r ϑ,
ϑ = ϑ b , at z = f b (x, y, t), -k r ∂ϑ ∂z = Q geoth , at z = f r = const.,    (3.37)
where f r is the lower position of the rock-bed (approximately 5 km below the bottom interface). This completes the formulation of the boundary conditions, which are summarized by the following statements:

(1) at the top interface: the kinematic relations (3.2), (3.3); the jump condition of fluid mass (3.13) with specifications (3.11), (3.12); the continuity requirement of the total traction (3.14), its distribution among the solid and fluid constituents according to (3.15); the jump condition of energy (3.17 (3) in the rock-bed: the kinematic equation (3.35) and the heat conduction problem (3.36), (3.37).

Scalings: relative orders of magnitude of the different dissipative phenomena

In this section estimates are given for the various different dissipative terms arising in the formulation of the theory outlined in § § 2 and 3; the intention is to see which of these numerous expressions are likely to contribute in a non-negligible fashion to the overall behaviour of the till-layer. Of the four items listed in § 1 we find the first two points are automatically taken into account by the formulation developed so far. In particular, the melting rate of ice and the production rate of sediment at the top and bottom interfaces are parametrized and, through their kinematic equations, the layer thickness evolves in conformity with the mass balances of the sediments and the water. More critical are the last two items listed in § 1, i.e. an estimation of the dominant dissipative mechanisms; how much ice is melted and pressed as water into the till-layer is determined by these terms.

Table 1 lists numerical values for physical constants, and table 2 gives estimates for process quantities arising in the formulation, which we shall now discuss in some greater detail. The moisture content in temperate ice is known to be only a few percent, typically 1% but no more than 2% [START_REF] Hutter | Thermomechanically coupled ice sheet response, cold, polythermal, temperate[END_REF]Greve 1997a, b). The solid volume fraction in the till deposition cannot exceed 0.85 and is probably seldom below 0.4 (G. K. C. Clarke, personal communication); thus, mass fractions ξ s are between 0.7 and 0.9. The volume fraction of the rock-bed is basically unknown, but can vary according to the amount of cracks present in the rock. Here, it will be assumed that the dynamical effects of the presence of water in the rock-bed can be ignored (ξ f = 0) except for the drainage motion which operates as a sink to the water content in the till-layer.

As for the melting temperature at the base of the ice sheet, note that it depends on pressure via the Clausius-Clapeyron gradients (β), which at 3000 m depth yields -2.6 • C for its value. The geothermal heat flow for pre Kambrian rock is 42 mW m -2 , but it is known that values as high as 55 mW m -2 yield better results in ice sheet modelling of Antarctica [START_REF] Greve | Polythermal three-dimensional modelling of the Greenland Ice Sheet with varied geothermal heat flux[END_REF][START_REF] Hansen | Application of a polythermal ice sheet model to the Antarctic Ice Sheet: steady-state solution and response to Milankovic cycles[END_REF][START_REF] Huybrechts | The Antartic Ice Sheet and environmental changes: a three-dimensional modelling study[END_REF]). On the other hand, the conductive heat flow from the top interface into the ice is practically zero when the temperate ice has finite thickness, but may be of the order of the geothermal heat flow when the ice is cold. The heat generated by the sliding of the basal ice over its bed depends on the sliding velocities and the frictional coefficient, C t , which itself may vary owing to the roughness of the bed. Taking it as given in table 1, Calov & Hutter (1996b) obtained for the Greenland ice sheet maximum values of the frictional heat Q fr more than ten times the geothermal heat flow, which via the relation a ⊥ = Q fr /(Lρ) corresponds to approximately 4 mm ice that is melted per year.

The stress state at the base can best be estimated by setting the ice pressure at the top interface equal to ρ i gH, the over-burden pressure of a column of ice with thickness H, and the corresponding shear traction equal to ρ i gH H , where H is the surface slope of the ice sheet. This yields, for H ≈ 3500 m, pressures as high as 3 × 10 7 Pa and shear tractions not larger than 10 6 Pa. Tractions up to this order must be transmitted through the till-layer to the bottom interface. They give us a means to estimate the coefficients m b s and m b f in the phenomenological laws for the abrasion rate and the drainage function. In fact, substituting in formula (3.32) 1 for the quantities σ and |τ | the above two estimates and assuming m s 1 = m s 2 = 2, as well as estimating the abrasion rate to be 1 mm a -1 and 1 m a -1 (two extrema) yields m b s = 10 -35 and 10 -31 (kg -3 m 2 s 7 ) as a lower and upper bound, respectively. Similarly, with m f 1 = 0, m f 2 = 4 and with m b f = 3 × 10 -31 to 3 × 10 -28 (kg -3 m 2 s 7 ) and the above stresses, drainage velocities are of order 10 -3 to 1 m a -1 .

As for the interaction forces, it is known from [START_REF] Wu | Schwerkraftgetriebene Scherströmungen in gestättigten Binärmischungen nicht-Newtonscher Fluide mit Anwendungen auf das Fließen von Sedimentverschmutztem Eis[END_REF] that β s should be negative and its absolute value small as compared to the saturation pressure. Setting the latter equal to the over-burden pressure yields β 1 s 10 7 Pa, justifying the estimate listed in table 1. On the other hand, the Darcy coefficient m d may be obtained from well-know values of soil permeabilities (see, for example, Vulliet & Hutter (1988) and table 1) yielding values from 1.0-5.0 × 10 10 (kg m -3 s -1 ). Finally, the constitutive relation for the sediment stress tensor must be specified. In a first approximation we shall assume that A s (ν) is constant; the work of Vulliet & Hutter (1988) then suggests the values in table 1.

We next turn our attention to the energy jump condition (3.15) at the top interface. The dominant term of [[ε I ]] in (3.17) is the latent heat of ice (ξ + f -1)L 330 × 10 3 J kg -1 with a comparable somewhat-smaller contribution from the term

ξ - s (ε - f (θ 0 ) -ε - s (θ 0 )).
The term in (3.17) involving the specific heats is about two orders of magnitude smaller, so that [[ε I ]] 5 × 10 5 J kg -1 . Estimating the melting rate of ice and the flow of moisture from the ice side to be of the order of Table 1. Physical parameters used in the model, taken in part from Greve (1997b) and [START_REF] Buntenbarth | Eine Einführung in die allgemeine und angewandte Wärmelehre des Erdkörpers[END_REF], with our extensions quantity symbol value units gravity acceleration g 9.81 (m s -2 ) universal gas constant R 8.314 [35,55] mW m -2 heat flow from the ice into top layer [0, 100] mW m -2 heat generated by the ice sliding over its bed [0, 1000] mW m -2 basal ice melted away [0, 5] mm a -1 basal pressure, (∼ ρigH)

(J mol -1 K -1 ) density of ice ρi 910 (kg m -3 ) water ρf 1000 (kg m -3 ) sediment ρs 2700 (kg m -3 ) rock-bed ρr 2700 (kg m -3 ) asthenosphere ρa 3300 (kg m -3 ) specific heat of ice c i 146.3 + 7.253(K) (J kg -1 K -1 ) 2127.5 + 7.523( • C) water c f 4200 (J kg -1 K -1 ) sediment c s 700 ÷ 1000 (J kg -1 K -1 ) rock c r 750 (J kg -1 K -1 ) latent heat of ice L 335 × 10 3 (J kg -1 ) Clausius Clapeyron β 8.7 × 10 -4 (m -1 K) gradient of ice internal energy of water ε f (ϑ0) 750 × 10 3 (J kg -1 ) at 273.15 K internal energy of sediment εs(ϑ0) 150 × 10 3 (J kg -1 ) at 273.15 K heat conductivity of ice ki 9.828 exp(-0.0057ϑ (K)) (W m -1 K -1 ) water k f 0.552 (W m -1 K -1 ) rock k r 3.0(∼ 2.1 ÷ 5.5) (W m -1 K -1 ) sediment k s 2.5 ( W m -1 K -1 ) dynamic shear viscosity µ f (1.4 ÷ 2.0) × 10 -3 (kg m -1 s -1 ) of water (at 273.15 K) dynamic bulk viscosity λ f (neglected) (kg m -1 s -1 ) of water (at 273.15 K) representative apparent µs 1.0 × 10 7 (kg m -1 s -1 ) viscosity for sediment fluidity coefficient As 1.2 × 10 -23 ((K Nm -2 ) -9.6 ) for sediment k 0.7 × 10 -10 (K N -1 m 2 s) k 1 1.0 - n 9.6 - thermodynamic pressure coefficient |β 1 s | 10 5 ÷ 10 7 (Pa) soil permeability K 10 -7 ÷ 5 × 10 -6 (m s -1 ) Table 1. Continued quantity symbol value units abrasion coefficients m b s 3 × 10 -31 ÷ 3 × 10 -28 (kg -3 m 2 s 7 ) m s 1 2.0 - m s 2 2 . 0 - drainage coefficients m b f 3 × 10 -31 ÷ 3 × 10 -28 (kg -3 m 2 s 7 ) m f 1 0 - m f 2 4 - friction coefficient C t 2 × 10 -3 (s -1 ) (top interface) friction coefficient C b = λCt λ 5 ÷ 100 - (bottom
[ 0 , 3 × 10 7 ] kg m -1 s 2 basal shear traction [0, 10 6 ] kg m -1 s 2 abrasion rates [2.7, 2700] kg m -2 a -1 abrasion velocity [10 -3 , 1] m a -1 drainage rates [1, 1000] kg m -2 a -1 drainage velocity [10 -3 , 1] m a -1 2 × 10 -7 kg m -2 s -1 , the left-hand side of (3.15) yields [[ε I ]]M t = 5 × 10 -2 W m -2 .
If this term contributes to the energy jump condition, then [[q • n]], the conductive contribution, is of the order of the geothermal heat flow, 5 × 10 -2 W m -2 . This would be exact if no motion takes place, and steady state conditions prevail. The heat flow into the ice is negligible in comparison to this, because of the very small value of the Clausius-Clapeyron constant in the temperate ice. To estimate the contribution from the moisture flux, Lj • n, let the moisture diffusion velocity be 1 m a -1 = 0.32 × 10 -7 m s -1 ; this yields for Lj • n the order of magnitude of 10 -1 to 10 0 W m -2 which is not negligible in comparison to the conductive contribution.

A typical member of the stress or pressure dependent term in the expression (3.17) for the heat flow jump is νp(v s -v) • n. Using ν 0.5, p 2 × 10 7 Pa and an upper bound for the diffusion velocity u s • n 1 m a -1 , we obtain |νpu s • n| ∈ [10 -3 , 10 -1 ] W m -2 , indicating that the omission of these terms is not justified either. Using the same value for the diffusion normal speed u s • n and applying the values of ε a as listed in table 1, it may be concluded that

a ρ a ε a u a • n ∈ [10 -2 , 10 0 ] W m -2 ,
again a non-negligible contribution. It follows that none of the terms that define the energy flux vector q of the mixture ought to be ignored.

Next, to extimate the dissipation rate [[v • T n]], we restrict this expression to the sliding motion [[v ]]T n and substitute values 10 4 -10 6 Pa for the shear traction T n , and 1-10m a -1 for the velocity difference due to sliding, [[v ]]; with these,

[[v • T n]] ∈ [ 1 3 10 -3 , 1 3 ] W a -2
, making the stress power equally significant. It follows that all terms in the energy jump balance (3.15) contribute to this balance with comparable magnitudes; only the terms in (3.17) involving the specific heats and the jump [[1/ρ]]M t σ in (3.22) may safely be ignored.

A similar analysis of order of magnitude estimates can also be performed with the energy jump equation (3.27) that applies at the bottom interface. Because no phase changes occur, the contribution on the left-hand side is negligibly small in comparison to those terms arising on the right. Thus we have the approximate balance [[q

• n]] = [[v]] • T n.
The geothermal heat flow from below is approximately 5 × 10 -2 W m -2 . On the layer side, since the dominant gradients are across the interface, there is no evident reason to ignore any one of the non-conductive energy flux terms in the respective definition (2.19) 3 . Thus, the same stress powers on the diffusive motions and the diffusively convected energy transports are likely to contribute non-negligibly to the energy flow jump, as do those at the top surface. Furthermore, as abrasion occurs, its contribution to the generation of heat is non-negligible and produces heat of similar, if not higher, order of magnitude, as does the frictional sliding motion of the ice over its bed.

The above analysis shows that a detailed account of all the terms in the energy jump conditions is significant, if a proper thermal treatment is to be achieved. What about the local energy balance equation (2.5)? As the dominant heat flow component is transverse to the sediment-water layer, assume it to be between one and ten times the geothermal heat flow and take a layer thickness of five to ten metres to estimate its gradient; thus | div q| ∈ [10 -3 , 10 -1 ] W m -3 . This must compare with the dissipation rate tr(T D), with dominant term τ γ, where τ is a representative shear stress and γ the corresponding shearing. With τ 10 6 Pa and γ ∼ v/H, where v ∼ 100 m a -1 and H ∼ 5-10 m this yields τ γ ∈ [10 -2 , 10 0 ] W m -3 , which is comparable to the value of | div q|.

This completes the order of magnitude analysis of the energy equations. It demonstrates that an adequate treatment of the saturated sediment-water layer beneath an ice sheet must be a binary mixture of a saturated mixture in which proper energy balance must involve the local balance law plus the energy jump conditions at the top and bottom interfaces.

Steady shear flow

Consider this saturated soil-water mixture subject to steady loading conditions from above, assume the upper and lower interfaces to be plane and their distance h to be constant in space and in time, ignore for the moment abrasion of till from the lower interface and suppose the horizontal velocity of the upper edge to be imposed via the application of shear traction and to be unidirectional. Assume, moreover, processes to be uniform in the horizontal direction except for the constrained pressure p. If we then choose x, y as the horizontal axes and z vertical, we have ∂(•)/∂x = ∂(•)/∂y = 0 except for p. It suffices to consider plane flow in the direction e 1 (horizontal) and e 3 (vertical) so that the peculiar velocity fields take the forms (a = f, s) v a = u a e 1 + w a e 3 , u a = u a (z), w a = w a (z),

(5.1) from which we readily deduce

D * a = 1 2 u a (e 1 ⊗ e 3 + e 3 ⊗ e 1 ) -1 3 w a e 1 ⊗ e 1 -1 3 w a e 2 ⊗ e 2 + 2 3 w a e 3 ⊗ e 3 , II(D * a ) = 1 4 (u a ) 2 + 1 3 (w a ) 2 , (5.2)
in which primes on u a , w a denote differentiation with respect to z.

With the above suppositions, the balance laws of mass (2.20) 1,2 take the forms

(νw s ) = 0, ((1 -ν)w f ) = 0, (5.3)
or upon integration νw s = V so , (1 -ν)w f = V fo , (5.4) in which V so and V fo denote constants of integration, the vertical volume flow of the solid and the fluid, respectively. Because the top interface is stationary and ν = 0 there, equation (5.4) 1 implies V so = 0; it follows, as ν = 0 ∀z ∈ [0, h], that w s = 0 ∀z ∈ [0, h].

Given the specific form (5.2) 1 of the stretching deviator D * a the materially dependent stress contribution, since t * a is affine to D * a , can be written as t * a = τ a (e 1 ⊗ e 3 + e 3 ⊗ e 1 ) + σ xa e 1 ⊗ e 1 + σ ya e 2 ⊗ e 2 + σ za e 3 ⊗ e 3 , (5.5) and imposing the constitutive relations (2.12) yields

σ xf = σ yf = -2 3 µ f w f , σ zf = 4 3 µ f w f , σ xs = σ ys = σ zs = 0, τ a = µ s u a .
(5.6) Balances of momentum now imply the following.

In the horizontal direction: (5.7) in which µ a are effective viscosities and m d is given in terms of the permeability in (2.14). From these it follows that ∂p/∂x is only a function of z.

-ν∂p/∂x + (µ s u s ) + m d (u f -u s ) = 0, -(1 -ν)∂p/∂x + (µ f u f ) -m d (u f -u s ) = 0,
In the vertical direction:

(-ν(β s + p)) + (p + (1 -ξ s )β s )ν = -m d w f , (-(1 -ν)p) -(p + (1 -ξ s )β s )ν + 4 3 (µ f w f ) = m d w f .
(5.8)

Adding the two equations implies

(-p -νβ s + 4 3 (µ f w f )) = 0 ⇒ p = 4 3 µ f ν (1 -ν) 2 V fo -νβ s + p(x),
(5.9) where p = ∆p x + p 0 and ∆p is the assigned pressure gradient and p 0 is to be determined by the pressure boundary condition at the top interface. Substituting this result back into one of equations (5.8) allows elimination of the pressure and yields, after some manipulations, a second order ODE for the function ν(z), for which the following first integral can be found:

ν = 3(1 -ν) 2 4µ f V 0f [2β 1 s F (a, ν) -k + αV 0f z],
(5.10)

where a := ρf /ρ s , α := ρf g/K (5.11) and

F (a, ν) = (1 -ν) ν + a 2(1 -a) ν + a 1 -a dν = -1 2 ν 2 + 1 2 a -2 a -1 ν - 1 2(a -1) 2 a ln |νa -ν -a|;
(5.12) k is a constant of integration. The solution to the first-order differential equation (5.10) requires determination of an additional constant of integration: for instance ν 0 , the value of the solid volume fraction at z = 0: ν 0 = ν(0). Once this constant is assigned, the ODE (5.10) admits for every k a unique solution ν(z; ν(0), k). In particular, we may define the function Φ(ν(0), k) = ν(h; ν(0), k) whose unique value is the solid volume fraction at z = h: ν h = ν(h). The values ν 0 = ν(0), ν h = ν(h) and the pressure p 0 can be expressed in terms of boundary tractions by imposing the stress boundary conditions (3.14), in which the basal normal traction of the ice is distributed between the normal tractions of the solid and water stresses, governed by the parameter α and by imposing condition (3.29) 2 at the bottom. Considering equations (3.14) and (3.32) at z = h and z = 0, respectively, recalling that the constant k must be determined by a shooting technique to attain the required ν(h), and evaluating equation (5.10) for z = h, we obtain the following system of equations for ν (h), ν 0 , ν h and p 0

ρ i gH(ν h -αν 2/3 h ) = β 1 s ν h (ν h -1), p 0 = ρ i gH, (ν -1)ρ i gH + β 1 s βs ν(1 -ν) + 4 3 V fo µ f ν (1 -ν) 2 = 4 ρf V fo /m b f , Φ(ν 0 , k) = ν h ,     
(5.13) recalling that the parameter α is defined by (5.11), and is not related to α. The last of equations (5.13) requires determination of a general solution of the family of ODEs (5.10) parametrized by k. Such a solution does not seem to be available in closed analytical form, however we can determine a field ν(•) which very closely approximates the solution to the ODE (5.10) and verifies equations (5.13).

Indeed, we start by remarking that when a ∈ (0, 1) then the function F (a, ν) with ν ∈ [0, 1] can be very closely approximated by a linear function of ν. This simply follows by direct evaluation. From now on we will thus assume

F (a, ν) F (a, ν h ) + ∂F ∂ν (ν h )(ν -ν h ).
(5.14)

Using this approximate expression in (5.10), differentiating the result with respect to z, introducing the dimensionless variable ζ := z/h and the new dependent variable y = 1/(1 -ν) we obtain the ODE

1 b y (ζ) + C b y (ζ) y 2 (ζ) -1 = 0, (5.15) in which b = 2αh 2 µ f , C = 3hβ 1 s 2µ f |V 0f | .
(5.16)

The estimates obtained in § 4 clearly show that equation (5.15) defines a singular perturbation problem; indeed, using the numerical values from tables 1 and 2 we find that |C| ≈ 10 20 , b ≈ 10 16 . A simple scaling analysis of (5.15), employing matched asymptotic expansions, shows that ν is to a very high degree independent of z (i.e. constant) except in an extremely thin boundary layer which can be situated only at the top interface. The above equations allow for the determination of all the constants needed for determining the pressure and volume fraction in the till-layer.

Next, consider the horizontal momentum equations (5.7) and assume in a first qualitative analysis that the viscosities µ f and µ s are constant, i.e. independent of stretching. Then, multiplying (5.7) 1 with µ f and (5.7) 2 with µ s and subtracting the emerging equations from one another and, in turn, adding the two equations yields

(u f -u s ) -m d 1 µ s + 1 µ f (u f -u s ) + 1 -ν µ f + ν µ s ∆p = 0, -∆p + (µ s u s + µ f u f ) = 0,    (5.17)
of which the second possesses the solution (5.18) where u 0 and u 0 are constants of integration. Equation (5.17) 1 is a second-order ODE for the difference velocity. Through m d it depends on the solution for ν(z). Two further constants of integration will enter, so that equations (5.17) comprise four such constants to be determined by, e.g. the no-slip conditions at the bottom interface and the shear part of (3.12). Locally, (5.17) 1 exhibits exponential behaviour, while (5.17) 2 is of Couette-type. So both constituent horizontal velocities will consist of the superposition of Couette and exponential-type behaviour. More realistically, the soil viscosity is stretching dependent. Such a nonlinear influence can iteratively be accounted for by evaluating µ f and µ s for this linear solution and subsequently reintegrating (5.7). This application, as well a more detailed analysis of the behaviour of the layer, is reserved to a future paper [START_REF] Dell'isola | A qualitative analysis of the dynamics of a sheared and pressurized layer of saturated soil[END_REF]). Alternatively we could have solved equation (5.7) subject to viscous sliding conditions at the bottom interface. Strictly, for this analysis to be reasonable, abrasion must be ignored in this case, as otherwise w s = 0.

µ s u s + µ f u f = u 0 + u 0 z + 1 2 ∆pz 2 ,
The solid-volume fraction and the horizontal velocity components u s and u f that are determinable with the aid of (5.17) depend on the vertical fluid volume flow V f0 that has so far been assumed to be prescribed. This latter quantity is now obtainable from an energy analysis: more specifically, consider the balance of energy (2.20) 5 .

Lengthy algebraic computations in which the prerequisites of the considered steady flow are incorporated show that the energy equation takes the following form:

V f0 (1 -ν) E(ν, ϑ, ϑ ) = (k s + k f )ϑ + U(µ a , u a , ν) + V f0 (1 -ν) D(ν, ϑ, ϑ , p), (5.19)
where

E(ν, ϑ, ϑ ) := ρf (1 -ν)[ξ s ∆ε + ϑ (ξ s c s + ξ f c f )], ∆ε := ξ s [ε s (ϑ 0 ) -ε f (ϑ 0 )] + ξ s [(c s -c f )(ϑ -ϑ 0 )], U(µ a , u a , ν) := a,b {[µ a u a ( 1 2 u a -ξ b u b )] + (ξ b u b ) µ a u a }, D(ν, ϑ, ϑ , p) := 2ξ s (2νβ s + p) + (νβ s + ν 2 β s + ν p) , +ρ f (1 -ν)c f ϑ + ξ s (ρ f (1 -ν)ε f (ϑ, ϑ 0 ) + ρs νε s (ϑ, ϑ 0 )), -ξ f (ρ f (1 -ν)ε f (ϑ, ϑ 0 ) + ρs νε s (ϑ, ϑ 0 )) + ρf ν ε f (ϑ, ϑ 0 ).                        (5.20) Equation (5.19
) is a second-order nonlinear ODE for ϑ in which E, U and D are known functions of the arguments as indicated in equation (5.20). As with the original energy equation the left-hand side of (5.19) expresses the rate of change of internal energy; the first term of the right-hand side of (5.19) accounts for the conductive heat whereas the remaining two terms express the energy transport on the diffusive motion as well as the dissipation. The equation is inhomogeneous through the appearance of the known function U(µ a , u a , ν) but it contains the yet still unknown parameter V f0 . This parameter will be determined along with the solution of (5.19), as there are three boundary conditions at our disposal; in fact these are: at the bottom interface: for no-slip [[q • n]] = 0, (5.21) as given in (3.27) in which q + • n = -Q geoth ; at the top interface:

ϑ = ϑ m , -[[ε I ]]V f0 (1 -ν) -ρf = [[q • n]] + [[v • T n]] , (5.22) [1] [2] [3]
where ϑ m is the melting temperature of the ice at the pressure of the interface and in which [2] is given in (3.19), [3] is given in (3.22) and [1] equals -M t ; notice that [2],

[3] are functions of V f0 via the solution of the previous ordinary differential equations. We view the above as a two-point boundary-value problem with ϑ = ϑ m given at the top and Q geoth given at the bottom interface and V f0 obtained as an output of the solution using (5.22) 2 .

To close this section we state that the analysis of this problem has not only been presented to obtain a qualitative physical understanding of the deformation of the sediment layer, but equally to show that the mathematical problem for the simplest set of physical approximations retaining key features is well posed.

Concluding remarks

The focus of this paper is the derivation and justification of a continuum thermodynamic theory for the relatively thin layer of sediment and water below the temperate basis of a large ice sheet. In the glaciological jargon, it is thus a contribution to the physics of soft glacier beds and attempts to identify their dominant governing processes.

The layer is modelled as a binary mixture of a viscous heat-conducting fluid and an isotropic viscous heat-conducting solid with stresses and interaction force given, in equilibrium by the saturation pressure and contributions from a free energy, and in non-equilibrium by additional contributions representing Newtonian or polynomial as well as Darcy-type behaviour. The boundary conditions at the non-material top and bottom interfaces include the kinematic equations describing their evolution and the mass, momentum and energy jump conditions for the individual constituents and/or for the mixture as a whole. Melting of ice at the top, and abrasion of rock at the bottom interface are the two physical mechanisms that add fluid and solid mass to the till-layer and thus determine, in conjunction with mass balance and kinematic equations of the surfaces, the till-layer thickness. On the other hand, the sliding of temperate ice over the top interface, and of the sediment-water mixture over the rock-bed comprise dissipative mechanisms that contribute to the energy jump conditions: (i) at the bottom surface thereby enlarging the heat flow from the rock-bed into the till-layer; and (ii) at the top interface and thus contributing to the melting rate of ice. Within the layer, heat is conducted along its temperature gradient (primarily transversely to the layer) and convected by the specific internal energy and the stresses on the diffusive motion; the latter corresponds to the power expended by the Darcy interaction force on the constituent relative velocity. Within the layer, heat, however, is also generated by the dissipation due to the power of the stresses. All these processes have been shown to be likely significant in an adequate thermomechanical description of the layer physics.

Even though the proposed binary-mixture model with two interfaces at the bottom and the top is complicated, a steady-state analysis of a shear flow within a layer of constant thickness has been analytically amenable, thus giving at least a partial proof of the usefulness of the model. On the basis of the supposition that the peculiar traction and the heat flow into the ice are prescribed at the top interface and that the geothermal heat flow is prescribed, the across-layer solid-volume fraction profile and the profiles of the along-layer sediment and water velocities were determined, as was the across-layer water velocity and temperature profile, which are coupled through the Stefan-type boundary conditions. The above analysis and numerical computations show that ν is essentially constant while u s and u f are linear combinations of parabolic and exponential profiles. This is a significant result because in earlier and much simpler works of MacAyeal (1992), [START_REF] Boulton | Sediment deformation beneath glaciers: rheology and geological consequences[END_REF] and Alley et al. (1987a, b) a linear barycentric velocity profile was supposed ab initio.

The shear flow problem is believed to capture the essential physics at least as far as the geometry of the profiles of the solid volume fraction and the velocity components are concerned. They will give us guidelines on how to treat the full time-dependent problem. Nevertheless, despite its physical transparency, the above theory is rather formidable: a complete three-dimensional (numerical) procedure for solving the coupled ice-sheet-till-layer-rock-bed problem is not warranted, because the layer with 5-10 m thickness falls below the mesh size of usual numerical integration procedures of ice sheet modelling. So, we wish to think of the layer, together with its top and bottom boundaries that connect it to the ice sheet and the rock-bed as a single (ideally) infinitely thin interface, that is equipped with surface fields, which, in their entirety, replace the layer and boundary fields. Equations to be derived for these fields will

  Dedicated to Professor Peter Haupt on the occasion of his sixtieth birthday. Proc. R. Soc. Lond. A (1998) 454, 1169-1195 c 1998 The Royal Society Printed in Great Britain 1169 T E X Paper

Figure 2 .

 2 Figure 2. Close-up of the basal region of the ice sheet where a layer of temperate ice exists. Definition of a Cartesian coordinate system and identification of the top and bottom interfaces as z = ft(x, y, t) and z = f b (x, y, t), respectively.

  .11) Turning next to the momentum jump condition (3.4) 2 , we neglect the term [[v]]M t in comparison to [[T n]], implying that the total traction must be continuous across the interface [[T n]] = 0.

  ) with specifications (3.19), (3.20), (3.21) and (3.23), (3.24); (2) at the bottom interface: the kinematic relation (3.27); the jump conditions of solid and fluid mass (3.29) with the specifications (3.32)-(3.34); the energy jump condition (3.30) with the specifications (3.31), (3.32)-(3.34);

Table 2 .

 2 Process quantities and their approximate ranges of numerical values

	interface)
	time lag for bed adjustment	τv	9 × 10 10 ÷ 15 × 10 10	(s)
	(relaxation time)
	process quantity		range of values
	moisture content of ice, ξ f		[0, 0.02]
	volume fraction of sediment, ν	[0.5, 0.85]
	mass fraction of sediment, ξ s		[0.7, 0.9]
	water content in rock, ξ f		[0, 0.001] (negligible)
	melting temperature at base, ϑ	[273.15, 270] K
	geothermal heat flow, Q geoth

† This means that we assume no transition of a sediment water layer into a sediment-ice layer with a phase change surface separating the two. This more general situation can be treated, but for simplicity will not be dealt with here.Proc. R.Soc. Lond. A (1998)

Proc. R. Soc. Lond. A (1998)

† Note that for a phase interface of a compressible fluid[ [ε] ] = L -[ [1/ρ] ]p ⊥ , where [ [1/ρ]] is the jump in specific volume and p ⊥ the pressure perpendicular to the interface. Proc. R.Soc. Lond. A (1998)

F. dell'Isola and K. Hutter

form the new basal interface conditions that must connect the ice sheet variables with those of the rock-bed. It is hoped that such a two-dimensional description is simpler to implement in a numerical code, but still able to capture the essential physics outlined in this paper. Our own past experience with interface theories [START_REF] Alts | Towards a theory of temperate glacier. Dynamics and thermodynamics of phase boundaries between ice and water[END_REF], 1988a, b, 1989;[START_REF] Dell'isola | Deduction of thermodynamic balance laws for bidimensional nonmaterial directed continua modelling interphase layers[END_REF] gives us reason for optimism.

While performing this work K.H. was supported by a Grant from Università di Roma LA SAPIENZA. Both authors thank Professor C. Gavarini, Director of the Dipartimento di Ingegneria Strutturale e Geotecnica for securing that Grant.

We thank Professor L. W. Morland and Dr H. Blatter for their reviews of an earlier version of this paper and Dr Y. Wang for his help with T E Xing the manuscript.