Symbols Used
The following are the main symbols used in the paper:
1. Introduction
The formulation of a basal sliding law which realistically describes the physical processes at the base of glaciers and ice sheets has long been of interest. Reference WeertmanWeertman (1957) was the first to model the basic physical processes involved, which include the enhancement of flow round protuberances by stress-dependent viscosity, and lubrication of the bed by a very thin water Film (the mechanism of regelation). Various authors have since tackled different aspects of this problem. Notably, Reference NyeNye (1969, Reference Nye1970) and Reference KambKamb (1970) considered the problem of slow flow of a constant-viscosity fluid over a bedrock of small slope, and were thus able to linearize the problem, and establish a sliding law for very general bed shapes. Kamb extended his solution to nonlinear viscosity (Glen’s law) in an ad hoc manner, Reference MorlandMorland (1976[a], Reference Morland[b]) solved the same problem using the methods of complex variable theory, and included, in the latter paper, the frictional effects of basal debris.
Two major effects, those of a non-linear flow law, and of subglacial cavitation, have been seriously addressed by Reference LliboutryLliboutry (1968, Reference Lliboutry1975, Reference Lliboutry1978, Reference Lliboutry1979) and by Reference WeertmanWeertman (1964). Reference LliboutryLliboutry (1979) reviewed his work, and Reference WeertmanWeertman (1979) gave an overall view of the subject at that time. The non-linearity of the flow law incapacitates exact solutions of the Nye‒Kamb type, but Weertman’s notion of a controlling obstacle length is very useful, and the general “Weertman‒type” sliding laws
where τ is stress, u is basal velocity, and c is constant, have proved very popular. Here a is an exponent which ranges from one for a Newtonian fluid, to 1/n (the exponent in Glen’s law) for a non-linear fluid, the latter in the particular case that roughness is absent at low wavelengths. More generally, 1 ≥ a ≥ 1/n. Reference FowlerFowler (1981) introduced complementary variational principles for the non-linear problem. He showed that, in the absence of cavitation, an important dimensionless parameter was σ/ν n+1, where σ = x 0/d, ν = y 0/x 0, and x 0, y 0 represent a typical bedrock wavelength and bedrock amplitude, respectively. The length scale d represents the length scale of velocity variations in the “outer” flow, away from the bed. It may be thought of as a typical glacier depth, although it can be defined purely in terms of the prescribed bed shape, as is done below, just before Equation (2.5). Thus ν is the typical magnitude of the bedrock slope relative to the mean. The parameter µ 1 = σ/ν n+1 can easily be either large or small for reasonable values of ν and σ, so both large or small sliding velocities can be understood in terms of variable bedrock roughness. However, Fowler neglected the extra functional effects of regelation, which is tantamount to neglecting roughness at short wavelengths (< 1 cm).
The importance of cavitation has been stressed by Reference LliboutryLliboutry (1968, Reference Lliboutry1979). The latter paper affords the most complete analysis of this subject. Again, the presence of cavities renders analytic solutions intractable, but Lliboutry was able to derive functional representations for the sliding law by making judicious assumptions about the form of the flow. In particular, he introduced a shadowing function s(t), where if t is the angle of declination at which a light were to be shone at the bedrock, then 1 − s(t) is the fraction of the bed in shadow. Thus, s(t) is purely a geometrical property of the bed. Next, he identified t with the (assumed constant) slope of the cavity roofs, and assumed that t is inversely proportional to the velocity (and depends on other things as well), as is at least qualitatively reasonable.
From this, he could determine the sliding law for different bedrocks, characterized by their shadowing function. In particular, he made the important observation that the sliding law over a sinusoidal bedrock, and that over a general (Gaussian) bedrock, are radically different. To quote, for a sinusoidal bedrock, “all the bumps are drowned simultaneously [as t → 0] and the friction drops to zero”. Thus Lliboutry found, for a sine profile, that τ increases to a maximum and then decreases to zero as u tends to infinity (at fixed effective pressure N). However, for a “Gaussian” bedrock, he found that τ/N is essentially constant as u → ∞, giving the “solid friction law”, which has been used, for example, by Reference ReynaudReynaud (1973).
The present paper is inspired by the above results of Reference LliboutryLliboutry (1979). We wish to examine the sliding law for a fairly general bed form, but without using the arbitrary assumption that cavity roofs all have (the same) constant slope, which will not in general be true. Rather, we base our method on another suggestion in the same paper, which is to consider a bed consisting of superimposed bumps of different scales. As Lliboutry pointed out for the case of no cavitation, this is strictly an asymptotic procedure, which assumes that the successive length scales λ n of the different undulations form an asymptotic sequence. However, we can hope that some useful results can still be obtained even if λ n+1/λ n is not too small. As such, this procedure is akin to methods of re-normalization group theory, which have been recently exploited in understanding a wide variety of physical phenomena. For example, see Reference Smalley, Smalley, Turcotte and SollaSmalley and others (1985).
The basic building block for this model is that of flow over a periodic “humped” bedrock. Results for the Nye‒Kamb version of this problem (small slopes, constant viscosity) have recently been given (Reference FowlerFowler, 1986), and they corroborate Lliboutry’s results for a sinusoidal bedrock. The model of the present paper essentially corroborates Lliboutry’s results for a general bedrock, although different cases present themselves.
One important purpose of analysing sliding with cavitation is to examine whether a multiple-valued sliding law can realistically occur. If it could, then one could explain surges almost immediately (Reference HutterHutter, 1982). Reference LliboutryLliboutry’s (1979) “Gaussian” law is marginal in this respect. However, another point arises, in that the effective pressure N appears in the cavitational sliding law. One can then enquire as to whether N itself depends on u. The answer has been persuasively shown to be yes (Reference KambKamb and others, 1985), using both detailed field measurements and theoretical analysis. In this paper we will take the simplest model, both of Reference RöthlisbergerRöthlisberger’s (1972) tunnel drainage, and Kamb’s “linked cavity” drainage, to show that N drops catastrophically at the onset of cavity drainage. As a result, we show that a multi-valued sliding law is a realistic possibility, and we offer some comments on the ramifications of this for surges, seasonal waves, hydraulic pressure fluctuations, and also to sliding of ice sheets.
2. Solutions of the Nye–Kamb Problem with Cavitation, Periodic “Humped” Bed Forms
Reference FowlerFowler (1986) has presented the solution of the following version of the Nye‒Kamb sliding problem. Given a rough bed of amplitude y 0 and wavelength x 0 written in the form
(so h, h′ ~ 0(1)), find the relation between shear stress τ and velocity u far from the bed, when ice, considered to be of constant viscosity, flows over it. The boundary conditions at the bed are that there is no shear stress, and that the normal velocity is zero. This latter condition, in particular, is equivalent to neglecting the drag due to regelation (see Reference FowlerFowler, 1981, p. 672, for discussion), which can be shown to be a good approximation unless significant roughness is present at vertical scales of less (perhaps a good deal less) than 1cm (Reference FowlerFowler, 1981). Even in that case, the regelative drag can then be added separately (Reference LliboutryLliboutry, 1979); it is not germane to the problem of cavitation, which involves larger length scales.
By the nature of the analysis, it was convenient to restrict the results to periodic bedrocks with one cavity per period. Roughly, this means one hump per period. A typical result is shown in Figure 1. Sliding laws for other bed forms also exhibit a single maximum of stress, followed by a decrease to zero as u → ∞.
Basal stress versus basal velocity for the skewed, periodic bed, shown in the inset.
Figure 1 shows a curve of stress versus velocity, which also involves the effective pressure N, defined by
where p i is the ice-overburden pressure and p w is the local hydraulic drainage pressure in the water. Precisely, Figure 1 depicts the relation between stress, velocity, and effective pressure in the following dimensionless way:
where τ, u, p c are dimensionless stress, velocity, and effective pressure, respectively. To convert Equation (2.3) to dimensional units, we use the definitions of Reference FowlerFowler (1986):
where a suffix d denotes a dimensional unit. The scales [τ] and U represent typical values of these quantities, e.g. [τ] = 1 bar, U = 100 m year−1; these may vary according to the context. Introducing the depth scale d via [τ] = ηU/d, where η is the viscosity (assumed constant in the calculation), and using the definition σ = x 0/d = ν 2 to define d (Reference FowlerFowler, 1986), we finally rewrite Equation (2.3) in the dimensional form
The relation (2.5) is independent of a particular glacial context and could be compared with (for example) the experiments of Reference Budd, Budd, Keage and BlundyBudd and others (1979). Equation (2.3) is more useful to see how sliding velocities are likely to compare with surface velocities in a particular glacier or ice sheet. Notice f(β) (where β is the argument of f in Equation (2.5)) is a dimensionless function, having a maximum of f ≈ 0.4 at a value β ≈ 4. Other bedrocks (sine wave, symmetric humps) betray similar features (Reference FowlerFowler, 1986). The variation of the point of maximum τ is somewhat lessened if ν is defined in terms of up-stream slopes.
The point of onset of cavitation is not marked in this figure. For small β (no cavitation), f(β) is linear (for this Newtonian viscosity), and then τ ∝ u. Cavitation sets in for β > β c, where in Figure 1, β c ≈ 0.09. For other more symmetric bedrocks, the values are a little higher. To gain an idea of the size of β, we write
where ϵ would be the surface slope of a glacier. If [τ] is a (typical) shear stress, then N* ≈ 1 if p w is atmospheric (or zero) (free drainage or no water at the bed), whereas N* « 1 for high water pressure. Then
where
measures the bedrock-slope variations relative to mean surface slope. For glaciers, we might have ν ≳ 0.1, ϵ ~ 0.1, so that µ 2 ≳ 1 is typical. Then cavitation occurs if u > β c µ 2 N* ~ 0(1) for a low water-pressure case. As is well known, cavitation does often occur in glaciers. For large ice sheets, ν tends to be larger (corresponding to much larger irregularities), whereas ϵ is smaller, typically 10−3. Then µ 2 ~ 103, and cavitation would occur if u ≳ 0(103), where as u ~ 1 for τ ~ 1. Cavitation in this sense is ruled out for large ice sheets, though some lubrication is afforded by basal lakes, for example, in Antarctica. However, “temperate” sliding theories, such as developed here, are likely to be irrelevant for cold ice masses which only reach melting temperature at their beds.
“Shadowing” Function
For later purposes, we need to define a function which measures how much of the bedrock is overlain by cavities. Such a function was introduced by Reference LliboutryLliboutry (1968), and is called the shadowing function, and is denoted by s(β), since it will be a function of the (dimensionless) velocity. Specifically, s(β) is the fraction of the bedrock which is free of cavities. For the constant-velocity Nye–Kamb problem solved by Reference FowlerFowler (1986), we can calculate s(β) as part of the solution. We find that β < β c, s = 1, and for β > β c, it decreases monotonically to an asymptotic limit of zero, as β → ∞. We portray this in Figure 2, where we have defined an auxiliary function
It may be pointed out that Reference LliboutryLliboutry’s (1968, Reference Lliboutry1979) estimation of the shadowing function appropriate to the flow is heuristic, whereas, the function in Figure 2 is derived from the solution of the flow problem. Nevertheless, Lliboutry’s estimation of the form of this function is essentially correct.
The function S(β) for the same cases as in Figure 1
The details of S(β) near β = 0 are obscured by carrying out a reasonable number (e.g. 20) of computations for different values of β (this is why Figure 2 has kinks ‒ it is not smoothed as Figure 1 is), the approach of s to zero as β tends to β c ≈ 0.09 has been omitted. In fact, S approaches zero with infinite derivative as β → β c, and it is zero for β < β c. Also, S tends to infinity as β tends to infinity. Thus, the important feature is that S = 0 for β < β c(no cavitation), and increases rapidly, and monotonically, for β > β c.
Non-linear viscosity
We are not able to legitimately generalize Equation (2.5) for a non-linear flow law. Nevertheless, it is tempting to do so in an ad hoc manner, and we shall follow this temptation here. If we consider Glen’s law
then the viscosity η is
Thus the far-field viscosity in the sliding problem is η ∞ = 1/Aτ D n−1 ; however, stresses are higher by 0(1/ν) near the bedrock, since the drag in the far field is balanced by the integral of the normal stress along the bedrock; this must be higher by 0(1/ν), since the bedrock amplitude is only 0(ν), We will assume that deviatoric stresses are also 0(1/ν) near the bedrock, (This is true when n = 1, and can be deduced when n > 1 from the scaling of the problem (Reference FowlerFowler, 1981). Thus we choose
to generalize Equation (2.5). There follows
which is in fact consistent with Reference LliboutryLliboutry’s (1979) nonlinear result. The dimensionless version of Equation (2.13) is, using Equations (2.4) and (2.6),
where
and
where d is the depth scale, and the velocity scale U now satisfies
Equation (2.14) generalizes Equation (2.3), both for n = 1 and for d defined independently of x 0 and y 0, as might be appropriate for glaciers and ice sheets, where d is determined by the large-scale flow.
For small β < β c, we assume
where c ~ 0(1) is a roughness coefficient. This holds for n = 1, and extends to β ~ β c. Adopting Equation (2.18), Equation (2.14) reduces to
This is, in fact, exactly the law found by Reference FowlerFowler (1981) for non-linear sliding without cavitation, and therefore gives some confidence in its applicability. Cavitation sets in if µ 1 µ 2 β c N*τ n−1 ~ µ 1 µ 2 (if N* ~ 1). For glaciers, one can have µ 1 µ 2 ~ 1, but for ice sheets, µ 2 ~ 103 as before, and µ 1 ~ σ /ν n+1 ~ 1/20 for d ~ 2 km, [x] ~ 100 m, ν ~ 1. Again µ 1 µ 2 ≫ 1, and cavitation is unlikely. (Both cavitation conditions, for n = 1 and n > 1, can be written as τ > β c µ 2 N*, which also shows the difficulty for ice sheets.)
In summary, sliding with cavitation over single-humped bedrocks is controlled by two dimensionless groups, µ 1 and µ 2, which may be thought of as between them controlling the extent of cavitation, and the relative size of the resultant velocity. An important point which emerges from inspection of the S(β) graphs, and more dramatically from the plots of the cavities themselves (Reference FowlerFowler, 1986), is that the stress only starts to decrease when the cavity from one bump reaches the stoss face of the next bump. The drag starts to decrease, since it is all being contributed by the high pressure on the stoss face, and the extent of influence of the high pressure wanes due to the invasion of the stoss face by the up-stream cavity. As Lliboutry put it, the bumps get drowned. It is because of the simultaneous drowning that the stress decreases; but in reality, when bumps of one amplitude and wavelength become drowned, larger bumps of longer wavelength, with a correspondingly smaller β value, will not be drowned, and thus will take up the excess stress from the smaller bumps. Consequently, we can realistically expect stress to continue to increase, so long as the entire bed is not almost drowned, as indicated by Lliboutry. In the following section, we use a re-normalization approach, based on the above discussion, to get some idea of how a sliding law for a more realistic bedrock will behave.
3. Re-Normalization Approach for a General Bed: the Sliding Law τ = kN a1 u a2
Let us now consider a more general bedrock, which we consider to be a superimposition of similar bedrocks of different scales – for example, sine waves, corresponding to Fourier decomposition. In order to be able to superimpose results for the beds of the preceding section, it is necessary that the length and height scales of the different components be asymptotically distinct. In that case, the flow over a bump of a certain scale is independent of the flow over larger-scale bumps. The smaller-scale bumps exert a traction on the flow, which, however, can be combined additively (cf. Reference FowlerFowler, 1981). In this sense, the different components respond independently to the flow.
To scale the bed, we suppose there is some preferred wavelength x 0 and associated amplitude y 0; for example, x 0 might be the scale at which some measure of roughness (e.g. the spectral power density of the slopes) was maximum. Then each component is of the form
where the dimensionless wavelength λ ranges from 0 to ∞, and γ(λ) is a prescribed roughness spectrum. Following Equation (2.13), the dimensional drag offered to the farfield flow is
where we have defined
and
As before, we non-dimensionalize τ D, µ D, and N using
and
then
is the dimensionless drag offered by the component of wavelength λ.
The above is for the case that all the bumps are in contact with ice, and that bumps of smaller scales are absent. In reality, for a decreasing sequence {λr} (such that λ r » λr+1), some of the bumps will be ubmerged by cavities, and the cavity-free parts of unsubmerged bumps will have smaller bumps superimposed on them. Let ρ r be the fraction of λr bumps which are unsubmerged by larger cavities, and let s r = s(β r) be the fraction of uncavitated bedrock for flow over λr bumps, where β r is the argument of f in Equation (3.7), that is
where N
r (= N
λr
) is the effective pressure experienced by the λr bumps. Further, let
τ
r be the overall stress actually experienced by the λ
r
bumps (including the additional stress due to smaller bumps), and let 
be the actual stress experienced by one of the unsubmerged bumps. Then we have
In order to calculate 
, the problem studied by Reference FowlerFowler (1986) must be amended to allow for a non-zero stress at the attached parts of the ice flow. For the case of a Newtonian rheology, this can be done by subtracting off a linear shear flow, and the result is simply that the extra stress at the interface must be added to that computed from Equation (3.7). For the case of a non-Newtonian rheology, Reference FowlerFowler (1981, p. 668) derived the same result assuming small slopes and no cavities. Allowing for cavities, the result is then
Now, those λr+1 bumps which are unsubmerged exist on the uncavitated parts of unsubmerged λr bumps. It follows that
Combining Equations (3.9) to (3.11), we have
which may be considered to represent the overall stress contributed by the λr bumps alone. (This means, below, that we compute the total stress by summing the “individual” contributions.)
Let p r be the mean ice pressure on the uncavitated parts of the unsubmerged λr bumps, and let p be the cavity pressure (assumed the same in all cavities). Then, as pointed out by Reference LliboutryLliboutry (1979), we must have, by an approximate averaging,
that is
This follows because far from the λr+1 bumps, the pressure (on unsubmerged λr bumps) is p r whereas this same pressure is transmitted to the bed by a pressure p r+1 over the uncavitated part, and p on the cavitated part, of λr+1 bumps. Then the total stress is obtained by summing τ r over r. As r → — ∞, (λ→∞), we expect m → 0, and N r→ N*, the (dimensionless) effective pressure, s r → 1, and ρ r → 1 (no cavitation, since β r → 0).
Notice immediately that Equations (3.14) and (3.11) imply that
is constant, and thus τ is defined by
it remains to solve for β as a function of λ, and then sum Equation (3.16).
First, we note that Equation (3.7) may be rewritten as
If we put
and define
then Equation (3.17) takes the form
which serves as an implicit definition of g(b r). Then Equation (3.18) gives τ λ as an explicit function of b r. Further, differentiation of Equation (3.20) with respect to b r yields, using the chain rule,
(note βr = b r/g n-1). Thus g′/f′ is positive providing (n − 1)βr f′ + f > 0, and g has the same shape as f in Figure 1. This condition on f states that (if n ≥ 1) β 1/(n−1) f is increasing, which is certainly true for small enough β. If this condition is not satisfied, then g(b r) will increase for large enough b r, as was found by Reference LliboutryLliboutry (1968, fig. 12), giving a multi-valued sliding law. Below, we suggest that such subtleties are irrelevant in reality, as the values of b r are restricted to values for which f′ > 0. We now write Equation (3.16) in the form
We define
then Equations (3.19) and (3.14) imply that b r satisfies the difference relation
where
and we consider S as a function of b r+1; it will be of the same general form as in Figure 2. Then the graph of b r+1 versus b r is shown schematically in Figure 3, for q r+1/q r less than or greater than one. In this figure b c is the value of b at the onset of cavitation. Thus S = 0 for b < b C, and increases rapidly thereafter.
Schematic graph of br+1 versus br.
Now the behaviour of maps such as Equation (3.24) is well known (e.g. Reference Collet and EckmanCollet and Eckmann, 1980) and is illustrated on Figure 3. For fixed q r+1/q r, b r → 0 when r → ∞ (smaller wavelength) if q r+1 < q r, whereas if q r+1 > q r, b r will tend towards b*, the intersection of Equation (3.24) with the line b r+1 = b r, which is the non-zero fixed point of Equation (3.24), viz.
Thus, it is reasonable to consider these two states as approximate alternatives. As q r+1/q r increases through one, there will be a transition of b r from 0 towards b*.
The picture which then emerges of the degrees of cavitation is that if q decreases with r, there is essentially no sliding. For q increasing, b increases from 0 to b c: these bumps are not cavitated and contribute the usual viscous drag. For b > b c, the bumps are cavitated, but as further wavelengths are considered (larger r), the degree of cavitation (as measured by b r) rapidly tends towards a (relatively) fixed value. Drowning does not occur, even of the smallest bumps. For S increasing quite rapidly, and particularly for n > 1, the value of b* is quite close to b c.
Our approximation to the sliding law now follows from Equation (3.22), using the fact that s = 1 for b < b c:
where b* is given by Equation (3.26), and s* = s(b*). For the uncavitated bumps, b < b c, we suppose, following Equations (2.18) and (2.19), that f(β) = cβ, and hence, from Equation (3.20)
where R = c 1/n . For a sine wave, R = 1 for n = 1. R ≈ 1.46 for n = 3 (Meysonnier, 1983). Noting that N r = N* for b < b c, we finally have
where g* = g(b*)/s*, and can be approximately evaluated using Equation (3.28).
There remains the task of evaluating (3.29) for various m(λ). Before doing so, notice that (3.29) is quite similar to the form of Reference LliboutryLliboutry’s (1979) result. The first term in (3.29) is a viscous drag over non-cavitating bumps, the second due to cavitation. Lliboutry has the first term split in two, one for large wavelengths (τ ∝ u 1/n ), the other for very small wavelengths which do not cavitate, for which the Weertman law τ ∝ u 2/(n+1) includes the effects of regelation. Whether we find this lack of cavitation at small wavelengths depends on our choice of q(λ).
Let us now consider power-law bedrocks. These are bedrocks for which the amplitude function
Typically, one considers superimposed sine waves but we have the advantage that we are not restricted to particular wave forms of this type, providing we adjust f(β) accordingly. The case α = 3 corresponds to the so–called “white roughness”, because it contributes equally to the viscous drag at all wavelengths. Measured values are typically in the range 2–4 (Melosh and Kamb, unpublished), with (perhaps) preference towards the lower end α ≈ 2.4 (Reference BenoistBenoist, 1979; Reference HalletHallet, 1979). Self-similar profiles obeying (3.30) have
Since we suppose λ decreases as r increases, this implies q increases for all r if α < 4, and thus that cavitation extends to all (small) wavelengths. These profiles also have m tending to zero at low wavelength (cf. Reference Budd, Budd, Keage and BlundyBudd and others 1979) if α > 2.
If q increases with r for all r, there is a unique critical λ such that b = b c, given by (from Equation (3.19))
Cavities exist over bumps of wavelength 
.
To estimate the sums in (3.29), we suppose λ ∼ exp(−r/r 0), and then approximate the sums by integrals, using −dr = r 0dλ/λ. The reason for choosing a logarithmic dependence was given by Reference NyeNye (1970); see also Reference WeertmanWeertman (1979). Nye recommended a value of r 0 slightly less than one, depending on the method used to extract the rough bed from its smooth average. In view of other approximations involved, we will put r 0 = 1 below. Thus in (3.29), summation
and we have m = λq. Then (3.29) can be approximated as
The sums converge or diverge with the integrals, and for convergence we require, if q ~ λ−Σ(Σ = 2 ‒ α/2), Σ > 1, and (1 + 1 /n) Σ < 1, that is
For n = 1, this is 2 < α < 3, for n = 3, it is 2 < α < 2.5. Assuming (3.34), (3.33) is
where from Equations (3.32) and (3.31) (q = λ−Σ, with co-efficient 1 by choice of x 0 and y 0)
After simplification using Equation (3.36), (3.35) can be written in the form
where
and
is an 0(1) constant. Equation (3.37) is a very simple generalization of the functional form of Weertman’s law,
. Notice some important features.
-
(i) τ increases with u (at fixed N*) for all n, and levelling off does not occur (in contrast to Reference LliboutryLliboutry’s (1979) result). This assumes that q → 0 and λ →∞, so that b in Equation (3.17) tends to zero (no cavitation). This seems very reasonable.
-
(ii) we have neglected regelation in deriving an exponent different to 1/n for u in Equation (3.37). Regelation will be important at low wavelengths, if the ice is otherwise blocked. In the classical theory without cavitation (Reference NyeNye, 1969, Reference Nye1970; Reference KambKamb, 1970), regelation is important for length scales less than perhaps 50 cm. Non-linear ice rheology seems to lower this value (Reference FowlerFowler, 1981). The importance of regelation is registered by how fast the bedrock amplitude spectrum γ(λ) decreases to zero as λ → 0. If it decreases rapidly enough, then the Nye–Kamb drag without regelation converges as λ → 0, which is to say that the effect of regelation in this case will be small and may be neglected. The condition that this should apply (when there is no cavitation) is that γ << 0(λ2) as λ → 0. Alternatively, this requires q → 0 as λ → 0, which is only satisfied by (3.31) if α > 4.
Conversely, if α < 4 in (3.31), then the small-scale roughness is enough to stop effectively viscous deformation, and it is only because of regelation that there can be any flow at all. Thus, if α < 4, the classical theory cannot neglect regelation, and the drag converges for (3.31), providing now (when n = 1) α >2.
Now, if we suppose α < 4, and consider in addition cavitation, then we find that q increases as λ decreases, so that cavitation occurs at a length scale much larger than the regelation length scale. As the length scale decreases, cavitation facilitates sliding, and the drag converges at small wavelengths for α < 2. Therefore, it seems that when cavitation is added to the classical sliding problem, we can consistently ignore the effect of regelation as a small one on the sliding law, providing now α > 2. Whether regelation as well as cavitation will extend this range of convergence downwards is unknown, but we have at least shown that for α > 2 regelation may be neglected.
-
(iii) Even for the general bedrock considered here, the role of the two parameters μ 1 and μ 2 is crucial, μ 1 governs the actual magnitude of the velocity, whereas μ 2 measures the role of the effective pressure. The 0(1) parameter k can only be estimated in practice by fitting experimental data; in fact, the best one could do would be to fit the constants a 1, a 2, a 3 in a law of the form
to the data. This has been done to some extent (Reference Budd, Budd, Keage and BlundyBudd and others, 1979; Reference BindschadlerBindschadler, 1983).
-
(iv) It is obvious that more general beds will produce sliding laws which, computed from (3.33), can be much more complex. For example, a more general bed will have a distribution ρ(m, λ) of wave forms of “slope” m and wavelength λ; it is not then entirely obvious how to generalize (3.33), in particular, how to transform sums over m to integrals over m; and it may be that the model is too coarse to allow useful indulgence in such generalizations.
Finally, notice that (3.33) can be written in the form
where, using Equation (3.32),
and c 1, c 2 are 0(1) constants, F(T) decreases monotonically to 0 as T increases to ∞, G(T) increases monotonically in the opposite manner. The sliding law in Equation (3.37) follows from the approximations
which may usefully serve except (perhaps) at the extreme ranges of T.
Finally, it may be useful to rewrite Equation (3.40) in dimensional terms, using Equations (3.5) and (3.6). We find
and then
which may be compared with Reference LliboutryLliboutry’s (1979) equations (93), (103), and (108), or (110). The expression of the coefficient of u D 1/n seems to be new.
4. A Simplified Theory of Basal Drainage: the Multi-Valued Relation N(u)
Since the sliding law depends on the effective pressure N, it is necessary to have some idea as to how N will change in response to a change in the hydraulic regime. That is the aim of the present section. First, we consider the drainage through a central conduit (Reference RöthlisbergerRöthlisberger, 1972), as appropriate to non-surging glaciers. We first assume that the conduit is flooded. Let S R denote the cross-sectional area of the channel, m the mass of ice melted per unit length per unit time, P i the ice-overburden pressure, p the water pressure, M the influx of water from the glacier to the channel per unit length per unit time, and Q R the volume flux of water. The equations governing the variables have been given by Reference NyeNye (1976), following Reference RöthlisbergerRöthlisberger (1972):
in which X is the mean bedrock slope, K is a constant related to the viscosity of ice, ρ w is the density of water, g s is the component of gravity down the bedrock slope, N is an empirical constant related to turbulent channel flow, and L is the latent heat. These equations represent a simplified version of the more complete model given by Spring and Hutter (1981, 1982), and in particular they neglect terms representing thermal inertia, cf. Reference ClarkeClarke (1982). The “momentum” Equation (4.1) 3 assumes the Manning formula; other choices of friction factor are possible, leading to slightly different exponents, cf. Spring and Reference HutterHutter (1981, 1982)
If d is a typical glacier depth and ℓ is a typical length, and X is the mean bedrock slope, then p i ~ p i gd, so that p ≤ ρi ~ p i gd and hence ∂p/∂x ~ p i gd/ℓ; consequently, (∂ p/∂ x)/ρ w gX ~ d/ℓX. For typical values d ~ 100 m, ℓ ~ 10 km, X ~ 0.1, d/ℓX is small and so we can neglect ∂p/∂x in Equations (4.1). This is a singular approximation which will cause a boundary layer to form at the glacier snout, where p is prescribed. Away from this layer, and neglecting ∂p/∂x. Equations (4.1) imply successively
(assuming a steady state),
and finally
gives the effective pressure. Thus, if M is known, we solve (4.3) to find Q R as a function of x, the distance downstream of the head with a boundary condition, e.g. Q R = 0 at x = 0 Notice that if Q R is constant, N R is also, as noted by Reference RöthlisbergerRöthlisberger (1972). However, in general, Q R will vary with x. We observe that the length L/gX, with L ≈ 3.3 × 105 J kg−1 and gX ≈ 1 m s−2 (X ≈ 0.1), is 330 km, so that Q R will not vary much due to frictional melting, and most of the variation of Q R is due to external sources — surface melt water, for example. To determine M in any reasonable manner would be prohibitive, but luckily we see that N R is proportional to Q R 1/4n which for n = 3 is very insensitive. This suggests that we simply treat Q R as constant, and equal to its outlet value, in (4.4) and Equation (4.5). Then N R is essentially constant.
With ρ i ≈ 900 kg m−3, ρw ≈ 1000 kg m−3, gX ≈ 1m s−2, K ~ 0.02 bar–3 year−1 (Röthlisberger, 1972), L as above, we find
Thus, if Q R/S R = ū m s−1, then, with n = 3,
Further, with N = 700 m−8/3 kg (Reference NyeNye, 1976), we have
so that, if QR = Ǭ m3 S−1, (4.4) implies
and thus
Thus, with these parameter values, if Q R = 10 m3s−1 we would have S R ~ 2m2, u ~ 4.8 m s−1; then (4.7) would imply N R ≈ 30 bar. By comparison, a glacier of depth 100 m has p i = ρi gd = 9 bar. It should be said that indeterminacy of some of the above parameters (notably N and K) will cause different precise values, but it is clear that it is plausible that many glaciers should have N R > p i, i.e. p w < 0! The resolution of this paradox is that if Equation (4.5) implies N R > p i, then the appropriate condition is that p = p A, the atmospheric pressure, the channel is not filled as previously assumed, and consequently (4.4) must be discarded. Additionally, the form of the channel will be modified, and the Equations (4.1) 1,4 will have to be changed. The momentum Equation (4.1) 3 will involve the water depth, which will be determined in place of p. We need not pursue this, as these details do not concern us. The fact that we can commonly expect unflooded conduits has been pointed out by Reference HookeHooke (1984) and by Reference LliboutryLliboutry (1983). The latter paper suggests that determination of N is inherently a time-dependent (seasonal) problem. It remains to be seen whether an annual average can satisfactorily describe long-term processes such as surges.
Let us now turn to the situation as it pertains during a surge (Reference KambKamb and others, 1985). Observations of Kamb and others suggest that, when cavities are widely present, most of the drainage occurs through the cavities, which are linked by small conduits, or by thin water sheets (Reference WeertmanWeertman, 1972). One difference this generates is that flow is much slower, since ū, the mean velocity, is much lower than in the Röthlisberger case. Drainage now occurs between cavities and one can in a simplistic manner apply the Röthlisberger theory with a suitable change of the important variables. It is important to realize that the actual situation will be highly time-dependent, and again we hope that our steady-state model reasonably represents the average state of the system. Suppose there are n K cavities of largest dimension across the glacier width. We can approximately identify this largest dimension with the critical wavelength λc x 0. As an approximation
in the notation of section 3, and where W is the glacier width. We again assume a steady state (thereby excluding discussion of transient effects due to rainfall, for example), so that the volume flux Q K of water is prescribed; we take it as constant. If the largest cavities have one main drainage conduit linking them to each other, then the flux through each conduit is ≈ Q K/n K, and the Reynolds number is ≈ Q K/S K 1/2 v w n K where v w is the kinematic viscosity of water ≈ 2 × 10−6 m2 s−1, and S K is the cross-section of each inter-cavity conduit. Let us assume that the Reynolds number >103, so the flow is turbulent. For example, if the “Kamb orifices” have cross-section S K ~ 10−2 m2, Q K ≥ 2 m3s−1, nK ~ 103, then Re ~ 104.
The determination of the effective pressure N K is much as before. However, it is no longer obvious that the average water-pressure gradient is the slow drift due to gravitation, since the ice pressure now varies on the cavity-length scale. Nevertheless, we take ρ w gX – ∂p/∂x as being approximately ρ w gX, since if |∂p/∂x| » ρ w gx, then −∂p/∂x is positive, but also we can suppose p equal at the end of orifices. This is inconsistent unless p actually varies little. Nevertheless, it must be said that a more thorough description of orifice dynamics in a three-dimensional cavity field would be welcome. A recent paper by Reference WalderWalder (1986) contains further discussion.
As before,
where m k is the melting rate, and
where P is the local ice pressure, which is given by (Reference LliboutryLliboutry, 1979, compare Equation (3.13) and subsequent discussion above)
where s is the same “shadowing” function as in Equation (2.9); thus
replaces Equation (4.5). From Equations (4.16), (4.5), (4.12), and (4.4)), we find
where we write s(N K) to emphasize the dependence of s on the effective pressure, and where we have taken Q K ≈ Q R. A typical value of n K is perhaps 103 (width 1 km, bedrock wavelength 1 m), and then N k is reduced at the onset of cavitational drainage by a factor (for n = 3) of at least two. Strictly speaking, n K varies with λc, but such variation is largely irrelevant compared to its magnitude, and a reasonable guess is that N K is approximately constant given by Equation (4.16), with n K ≈ W/x 0. Notice that N R is the value given by Equation (4.5), which may be greater than the overburden pressure. Consequently, in spite of the reduction, N K may also be greater than overburden. In that case, presumably p w = p A, the atmospheric pressure.
Stability
The above estimates give approximations for the effective pressure in the cases where either Röthlisberger’s channels or Kamb’s linked cavities provide the basic drainage system. However, one must also examine the stability of the system when both cavities and a channel are present. Our analysis follows that of Röthlisberger and Kamb, who showed that a drainage system which has p′ (Q) negative (channels) will form a single drainage conduit. On the other hand (Reference KambKamb and others, 1985), cavities have p′ (Q) positive, and so can exist stably with many different cavities interlinked. This follows from (4.17) providing ∂N K/∂N R < 0. By considering the graph of ln(1/s) in Figure 2, it is not difficult to see that this will typically be true, since near the critical value of u/N n for which cavitation is initiated, s decreases very rapidly as u/N n increases, that is ∂s/∂N is large and positive. Therefore, also ∂[N/s]/∂N will typically be negative. Nevertheless, this stability characteristic is rather different from that of Reference WalderWalder (1986) and Reference KambKamb and others (1985), which is based on a Röthlisberger-type balance of flow within the cavities themselves. This involves the assumption that melting at cavity walls can be significant. To estimate its significance here, we do a Röthlisberger-type calculation for a cavity. The melting rate m ~ ρ w gXq/L, where q is the flux through the cavity. The normal velocity corresponding to this is v n ~ m/ρ iℓ. The ratio of v n to νU, where ν is the bump-aspect ratio and U is the sliding velocity, measures the importance of melting. Thus, v n * = v n/νU ~ ρw gXq/ρ iℓvUL, where ℓ is cavity length. With values q = 0.02 m3 s−1, ℓ = 1 m, ν = 0.2, U = 100 m year−1, we get v n * ~ 10−1. On this basis, melting is small, and we thus visualize the break-down of channel flow as being due to the enlargement of Kamb orifices rather than of the cavities themselves (on this latter topic, a paper by B. Kamb is in preparation). There remains the analysis of the situation when both systems are present.
In this case, let p K be the cavity/conduit water pressure in the cavity system, and let p R be the pressure in the Röthlisberger channel. If these are not equal, we suppose flow occurs from one to the other at a rate proportional to p R ‒ p K. Let the cross flow per unit length down-glacier towards the Röthlisberger channel be k(p K ‒ P R). Then, conservation of water volume yields the following equations:
where ℓ is the mean cavity spacing, and V K is the total cavity volume per lengtht ℓ. We take
The latter relation reflects the fact that cavities (of λ-bumps) are of height γy 0, length (1 − s)ℓ, and of typical lateral dimension ℓ. We revert to summation over λ, as in section 3, as the total volume is important. We did not really need this for calculating N K, since n K really measures the number of main inter-cavity channels, without reference to individual cavities. Except for very, very small cavities, we may take V K ≫ ℓsS K, and then Equations (4.18) and (4.19) simplify to
where we have put k = 1 without loss of generality (it just measures the time unit), and
As in section 3, we approximate the sum (4.21) by an integral over λ. Following section 3, we suppose b ≈ b* for λ < λc, so that s ≈ s(b*) = s*, and then
For a power-law bedrock, λ c is given by (3.36); with m = λ1−Σ (cf. (3.31)), we have
Where
and ∧ is given by
(hereN = N K). Then
and, using Equations (4.27) and (4.28), we find that (4.23) can be written as
where
A* is essentially of the order of the total cavity cross-sectional area.
By eliminating Q R from (4.4) and (4.5), we find
where
and thus Equation (4.22) becomes
The stability of the co-existing state N R = N K = N is easily examined (see Appendix). It is stable (i.e. tunnel drainage is stable) if
For illustrative purposes, note that
For typical values n = 3, Σ = 1, S R ≈ 5 m2, A* ≈ 200 m2, we have ᴧc ≈ 1/4, and this value can thus be realistically exceeded. For values of ᴧ less than ᴧc, cavities and a tunnel can co-exist, stably. It should be emphasized that cavities continue to exist when tunnel drainage is operative.
Thus, our conclusion concerning the effective pressure is that
where N R is given by (4.4) and (4.5) if N R < p i and N = p i if N R > P i. When cavitational drainage occurs we have
with a value
for example, ifs ≈ 1/2, W = 1000 m, x 0 = 1 m. Notice that s in (4.39) is really a function of ∧. The effect of this abrupt transition on the sliding law is examined in the following section.
5. A Multi-Valued Sliding Law Incorporating Drainage Characteristics
Assuming the ad hoc generalization of the Newtonian result, (2.13), is appropriate, the sliding law for a general bed can be written in the (dimensional) form
and we can assume φ is realistically monotone increasing, e.g. φ ~ ᴧγ, γ > 0, as in section 3. The simple theory of section 4 suggests that a qualitively reasonable approximation to the effective pressure N is given by Equations (4.34) and (4.35). Thus N ≈ N R for Λ < Λc, N ≈ δN R for ᴧ > ᴧc as shown in Figure 4, which shows that N versus u lies on a z-shaped curve. For values of u < u 1, N ≈ N R, for u > u 2, N ≈ δ N R, but for u 1 < u < u 2, N can take any of three values, of which we can expect the middle one to be unstable. However, in this range, two different values of N are possible.
Typical representation of effective pressure N versus velocity u when a sudden increase in N occurs at a critical value of ᴧ = u/Nn.
The effect of Figure 4 on the sliding law is shown in Figure 5. For simplicity, we use the generalized Weertman’s law
given in (3.37). For u < u 2, stable tunnel drainage is possible, and the stress follows the upper curve. For u > u 1, stable cavity drainage is possible, and the stress follows the lower curve. Between u 1 and u 2, an unstable region is possible on which ᴧ = ᴧc, and consequently, the variation of N with u implies that the stress τ varies non-monotonically with u. Assuming the derivation of Equation (5.2) is reasonable, this is a general conclusion, provided only that the cavity-drainage pressure N K ≈ δN R is less than the overburden pressure. With the estimates in section 4, we might have N K ≈ 10bar, which would imply multi-valuedness is only possible for sufficiently deep glaciers ≥100m deep. This estimate is of course not precise, and will vary with drainage characteristics, etc. Whether surging behaviour can actually occur or not will depend on whether the stress can reach high enough levels to trigger the transition to cavity drainage. This again depends importantly on the depth and is studied more closely in the following section.
A multi-valued sliding law for τ versus u. when the effect of the hydraulic switch portrayed in Figure 4 is taken into account.
6. Discussion
Much of the preceding mathematical development will be hard to follow, so let us attempt to summarize the main points of the argument.
First, the Nye–Kamb sliding problem with cavities can be solved for periodic bedrocks with isolated bumps. For this case, one finds τ/N = g(u/N n), where τ is shear stress, u is velocity, and N is effective pressure. The function g(b) has a single maximum, and tends to zero as b → ∞, when all the bumps become simultaneously drowned. As Reference LliboutryLliboutry (1979) emphasized, such a result is not realistic for bedrocks with many different obstacle sizes. Then, when small obstacles become drowned, larger obstacles take up the extra stress, and thus basal stress continues to increase with velocity. A general way of deriving a sliding law for a non-uniformly bumpy bed is the re-normalization principle enunciated by Reference LliboutryLliboutry (1979). Application of this yields a sliding law which may be generally written in the form
where τ, N*, u are dimensionless variables,
and v = y 0/x 0 is a typical bedrock roughness slope, ϵ is the (typical) ice-surface slope, σ is a measure of the coarseness of the bed, σ = x 0/d, and d is the depth. Using Equations (2.4), (2.6), and (2.17), Equations (6.1) can be written in the dimensional form
but it is actually more useful to stay with Equations (6.1). For a power-law bed of spectral exponent α (2 < α < 4), one has
so that
where
Generalized power laws of the form in Equation (6.4) have in fact been examined experimentally by Reference Budd, Budd, Keage and BlundyBudd and others (1979), and in field work by Reference BindschadlerBindschadler (1983). Budd and others’ experiments were fitted by them to the form of Equations (6.1), with n = 2, and φ ~ ᴧ1/3, corresponding to α = 2.5 (with n = 2). For lower effective pressures, they found τ ∝ Nu, which is not really compatible with Equations (6.1). Reference BindschadlerBindschadler (1983) found field data to be reasonably consistent with the same power law, φ ~ ∧1/3. It must be said, however, that these results are hardly corroborations of Equation (6.5). Budd and others’ experimental results are so far removed in scale from glacier sliding that it may be doubted if even the same physics is responsible. Reference BindschadlerBindschadler (1983) did not actually present his data, and it may be that a power law is not the best-fit law.
For beds which have a cut-off in “power” at high and low wavelengths, a more general form φ(∧) may be obtained from Equation (3.40), φ(ᴧ) = c 1 F(1/ᴧ) + c 2ᴧ1/n G(1/ᴧ). Typically, this is a monotone increasing function from φ(0) = 0 to φ(∞) = F < ∞. The behaviour is thus broadly similar to the power law in Equation (6.4), except that φ saturates as u → ∞, leading to the solid friction law (Reference ReynaudReynaud, 1973)
However, Equation (6.7) is not sufficient on its own for a complete analysis, since τ = 0 as u → 0. For constant N, this sliding law is typically as shown in Figure 6. For many applications, where the drainage pressure is atmospheric, then N* ≈ H, the (dimensionless) depth, and also τ ≈ H (neglecting variations of surface slope due to bedrock topography), so that Equations (6.1) become
In the particular case that φ saturates, φ(∞) = F, then (6.8) indicates that no finite sliding velocity exists if
i.e. the bed is too smooth. This is a conceptual draw-back to a solid friction law, and would lead to unlimited acceleration. In reality, corrective terms involving surface-slope variation and longitudinal stress would yield a meaningful law, but then the velocity would increase rapidly with distance down-stream, and it seems unlikely that such a situation could occur in a steady state. Since the approach to a saturation limit as ᴧ → ∞ assumes that cavitation extends to all wavelengths, and floods the bed, it may be less realistic on the whole than the assumption of a continually increasing φ. Notice that the condition (6.9) is that v ≤ ϵ, which essentially says that locally the whole bed slopes down-stream. Now one can imagine that this condition may hold locally, but not globally, due to the common occurrence of pinning points (controlling zones), where the flow is effectively resisted, e.g. bedrock rises or riegels. We thus reiterate the point of view of Reference Lliboutry and ReynaudLliboutry and Reynaud (1981). Such pinning points are the large-scale cavity-free bumps of the theory, and it is simply a matter of what length scale they occur on which decides whether a “local” friction law such as Equation (6.5) is globally appropriate with the usual τ ≈ H, or whether the large-scale flow with variable topography must be solved as well to determine the basal shear stress.
A possible sliding law with constant N, leading to “solid” friction at high velocities.
There has been much interest recently in the effect of water pressure on the sliding velocity (Reference IkenIken, 1981; Reference Iken and BindschadlerIken and Bindschadler, 1986). The latter authors found quasi-steady velocities to vary with effective pressure in a way consistent with Equation (6.5) at fixed stress τ. However, Iken reported, at the recent Interlaken meeting on “Hydraulic Effects at the Glacier Bed”, that later measurements in 1984 and 1985 seem rather different.
In order to understand the surging process in Variegated Glacier, Kamb has developed a model of the subglacial drainage system rather similar in concept to that presented here (see Reference KambKamb and others, 1985). Similar work has been done by Raymond and Harrison, and was reported at the Interlaken meeting. In concept, it seems closer to the simple model adopted here than Kamb’s detailed treatment. The following picture is suggested, which gives a physical description of the analysis of this work. Cavities are always present, over sufficiently small bumps, or steps. Joints in the bed ensure these are interconnected but at low velocities drainage is primarily through a central Röthlisberger channel. Drainage through cavities is very small, because the joints, etc. are small and tortuous, and viscous heating in the cavity system is negligible.
The stability of this drainage system depends on the flow rate, as follows (I am grateful to C. Raymond for this observation). We examine the stability of the system considered as two reservoirs (cavities and tunnels). These are connected by the joints and striae (Kamb’s “step cavities”). Let p K be the cavity-system pressure, and p R be the tunnel pressure. If p K and p R are different, water flows from one reservoir to the other. Suppose V K and V R are cavity volume and tunnel volume per unit length. Presuming a constant flux through the system, one finds dp K/dV K > 0, dp R/dV R < 0. For this reason, a single tunnel on its own is unstable, and many linked cavities are stable. Now suppose a steady-state p K = p R is perturbed, say by increasing p R by Δp R > 0. Then water flows from tunnel to cavity, leading to an increment of cavity volume, ΔV K > 0, and thus an increment of cavity pressure Δp K > 0. If Δp K > Δp R, the system is stable, and will return to its steady state. Specifically, stability is ensured if dp K/dV K > |dp R/dV R|. The right-hand side of this inequality is independent of sliding velocity, but of course V K depends on p K and u, in fact ∂V K/∂p K > 0, ∂V K/∂u > 0. For example if V K = p K u, then stability ensues if u < |dV R/dp R|. In this way, one finds that the stability of the combined tunnel/cavity system depends on sliding velocity.
When the combined system breaks down, tunnel(s) rapidly collapse, and the drainage is forced through the joints at the bed. These then must enlarge by viscous heating and eventually a system of interconnecting Röthlisberger-type channels will connect the cavities. The resultant drainage pressure can be calculated, and is much higher. We have implicitly assumed that the cavity system becomes unstable when the flow velocity is sufficiently reduced, though this is not obvious. A detailed examination of such instability has been carried out by Kamb, who found that cavities were unstable if the sliding velocity u is low enough. It seems likely that this mechanism of instability is not essentially the same as that considered here.
The analysis of the Variegated Glacier surge by Reference KambKamb and others (1985) focuses on the detailed day-to-day variation of velocity and water pressure. In order to study surges from the point of view of long time-scale, nonlinear oscillations, one needs a description of basal sliding which can incorporate quasi-steady fluctuations in water pressure. This is provided by the incorporation of the drainage pressure as a function of sliding velocity, and when this is done, one finds that the monotonic τ versus u curve of Equation (6.5) becomes multi-valued, due to the drop in effective pressure when the tunnel system disappears.This multi-valuedness alone is sufficient to understand completely the large-scale properties of surges (Reference HutterHutter, 1982), that is, why they recur periodically, why they have two distinct time-scales of behaviour, and so on. The issues of how this can be shown, and what parameters may control surges, will be developed in a separate paper.
For different glaciers or ice sheets, the important parameters which control sliding are µ 2 and µ 3, given by Equations (6.2). In addition, the subglacial drainage system will switch to a linked-cavity system if ∧ in Equation (6.1) is large enough (see Equation (4.35)). For glaciers, we might typically expect µ 2 ≥ 0(1), so that ᴧ ~ 0(1). Then u ~ 0(µ 3), and if x 0 ~ 1 m, d ~ 100 m, y 0 ~ 0.2 m, ϵ ~ 0.1, we have µ 3 ~ 50, so that one expects large sliding velocities (compared with shearing motion). For such beds, the bulk of the resistance may then in fact be due to pinning points (e.g.riegels) as mentioned earlier (personal communication from I. Whillans). However, when sliding velocities are predominant, it is no longer appropriate to choose the velocity scale U from the flow law, as in Equation (3.6). The correct method has been given by Reference FowlerFowler (1982, appendix A) and is tantamount (when µ 3 ≥ 0(1)) to putting τ = (1/µ 3)1/n τ, N * = (1/µ 3)1/n N *. Then µ 3 measures the importance of sliding to shearing and the Sliding law is given by Equation (6.1), ᴧ = u/N *n , and u, N * ~ 0(1). Instability still occurs if ᴧ > ᴧc given by Equation (4.36), since S R is determined directly by the flux Q R, see Equation (4.4). We thus have he following picture. For temperate glaciers with smoothlarge-scale bedrock topography, we expect sliding to be predominant, and of 0(µ 3) times the differential shearing motion. Tunnel instability leading to reduced effective pressure and thus surges, will occur if ᴧ > ᴧc, that is for small water outflux, wide glaciers, steep slopes, and smooth beds. In addition, surges require that the linked cavity water pressure δN R is less than overburden. A quantitative compendium of these predictions and comparison with actual characteristics of surging glaciers is deferred to a separate paper.
Lastly, let us comment on the assumed morphology of the bed. Our analysis has assumed a hard bed, without basal debris. This is the classical assumption, and is not generally realistic. Rock fragments and erosive debris may provide an active (deformable) sub-sole drift, or may be incorporated into the basal ice. In either event, effective water pressure will affect sliding, and it is reasonable that any sliding law will still have ∂τ/∂u > 0, ∂u/∂N < 0 as for a hard bed. Thus, one does not expect soft beds per se to affect the qualitative shape of the sliding law. For surging sub-polar glaciers, e.g. Trapridge Glacier (Reference Clarke, Clarke, Collins and ThompsonClarke and others, 1984) the important difference is in the drainage mechanism which may be blocked off by a cold ice front. Thus a tunnel system is inappropriate, and drainage is likely to be through a permeable bedrock.
In this case, explanation of surges by the kind of switching mechanism between drainage systems, described byReference KambKamb and others (1985), requires some kind of mechanism for blocking the escape route. For a given water flux through the glacier, p w (the pore-water pressure) remains constant, determined by the subglacial bedrock permeability. Thus, as the glacier thickens, N increases, and therefore the sub-sole drift will begin to consolidate, thus decreasing its permeability. One may speculate that eventually the drift permeability becomes so low that water is unable to escape and so the bed is flooded, p w increases rapidly, and a surge is initiated. These ideas follow those of Clarke (Reference Clarke, Clarke, Collins and ThompsonClarke and others, 1984), as reported at the Interlaken meeting.
In conclusion, it seems that a (relatively) unified mathematical treatment of sliding and surges may be feasible, but the detailed physics is likely to depend on the particular kind of bed considered.
Appendix
Consider the system in Equations (4.30) and (4.34):
where
There is a steady-state solution
to study its stability, we put
and linearize for smalln 1, n 2 to obtain
Solutions proportional to exp(λt) exist, where λ is an eigenvalue of the matrix
thus
whence
It follows that the co-existing system is unstable if λ2 > 0, that is if g1 > g 2 N, i.e.
simplifying, we obtain Equation (4.35).
Acknowledgements
This work was reported at the Interlaken meeting on “Hydraulic Effects at the Glacier Bed”, September 1985. Many thanks to A. Iken and H. Röthlisberger for a splendid time, and to C. Raymond for interesting discussion.
