1. Introduction

In seeking the simplest realistic model of glacier flow, an attractive procedure is to analyse a given set of equations and boundary conditions in a mathematically consistent fashion, rather than make physically plausible, but nevertheless ad hoc, assumptions and approximations. The standard method of doing this consists of first non-dimensionalizing the variables, using the given parametric “inputs” to the problem, and then making formal asymptotic approximations to the resultant model using the appropriate parametric limits that suggest themselves from the numerical order of the various dimensionless parameters which naturally arise. It can be checked a posteriori whether terms neglected in this fashion are in fact small, and in this sense the procedure is a rational one: no a priori assumptions are made.

A rational derivation of a “reduced” model of a two-dimensional, incompressible glacial ice flow has been discussed at length elsewhere (Reference Fowler and LarsonFowler and Larson, 1978), and so we content ourselves here with merely presenting the resulting equations and boundary conditions, and commenting on some of the novel features which they introduce: this is done in Section 2.

In Sections 3–5, we consider how the proposed model may be used to examine various phenomena of dynamic interest. Specifically, we treat the kinematic waves studied by Reference NyeNye (1960) from a non-linear point of view; we consider the important dynamic effect of using as a basal boundary condition a sliding law which is continuously dependent on the temperature below the pressure-melting point; and we show how consideration of a kinematic wave equation can explain the essential features of the observed seasonal waves (Reference HodgeHodge, 1974), on the basis that cavitation has a major effect on the sliding law, as is suggested by previous work (Reference FowlerFowler, unpublished).

2. Mathematical Model

We consider a two-dimensional incompressible ice flow whose geometry is typically as shown in Figure 1. The head and snout of the glacier are identified by the points x ₀ and x_S; x _Q, X_Z ,and x _M are identified below. In general, the glacier will consist of distinct cold and temperate regions, divided by a melting surface; it should be emphasized that we are not necessarily assuming that Figure 1 gives an accurate representation of every glacier’s temperature profile, but it serves to portray the fact that we must normally expect regions of both cold and temperate ice to exist. The shape of the melting surface is to be determined in the solution of the problem, and will not necessarily be as in the figure.

Fig. 1.

Glacier geometry.

If we let the two-dimensional velocity field q have components u and v parallel to the x and y axes respectively, then, using suffixes to denote differentiation, the continuity equation expressing conservation of mass,

is automatically satisfied by defining a stream function ψ such that

In the cold zone, one then obtains coupled equations for the stream function ψ and the temperature T. It is convenient to introduce a change of coordinate

where y ═ η(x, t) is the equation of the top surface of the glacier, so ξ measures the distance downwards from the top surface, and an associated change of variable

where H(x, t) is the depth (measured perpendicular to the x-axis). In a shallow, shearing flow such as in a glacier, the constitutive relation between the strain-rate tensor e_ij and the stress tensor T_ij is approximated by the same relation between the shear strain-rate u_y ═ Ψ_yy ═ Ψ_ξξ, and the shear stress T ₂. Since the surface shear stress is approximately zero, the balance between T _2y and gravity implies that T₂ is proportional to the vertical coordinate ξ; when the variables are non-dimensionalized, scaled, and reduced (by setting certain small parameters equal to zero), the usual Glen’s law for ice with exponent n may be written in the approximate form

The energy equation is similarly found to be

where the left-hand side is the convective derivative (═ T_t+uT_x+uT_y), and the two terms on the right represent the viscous dissipation (strain-rate times stress) and heat conduction (the vertical component of which dominates the horizontal component due to the shallowness of the flow).

The conditions that we impose for these equations on the boundary of the cold zone are:

• (i) on the top surface, the surface accumulation/ablation rate and temperature are prescribed;

• (ii) on the melting surface, the temperature is equal to the pressure-melting point, and the heat flux is continuous;

• (iii) on the (unknown) base, a geothermal heat flux is specified, the normal velocity is zero, and the tangential velocity is a prescribed function of the basal stress and the temperature.

If we define the surface-flux function s(x, t) by

where a is the accumulation/ablation rate (a > 0 in the accumulation area, a < 0 in the ablation zone), then the kinematic boundary condition at the top surface requires

When this is transformed using Equations (2), (3), and (4) we find

where we may define Ψ to be zero at x₀. We also require that

On the melting surface, scaling shows that the melting temperature is effectively constant (equal to zero when normalized) and if we assume there is no heat flux into the temperate ice by moisture transport away from the melting surface, suitable conditions there are

On the bedrock ξ ═ H(x, t), the no-flow-through condition implies that ξ ═ H is a streamline, ψ═ 0 there, that is, using Equation (4),

A sliding law for the horizontal velocity u ═ ψ_y ═ —Ψ_ξ which is dependent on temperature T and basal stress ≈ H (the depth) is

Lastly, with zero geothermal heat flux (from scaling considerations), the heat flux into the cold ice above is equal to the viscous heating generated by basal sliding: this implies

which is valid until xz where the basal temperature reaches melting point, beyond which the appropriate condition is

The approximations (and their physical interpretation) on which these equations are based are discussed below.

The above equations and boundary conditions are not relevant in zones of temperate ice: in these the temperature is effectively constant (T ≈ o°C) and the role of enthalpy variable is taken on by the moisture content (Reference LliboutryLliboutry, 1976). It appears that the flow law of temperate ice also depends on the moisture content, and so the equations for ψ and the moisture content w are once again coupled. In this case it becomes necessary to specify a term in the energy equation which describes the hydrological flow of moisture through the ice (this plays a similar role to that of heat conduction in cold ice). It is clear that a proper description of such a term is necessary before the dynamics of temperate ice can be usefully studied. In the earlier paper (Reference Fowler and LarsonFowler and Larson, 1978), it was assumed, for want of any better information, that the transport by this means was negligible; in this case the cold zone effectively uncouples from the solution in the temperate zone, and the latter is of little further interest: subsequent analysis is then largely concerned with the cold zone—or with completely polar glaciers.

It will be noticed in Equation (13) that the sliding law is taken to be a function not only of the basal stress (≈H) but also of the temperature T: this is to accommodate the realistic physics (Reference FowlerFowler, 1979) which demands that the basal velocity should increase from zero to its full temperate value over a small (typically 10⁻¹ deg) but crucially finite temperature range just below the pressure-melting point: the important consequence of this novel assumption is discussed below in Section 4. The points xQ, x_Z, x_M in Figure 1 and in Equation(14) may now be interpreted as follows: x_Q is the point where the basal ice first begins to slide; the ice is wholly frozen to the bedrock in (x₀, x_Q). The point x_z is where the temperature reaches melting point and the sliding law is the fully temperate one; thus (x_Q, x_Z) is the region of “sub-temperate” sliding. Lastly, x_M is the point where the melting surface on which T ═ 0 “breaks away” from the bedrock; in (x_Z, x_M) the basal temperature is at the melting point, but the viscous source heating at the boundary is used up in warming the cold ice above until T_ξ ═ 0 (i.e. the heat flux into the ice decreases to zero) which is precisely when x ═x_M. The melting surface ξm dividing cold and temperate ice is a necessary constituent of the model, since in general cold ice will be heated by the frictional dissipation source term in the energy equation. It is however quite realistic to consider particular limits represented by the so-called “polar” and “temperate” glaciers. In view of the comments above on moisture transport, realistic models for the latter cannot be said to exist.

The assumptions made in deriving the reduced model presented above are the following: the Reynolds number (Re) ═ 0; δ (typical depth/length) ═ 0; μ ═δ/(mean bedrock slope) ═ 0; λ (dimensionless geothermal heat flux) ═ 0; and that Z, which measures the error in assuming an exponential form exp [kT] for the temperature dependence of the flow law rather than the more usual Arrhenius term, is zero. These approximations are based on the numerical estimates

The Reynolds number (Re) measures the slowness of the flow, and its neglect entails ignoring inertial terms in the momentum equations; δ is a measure of the shallowness of the flow; neglecting μ means that variations of the surface slope from the mean bedrock slope are ignored. (Note that this neglect of μ. is inapplicable to a flow with a horizontal base, e.g. an ice sheet, to which the definition of μ is not appropriate). The above approximations are expected to be good ones, in the sense that they are mostly of non-singular type; this is not strictly true of the neglect of δ and μ, but neglect of δ only appears to be associated with a degeneracy of the equations near the head x₀ and snout xs, and does not affect the bulk of the flow. The parameter μ represents a diffusional term for kinematic waves, and thus may become prominent if its neglect leads to the prediction of surface shock waves: for further discussion of this point, see Reference Fowler and LarsonFowler and Larson (in press).

The important parameters arising in the model are k, β1, and β ₂. These measure respectively the strength of the dependence of viscosity on temperature, the magnitude of the viscous source heating, and the magnitude of heat conduction. With typical values of the physical input parameters, it is found that these parameters are of numerical order one, although with a realistic variation of such quantities as mean accumulation-rate, mean bedrock slope, etc., larger or smaller values may easily be obtained: in particular, β1, and β ₂ may be small. Furthermore, when the analogue of Equation ((6) is studied for larger ice masses such as ice sheets, it is expected that realistically small values of β₂ may be obtained. For these (and purely mathematical) reasons, it is of interest to study the proposed model under various asymptotic limits of the given parameters. In general, this is not an easy task.

3. Kinematic Waves in the limit K → 0

If we suppose that the surface-flux function s(x, t) is independent of t, that is s(x, t)≡s(x), then on integrating Equation (5) twice and using the boundary conditions for Ψ, Equations(9), (12), and (13), we find that in x < x _m conservation of mass takes the form

where the flux Q is given by

and ub ═ F(H, T) is the basal sliding velocity. Equations of the form of Equation (16) were studied by Reference Lighthill and WhithamLighthill and Whitham (1955) and applied to glacier flow by Reference WeertmanWeertman (1958) and Reference NyeNye (1960, Reference Nye1963). There is a slight subtlety in the present case, since the basic (“datum”) steady state is of finite extent, and so linear analyses of the type proposed by Nye lead to non-uniformly bounded solutions at the snout, where H → 0. One can easily avoid this pitfall by applying the method of characteristics to the essentially non-linear system of Equation ((16)16). As an illustration, we suppose x < x_Q , that is ub ═ 0, and take the formal asymptotic limit k → 0. In this case we obtain

of which the characteristic equations are

Equation ((19) immediately gives the wave speed, which is approximately four times the surface speed (this being Hⁿ⁺¹/n+1). The solution of Equations (19) and ((20) may be written down in terms of a characteristic parameter σ ϵ [x₀ , xs(o)]. It is

and is valid for all times in which σ as defined by Equation (22) is a single-valued function of x and t: if σ becomes multi-valued, then one can make the solution single-valued but discontinuous by the insertion of appropriate shocks (which would physically be smoothed out by considering a non-zero μ). It is possible to show that such shocks will form if the initial data s₁ (σ) is such that s₁' (σ) > 0 for some value of σϵ[x₀, xs(0)]: if the initial data represents a sudden change in climate, this condition may be interpreted as representing an initial increase in the accumulation-rate (or decrease in the ablation-rate).

Whether shocks form or not, Equations(21) and ((22) imply that any initial disturbance reaches the snout in a finite time, and for small disturbances, all shocks which do form must do so near the snout (where in any case the equations are not strictly valid): over the remainder of the glacier one can obtain an explicit form for the travelling wave nature of the solutions by linearizing the characteristics rather than H. If we denote the steady state by H ═ H₀(x), these solutions may be written (Reference Fowler and LarsonFowler and Larson, in press) as

where ø ⪡ 1 is the initial perturbation from the datum profile.

The form of these travelling waves incidentally shows that the glacier profile is stable to small perturbations. In the same limit k → 0 with u _b ═ 0, one can also show that the steady-state temperature profile, although not explicitly known, must also be linearly stable: this is done by Reference Fowler and LarsonFowler and Larson (1978), who show that (at least for completely polar glaciers) consideration of non-zero k, i.e. temperature-dependent viscosity, is a necessary constituent of any realistic discussion of thermal instability, proposed by Reference RobinRobin (1955) as an initiating mechanism of glacial surges. In other cases where the temperature and flow fields do happen to be stable, we expect that the mechanism of kinematic waves proposed above will remain qualitatively valid even when u _b and k are non-zero.

4. Dynamic Effect of a Temperature-Dependent Sliding Law

The apparently minor refinement of prescribing a sliding law which is continuously dependent on the temperature near the melting point, rather than discontinuous there, leads to a major change in the nature of the temperature field and dynamic behaviour of a glacier. It should be emphasized that a continuous sliding law was included in the model, not because any such effect was envisaged, but simply because the physics seemed to demand it. This illustrates the need for formulating an appropriate model before taking any limits that suggest themselves, since in this case it is found that the limiting behaviour of a glacier as the sliding law becomes discontinuous is not the same as that of a glacier with a discontinuous sliding law. Such non-uniform limits are commonly encountered in other areas of applied mathematics; in the present case the above discrepancy occurs because the discontinuous sliding law has a finite jump at T ═ 0 and no intermediate basal velocity between 0 and F[H, 0] is possible, whereas in the limiting form of the continuous law, any such intermediate velocity is admissible at T ═ 0, and thus the discontinuous law does not really represent the correct limit of the continuous law. More specifically, if we admit a discontinuous sliding law of the form

then, since the basal velocity is discontinuous, there must be some kind of discontinuity in the surface behaviour, and if μ ═ 0, this takes the form of a finite jump in the depth at the point x Q (═ xz in this case). Now if we examine the particular asymptotic limit of the steady-state versions of the model given by Equations (5) and (6) in which β ₂ →∞ with β1/β2 finite (admittedly this is unlikely to be a realistic case), it is found that just as for k ═ 0, the equations uncouple (the convective terms vanish) and the solution for the entire problem in the cold zone may be constructed from a knowledge of the temperature, which may be found explicitly as the solution of a second-order differential equation in ξ which is parametrically dependent on x through the boundary data on ξ ═ 0 and ξ ═ H. It is then found from this explicit solution that if (—v, 0) is the dimensionless temperature range over which the sliding velocity increases from zero to its full temperate value, i.e. T ═ —v at xq and T ═ 0 at xz, then in the limit as v → o, xQ → xz, and so the glacier behaviour with v ═ 0 is not the same as that when v →o. A reasonable estimate for v is 10⁻², so we see that v is indeed small. From the above result, it is clear that rather than apply the boundary condition u _b ═ F[H, T] in (x_q, x _z), we should in the limit v → 0 specify that T ═ 0 there. Thus the effective boundary conditions to replace Equations (13) and (14) are

and we see that there are three separate sets of boundary conditions to be applied on ξ ═ H. Equation (25) is the correct limiting version of Equation (13); it obviates the need to specify F[H, T] for T < 0, which is just as well, and it is fairly clear that a numerical solution of the equations will be much simpler with these effective boundary conditions. Questions of stability are slightly more subtle. These are to be discussed in a paper in preparation by Reference Fowler and LarsonA. C. Fowler and D. A. Larson.

What happens in (x _q , xz) is clear: the basal velocity is essentially governed by the bulk ice flow rather than the basal “inner flow” (which really specifies the precise basal temperature in (—v, 0)), and u_b gradually increases in (x _Q, x _z) until it becomes equal to the full sliding velocity F(H, 0). Thus we must expect that, in general, glaciers may have substantial parts of their bases at temperatures within this “sub-temperate” range, and accordingly basal velocities within these regions will not appear to be functionally related to the basal stress (and hence the depth). Although this result is explicitly derived in the particular limit β ₂ —> ∞, it depends only on the continuity of the sliding law, and hence will be valid for any value of β₂ whatsoever, as long as v ⪡ 1. This result in itself justifies the study of the limit β ₂ → ∞, and answers negatively the question posed by Reference RobinRobin’s (1976) paper—at least in certain regions.

5. Seasonal Waves

Seasonal waves are waves of velocity disturbances which propagate at very high speeds down glaciers (typically 20 to 150 times the surface speed). They have been reported by, amongst others, Reference Deeley and ParrDeeley and Parr (1914) on Hintereisferner and Reference HodgeHodge (1974), who observed such waves on Nisqually Glacier with an approximate speed of 20 km a⁻¹.

To date, no satisfactory theoretical explanation of their behaviour has been forthcoming. Reference WeertmanWeertman (1962) considered the effect of kinematic waves in the subglacial water film, while Reference HodgeHodge (1974) considered that the waves might be due to a correlation between the sliding velocity and the liquid water stored in the glacier, since the latter varies on a seasonal basis. While this may be an important part of the phenomenon, it does not in itself answer the fundamental question: why are the wave propagation speeds so large? And why do the waves consist of velocity perturbations with no apparently related wave motions in the depth?

One possible answer to these questions is provided below. We will essentially concentrate on Hodge’s paper, and in particular seek to understand his figure 6, which is a contour map in (x, t) space of the surface velocity. The lines of constant velocity mostly oscillate in space, but some form closed loops: the contours are skewed slightly, which represents a seasonal wave with a transit time of about a month from the equilibrium line to the snout.

Our analysis is based on the supposition that the basal sliding velocity becomes a rapidly varying function of the basal stress, and hence the depth, when subglacial cavitation occurs. This was suggested as a theoretical possibility by Reference LliboutryLliboutry 1968, and it was shown by Reference FowlerFowler (unpublished, Reference Fowler1979) that in the particular case of Newtonian flow over a sinusoidal bedrock of small roughness, the leading-order solution for the basal flow leads to a sliding law in which the sliding velocity increases without limit at a constant stress when cavitation is present. It should be emphasized that this result may well be importantly affected by the consideration of a more realistic bedrock in which more than one obstacle size is present. Nevertheless, the important point is that (for this particular case) the velocity increases rapidly with stress when cavitation is present, and it seems reasonable to suppose that this qualitative feature of the sliding law may well be representative of the flow over more realistic bedrocks also. Based on this, we here suppose that the sliding law is of the form shown in Figure 2, that is the velocity increases rapidly with basal stress but not at an infinite rate. We suppose that μ ═ 0 as before so that the basal stress T_b is approximately equal to H, and we assume cavitation sets in when H ≈ H_c. We will formally assume that shearing within the ice mass is negligible, so that the flux

and Q increases rapidly near H_c. (This assumption is not necessary for the subsequent analysis.) The kinematic wave equation describing conservation of mass is then

if we suppose that the flux function s is independent of time. The steady-state solution of Equation (27) is

Fig. 2.

Possible form of the sliding law when cavitation is present.

If we denote the critical value of Q at H_c by Q ,c, then in regions where Qc < s(x), the glacier will be of effectively constant depth Hc. In such regions, a small perturbation in the depth will generally have a large effect on Q (and hence the surface velocity).

Now the characteristic equations (see, e.g. Reference WhithamWhitham, 1974) of Equation (27) may be written as

The first of these immediately states that finite disturbances in the velocity of the type described above propagate at a speed Q’(H). From what has been said, Q’(H) ⪢ 1 in regions where Q > Q c, and so these velocity waves will travel down glacier at very high rates without any apparently related waves in the depth profile (since perturbations in H are of small magnitude): this seems to explain immediately the basic phenomena observed. Of course, no analytic estimate of Q’(H) is available when cavitation is present, and to this extent the above theory remains hypothetical (as compared, for example, to Reference NyeNye’s (1960) specific prediction of a wave speed of three to five times the surface speed): to the wave speed of 20 km a⁻¹ described by Reference HodgeHodge (1974) corresponds a value of Q’(H)≈ 200. However, it is worth pointing out that there does not appear to be any other obvious dynamic mechanism which can supply the high rates of propagation observed (mere seasonal variation of the inputs cannot do this). The only other possibility seems to be that the motion becomes highly enhanced due to the presence of appreciable quantities of melt water (as suggested by Hodge): a quantitative description of such a proposed process must await the appropriate modelling of hydrological moisture transport.

If we write Q’(H) as a function of Q, say

then we see from Figure 2 that c(Q) will have the form shown in Figure 3, with critical behaviour at Q ═ Q c. Multiplying Equation (27) by Q’(H) ═c(Q), we obtain the wave equation for the flux in the form

Fig. 3.

Schematic form qf Q’(H) as afunction of Q.

For regions in which Q > Q c, we have by assumption

If we define a small parameter ϵ by considering c(Q) ≈ 1/ϵ, and let

then the correctly scaled version of Equation (31) is

valid for the description of seasonal waves over the fast time scale T. Equation (34) is applicable in regions where s(x) > Q c. Denote such a region by (x ₁, x ₂), and suppose s(x) < Q c outside this region. In the limit as ϵ→ o, we require that flux changes outside (x ₁, x ₂) should occur over the longer convective time scale t, and so appropriate conditions on Q are that Q ═ Q c at x ₁ and x ₂. We can satisfy the up-stream boundary condition at x ₁ , but not in general the second: this may require the reintroduction of the diffusion-like coefficient μ, which would correspondingly introduce a small term proportional to Q _xx in Equation ((34), and an associated local analysis in the vicinity of x _2; we shall not consider this subtlety further in the present instance.

For want of any better information, we now choose

to illustrate the sort of behaviour we expect. For Nisqually Glacier, a typical value of ϵ is 5X10⁻³: the bulk-flow time scale may be estimated at 20 years, so that this ϵ corresponds to a transit time down-glacier of c. 1 month (as observed). With Equation (35), Equation (34) is

with solution

where ϕ represents a travelling wave generated by an initial disturbance ϕ(x), ϕ(x ₁ ) ═ 0. Let us now seek the solution corresponding to a regular variation in the flux function. Like Reference Fowler and LarsonFowler and Larson (in press), we model such a seasonal variation by adding a sinusoidal term to s(x), and thus seek solutions of

where the real part of Q is to be taken.

With ϵ ≈ 5 × 10^{– 3} and t associated with a time scale of 20 years, T is associated with a scale of c. 1 month; the seasonal variation defined by Equation ((38) has a dimensionless period (in T ) of 2 π/Ω ≈ 6/Ω: since this corresponds to 1 year ═ 12 months, we have 6/Ω ≈12, i.e. apparently Ω ≈ ½ for Nisqually Glacier. Typically we may take Ω ≈ 1. (Contrast the effect of seasonal variation on Nye’s kinematic waves, where the seasonal frequency relative to the long time scale is so fast that its effect is “averaged out” and may be neglected (Reference Fowler and LarsonFowler and Larson, in press).)

The solution of Equation ((38) is (if we impose Q ═ s(x) at x ₁),

Since Ω≈ 1, Equation (39) represents a regular seasonal wave of finite amplitude (if s’ ₁ ═ 0(1)). To be more specific, let us choose

where a is constant. Equation (40) reflects a steady oscillation between extra accumulation and ablation: realistically we expect a≈1. We then find, on substituting Equation (40) into Equation (39), evaluating and taking the real part of the resultant expression, that

where for convenience we choose x 1 ═ 0. Thus Q may be thought of as a modulated wave of speed 2, or as a travelling wave of constant shape and speed 1 together with a superimposed oscillatory component. Reference HodgeHodge (1974) gives data on a relationship of this type in the form (his fig. 6) of a contour map of the surface velocity in the (t, x) plane. We may obtain the characteristic features of such a map in the present case by writing

and seeking intersections in the (T, X) plane of the two surfaces

and

for various different (constant) values of Q (corresponding to different values of the surface velocity). The surface described by Equation (44) is a smooth sheet which, for reasonable functions s bends concavely upwards in the X direction (if we suppose x ₁ lies in the ablation zone). On the other hand, sin X cos T represents an “egg-box” curve in the (T, X) plane consisting of a checkered formation of undulating peaks and troughs. Changing this to sin X cos (T—X) simply has the effect of skewing the undulations in the T direction, but does not materially alter the properties of the intersection of Equations (43) and (44). The intersection determines two types of curves depending on the value of Q. As T—X increases from zero on curve (1) in Figure 4a, the locus of the intersection which lies on the point b₁ on the peak (solid line) decreases in X; hence Z_z decreases until Z ₂ ═ O when the locus leaves the peak and enters the trough, the axis of which is the dotted line in the figure at T—X ═ π. X further decreases till the point a₁ is reached at T—X ═ π, after which the locus retraces its steps to B₁ at T—X ═ 2π. Thus the locus is an oscillating curve in (X, T—X) space.

Fig. 4.

a. Sectional view of intersection of the two curves representing Equations(43) and (44).

b. Possible form of contours Q ═ constant in the (T, X) plane.

The same is true of the locus bounded by the points A₂ B₂ on curve (2) in Figure 4a. However, that bounded by M and N consists of a series of closed loops, because these loci cannot leave the peaks on which they lie, since MN does not cross the X axis. The curves bounded by MN, A₂B₂, and A₁B₁ are shown from top to bottom in Figure 4b, and one can immediately note the resemblance to figure 6 of Reference HodgeHodge (1974): however, there are also striking differences, and one should not attempt to seek too great a quantitative comparison.

An alternative way of examining Equation (41) in relation to this figure is to consider the “height” Q to be the superimposition of the decreasing function s(x) and the skewed egg-box. The overall “topography” is then that of a slope downwards to the snout x ═ xs, with various hills and basins evident due to the seasonal variation: again this is essentially what is represented in Hodge’s figure. Of course the present model is a huge simplification of any realistic seasonal waves, and cannot be used for predictive purposes (at least not quantitatively). Even so, there are two particular facets of the above explanation that must be considered in terms of the qualitative applicability of the theory. One is that a much larger value of Ω than that given here (≈½) must be invoked in order to obtain the secondary “loops” evident in Hodge’s figure: this is probably not very important in view of the quantitative assumptions of the model. Of more importance is any possible discrepancy between the phase of the level curves in Equation (41) and in Hodge’s diagram, since as he says, “the acceleration of the glacier throughout the winter in the ablation zone is crucial to a correct interpretation of the velocity variations”. It is therefore important to try and relate the phase of Equation (41a) to that of the experimental data. We see that sin ½Ωx is a maximum when x ═ π/Ω. At this point

Now we see from Equations (38) and (40) that the variation in the seasonal accumulation-rate is a cos Ωτ, so that the origin τ ═ 0 corresponds to the time of maximum accumulation-rate, i.e. winter. On examining Hodge’s figure 6, it seems that in the above context, a suitable choice of “x ═ π/Ω” is about 300–400 m along the centre line, and here the Q profile in τ is indeed “sinusoidal” with the origin taken about February. It is hard to be more specific, but this seems to agree to a reasonable extent with the theory. (It is difficult for example to say where the origin x ═0 at which the waves are initiated should be put.)

More detailed analysis is necessary to determine whether such apparent phase agreement is adequate to explain the observed results: this is rather important in view of Hodge’s idea that liquid water storage, in particular at the glacier base, may have a major effect on the dynamics of temperate glaciers.

Lastly we consider Hodge’s figure 13, in which the surface speed averaged over 19 roughly equally spaced x-posilions is plotted against time. We compare this with taking the average in x of Q as defined by Equation (41b), with the assumption that x ═ xs is close to a period of sin Ωx. The result is (denoting averages by overbars)

which is very similar to the second graph in Hodge’s figure 13 for 1969, if we choose an origin τ ═0 at about March. This is again suggestive rather than predictive, but the agreement is nevertheless quite satisfactory.