Introduction

Deformation of unlithified sediment beneath wet-based ice masses may contribute significantly to the motion of modern ice streams (e.g. Reference Alley, Blankenship, Bentley and RooneyAlley and others, 1986) and surging glaciers (e.g. Reference Clarke, Collins and ThompsonClarke and others, 1984). Bed deformation may have also enabled fast flow of Pleistocene ice masses, perhaps contributing to the 100 000 year climate cycle, to extreme millennial climate variability and to the low aspect ratios of ice sheets during the Last Glacial Maximum (Reference Clark, Alley and PollardClark and others, 1999). There is therefore widespread interest in modeling accurately the flow of glaciers that rest on sediment.

In most models, the till that commonly underlies glaciers has been assumed to behave as a viscoplastic fluid, such that at shear stresses exceeding its Coulomb yield strength, it shears at rates proportional to the excess shear stress (Reference Boulton and HindmarshBoulton and Hindmarsh, 1987; Reference AlleyAlley, 1989; Reference Walder and FowlerWalder and Fowler, 1994; Reference BoultonBoulton, 1996; Reference Clark, Licciardi, MacAyeal and JensonClark and others, 1996; Reference Jenson, MacAyeal, Clark, Ho and VelaJenson and others, 1996; Reference HindmarshHindmarsh, 1998a–c; Reference Licciardi, Clark, Jenson and MacAyealLicciardi and others, 1998). Support for this assumption has come from the field measurements of Boulton and co-workers at Breidamerkurjökull, Iceland (Reference Boulton and HindmarshBoulton and Hindmarsh, 1987; Reference Boulton and DobbieBoulton and Dobbie, 1998). By inserting segmented rods and cable extensometers in till beneath the glacier margin, they measured displacement as a function of depth and found that most glacier motion resulted from pervasive shearing of the upper 0.5–1.0 m of till. Rates of shear strain increased upward toward the glacier sole. These displacement profiles have been interpreted as evidence of linear or slightly non-linear viscoplastic behavior (stress exponent < 2.5) (Reference AlleyAlley, 1989; Reference Boulton and DobbieBoulton and Dobbie, 1998).

As was first noted by Reference KambKamb (1991), however, most experimental data from soil mechanics (e.g. Reference MitchellMitchell, 1993, p.343) and geophysics (e.g. Reference Biegel, Sammis and DieterichBiegel and others, 1989) contradict such inferences. These data indicate that strengths of slowly shearing granular materials do not increase significantly with deformation rate, contrary to the viscoplastic behavior inferred from displacement profiles. Laboratory experiments with various testing devices on a variety of tills confirm that the strength of till is insensitive to its deformation rate (Reference Tika, Vaughan and LemosTika and others, 1996; Reference Iverson, Hooyer and BakerIverson and others, 1998; Reference Tulaczyk, Kamb and EngelhardtTulaczyk and others, 2000). These studies indicate that till is best idealized as a Coulomb material, such that its strength is independent of deformation rate and linearly dependent on effective normal stress. This is not a surprising conclusion: soils engineers have characterized soil deformation successfully with Coulomb models for many decades (e.g. Reference Terzaghi, Peck and MesriTerzaghi and others, 1996), and recent studies of till deformation beneath modern glaciers have yielded results that are difficult to explain with viscoplastic models (e.g. Reference Hooke, Hanson, Iverson, Jansson and FischerHooke and others, 1997; Reference Truffer, Harrison and EchelmeyerTruffer and others, 2000).

Despite the clear motivation for Coulomb models, viscous deformation resistance continues to be assumed in models of subglacial till deformation (Reference NgNg, 2000; Reference Thorsteinsson and RaymondThorsteinsson and Raymond, 2000). This, in part, may reflect the lack of an explicit and simple analytical description of how distributed till deformation can occur in subglacial till without viscous deformation resistance.

The objective of this investigation is to demonstrate that displacement profiles like those observed in till beneath Breidamerkurjökull may arise from deformation resisted only by Coulomb friction and that, therefore, a viscoplastic till rheology cannot be inferred confidently from such profiles. Our analysis also indicates that neither dilatant hardening (Reference Iverson, Hooyer and BakerIverson and others, 1998) nor pore-pressure diffusion into the bed from the ice–till interface (Reference Tulaczyk, Mickelson and AttigTulaczyk, 1999), two proposed mechanisms for distributing deformation in Coulomb models, is necessary to account for the observed distribution of displacement.

Background

Yielding will occur in a granular material when the shear traction applied to a planar element in the material equals the Coulomb yield strength, proportional to the effective normal stress on the element. Everywhere the material is near this condition, called limiting equilibrium, it will be on the verge of deforming irreversibly. As many others have noted (e.g. Reference Boulton and JonesBoulton and Jones, 1979; Reference Alley, Blankenship, Bentley and RooneyAlley and others, 1986), for a subglacial till layer to be weakened sufficiently to approach this limiting state, till pore-water pressure must be a large fraction of the total normal stress.

Laboratory experiments, motivated by diverse problems in physics, geophysics and engineering, provide insight into how deformation occurs as limiting equilibrium is approached. These studies highlight the stress heterogeneity that characterizes granular materials and the temporal variability of intergranular stresses during deformation. Experiments with particles in static and slowly shearing packs indicate that stresses are concentrated along “force chains” which coalesce and form complex networks that include many particles (Fig. 1a) (e.g. Reference Dantu and GlanvilleDantu, 1957; Reference Drescher and de Josselin de JongDrescher and de Jong, 1972; Reference Allersma, Vermeer and LugerAllersma, 1982; Reference LiuLiu and others, 1995; Reference Painter, Tennakoon, Behringer, Herrmann, Houi and LudingPainter and others, 1998; Reference Veje, Howell, Behringer, Schollann, Luding, Herrmann, Herrmann, Houi and LudingVeje and others, 1998; Reference Howell, Behringer and VejeHowell and others, 1999). During shearing, these chains reorganize continuously and cause intergranular stresses to fluctuate significantly (Fig. 1b) (e.g. Reference Painter, Tennakoon, Behringer, Herrmann, Houi and LudingPainter and others, 1998; Reference Howell, Behringer and VejeHowell and others, 1999). Such fluctuations have been measured normal to the shearing direction in ring-shear experiments with simulated and natural tills (Fig. 1c) (Reference Iverson, Hooyer and HookeIverson and others, 1996, Reference Iverson, Baker and Hooyer1997; Reference HookeHooke, 1998, fig. 7.16). Not surprisingly, these spatial and temporal variations in intergranular stresses give rise to instantaneous velocity fields that are highly heterogeneous, with groups of particles pushed forward sporadically along surfaces of concentrated particle slip or rolling (Reference Drescher and de Josselin de JongDrescher and de Jong, 1972; Reference Misra, Herrmann, Houi and LudingMisra, 1998).

Fig. 1. (a) Force chains in 2900 photoelastic (birefringent under strain) disks sheared within a Couette apparatus. The outside diameter of the sample chamber is 0.38 m. Disks are 7 and 9 mm in diameter. Lighter disks are under larger forces. The rate of shearing at the specimen center line is 2 mm s⁻¹ (adapted from Reference Howell, Behringer and VejeHowell and others, 1999). (b) G ², an optical measure of the force on disks proportional to the number of birefringent bands they display, as a function of displacement during an experiment like that shown in (a). Over the measured range of G ², its value is nearly linearly proportional to the force on disks. The value plotted is an average over one-tenth of the samplers circumference (adapted from Reference Howell, Behringer and VejeHowell and others, 1999). (c) Local stresses measured normal to the shearing direction at the upper boundary of a till specimen during deformation in a ring- shear device. Forces are measured over circular platens with diameters of 18 and 36 mm. The rate of shearing is 0.01 ^{mm s−1} (adapted from Reference Iverson, Baker and HooyerIverson and others, 1997).

Results of distinct-element numerical models, in which equations of motion (forces and moments) are solved simultaneously for every particle assuming appropriate particle-contact rules (e.g. Reference Cundall and StrackCundall and Strack, 1979), reinforce these laboratory results. Simulations of slow deformation with such models indicate that force chains always develop with the associated tendency for grains to move and rotate sporadically in coherent clusters (Reference Cundall, Drescher, Strack, Vermeer and LugerCundall and others, 1982; Reference Misra, Chang, Babic, Liang and MisraMisra, 1997; Reference Williams, Rege, Chang, Babic, Liang and MisraWilliams and Rege, 1997). In especially relevant simulations of fault-gouge deformation in a shear zone, Reference Morgan and BoettcherMorgan and Boettcher (1999) applied a distinct-element model to grain-size distributions similar to those of till and found that force chains formed and failed continuously, resulting in local stress fluctuations. Shearing was caused by the failure of force chains and focused along planes that became parallel to the shearing direction at low strains (0.1). Although strain accrued through sporadic slip along individual planes, the cumulative effect of these slip events (shear strain = 2.0) was a relatively uniform distribution of finite strain.

Hypothesis

We explore the hypothesis that local stress fluctuations like those observed in experiments cause bulk deformation of till subglacially (e.g. Reference Boulton and HindmarshBoulton and Hindmarsh, 1987). Mean shear and effective normal stresses on a till layer control whether it is near limiting equilibrium and on the verge of deforming. However, particle movement at a particular depth is driven by brief reductions in effective normal stress caused by the failure of force chains. Consequent force imbalances on groups of particles push these groups forward episodically along slip planes at various depths, resulting in bulk deformation of the bed over time.

The objective of our analysis is to calculate the distribution of motion with depth in the bed that results from these local departures from static equilibrium. Rather than consider a viscoplastic till rheology, we assume simple Coulomb friction, in which rate-independent slip between grains depends linearly on effective stress. Non-synchronous till displacements along hypothetical slip planes are calculated using conservation of momentum, which also dictates the maximum depth to which deformation extends in the bed. Integrating till displacements from this depth to the ice–bed interface yields the bed-normal displacement profile.

Our general approach is similar to that of distinct-element models in that, rather than apply a continuum model for flow, we solve the equation of motion for discrete portions of the till while applying an appropriate friction relation for slip. Although our simple analytical model falls well short of the realism achieved using numerically intensive alternatives like the distinct-element method, it illustrates with a minimum of complexity an important concept in interpreting the subglacial deformation of till: that fluid-like deformation profiles need not be a consequence of a fluid rheology.

Analysis

Consider a thick till bed that slopes at an angle β (Fig. 2). The depth, y, below the ice–till interface is measured normal to the bed surface. Motion in the bed is assumed to occur along a series of bed-parallel slip planes separated by a uniform distance, δ, much smaller than the bed thickness. Although local slip surfaces in granular media are sometimes inclined at an angle to the macroscopic shearing direction, experimental evidence (Reference Mandl, de Jong and MalthaMandl and others, 1977; Reference Logan, Dengo, Higgs, Wang, Evans and WongLogan and others, 1992; Reference Beeler, Tullis, Blanpied and WeeksBeeler and others 1996) indicates that at sufficiently large strains slip surfaces tend to align with the shearing direction, consistent with the results of distinct-element models (Reference Morgan and BoettcherMorgan and Boettcher, 1999). Thus, idealizing the bed as essentially a deck of cards, as we do here, is a reasonable abstraction. The till is assumed to be water-saturated with a constant bulk density, ρ _t. Possible ramifications of bulk density varying temporally are discussed later in the paper.

Fig 2. Model bed geometry.

As in distinct-element models of granular shearing, we calculate displacements using Newton’s second law. In the most general case, slip in the bed may occur along neighboring slip planes concurrently. Thus, a slice of till bounded by slip planes may have a frame of reference that is accelerated relative to that of the till slice beneath it. For this case, if V _r is the bed-parallel velocity of a till slice relative to the slice immediately underlying it and t is time, the equation of motion for slip is

where F is the net bed-parallel force acting on till and ice above the slip plane, m is the mass of till and ice above the slip plane, and V is the velocity of the till slice beneath the slip plane relative to a fixed reference frame. The second term on the righthand side of this equation results from the accelerated reference frame (Reference AharoniAharoni, 1972, p. 192). This term can be eliminated, however, if we assume that slip events at various depths occur at discrete times. Distinct-element modeling of fault-gouge deformation by Reference Morgan and BoettcherMorgan and Boettcher (1999) indicates that, although more than one slip surface can be active at any time in a shear zone, there is a general tendency for slip to alternate from one plane to another. Similar behavior has been assumed in models of subglacial till deformation (Reference Tulaczyk, Mickelson and AttigTulaczyk, 1999). For non-synchronous slip events, till below an active slip plane will always be at rest, and Equation (1) reduces to the standard form of Newton’s second law:

A net driving force F on till and ice above a particular slip plane exists when the Coulomb shear strength becomes less than the downslope weight per unit area of the overlying till and ice. If σ _e is the effective normal stress on planes at depth y, and φ is the Coulomb friction angle of the till, then the Coulomb shear strength is σ _e tan φ. The till is considered to be deformed sufficiently to be in its critical state, and thus φ is a constant and cohesion is negligible (Reference SkemptonSkempton, 1985). If W _d is the downslope component of the weight per unit area of the overlying till and ice at depth y, the net force per unit area that causes motion is W _d =− σ _e tan φ. If H is the ice thickness, the mass per unit area of the till and overlying ice at depth y is ρ _t(y + CH), where C = ρ/ρ _t and ρ is the density of ice. Incorporating these relations in Equation (2) yields

V _r can only be positive because up-glacier motion does not occur if W _d < σ _e tan φ.

To characterize σ _e and W _d as functions of depth, we first emphasize that σ _e has a fluctuating component due to the formation and failure of force chains. Motion is expected to occur when force chains fail (e.g. Reference HookeHooke, 1998, p. 109; Reference Morgan and BoettcherMorgan and Boettcher, 1999). Thus, we define σ _e such that it includes a transient perturbation, . Such perturbations are expected to occur randomly with depth and to have a wide range of magnitudes (Reference Painter, Tennakoon, Behringer, Herrmann, Houi and LudingPainter and others, 1998; Reference Veje, Howell, Behringer, Schollann, Luding, Herrmann, Herrmann, Houi and LudingVeje and others, 1998). Rather than attempt to characterize these perturbations statistically, which would add realism but sacrifice clarity, we make the simplest possible assumption: that perturbations, although not concurrent, are of equal magnitude at all depths. If is taken to be positive when the effective normal stress is increased, then the effective stress in the till on any plane normal to y during a perturbation is:

where g is the gravitational acceleration, σ ₀ is the total normal stress at the ice-till interface, and λ is the ambient pore-pressure ratio, which equals the fraction of the total normal stress supported by pore-water pressure in the absence of an effective-stress perturbation. The value of λ can vary from 0 to 1 and accounts in a general way for the full range of possible pore-water pressures, given that the steady-state pore pressure is not necessarily hydrostatic. The downslope component of the weight per unit area of ice and till at depth y is

where W ₀ and ρ _t gy sin β are the contributions due to the weights of the ice and till, respectively. Substitution of Equations (4) and (5) into Equation (3) yields

To condense Equation (6), we define the following parameters that describe the state of stress on potential slip planes:

where S ₀ is the difference between the downslope component of the weight of the ice per unit area and the available shear strength of the bed at y = 0 in the absence of an effective-stress perturbation, S′ is the shear-strength perturbation caused by the effective-stress perturbation, and a is a constant (0 < α < 1) that characterizes the steady-state stress due to the weight of the till. Since we are assuming that the bed is at rest unless perturbed sufficiently, the shear strength of the bed exceeds the downslope weight per unit area of ice, and thus S ₀ ≤ 0. Only a reduction in effective normal stress, corresponding to a negative value of , can cause slip, and thus we consider only the case in which S′ ≥ 0. Using Equations (7a–c) allows the equation of motion (Equation (6)) to be written in a compact form:

If the shear-strength perturbation, S″, is sufficiently large (S″ + S ₀ > ρ _t gαy)), W _d (Equation (5)) exceeds resisting stresses at depth y, slip occurs there in the till, and deformation will extend to a depth in the bed, y ₀, that depends on the magnitude of the numerator of Equation (8). Note that at this depth V _r will always equal zero, so this depth can be determined by setting the lefthand side of Equation (8) equal to zero and solving for y to obtain

Deformation decreases to zero at depth because the bed shear strength, controlled by intergranular friction, increases more rapidly with depth than W _d. At depth y ₀, friction is sufficiently large that W _d is balanced exactly by the till shear strength, and there is no motion.

Next we determine the relative velocity between neighboring till slices for a given shear-strength perturbation, S′, by solving Equation (8). Two time periods are relevant: the period during the shear-strength perturbation when the till slice (and overlying till and ice) are accelerating and the period immediately after the shear-strength perturbation when the till slice is decelerating but still moving owing to the inertia of the till and ice. For the duration of the perturbation, V _r is obtained by integrating Equation (8) and applying an initial condition which states that motion commences from a state of rest. The integration yields

where T is the duration of the perturbation. During the second stage of motion, S′ = 0 because the perturbation has ceased. Integrating Equation (8) with S′ = 0 and applying Equation (10) as an initial condition at t = T yields

where t _stop is the time required for till at a particular depth y to come to rest. This time is determined by setting V _r = 0 in Equation (11) and solving for t to obtain

Note that values of t _stop diminish with depth, indicating that motion events caused by effective-stress perturbations near the glacier sole persist longest. This persistence indicates that effective-stress perturbations of a given magnitude induce larger force imbalances near the glacier sole.

The bed-normal profile of displacement can be obtained by first determining the relative displacement, X _r(y), between till slices. To do so, we integrate Equation (10) from 0 to T and Equation (11) from T to t _stop, and add the two results to obtain

Writing this equation in terms of t _stop (Equation (12)) yields a simpler expression,

which illustrates that at large depths where t _stop approaches T, X _r(y) approaches zero. At depth y ₀, T = t _stop, and X _r(y) = 0.

To calculate the total displacement, X(y), the relative displacements on slip surfaces between y ₀ and y must be summed. The number of slip surfaces is given by the range of depth over which the summation is performed, y ₀ − y, divided by the distance between slip surfaces, δ. Thus, keeping in mind that δ is assumed to be very small, the summation is expressed by the integral,

where the negative sign accounts for the negative (upward) direction of the integration. Substituting Equation (13) into this equation and integrating results in a lengthy expression for the desired displacement profile (see Appendix).

For the important special case in which the downslope component of the weight per unit area of the glacier (W ₀) and the bed shear strength are balanced sufficiently at the ice–till interface to assume S ₀ = 0, Equation (13) can be written conveniently in terms of y ₀ (Equation (9) ) to obtain

Substituting this equation into Equation (15) and integrating yields the displacement profile for this special case:

Equation (17) is similar to the general expression presented in the Appendix, but that expression contains a more complex logarithmic term and cannot be written as concisely in terms of y ₀. Note that displacements indicated by Equations (16) and (17) become infinite as y → 0, so the model results are meaningful only at depths y ≥ δ/2. Displacements at y = 0 depend on interactions between basal ice and till particles that are not included in the model.

Equation (17) indicates that total displacements depend on the square of the duration of the shear-strength perturbation, T ², and depend non-linearly on the magnitude of the shear-strength perturbation, S′, through dependence on y ₀. Displacements are largest if the slip-plane spacing δ is small because δ is inversely related to the number of slip planes, and the total displacement of a given till slice increases with the number of slip planes beneath it. Note also that the normalized total displacement, X(y)/y ₀, scales with the dimension-less parameter αgT ²/2δ.

Parameter Choices and Results

The objective now is to explore whether profiles generated by the model are consistent with those observed at Breidamerkurjökull, Iceland. Reference Boulton and DobbieBoulton and Dobbie (1998) presented the most complete data on the distribution of till motion beneath that glacier (Fig. 3a). Through a borehole about 1 km from the glacier margin, these authors placed a displacement marker at the level of the glacier sole, and four other markers at various depths in the underlying till. The markers were connected with wires to potentiometer-gauged spools that were held in the borehole at the base of the ice. As the bed sheared, the spools were carried away from the markers, and the length of wire paid out was monitored for 17 days to obtain the distribution of displacement in the bed over that period. Unfortunately, certain details of the experiment were not discussed, such as how the spools were held fixed at the bottom of the borehole and how the depth of marker emplacement in the bed was measured. In addition, the trajectories of the wires between the markers and the spools were unknown as a function of time. Owing to this significant uncertainty, others have restricted use of this technique to measurements in the bed only very near the glacier sole (Reference Blake, Fischer and ClarkeBlake and others, 1994; Reference Engelhardt and KambEngelhardt and Kamb, 1998). Despite these uncertainties, however, inferred profiles were qualitatively similar to those measured directly by excavation nearer to the glacier margin (Reference Boulton and HindmarshBoulton and Hindmarsh, 1987), with shear strain increasing upward toward the glacier sole. For the purposes of this paper, therefore, the displacement profiles of Figure 3a will be taken at face value.

Fig. 3. (a) Daily displacement of anchors in the bed of Breidamerkurjökull over a 17 day period, as reconstructed by Reference Boulton and DobbieBoulton and Dobbie (1998, fig 13). (b) Daily displacement indicated by the model over 17 days, assuming that limiting equilibrium was approached each day during one period of sufficient duration for a perturbation in effective normal stress to occur once at each depth in the bed. Parameter values used in the model are listed in Table 1.

Many of the parameters of the model are known from field observations at Breidamerkurjökull and from ancillary laboratory measurements (Reference Boulton and HindmarshBoulton and Hindmarsh, 1987; Reference Boulton and DobbieBoulton and Dobbie, 1998). These are H, σ ₀ , ρ _t , C, β and φ (Table 1). Other parameters, such as λ and W ₀, can be estimated. Although pore-water pressure, recorded by sensors in the displacement markers, fluctuated through a range of about 200 kPa (Reference Boulton and DobbieBoulton and Dobbie, 1998), averaging over the period of record indicates that λ = 0.9 is a good approximation. The value of W ₀ is specified such that it equals the resisting stress at the bed surface in the absence of effective-stress perturbations. This is a convenient but not essential assumption that satisfies the special case noted earlier in which S ₀ = 0, and allows use of Equation (17) in computing the velocity profile. Thus, using values from Table 1 and Equation (7a), W ₀ = σ ₀(1 − λ) tan φ = 58.7 kPa, a reasonable value. An alternative would be to calculate W ₀ from the glacier thickness and surface slope, but that approach would be suspect owing to unknown longitudinal stress gradients near the glacier margin.

Other important parameters in the model — the spacing of slip surfaces, δ, the shear-strength perturbation, S′, and its duration, T — can be estimated only within broad limits. The value of δ must be sufficiently large to approximate the movement of many grains in coherent groups (e.g. Reference Drescher and de Josselin de JongDrescher and de Jong, 1972; Reference Misra, Chang, Babic, Liang and MisraMisra, 1997) but sufficiently small to result in time-integrated strain that is pervasive macroscopically. Values of S′, driven by perturbations in effective normal stress as given by Equation (7b), are expected to be small for slow deformation. Moreover, distinct-element modeling of fault-gouge deformation indicates that cushioning of large particles by adjacent smaller particles results in short, discontinuous force chains that keep stress fluctuations small (Reference Morgan and BoettcherMorgan and Boettcher, 1999). A reasonable assumption, therefore, is that the value of S′ is a small fraction of the mean shear strength of the till layer. Values of T are likely to be small since individual force chains remain intact for only small increments of strain (Reference Morgan and BoettcherMorgan and Boettcher, 1999).

Given that values of δ, S′ and T are not well known, we now explore whether displacement profiles, calculated with Equation (17), can be fitted successfully to those from Breidamerkurjökull (Fig. 3a) by adjusting the values of these parameters through reasonable ranges. The value of S′ controls the depth to which deformation will extend after sufficient time has elapsed for perturbations in effective normal stress to occur at all depths in the bed (Equation (9)). Thus, the value of S′ is adjusted to match the observed depth of deformation. Then T is adjusted until the total displacement of the uppermost till slice agrees with the observed displacement there. The value of δ is taken to be 0.01 m.

Figure 3b depicts displacement in the bed indicated by the model if this fitting procedure is used. We assume that limiting equilibrium was approached for one period each day at Breidamerkurjökull, consistent with the roughly diurnal variation in pore-water pressure measured there (Reference Boulton and DobbieBoulton and Dobbie, 1998, fig. 11), and that during each period a perturbation in effective normal stress occurred once at every depth in the bed. Fitted values of S’ and T are 3.1 kPa and 0.16 s, respectively. Somewhat larger values of T are obtained if a larger value of δ is chosen, but since displacements vary with T ², T cannot be much larger without producing displacements that are unrealistically large. Conversely, the fitted value of T is smaller if a smaller value of δ is used.

The model results agree with the observed convex-upward displacement profile. This displacement pattern is independent of the fitted values for S′, T and δ and is therefore a fundamental feature of the model. The profile is convex upward because friction along slip surfaces increases with depth in the bed more rapidly than the downslope component of the weight per unit area of till and ice. Thus, a given perturbation in effective normal stress applied to each depth in the bed causes larger force imbalances and relative displacements between till slices near the glacier sole than at depth.

Discussion

The model presented here should be viewed as illustrative but not predictive. It is based on a fundamental principle, conservation of momentum, and a well-established empirical rule, Coulomb friction. It illustrates that the observed distribution of strain at Breidamerkurjökull arises naturally from these principles if combined with the fact that brief, local departures from static equilibrium instigate motion in the bed. On the other hand, the model is not a useful predictive tool at this time. It contains several parameters that are difficult to measure or estimate and that we have adjusted to bring the model results into quantitative agreement with observations. In addition, as discussed hereinafter, certain subglacial complexities are neglected.

An important tenet of the model is that over brief intervals and at specific depths in the bed, driving stresses exceed resisting stresses as force chains rearrange and local normal stresses fluctuate. Although such departures from quasi-static stress equilibrium are usually not considered in glaciology, they drive the slip events that cause the slow deformation of granular materials. Indeed, this concept is essential in distinct-element models of slow granular flow, in which net forces and torques on individual particles are computed and used to calculate particle displacements and rotations from Newton’s second law (e.g. Reference Cundall, Drescher, Strack, Vermeer and LugerCundall and others, 1982; Reference Morgan and BoettcherMorgan and Boettcher, 1999).

An idealization of the model is that perturbations in effective normal stress are assumed to occur at discrete times with equal frequency and magnitude at all depths. A more realistic model, which would require a numerical approach, would apply perturbations randomly with depth and time and a distribution of perturbation magnitudes consistent with that observed experimentally (e.g. Reference Veje, Howell, Behringer, Schollann, Luding, Herrmann, Herrmann, Houi and LudingVeje and others, 1998). There is no obvious reason to believe, however, that the displacement profile in that case would differ significantly from the concave-down profile obtained with our simpler model. If there were independent evidence that the character of force chains varied systematically with depth, causing associated variations in the magnitude or frequency of stress fluctuations, a more complicated model would be justified. Although such a situation is possible and might be caused, for example, by a progressive change in grain-size distribution with depth, it is probably not typical.

Another neglected complexity is that during a shearing episode, till porosity may vary, contrary to the constant bulk density assumed in the model. Local variations in till porosity are likely as force chains form and fail causing particle rearrangement. In addition to affecting the local friction angle of the till, a change in porosity may cause a change in pore-water pressure, owing to the generally small hydraulic diffusivity of till. If, for example, porosity increases during slip, pore-water pressure may be reduced, thereby strengthening slip surfaces and reducing the displacement expected from a given perturbation in effective normal stress. Reference Iverson, Hooyer and BakerIverson and others (1998) suggested that this process, called dilatant hardening (e.g. Reference ReynoldsReynolds, 1885), might distribute strain in till due to strengthening of active slip surfaces and their consequent migration. Although this process may well occur and help distribute deformation, our model illustrates that even in the absence of this effect, distributed deformation will occur in a Coulomb till.

Reference Tulaczyk, Mickelson and AttigTulaczyk (1999) and Reference Tulaczyk, Kamb and EngelhardtTulaczyk and others (2000) emphasized that water-pressure fluctuations at the ice–till interface and consequent pore-pressure diffusion into the bed will help distribute strain in a basal till layer. Pore-water pressure variations will decrease in amplitude with depth and arrive later at depth than at shallower levels. These pore-pressure variations cause temporal variations in effective normal stress with depth that are hypothesized to drive vertical migration of bed-parallel shear zones in the bed. Like dilatant hardening, this process may be important, but our model illustrates that it is not necessary to distribute deformation.

Conclusion

Despite its limitations, our model of subglacial till deformation demonstrates that displacement profiles like those observed at Breidamerkurjökull do not contradict the results of laboratory tests that indicate till behaves as a Coulomb material. This is an important conclusion since Coulomb behavior has widespread empirical support and has for many decades been used successfully to explain diverse phenomena associated with soil deformation. Coulomb behavior is, indeed, the appropriate null hypothesis for analyzing till deformation, so a viscoplastic till rheology should not be invoked if Coulomb models can account for observations. Our Coulomb model of the bed can account for observations at Breidamerkurjökull. Thus, fitting viscoplastic flow relations to displacement profiles measured in till there or in till beneath other glaciers may yield little insight.