Skip to main content Accessibility help


  • Access


      • Send article to Kindle

        To send this article to your Kindle, first ensure is added to your Approved Personal Document E-mail List under your Personal Document Settings on the Manage Your Content and Devices page of your Amazon account. Then enter the ‘name’ part of your Kindle email address below. Find out more about sending to your Kindle. Find out more about sending to your Kindle.

        Note you can select to send to either the or variations. ‘’ emails are free but can only be sent to your device when it is connected to wi-fi. ‘’ emails can be delivered even when you are not connected to wi-fi, but note that service fees apply.

        Find out more about the Kindle Personal Document Service.

        Channelized subglacial drainage over a deformable bed
        Available formats

        Send article to Dropbox

        To send this article to your Dropbox account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you use this feature, you will be asked to authorise Cambridge Core to connect with your <service> account. Find out more about sending content to Dropbox.

        Channelized subglacial drainage over a deformable bed
        Available formats

        Send article to Google Drive

        To send this article to your Google Drive account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you use this feature, you will be asked to authorise Cambridge Core to connect with your <service> account. Find out more about sending content to Google Drive.

        Channelized subglacial drainage over a deformable bed
        Available formats
Export citation


We develop theoretically a description of a possible subglacial drainage mechanism for glaciers and ice sheets moving over saturated, deformable till. The model is based on the plausible assumption that flow of water in a thin film at the ice-till interface is unstable to the formation of a channelized drainage system, and is restricted to the case in which meltwater cannot escape through the till to an underlying aquifer. In describing the physics of such channelized drainage, we have generalized and extended Röthlisberger’s model of channels cut into basal ice to include “canals” cut into the till, paying particular attention to the role of sediment properties and the mechanics of sediment transport. We show that sediment-floored Röthlisberger (R) channels can exist for high effective pressures, and wide, shallow, ice-roofed canals cut into the till for low effective pressures. Canals should form a distributed, non-arborescent system, unlike R channels. For steep slopes typical of alpine glaciers, both drainage systems can exist, but with the water pressure lower in the R channels than in the canals; the canal drainage should therefore be unstable in the presence of channels. For small slopes typical of ice sheets, only canals can exist and we therefore predict that, if channelized meltwater flow occurs under ice sheets moving over deformable till, it takes the form of shallow, distributed canals at low effective pressure, similar to that measured at Ice Stream B in West Antarctica. Geologic evidence derived from land forms and deposits left by the Pleistocene ice sheets in North America and Europe is also consistent with predictions of the model.

Mathematical Symbols

  • α Exponent in till flow law

  • A i Ice flow-law parameter

  • A s Till flow-law parameter

  • b Exponent in till flow law

  • b 1,b 2,b 3 Numerical constants

  • c Average suspended-sediment concentration as mass fraction

  • Cp Specific heat of water

  • C t Rate of change of melting temperature with pressure

  • d Ice thickness

  • D B “Characteristic” size of suspended sediment

  • D 50 Median grain-size of sediment

  • D Deposition rate of suspended sediment

  • Ē Rate of suspended-sediment erosion and entrainment

  • f R Dimensionless stream-bed roughness factor

  • g Acceleration due to gravity

  • h Stream depth

  • [h] Magnitude of stream depth in metres

  • h T Heat-transfer coefficient

  • H Conduit spacing transverse to ice-flow direction

  • H t Thickness of till

  • k w Thermal conductivity of water

  • K Dimensionless constant in bed load-transport equation

  • K i,K S Dimensionless constants related to conduit closure

  • K 1,K 2 Dimensionless constants related to suspended-sediment transport

  • l Perimeter of subglacial drainage conduit

  • l f Characteristic distance along glacier

  • l i Part of l bounded by ice

  • l s Part of l bounded by till

  • L Heat of fusion of ice

  • Rate of meltwater input to basal till layer

  • Rate of melting of ice at tunnel wall

  • Rate of sediment erosion on conduit floor

  • Rate of sediment input from tributaries

  • Rate of meltwater input from tributaries

  • n Exponent in flow law for ice

  • Ν Pi - Pw

  • N c Pi - Pc

  • ΔN C Variation in N c

  • [N c] Magnitude of N c in bar

  • Nu Nusselt number

  • P c Water pressure in subglacial tunnel

  • P i Ice-overburden pressure

  • p w Water pressure in till far from a conduit

  • Δp w P w - P c

  • Characteristic value of N c separating channel and canal solutions

  • Q Water flux

  • [Q] Magnitude of Q in m3 s −1

  • Q b Bed-load flux

  • Q S Suspended-load flux

  • R Conduit radius for idealized semi-circular geometry

  • Re Reynolds number for flow in subglacial conduit

  • R f Dimensionless constant related to suspended-sediment transport

  • s Curvilinear co-ordinate measuring distance along conduit

  • S Conduit cross-sectional area, equal to S i+S s

  • S i Conduit area bounded by ice

  • S s Conduit area bounded by till

  • t Time

  • T Celsius temperature

  • ΔT Temperature drop across boundary layer between water and ice at wall of conduit

  • u Mean speed of water flow in conduit

  • u b, Glacier-sliding speed

  • v S Particle-settling speed

  • V i Rate at which ice flows into conduit

  • V S Rate at which sediment flows into conduit

  • α Slope of glacier surface

  • αB Suspended-sediment concentration as volume fraction

  • γ ρwCpCt = 0.32

  • έ Shear strain rate

  • θ Slope of glacier bed or stream bed

  • θ is Angles characterizing conduit shape

  • K W Thermal diffusivity of water

  • ^ Dimensionless measure of till permeability

  • Dimensionless parameter related to critical stress for bed-load transport

  • ρ i Density of ice

  • ρ s Density of sediment

  • ρ w Density of water

  • Δρ ρs - ρw

  • T Average shear stress exerted by stream on boundary of conduit

  • Tb hear stress exerted by glacier on its bed

  • Tc Critical value of T to mobilize bed-load

  • 〈〉 Designates average over tunnel cross-section

1. Introduction

When seeking to understand the dynamics of ice sheets and glaciers, it is necessary to pose a sliding condition for the basal velocity when the basal ice reaches the pressure-melting point. This relates the velocity ub to the basal stress T b via a sliding law. It is now well recognized that this sliding law, which gives u b as a function of T b, also depends on the basal water pressure through the effective pressure N, defined as the difference between the overburden pressure p i and the water pressure p w

. (1.1)

This dependence has been elucidated both theoretically (Lliboutry, 1968) and observationally (Iken and Bindschadler, 1986). On the assumption that such a sliding law is available, it is then necessary to provide a prescription for N and this requires a theoretical model for the subglacial drainage pressure.

The classical theory of subglacial drainage is due to Röthlisberger (1972) and Shreve (1972), and envisages channels cut into the ice, and maintained open (with positive N, which promotes channel closure through ice creep) by melt-back due to the frictional heating by the turbulent water flow. This theory applies to hard beds, that is to say, the substrate over which the ice moves is considered rigid and impermeable. The same assumption is made in the classical theory of sliding. In the last several years, however, a radically different viewpoint has arisen. Observations of subglacial till deformation by Boulton and Hindmarsh (1987), Blake and Clarke (1989), Kohier and Proksch (1991), Blake (1992) and Humphrey and others (1993), along with the recognition that the rapid movement of Ice Stream B in West Antarctica may be due to shearing of weak subglacial sediment (Blankenship and others, 1986, 1987; Alley and others, 1986; Engelhardt and others, 1990), have caused a realization that the traditional theory of glacier sliding over an undeformable bed is of questionable relevance in some circumstances.

In this paper, we develop a model for subglacial drainage for the case in which the ice is underlain by a layer of sediment, which for convenience we will term till. We modify the Röthlisberger analysis to show that, in addition to channels carved into ice, “canals” can exist which are cut into the till, and which are maintained open by fluvial removal of inwardly creeping till. In the following section, we discuss the physical assumptions and ingredients which such a model requires. In section 3, we present a mathematical model for channelized flow, then show how it may be approximately “solved” in section 4. Our conclusions are presented in section 5. A list of mathematical symbols used is also provided at the beginning of this paper.

2. Physics of Drainage Through Till

Before addressing a quantitative model of subglacial drainage, we must consider the form that this drainage takes. A number of different possibilities have been considered in the literature. Most simply, drainage could take place by evacuation through a subglacial aquifer or by channels incised into the ice or the bed (here the till), through a system of linked cavities, or through a water film of variable thickness.

In the present analysis, we assume that there is no aquifer beneath the till and that all drainage is through the till or along the bed. The case in which there is a sub-till aquifer has been considered by Boulton and Hindmarsh (1987), Clarke (1987) and Shoemaker and Leung (1987).

A number of earlier workers (e.g. Boulton and Jones, 1979; Alley and others, 1986; Brown and others, 1987) have shown that for plausible permeabilities and thicknesses of subglacial till, Darcian flow through the pore space of a subglacial till layer israrely sufficient to evacuate all the meltwater flowing at the glacier bed. In certain special cases, meltwater may flow englacially (Hooke and others, 1988) or, if the subglacial till is being sheared rapidly enough, be simply adverted along with the till (Alley and others, 1986, 1987); in general, however, meltwater must escape via some sort of localized drainage at the glacier sole.

We consider that there are three main likely types of drainage system for sediment-floored glaciers: conduits incised in the ice (R channels); conduits incised into the till, which we henceforth call “canals”; and a water film of variable thickness (Weertman, 1972; Alley, 1989).

We will not offer the last word on films versus conduits in this paper but a few comments are in order. It is possible to envisage a water film of variable thickness at the ice/till interface, and unpublished data of H. Engelhardt and B. Kamb (cited by Kamb (1991, p. 16,589)) are consistent with (but do not require) existence of such a film beneath Ice Stream B in West Antarctica. However, the existence of a variable-thickness film is threatened by an instability due to viscous dissipation in the water. Walder (1982) demonstrated that for a rigid bed this instability leads to incipient initiation of R channels. There is no reason to suppose the same instability would not exist in the case of deformable bed. Various heuristic arguments can be advanced to suggest that such an instability will not in fact be manifested until the film thickness exceeds several millimetres; the point is debatable. Weertman and Birchfield (1983) argued that, despite the film instability pointed out by Walder (1982), channelized flow itself is unstable if all meltwater is subglacially derived, due to the supposed inability of channels to collect water from large lateral distances. We are skeptical of this conclusion, because it depends on the assumption of straight and parallel channels, and it is more natural to suppose that a network channel will not be parallel. Moreover, Weertman and Birchfield’s argument, which relies on peculiarities of the stress field around a channel, applies only for rigid, impermeable beds. If the bed consisted of permeable sediment, there would be less impediment to lateral water movement to channels.

For a sediment bed, another instability should come into play: the tendency for rill formation, just as in subaerial overland flow (Smith and Bretherton, 1972; Loewenherz, 1991). Thicker parts of a film will exert a greater shear stress on the sediment bed than thinner parts and will thus be relatively more efficient at eroding the bed and transporting sediment. Again, flow localization occurs. Thus, there are two instabilities tending to produce channelized flow over a sediment bed. Accordingly,we focus in this paper on the characteristics of a channelized drainage system over a deformable bed.

3. A Mathematical Model

We consider the situation shown in Figure 1, where a conduit exists at the ice/till interface, in general with parts cut into both the ice and the till. The ice-overburden pressure (and thus also the far-field till pressure) is pi, while the water pressure in the channel is pc. In general, we expect pi > p c, and our aim is to find an expression relating the effective channel pressure N c = p ip c to the water flux Q through the channel, as well as to the effective pressure N = p ip w, where P w is now understood as the pore-water pressure in the till far from the conduit.

Fig. 1. Idealized subglacial conduit at the ice-till interface. Ice and till move at speeds Vi and Vs, respectively.

When N c < 0, the conduit tends to close by creep but is maintained open by melt-back due to viscous dissipation. In a similar way, it is recognized that wet till also creeps, with a yield stress on the order of 50 kPa. In many circumstances, the shear stress at the bed will far exceed this, so that the till also will tend to creep into the conduit. In this case, the till floor may be maintained by erosion and fluvial transport of sediment, either as suspended load or bed load. Our model is predicated on the effectiveness of these processes in maintaining a steady channel such as that in Figure 1. We now describe the several constituents of our model.

Creep of ice

The creep-closure rate of a tunnel in ice is given by


where R is the tunnel radius, and A i and n are the constants in Glen’s flow for ice

and T are the second invariants of the strain-rate tensor
and the stress deviator tensor Tij). This result applies directly to the closure of a cylindrical borehole (Nye, 1953) and is reasonably presumed to apply approximately in other geometric configurations.

Creep of till

Boulton and Hindmarsh (1987) proposed a flow law for till of the form



is the shear rate, T is the shear stress, N is the effective pressure P ip w and A s,a,b are positive constants. This law, based on a fit to seven data points, has been criticized by other authors (e.g. Kamb, 1991; Blake, 1992) but, in the present context, we consider it the simplest law which captures the idea that Equation 10 should increase with N, although the singularity at N = 0 is unphysical.

Boulton and Hindmarsh (1987) conjectured, following Nye’s (1953) solution for closure of a borehole in ice, that the consequent closure rate of a canal due to creep of till would be

where implicitly Nc = N, i.e. P c = p w With a = 1.33<b = 1.8 (as measured), sediment creep is then most effective at low values of N c, which has an important bearing on the stability of the canal system, as we shall discuss later.

More recendy, Fowler and Walder (1993) have solved the Nye borehole-closure problem for the case of saturated, deformable till where in general p cP w. They showed that the closure rate depends on the parameter ^ a dimensionless measure of permeability. Estimates of ^ are in the range 10-;6 < A < 1, for a permeability range 10-;19–10-;13m2 (Freeze and Cherry, 1979, p. 29) and an apparent till viscosity of 109Pas (cf. Blake and Clarke, 1989; Blake, 1992; Humphrey and others, 1993). For most tills, then, the asymptotic case

applies and the correct generalization of Equation (3.2). is as long as

, (3.4)

as long as

. (3.5)

If the condition Δpw is not satisfied, the analysis predicts the existence adjacent to the hole of a “failed” zone in which piping probably occurs and permeability is locally enhanced. An analogous situation arises subaerially where streams are recharged by ground water: as stage falls following a flood, the hydraulic gradient driving seepage through the stream banks may become great enough that piping and bank failure occur (e.g. Bradford and Piest, 1980). Fowler and Walder (1993) considered how the solution to the borehole-closure problem must be modified to account for this effect; they showed the “corrected” closure rate is not greatly different from Equation (3.4)., so for simplicity we assume Equation (3.4). is generally valid in the rest of our analysis.

Pore-water pressure

In general, we expect p w > p c: meltwater generated at the ice/till interface makes its way towards conduits (although the opposite situation could occur transiently if, for example, flooding occurs via an englacial waterway). If the till is very permeable, or the melt supply is low, then the water will be evacuated towards conduits through the till. If these are spaced a distance H apart and the till is of thickness Ht, then the excess pore-water pressure necessary to drain water supplied at a rate

per unit area of ice-till contact is
where k is permeability and μw is viscosity (Shoemaker, 1986). Taking H t = 10 m,
, μw = 2 x 10-;3Pas, k = 10-;16m2 (a median value), then μpw = 1 bar; for H = 100 m, μpw = 100 bar. In general, μpw/N may be significant; for example, under Ice Stream B, N is only about 1 bar (Engelhardt and others, 1990).

Various criteria may be offered for the initiation of conduits over a sediment bed (cf. Shoemaker, 1986). A physically reasonable criterion is that conduit spacing adjusts itself so that sediments at the edge of conduits are marginally stable against piping; from Equation (3.5)., this means Nc = aN. If conduits were closer, piping would occur, increasing the permeability of sediment near the conduit and thereby drawing down pore pressure pw at the drainage divides until the condition Nc = aN was satisfied. The sediment creep law would then become

. (3.6)

In subsequent analysis, we will adopt the form of Equation (3.6)..

Sediment erosion and transport

The inward creep of till into a conduit is balanced by erosion of the sediment and its subsequent transport along the channel. We base our discussion on that of alluvial channels in non-cohesive sediment and give only the basic details: in particular, the formulae we use are very rough approximations.

Erosion occurs only if the shear stress T exerted by the water flow is larger than some critical stress Tc, which for sorted sediment depends on grain-size Ds. Beyond this critical stress, particles can be dislodged and they are transported downstream either as bed load if they are large enough (Ds > 0.1 mm) or as suspended sediment.

Various empirical expressions are available for bed-load transport; we take, for simplicity, the common Meyer-Peter and Müller (1948) equation, which we write as

, (3.7)

where (a;)(x)+ = max(x, 0), Qb is the bed-load-transport rate (mass per unit time), K is a constant approximately equal to 10, ls is the width of the stream bed, g is acceleration due to gravity, ρw is the density of water and μρ = ρs - ρw.

Theoretical justification for expressions such as Equation (3.7). has been developed by Wiberg and Smith (1989). For stream beds with a mixture of sand-and gravel-sized clasts — probably a reasonable description for many subglacial streams (e.g. Engelhardt and others, 1990; Humphrey and others, 1993) — a condition of “equal mobility” holds, wherein all grain-sizes begin to move at about the same stress, namely,



is an empirically found, dimensionless constant, and D 50 is the median grain-size (Parker and others, 1982).

Expressions for the rate of erosion of suspendable sediment are difficult to come by. Parker (1978a) emphasized that this is pardy because the usual theory of suspended-sediment concentration (e.g. Richards, 1982, p. 101) is a local, steady-state theory in which erosion rate Ė and deposition rate Ḋ balance exactly. In principle, however, Ė may be prescribed as boundary conditions on the concentration profile. Parker (1978a) proposed that the erosion rate is given by


where K 1 ≈ 0.1 and v s is the particle-settling velocity, and that the deposition rate is


where c is the depth-averaged suspended-sediment concentration and K 2 ≈ 6. The net-erosion rate of suspendable sediment is then

. (3.11)

An important feature of subaerial stream beds is that they evolve their own shape through the processes discussed above. We can expect this to be true also for subglacial conduits cut into sediment. For non-cohesive sand/silt stream beds, Parker (1978a) argued that a stable bed shape could be maintained by a secondary lateral circulation of sediment, as follows. In the central parts of the stream, the shear stress acting on the bed is greatest and consequently greater erosion there leads to a larger suspended load than at the banks. This causes a lateral eddy diffusion of sediment towards the banks, which balances the bed-load transport that tends to bring sediment back down-slope toward the centre. Thus, even though the banks are being actively eroded, the river can maintain a constant cross-section via this dynamic equilibrium. In alluvial rivers, this leads to channel cross-sections with high aspect ratio (i.e. width l much greater than depth h). In gravel-bed rivers, the adjustment mechanism is more complicated. Since sediment is not being carried in suspension toward the banks, a stable channel cannot exist unless no erosion (and thus no bed-load transport) occurs on the banks. Accordingly, the shear stress must be great enough near the channel centre to transport the imposed sediment load but below critical on the banks. Parker (1978b) resolved this by showing that lateral turbulent flux of downstream momentum can produce the required shear-stress distribution. However, it is generally the case that T ≈ Tc for alluvial gravel-bed rivers. Again, large aspect ratios are typical.

Channels eroded into deforming subglacial till will differ from alluvial channels in one important aspect: whereas for an alluvial channel a steady-state cross-section exists only if there is no net erosion, for the subglacial channel a steady state will exist only if the net erosion balances the creep of till toward the channel. Thus, for a sand/silt-floored subglacial channel, there will not be a macroscopic balance between Ė and Ḋ, but rather


where the brackets

denote cross-sectional averages. For a gravel-floored subglacial channel, t will not drop below Tc on the banks, but rather T will exceed Tc by just enough to erode the banks and transport the inwardly creeping till at the appropriate rate. Accordingly, (Parker’s 1978a, b) analyses of the shapes of equilibrium stream channels cannot strictly apply to the subglacial case; in fact, his analyses should give lower bounds on the magnitude of bed shear stress. Qualitatively, however, the processes described by Parker for subaerial streams should apply to subglacial ones,hence we fully expect subglacial streams over non-cohesive sediment to tend towards high-aspect ratio cross-sections.

Governing equations

Having described the constituent ingredients of the process, we now write down governing equations for the flow. We consider a conduit at the ice–till interface (Fig.1). In general, we expect the conduit to extend both upwards into the basal ice and downwards into the till. We can thus define cross-sectional areas S i and S B of ice-or till-bounded conduits, respectively, and the total cross-sectional area S by

. (3.13)

The differential equations governing the physics are straightforward and are based on those elucidated by Röthlisberger (1972), Nye (1976), Spring and Hutter (1981) and Lliboutry (1983) but, here augmented by relations describing sediment creep, erosion and trans port. They must properly be framed in the context of a two-phase flow and we assume the simplest such model, that of homogeneous flow; in other words, we assume the sediment moves with the same velocity as the fluid. Let s be the downstream longitudinal coordinate and t be time. Equations of conservation of mass are, for the water,


and for the sediment



is the ice-erosion (i.e. melting) rate per unit downstream distance;
are, respectively, the water- and sediment-supply rates due to flow from tributaries, moulins or veins;
is the volume fraction of suspended sediment and Qa is the suspended-sediment flux, defined by


where u is the mean flow rate; other symbols are as defined previously. We include the bed-load transport rate Q b in Equation (3.15)., because we envisage the bed load as existing as a thin skin at the bed, occupying negligible volume.

Using Equations (3.2) and (3.6), the kinematics of conduit closure are described by

, (3.17)



where l i and l s are defined in Figure 1 and K i and K s and O(1) shape factors, partly representing the geometric difference between the conduit in Figure 1 and the idealized cylindrical shape assumed in the closure-rate calculations of Nye (1953) and Fowler and Walder (1993). Equations (3.17) and (3.18) are appropriate generalizations for non-cylindrical conduits. If the conduit is actually cylindrical, then

and the creep-closure rate
as expected. For a wide, flat channel in ice (cf. Hooke and others, 1990) with height h and width
(thus ll i), we expect
, because the average closure rate depends only on the width of the “hole” into which ice is flowing. A similar argument applies to a wide, shallow conduit incised into till.

Neglecting inertial effects, the momentum balance for the flow can be replaced by a force balance between pressure forces and drag on the channel walls,


where sin θ is the bed slope and l is the channel perimeter, defined by

. (3.20)

Assuming turbulent flow, which is valid except for very small tunnels (cf. Weertman, 1972), the average wall shear stress is given by


where the dimensionless friction factor fR is weakly dependent on the Reynolds number and is typically ~0.1 for subaerial streams. In general Tw (the average shear stress over the sediment bed only) and T may differ but we assume for simplicity that they are the same. This seems reasonable for pressurized flow but may be less realistic for open-channel flow.

We write the energy equation in the form


where Cp is specific heat, T is temperature, ΔT is the temperature drop across a thermal-boundary layer adjoining the ice wall (but not the stream bed, the temperature of which is unconstrained) and h T is a heat-transfer coefficient. We neglect work done in eroding and transporting sediment (Bagnold, 1961). The temperature T will be related to the temperature of the ice, which we expect to be at the pressure-melting point, by


where c t is the rate of change of melting temperature with pressure and ΔT is the temperature step across the thermal boundary layer at the ice/water interface. The rate of heat transfer across the thermal boundary layer is given by


where L is the latent heat and we neglect any small temperature gradient in the ice (cf. Weertman, 1972).

The Equations (3.14), (3.15), (3.17)(3.19) and (3.22) are six differential equations for the variables

and T; in addition, the seven algebraic relations in Equations (3.7), (3.10), (3.13), (3.16), (3.21), (3.23) and (3.24) determine the subsidiary variables
. Conduit perimeter I is determined by Equation (3.20)., so provided that l i and l s are prescribed functions of S i and S s respectively, the model is complete.

4. Analysis

Our aim is to solve the model posed in section 3 to determine the effective channel pressure, N c = p ip c. This depends on the inputs and most notably on the water supply. Specifically, we can, in principle, determine the water flux Q by integrating Equation (3.11). and then Nc is determined in terms of Q. In what follows, we consider only steady-state solutions. We hope this provides some average description of the flow, although water pressure often fluctuates very irregularly, at least beneath alpine glaciers.

Critical effective pressure

From (3.17) and (3.18), we have at steady state

, (4.1)


whence the relative rate of ice and till creep is given by


where we define the critical effective pressure to be


Taking K iK s, l il s, and with ρi = 9 x 102 kg m-;3, ρs = 2.65 x 103 kg m-;3, A i = 7.36 × 10-;24 Pa-;3 s-;1 (cf. Lliboutry, 1983), n = 3 and A s = 3 × 10-;5 Pab-as-;1 (after Boulton and Hindmarsh, 1987), a = 1.33, b = 1.8 we find

. (4.5)

Thus, for

, till creep is larger, while for
, ice creep is. This leads us to consider two end cases, when either one or the other creep process dominates conduit closure.

Conduit geometry

In order to proceed, we need to describe the “constitution” of the conduit. Since the conduit in Figure 1 has a boundary partly of ice and partly of till, the respective closure rates will generally be different, i.e.

. Since both till and ice are viscous, a boundary condition applies at their interface. In what follows, we assume a no-slip condition. In this case, the ice and sediment speeds towards the conduit at the interface, V i and V s, must be equal. In a slightly crude way, we identify


where k = i or s and θ(θs) is the angle of the ice roof (till flow) to the horizontal. It then follows that


and thus the two extreme cases

correspond to values
, i.e. level-floored ice channels and flat-roofed till canals, respectively. That is to say, if
, the conduit is essentially a Röthlisberger channel in ice, floored by a relatively stiff till, while if
it is essentially a canal cut into soft till and roofed by relatively stiff ice.

Note that the determination of

depends critically on the assumed till rheology. However, the identification of two asymptotically distinct conduit shapes is robust.

Determination of water pressure–flux relationship

The way to solve the equations is now as follows. First, we recognize that, based on observations c

, and the volume flux Q = uS is determined from Equation(3.14).. From Equation (4.1)., it is easy to calculate that melting due to frictional dissipation is likely to be eligible in most circumstances, and therefore the water flux Q(s) in a tunnel is essentially determined by the input from surface-water and basal tributaries. In effect, we can consider Q(s) to be a prescribed function. Equation(3.15). determines the suspended-sediment load but this is not needed for the solution of the rest of the problem.

Using Equations (3.23) and (3.24), the energy Equation (3.22). can be rewritten in the form



(Röthlisberger, 1972), K w = k w/pwCp is the thermal conductivity. The Nusselt number Nu, a dimensionless measure of heat transfer at the ice–water interface, is defined by
. Equation (4.8). determines
in terms of T, u and p c, and in addition we have Equations (4.1), (3.19), (3.21) and


to determine the five unknowns S, u, T, pc and m, (given l i and ls ). It is straight-forward to combine these to form a single, non-linear, second-order differential equation for Pc but it is probably more illuminating, and certainly easier, to simplify the system using two judicious approximations. For the first of these, note that, since Nc = p ip c and P i = ρigd, where d is the ice depth, then


where α is the ice-surface slope (Weertman, 1972). For bedrock, that varies sufficiently slowly, the term

will be negligible in this expression. For example, if the bedrock slope changes significantly over a length l f along the flowline, then


where ΔN C is the variation of N c. For an alpine glacier, we can expect ΔN C <

is negligible if
, typically about 10 ice thicknesses. On an ice sheet, we might expect
, but for ice streams, at least, geophysical evidence from (Blankenship and others 1986, 1987) suggests strongly that
. In either case, it is therefore useful as a rough estimate to neglect
in Equation (4.10).. A similarly small error is made if we then approximate the hydraulic gradient using the surface slope, thus

. (4.12)

This approximation breaks down near the glacier terminus, where variations in p c are substantial (Fowler, 1987).

The second simplification results from Equation (4.8).. We note that Q/k wNu is dimensionally a length, and is in fact the length over which

relaxes to its equilibrium value. Taking values Q = 1 m3 s-;1, K w = 10-;6 m2 s-;1, Nu ≈ 0.1 Re08 (Spring and Hutter, 1981), Re ≈ 106 (based on u = l ms-;1, l = 1 m), so that Nu ≈ 104, we have Q/k wNu ≈ 102 m, so the derivative of
in Equation (4.8). is negligible almost everywhere. The other derivative term in Equation (4.8). is

. (4.13)

It is consistent with our previous approximation that

and thus using p c = p iN c,

. (4.14)

Finally, then Equation(4.8). can be approximated by


and we have only to solve algebraic equations to find Nc, since p i and Q are prescribed functions of s. The square-bracketed term in Equation (4.16). represents that part of the viscously dissipated energy that does not go into melting the ice roof. Röthlisberger (1972) showed that for the simplest case of a horizontal tunnel and constant discharge Q, this term would be a fraction γ of the total viscous dissipation. More generally, simple estimates of the square-bracketed term in Equation (4.16). suggest it is typically about 0.1–0.3 times the term Tul/L, so for algebraic convenience we neglect it.

Model summary

We now simply have five algebraic equations:


we take l il and

to find




Thus, in order finally to determine Nc in terms of Q, we must prescribe l in terms of S. It is at this point that the conduit shape becomes critical. Specifically, we now consider the two distinct cases

(channel) and
(canal), and seek to prescribe l(5) appropriately in each case.

Channel drainage

In the case

, the tunnel protrudes into the ice, forming a Röthlisberger-like channel. We then expect a more or less semi-circular cross-section, and accordingly we choose

. (4.19)

Equation (4.17). then yields

, (4.20)

which is essentially Lliboutry’s (1983) form of the standard Röthlisberger result. Note that as Q increases, Nc increases and thus pc decreases, which is consistent with formation of an arborescent network of channels, as we discuss below. Because Nc depends weakly on Q, we can reasonably take Q ≈ constant to evaluate Nc.

Canal drainage

Suppose now that

. In this case, drainage is through a canal incised into the sediment. The mechanics of sediment erosion and transport now become very important in determining the cross-sectional shape. In particular, canals are typically wide and shallow, in which case Equation (4.19). no longer applies. If h is the canal depth, then

. (4.21)

and the determination of h depends on the precise nature of the stream bed, gravel, sand/silt and clay beds having different prescriptions. We discuss these subsequently. For now, suppose that h is prescribed in terms of the bed composition. Then, using Equation (4.17)., we have

, (4.22)

Note the crucial result that, in contrast to Röthlisberger channels, canals incised into sediment lead to a water pressure that increases with increasing water flux. This is consistent with a distributed rather than an arborescent drainage network as for flow in cavities (Walder, 1986; Fowler, 1987; Kamb, 1987). We emphasize that the physics leading to a non-arborescent subglacial canal system is distinct from what is involved in the instability causing braiding of subaerial streams. That instability involves details of the sediment-transport process, includ ing non-steady phenomena such as growth and destruction of sediment bars and passage of bed-load pulses (e.g. Parker, 1976; Ashmore, 1991). The predicted non-arborescent nature of a subglacial canal system, on the other hand, is simply the result of the relation between Q and Nc.

Numerical estimates

The coefficients b 1, b 2 and b 3 depend only on sin α, and taking n = 3 they are given by



, Röthlisberger-like channels exist, and the relation in Equation (4.20). between Nc and Q may be expressed numerically as


where Q = [Q] m3 s-;1 and N c = [N C] bar. Inverting Equation (4.24). gives

. (4.25)

Therefore, for sin α = 0.1, [Q] = 1, we have [Nc] ≈ 38; for sin α = 10-;3, [Nc] ≈ 4. Thus, we can expect that in fact

for valley glaciers, although open-channel flow may be common, as pointed out by Hooke (1984). The calculated value of Nc is, however, probably less than
for very small slopes, contrary to assumption, so R channels are unlikely to exist beneath ice streams and ice sheets, where these are underlain by deforming sedimentsor till.

If, on the other hand, we assume

, then the expression in Equation (4.22). for water flux is, numerically,


where h = [h] m. For [Q] = 1, we find [N c] = 85 [h] for sin α = 0.1, [N C] = 4 [h] for sin α = 0.001. Thus, for steep valley glaciers, canals are only viable if h < 10 cm, while for ice sheets, canals can exist for stream depths up to about 2 m. In order to estimate what h is, we must consider further the mechanisms which determine it.

Determination of channel depth


As pointed out previously, alluvial gravel-bedded streams are wide and shallow, and adjust their shapes so that the average shear stress T exerted on the stream bed is approximately equal to the critical stress Tc for bed-load transport, given by Equation (3.8). for mixed grain-size beds. Thus, from Equations (3.8), (4.17)2 and (4.21), we have


for a subglacial gravel-bedded stream. The net input of sediment by creep actually requires a slightly higher average shear stress, as argued in section 3. A crude estimate using the Meyer-Peter and Müller bed-load transport Equation (3.7)., together with measured values of

on alpine glaciers (Gurnell, 1987), suggests that expression (4.27) need not be significantly modified. Thus, we have h/D50 = O(1) for alpine glaciers (sin α~0.1)— an estimate consistent with observations on steep alluvial streams — and h/D50~O(102) for ice streams (sin α~10-3).


Here the sediment is non-cohesive, but suspended-sediment transport is important. We again expect broad, shallow channels with S ≈ lh, much as in the alluvial case. An estimate of the flow depth comes from Parker’s (1978a) analysis of sand/silt bedded streams:


where D s is, after Parker (1978a), a “characteristic” size of suspended grains. (In his comparison of theory to field data, Parker chose Ds = D 15, i.e. the grain-size such that 15% of all grains were smaller.) Rf is a dimensionless parameter that depends on grain-size and is typically O(10-;1) for fine sand and coarse silt; hence

for sin α = 10-;3.


In contrast to a non-cohesive sediment, canals cut into a cohesive bed can reasonably be expected to maintain steep banks. In this case, the cross-section may be more or less semi-circular, and we obtain the Röthlisberger result (S ≈ l2) once again. However, if the sediment layer is thin and the canal erodes through it to bedrock then h is given by the thickness of the till layer and we would again have S≈lh.

Subglacial till

Till typically has a bimodal grain-size distribution (e.g. Drewry, 1986) with a matrix rich in clay and silt-sized grains, presumably derived from comminution processes, supporting coarse clasts. The cohesive nature of the matrix may allow steep banks to occur, but since till creep and erosion also occur, we can plausibly expect that a lag will tend to form at the stream bed, comprising sandy and gravelly material. The presence of such a layer will most likely cause the canal to evolve towards a broad, shallow shape.

At the risk of generalization in the absence of much observational verification, we therefore suggest that, on the whole, canals incised into subglacial till will be broad and shallow, with a prescribed depth h that is a function of slope and grain-size D s (and, probably, discharge), with

. Although we have applied the best present understanding of stream mechanics and sediment transport, our estimates of h are poorly constrained, with likely values in the range 10 mm-1 m, the value increasing as sin α decreases.Specifically, for gravelly till with D s~10 mm, 1 cm < h < 1 m, while for sandy/silty till with D s~0.1 mm, 1 cm < h < 10 cm.

Canal effective pressures

Since h will depend on sin a as well as Ds and bed type, it is impossible adequately to summarize the resulting variations in N c. Nevertheless, the gist of the above discussion can be summarized in Table 1, where we take Q = 1 m3s-;1, We see that

for both ice sheets and glaciers, so that we could expect a canal system beneath ice sheets, whereas for valley glaciers it appears at first sight that either R channels or canals might exist. We will have more to say on this matter below.

Table 1.

5. Discussion

The results of our analysis so far are summarized in Figure 2, which indicates the two limiting cases of canals

or channels
. In particular, both types of solution can exist for large enough slopes ~10-;1), while only the canal solution will exist for very small slopes (α~10-;3). It is particularly important to under stands how the rest of Figure 2 can be completed, since for given Q there are apparently two stable solutions when α is large.

Fig. 2. Plots of effective pressure vs discharge for Röthlisberger-like channels and canals, for several values of surface slope (sin α = 10-;3, 10-;2, 10-;1). The curves for R channels intersect the line

at lower values of discharge than do curves for canals.


In discussing stability, we distinguish between the stability of an individual solution branch (here used in a mathematical sense), i.e. the stability of a particular conduit to perturbations in that conduit, and what we will term collective stability, which is the stability of a particular drainage system against perturbation to another system, e.g. two (or more) channels, or two or more canals. First we discuss branch stability.

We suppose that each end-member solution (channel or canal) derived above is stable; we have no reason to suppose otherwise. Now consider what happens as Q increases from zero. For very small Q, Darcy flow through the till will suffice to drain all melt water. As Q increases, the till pore-water pressure will increase until N= 0, at which point flotation is initiated and a water film forms at the ice/till interface. As shown by Walder (1982), this sheet flow is always unstable, and we therefore expect that a bifurcating branch representing interfacial conduits exists with Nc > 0. This situation is portrayed in Figure 3.

Fig. 3. Probable conditions for very low values of discharge. Darcy flow can drain all meltwater for sufficiently low discharge. For Q = Q c, sheet-water flow commences at zero effective pressure (broken line) but this drainage configuration is always unstable relative to channelized drainage.

At larger Q, if canals and channels can exist, then we surmise that these solutions arise as bifurcating branches from an (unstable) solution in which

, with conduits more or less equally incised into ice and till. This branch will be the continuation of the bifurcating branch at low Q. The situation is portrayed in Figure 4. We have explored other possibilities, but this one seems the likeliest. There are certain immediate implications. In the case that channels and canals can coexist (as in Figure 4 ), then because the channel-solution curve intersects the line
at a lower value of Q than the canal curve, channels will be the preferred solution as Q is increased. Moreover, if canals and channels were to coexist at larger Q (each system being in itself stable), then because the effective pressure is higher for channels (and hence the water pressure is lower), leakage through the till will drain water from canals towards channels. (This is an example of collective stability of channels.) The net effect of these factors is that if channels can exist, they will be the preferred stable drainage style, and this is likely to be the case for valley glaciers. If the surface slope is sufficiently small, as for ice sheets or ice streams, canals will form the stable drainage network.

Fig. 4. Qualitative picture of the entire pressure vs discharge plot. The line

probably coincides (for large enough values of Q) with an unstable solution branch representing conduits inched equally into ice and till, and joins with the incipient-channel solution branch beginning at Q = Qc. For small slopes (sin α~10-;3), the channel solution does not exist, whereas for large slopes (sin α~10-;1), the canal solution does not exist. At intermediate slopes (sin α~10-;2), both solutions exist but canals are always unstable relative to channels.

Hydrological and geomorphological implications

Our analysis predicts that drainage at the base of ice sheets and ice streams should take place via canals of depth ~10 cm (plus, presumably, Darcy flow) at effective pressures of only a fraction of a bar, for fine-grained sediments. This seems consistent with basal water-pressure values at Ice Stream B, West Antarctica, inferred on the basis of geophysical measurements (Blankenship and others, 1986, 1987) and subsequently confirmed by borehole measurements (Engelhardt and others, 1990), although we should emphasize that actual data on the geometry of the basal drainage system is practically nil. Our qualitative prediction is independent of the exact rheological properties of till; as long as there is some value of effective pressure below which till can creep more quickly into conduits than can ice, distributed canals as we have described will be the stable drainage type. Strictly speaking, this requires

. The results of Fowler and Walder (1993) give
, but this particular value is by no means critical.

There is indirect sedimentological evidence for the existence of broad, shallow canals beneath ice sheets. Brown and others (1987) studied deposits left by the Puget lobe of the Pleistocene Cordilleran ice sheet. Their detailed stratigraphie work in the Vashon Till of the Puget Sound region, Washington State, U.S.A., revealed the existence of lenses of water-laid, sorted sediment within the till. Brown and others interpreted these as having been deposited in subglacial conduits and estimate the width of such conduits as ranging from 0.1 to 8 m. They further inferred, from several lines of evidence, that water pressure beneath the Puget lobe was very near the ice-overburden pressure. Such low effective pressures in channelized flow would be consistent with canals but not with R channels in the ice. We further note that the discontinuous character of the water-laid subglacial deposits described by Brown and others (1987) is consistent with a picture of broad, braided canals (similar to alluvial gravel-bedded streams), reflecting the constantly shifting geometry of the braids.

In northeastern England, an area of fine-grained till, Eyles and others (1982) described “shoestring” deposits, which they interpreted as channel fills, and which may also have been deposited within a canal system. These deposits are elongated sub-parallel to the former direction of ice flow. They have nearly flat upper surfaces and irregular, concave-up lower surfaces, and are much wider than thick. Gravel lags at the deposit base are overlain mostly by sand but also include laminated clay. Eyles and others (1982, p. 323) proposed that these deposits represent cut-and-fill structures formed by subglacial drainage, with grain-size variations within the deposit indicating variable discharge, including “ponding episodes possible through closure of the channel”. Diapiric deformation at the channel bases indicates low effective pressure in the till at the time of channel formation. Lenses of till that occur within some of these deposits may represent sediment masses that dropped into the channels from overlying ice.

Another line of geologic evidence supporting our picture of drainage beneath ice sheets comes from the distribution of eskers over the region of North America formerly covered by the Laurentide ice sheet. Eskers are very rare in areas where the ice flowed over thick till but ubiquitous where the ice-sheet bed comprised bedrock with only thin, discontinuous till cover. We would predict that, in the former instance, the stable drainage system would have comprised widespread, shallow canals incised downward into the till, whereas in the latter instance, with an essentially rigid bed, we would expect an arborescent R channel network, as in the classical view of subglacial drainage. Based on observed morphological and sedimentological aspects of eskers, we believe their formation would be very unlikely in the distributed-canal system, but definitely possible with the largest channels of an arborescent network, in accord with the observed distribution. This argument has been elaborated by Clark and Walder (in press).

Shoemaker (1986) and (Alley 1989, 1992) have previously presented analyses of channelized drainage for sediment-floored glaciers but with different results. Shoemaker (1986) never considered creep of till into sediment-floored tunnels; furthermore, he did not consider the way in which sediment-transport mechanics affect conduit shape. Consistent with the neglect of these factors, Shoemaker further assumed that the relation between Q and N c for sediment-floored conduits is simply Röthlisberger’s (1972) result. We have demonstrated that this is not necessarily correct. (Alley 1989, 1992) dealt more carefully with the issues of till creep and sediment transport, but did not take into account the fact that sediment-transport mechanics may strongly affect conduit shape.

Finally, we remark that the present analysis does not point at any sort of glacier-surge mechanism associated with a switch between R channel and canal drainage. Unlike the case of R channels and cavities over a rigid bed (Fowler, 1987), we do not expect the two different types of drainage systems to coexist stably. Additionally, nothing in the present analysis indicates a dependence of the drainage-system characteristics on the ice-sliding speed (as for cavities). A surge mechanism that involves a switch between R channels and cavities is still conceivable for a glacier resting on deformable sediment. The glacier bed will not be flat and cavities can still exist to the lee of relatively large rock clasts that protrude into the ice. Thus, the potential for a linked-cavity system still exists for a deformable sediment bed, as pointed out by Kamb (1987).


J. Walder thanks former USGS colleagues C. Paola and J. Pitlick for patiently explaining principles of sediment transport and for an introduction to the pertinent literature. N. Humphrey, A. Fountain, G. Clarke, R. Alley and D. MacAyeal made helpful comments on earlier versions of this paper. P. Clark brought to our attention the issue of esker distribution as related to bed type.

This work was begun while J. Walder was affiliated with the Department of Geological Sciences and Quaternary Research Center, University of Washington, Seattle, Washington. His trips to Oxford were supported in part by U.S. National Science Foundation grant INT-8732253.

The accuracy of references in the text and in this list is the responsibility of the authors to whom queries should be addressed

Note added at proof: there is no Equation (3.1).


Alley, R. B. 1989. Water-pressure coupling of sliding and bed deformation: I. Water system. J. Glaciol., 35(119), 108118.
Alley, R. B. 1992. How can low-pressure channels and deforming tills coexist sub-glacially. J. Gtariol., 38(128), 200207.
Alley, R. B., Blankenship, D. D. Bentley, C. R. and Rooney, S. T. 1986. Deformation of till beneath Ice Stream B, West Antarctica. Nature, 322(6074), 5759.
Alley, R. B., Blankenship, D. D. Bentley, C. R. and Rooney, S. T. 1987. Till beneath Ice Stream B. 3. Till deformation: evidence and implications. J. Geophys. Res., 92(B9), 89218929.
Ashmore, P. 1991. Channel morphology and bed load pulsee in braided, gravel-bed streams. Geogr. Arm., 73A(1), 3752.
Bagnold, R. A. 1961. An approach to the sediment transport problem from general physics. U.S. Geol. Surv. Prof. Pap. 4221.
Blake, E. W. 1992. The deforming bed beneath a surge-type glacier: measurement of mechanical and electrical properties. (Ph. D. thesis, University of British Columbia.)
Blake, E. W. and Clarke, G. K. C. 1989. In situ bed strain measurements beneath a surge-type glacier. Eos, 70(43), 1084.
Blankenship, D. D., Bentley, C. R. Rooney, S. T. and Alley, R.B. 1986. Seismic measurements reveal a saturated porous layer beneath an active Antarctic ice stream. Nature, 322(6074), 5457.
Blankenship, D. D., Bentley, C. R. Rooney, S. T. and Alley, R. B. 1987. Till beneath Ice Stream B. 1. Properties derived from seismic travel times. J. Geophys. Res., 92(B9), 89038911.
Boulton, G. S. and Hindmarsh, R. C. A. 1987. Sediment deformation beneath glaciers: rheology and geological consequences. J. Geophys. Res., 24(B9), 90599082.
Boulton, G. S. and Jones, A. S. 1979. Stability of temperate ice caps and ice sheets resting on beds of deformable sediment. J. Geophys. Res., 24(90), 29–43.
Bradford, J. M. and Piest, R. F. 1980. Erosional development of valley-bottom gullies in the upper midwestern United States. In Coates, D. R. and Vitek, J. D. eds. Thresholds in geomorphology. London, George Allen and Unwin, 75101.
Brown, N. E., Hallet, B. and Booth, D. B. 1987. Rapid soft bed sliding of the Puget glacial lobe. J. Geophys. Res., 92(B9), 89858997.
Clark, P. U. and Walder, J. S. In press. Subglacial drainage, eskers and deforming beds beneath the Laurentide and Eurasian ice sheets. Geol. Soe. Am. Bull.
Clarke, G. K. C. 1987. Subglacial till: a physical framework for its properties and processes. J. Geophys. Res., 92(B9), 90239036.
Drewry, D. 1986. Glacial geologic processes. London, Edward Arnold.
Engelhardt, H., Humphrey, N. Kamb, S. and Fahnestock, M. 1990. Physical conditions at the base of a fast moving Antarctic ice stream. Science, 248(4951), 5759.
Eyles, N., Salden, J. A. and Gilroy, S. 1982. A depositional model for stratigraphie complexes and facies superimposition in lodgement till. Boreas, 11(4), 317333.
Fowler, A. C. 1987. Sliding with cavity formation. J. Glaciol., 33(115), 255267.
Fowler, A. C. and Walder, J. S. 1993. Creep closure of channels in deforming subglacial till. Proc. R. Soe. London, Ser. A., 441, 1731.
Freeze, R. A. and Cherry, J. A. 1979. Groundwater. Engelwood Clifis, NJ, Prentice-Hall.
Gurnell, A. M. 1987. Fluvial sediment yield from alpine, glacierized catchments. In Gurnell, A. M. and Clark, M.J. eds. Glacio-fluvial sediment transfer; an alpine perspective. Chichester, etc., John Wiley and Sons, 415120.
Hooke, R. LeB 1984. On the role of mechanical energy in maintaining subglacial water conduits at atmospheric pressure. J. Glaciol., 30(105), 180187.
Hooke, R. LeB, Miller, S.B. and Kohier, J. 1988. Character of the englacial and subglacial drainage system in the upper part of the ablation area of Storglaciren, Sweden. J. Glaciol., 34(117), 228231.
Hooke, R. LeB, Laumann, T. and Kohier, J. 1990. Subglacial water pressures and the shape of subglacial conduits. J. Glaciol., 36(122), 6771.
Humphrey, N., Kamb, B. Fahnestock, M. and Engelhardt, H. 1993. Characteristics of the bed of the lower Columbia Glacier, Alaska. J. Geophys. Res., 98(B1), 837846.
Iken, A. and Bindschadler, R.A. 1986. Combined measurements of subglacial water pressure and surface velocity of Findelengletscher, Switzerland: conclusions about drainage system and sliding mechanism. J. Glaciol., 32(110), 101119.
Kamb, B. 1987. Glacier surge mechanism based on linked cavity configuration of the basal water conduit system. J. Geophys. Res., 92(B9), 90839100.
Kamb, B. 1991. Rheological nonlinearity and flow instability in the deforming bed mechanism of ice stream motion. J. Geophys. Res., 96(B10), 16,58516,595.
Kohler, J. and Proksch, R. 1991. In-situ measurement of subglacial till deformation beneath StorglaciÃren, N. Sweden. Eos, 72(44), 158.
Lliboutry, L. 1968. General theory of subglacial cavitation and sliding of temperate glaciers. J. Glaciol., 7(49), 2158.
Lliboutry, L. 1983. Modifications to the theory of intraglacial waterways for the case of subglacial ones. J. Glaciol., 29(102), 216226.
Loewenherz, D.S. 1991. Stability and the initiation of channelized surface drainage: a reassessment of the short wavelength limit. J. Geophys. Res., 96(B5), 84538464.
Meyer-Peter, E. and Müller, R. 1948. Formulas for bed-load transport. In Proceedings of the 3rd Conference. International Association of Hydraulic Research, Stockholm, Sweden, 3964.
Nye, J. F. 1953. The flow law of ice from measurements in glacier tunnels, laboratory experiments and the Jungfraufirn borehole experiment. Prot. R. Soc. London, Ser. A., 219(1139), 477489.
Nye, J. F. 1976. Water flow in glaciers: jökulhlaups, tunnels, and veins. J. Glaciol., 17(76), 181207.
Parker, G. 1976. On the cause and characteristic scales of meandering and braiding in rivers. J. Fluid Mech., 76(3), 457480.
Parker, G. 1978a. Self-formed straight rivers with equilibrium banks and mobile bed. Part 1, The sand-silt river. J. Fluid Mech., 89(1), 109125.
Parker, G. 1978b. Self-formed straight rivers with equilibrium banks and mobile bed. Part 2. The gravel river. J. Fluid Mech., 89(1), 127146.
Parker, G., Klingeman, P. C. and McLean, D. G. 1982. Bedload and size distribution in paved gravel-bed streams. J. Hydraul. Div., Am. Soc. Civ. Eng., 108(HY4), 544571.
Richards, K. 1982. Rivers . London, Methuen.
Röthlisberger, H. 1972. Water pressure in intra- and subglacial channels. J. Glaciol., 11(62), 177203.
Shoemaker, E. M. 1986. Subglacial hydrology for an ice sheet resting on a deformable aquifer. J. Glaciol., 32(110), 2030.
Shoemaker, E. M. and Leung, H. K. N. 1987. Subglacial drainage for an ice sheet resting upon a layered deformable bed. J. Geophys. Res., 92(B6), 49354946.
Shreve, R. L. 1972. Movement of water in glaciers. J. Glaciol., 11(62), 205214.
Smith, T. R. and Bretherton, F. P. 1972. Stability and the conservation of mass in drainage basin evolution. Water Resour. Res., 8(6), 15061529.
Spring, U. and Hutter, K. 1981. Numerical studies of jökulhlaups. Cold Reg. Sci. Technol., 4(3), 227244.
Walder, J. S. 1982. Stability of sheet flow of water beneath temperate glaciers and implications for glacier surging. J. Glaciol., 28(99) 273293.
Walder, J. S. 1986. Hydraulics of subglacial cavities. J. Glaciol., 32(112), 439445.
Weertman, J. 1972. General theory of water flow at the base of a glacier or ice sheet. Rev. Geophys. Space Phys., 10(1), 287333.
Weertman, J. and Birchfield, G. E. 1983. Stability of sheet water flow under a glacier. J. Glaciol., 29(103), 374382.
Wiberg, P. L. and Smith, J. D. 1989. Model for calculating bed load transport of sediment. J. Hydraul. Eng., 115(1), 101123