McCall Glacier

McCall Glacier is located at 69° 17′ N, 143° 50′ W in the Romanzof Mountains of the northeastern Brooks Range, Alaska. Il is situated on the northern front of the mountain range, about 100 km south of the Arctic Ocean across the coastal plain (Fig. 1. inset). The glacier is about 8 km long and has an area of 7.4 km². Ice originates from three cirques referred to, from northeast to southwest, as upper (UC), middle (MC) and lower cirque (LC) and it flows from an elevation of more than 2700 m to the terminus at 1350 m (Fig. 1). The surface is relatively steep, with a mean slope of 8°: there is a series of bulges along its length with slopes up to 15 separated by treads of 3- 5°. Seasonal mass exchange is quite small: winter balance is around 0.2 ma⁻¹ and summer balance is around −0.6 m a⁻¹. Because of mountain shading and wind- deposition patterns, the concept of an equilibrium- line altitude (ELA) is ill- defined on this glacier: the equilibrium line spans an elevation range of 350 m on average (Reference Trabant and Benson,Trabant and Benson, 1986; Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995). The mean annual air temperature is about −12° C at 1700 m elevation.

Fig. 1. Bed topography of McCall Glacier, interpolated from radio echo- soundings measurements made at the circles. The heavy line shows the glacier margin in 1958. Elevations are in meters above m.s.l. The ground traks of two elevation profiles along the glacier are shown.

Ice temperatures in the accumulation area arc between −1° and −15° C, both near the surface and at the base (Reference Orvig, and Mason,Orvig and Mason, 1963). These warm ice temperatures (more than 10 K warmer than the mean annual air temperature) are caused by latent- heat release of surface water that refreezes in the permeable firn. In the ablation area, surface- ice temperature is less than −10° C but the basal ice is temperate, at least in the vicinity of the center line. By extrapolation of a deep temperature profile at cross- section T6 (Fig. 2), Reference Trabant,, Harrison,, Benson,, Weller, and Bowling,Trabant and others (1975) predicted the 0° isotherm to be at 120 in. The present ice depth at T6, corrected for the known thinning of the glacier since 1972. suggests a basal layer of temperate ice that is about 30m thick at this location (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995). Surface temperatures suggest that this temperate layer exists beneath much of the ablation area, at least near the center line. Such an extensive temperate layer has been measured in the ablation area of poly- thermal White Glacier by Reference Blatter,Blatter (1987).

The drainage pattern observed in the 1990s is similar to that described by Reference Sater,Sater (1959). The drainage in the upper two- thirds of the glacier is almost entirely supraglacial. Meltwater drains by means of two major surface streams, one located near the center line of the glacier and the other on the east side. Both streams originaicai about 2200 m elevation in the lower and upper cirques, respectively. Below the confluence of the cirques, the central stream becomes strongly meandering and deeply incised (>10 m); it becomes englacial at about 1950 m elevalion, Qrabout 3.5 km from the snout. Crevasses in the upper ablation area, mainly formed during the cold season. fill with meltwater during summer. NVhiIe lilis process contributes to accumulation within the glacier, the dimensions of most crevasses (<0.5 m) seem to be too small to allow meltwater input to the bed.

Previous studies of the surface- velocity field of McCall Glacier, along with those of the surface geometry and mass, and energy balance, were made during the International Geophysical Year (IGY) in 1957- 58 (e.g. Reference Sater,Sater, 1959) and during the period 1969- 75 (e.g. Reference Wendler,, Fahl, and Corbin,Wendler and others, 1972: Reference Wendler, and Ishikawa,Wendler and Ishikawa, 1974). Our continuing investigations began in 1993 (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995) with measurements of surface and bed geometry, mass balance, meteorological variables, ice temperature and ice velocity.

Ice thickness

The thickness of McCall Glacier was determined in June 1993 by radio echo- sounding. A portable monopulse radar system operating at 10 MHz provided unambiguous identification of the ice bedrock interface. For a given bedrock echo, the locus of possible reflectors was taken to be a semi- ellipse, with transmitter and receiver at the focal points (Reference Echelmeyer,Eehelmeyer, 1983), The bed along a given profile is then the envelope of these intersecting ellipses. Accuracy of the ice- thickness measurements was conservatively estimated as ±10 m: this error is mainly a result of timing errors.

Cross- sections and longitudinal profiles of the glacier obtained from the radio echo- sounding arc- shown in Figure (2). The contour map of bed geometry (Fig. 1) was constructed by kriging the combined data set of the radio echo- sounding results and the mapped topography of the unglacierized valley (U.S. Geological Survey 1:63360 map). The maximum ice thickness is 250 m in the pronounced ovedeepening of the lower cirque. The average ice thickness along the center line is about 140 m. There are marked asymmetries in the valley shape near the confluence of the shallow eastern cirques With the deeper lower Cirque (cross- sectionT10): these are reflected in the velocity profiles described later. The deep trough situated on the west side of the glacier above the confluence gradually shifts towards the center below the confluence.

In a later section, we describe changes in annual ice velocity over time. As part of the interpretation of these changes, we need to know the changes in ice thickness since the 1970s. These changes were determined at about 55 locations on the glacier by comparing the surface elevation as surveyed in the 1970s by Reference Trabant and Benson,Trabant and Benson (1986) with that surveyed in the 1990s at the same horizontal position in space (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995). These elevation changes and, therefore, changes in ice thickness (assuming erosion of the bed to be small), are known to an accuracy of about ±0.3 m. The observed thinning of the glacier increases almost exponentially with decreasing elevation, from about 4 m on the upper glacier to 42 m m at the terminus (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995).

Fig. 2. Cross- sections (T−) and longitudinal (L−) profiles of McCall Claeier in 1993. The loci of potential radio- echo reflectors are shown as the dotted ellipses. Vertical exaggeration: cross- sections 2 x. longitudinal profiles 4 x.

Surface slope

Knowledge of present and past surface slopes, in addition to ice thickness, is a prerequisite for an interpretation of the velocity field and its temporal changes. The present surface slopes were determined from two elevation profiles of McCall Glacier acquired in 1993 using an airborne laser- profiling system (Reference Echelmeyer,Echelmeyer and others, 1996). These profiles approximately follow the center line, one descending through the Upper cirque to the terminus and the second through the lower cirque. Surface elevation along these profiles is accurate to about 0.3 m and measurements are spaced at approximately 1.5 m intervals along the flight path. For application to the flow- modeling described later. we have calculated the center- line surface slope over an averaging length of 150 m. corresponding to approximately one ice thickness. This averaging scale also corresponds approximately to the horizontal resolution of the depth measurements presented previously.

No continuous elevation profiles or topographic maps of this glacier were made in the 1970s and the density of surveyed markers on the surface was insufficient to determine the surface slope every where over an averaging length of 150 m. However, a 1:10 000 topographic map with a 5 m contour interval was constructed from aerial photographs taken in 1958 (Reference Bradenberger,Brandenberger, 1959). The overall vertical accuracy of this map is probably on die order of ±10 m (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995). For most areas of the glacier, this error is a slowly varying “offset” that should have little influence on the Calculation of surface slopes over a distance of 150 m. Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others (1995) have shown than more than 80% of the elevation change between 1958 and 1993 occurred between 1972 and 1993. Therefore, except in the terminal region, the 1958 surface slope should be a reasonable approximation to the surface slope in 1972. To check this approximation, we compared the surface slope between markers spaced 300- 500 m apart, which were surveyed in 1972. with the surface slope measured on the 1958 map at about the same locations. Except for the lowermost 1.5 km of the glacier, differences are random and generally smaller than ±0.5°. We therefore expect an accuracy of about ±0.5° in the 1972 slopes over much of the glacier. For the 1993 slopes, the error is ±0.3°.

The center- line surface slopes for the two time periods are shown in Figure (3). (In this figure and throughout the paper, .s is a curvilinear coordinate that follows the center line of the glacier, extending from the upper or lower cirque to the terminus, as indicated in Figure (1).) Like many other valley glaciers in the region. McCall Glacier shows an interesting Periodic succession of maxima and minima of surface slope. These bulges are accompanied by an almost linear increase in average slope towards the snout. There is a similar, but proportionally smaller, variation in ice thickness along the glacier. The wavelength of these oscillations is about 1500 m, or approximately 11 ice thicknesses. While some of the sleeper and thinner reaches of the glacier may he caused by resistant bedrock ridges along the valley walls (Fig. 1), there is not a one- to- one correspondence between die waxes of the glacier and these ridges.

Fig. 3. Surface slopes in 1972 and 1993 along the center line and transversely averaged slopes in 1993.

According to Reference Hooke,Hooke, 1991. waves in glacier beds may initially stem from existing perturbations in the bed profile. The series of ridges and overdeepenings then lends to be enhanced by a positive feed- back of erosion with variations in siibglacial water pressure. Abrupt variations in water pressure at the headwalls can cause erosion by quarrying, while on the down- glacier sides of overdeepenings the reversed bed slope causes a constant, high water pressure that minimizes erosion there. This process then causes amplification and broadening of the ovcrdeepenings.

An alternative theory that explains regular oscillations in valley- glacier geometry is dial proposed by Reference Mazo,Mazo (1989). He studied the coupled system between the flow of a glacier and the erosion of its bed. Bed geometry determines glacier geometry, which in turn controls the erosion of the bed through the distribution of basal shear stress. For a glacier with an initially uniform thickness h and a constant bed slope θ. a small sinusoidal perturbation in θ produces a corresponding oscillation in the oscillation in ice thickness then translates into a variation in basal erosion at the same wavelength but generally out of phase with respect to the initial disturbance. As a result, the disturbance in the bed migrates up- glacier with a speed that depends on its wavelength and the amplitude of the disturbance generally remains small. However, there is one wavelength. λ= 2h / tanθ that is stationary and whose amplitude grows by “constructive erosion For McClall Glacier, the mean slope is

and the mean thickness (h) is 140 m. This leads to λ ≈ 15 h or about 2100 m. Considering the approximate nature of Reference Mazo,Mazos’ (1989) theory, this is reasonably close to the observed wavelength of the topographic wave.

Transversely averaged surface slopes in 1993 are also shown in Figure (3). Transverse variations in surface slope lead to transverse shear stress gradients. Consequently, a transverse average of surface slope should enter into a one- dimensional Longitudinal flow analysis, such as that presented later. These slopes were obtained by averaging the longitudinal surface slope across a 500 m wide swath about the center line. We assumed that the ratio of the trans- Versely averaged slope to be center- line slope has not changed since 1958. The transversely averaged slope is not significantly different from the center- line value except at two local maximalocated at center- line positions of 2500 and 5200 m. where the two measures of slope differ by about 0.5- 1°.

As seen in Figure (3), there was a general steepening of Up to 2° along the lower glacier from 1972 to 1993. This is due to the observed thinning (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995). Such steepening should, to some extent, counterbalance the effects of the observed decrease in thickness on the ice velocity in this region. The 1993 terminus (s = 7000 m; Fig. 1) was considerably less sleep than the terminus show n on the 1958 map (at s = 7370 m). The terminus slope in 1972 is clearly not well represented on the 1958 map, because of the 8 m thinning and 70 m horizontal retreat from 1958 to 1972. On the upper glacier (above 1700 m center- line distance), the large differences between the 1972 and 1993 slopes are artifacts of known errors in the 1958 topographic map (Reference Rabus, Echelmeyer,, Trabant and Benson,Rabus and others, 1995).

Thrust faults

Conspicuous features in the mid- to lower- ablation area of McCall Glacier are discontinuities in the ice that extend over much of the glacier’s width (Reference Wakahama, and Tusima.Wakahama and Tusima, 1981: personal communication from L. Shapiro). These discontinuities form down- glacier ares and are marked by dirtladen ice bands. They dip up- glacier at shallow angles of about 30 . Topographic steps at these bauds sometimes occur, with the up- glacier sides being up to 0.3 m higher than the down- glacier sides. Also observed are dirt- laden lateral discontinuities, about 10 m from the ice edge, that parallel the eastern margin of the glacier over several kilometers. Both types of discontinuity resemble the surface expressions of the spoon- shaped thrust faults observed by Reference McCall,McCall (1952) on cirque glaciers. Several alternative explanations have been proposed for these features: (1) movement across a thrust fault in the ice, perhaps facilitated in differing rheological properties of the dirt- laden ice (e.g. Reference Nye,Nye, 1951); (2) differential ablation caused by albedo differences of dirty and clean ice. and the accumulation of heat- absorbing dust on die down- glacier side of the discontinuity; and (3) former crevasses that have been closed and are deformed by viscous flow. The latter explanation can be ruled out on McCall Glacier, as there are no large crevasses that Could be their predecessors.

About twenty to thirty of these features are found in a region of the glacier (s = 4000- 6000 m) where modeling suggests a high, non- plastic contribution to the How. This contribution is probably due to year- round basal sliding but alternatively might be caused by displacement across such faults.To investigate the latter idea, we made measurements on three such “faults” in This region, in order to determine whether significant differential motion exists across them. Pairs of wooden poles were drilled vertically into the ice, about 1 m apart, above and below the discontinuity. Kach discontinuity had a well- developed vertical offset. Relative motion between the pairs of poles was carefully determined over a 1 year period. At two of the sites, there was no differential motion and at the third there was a convergence of 0.03 ±0.01 m. Total surface displacement at each of these sites during the same time period was on the order of 10- 15 ma⁻¹, The average longitudinal strain rate in this region was about 0.005 a⁻¹ , corresponding to a not convergence of 10 m a₋₁ over a 2 km stretch. These measurements therefore show that, if there were 30 such faults, a maximum of 30 x 0.03 m a⁻¹ = 0.9 m a⁻¹ . or less than 10% of the net convergence in this region, is due to strain on the “thrusi- faults”. This lack of Significant strain agrees with measurements made on similar features found on Barnes Ice Clap (Baker, 1986). Either such features are entirely passive or they form by shear displacement over a relatively short time interval and then become inactive. In support of the latter mechanism, a new fault was observed to form on Hubley Glacier (3 km southeast of McCall Glacier) over a 4 week period in an area that was previously undisturbed (Fig. 4a and b) The fault showed an offset of about 0.1 m and it extended laterally over several hundred meters in clean ice. Ii had an up- glacier dip of about 20- 30°. Along this particular fault, there appeared to be little no dirty ice concentrated on the actual fault surface. No obvious dirt band had been observed before the offset formed and thus differential ablation could not have caused it.

Velocity measurements

Ice velocity was first measured on McCall Glacier by Reference Sater,Sater (1958). An unfortunate choice of survey geometry caused substantial errors in his data, and the 1957 reference frame is not given, so that these results could not be compared with later measurements. From 1969 to 1972. an extensive mass- balance network was maintained on McCall Glacier (Reference Wendler,, Fahl, and Corbin,Wendler and others, 1972: Reference Trabant and Benson,Trabant and Benson, 1986). The marker poles in this network were frequently surveyed by theodolite triangulation (Trabant, personal communication). We used these raw survey data to calculate mean annual surface velocities between 1970 and 1972. In some cases, we were also able to calculate seasonal velocities for the 2 month period from mid- August 1971. Our analysis suggests an uncertainly of about 0.3 m a⁻¹ in the 1970s annual velocities.

In 1993, a number of poles were emplaced in the ice at a sub- set of the 1970 locations. These were carefully surveyed 2- 5 times a year during early june to mid- August until 1995, and both annual and seasonal velocities were derived. The surveys were carried out using a 1 s theodolite and electronic distance meter. Based on an error of 5 mm in distance and 2 arcseconds in angle we found typical errors of 0.02 ma⁻¹ for annual velocity and 0.45 m a ⁻¹ for velocities determined over a 10 day period. Corresponding azimuthal accuracies are 1° and 10° for the animal and short- term velocity vectors, respectively. We used the same network control as in the 1970s; this was resurvcyed with Cilobal Positioning System CiPS methods in 1993. This allowed justification of the 1970s and 1990s surveys in an absolute coordinate frame, consisting of UTM northing and easting and height above- the WGS84 ellipsoid.

Fig. 4. (a) Arcuate “thrust fault” on South Hubley Glacier, newly formed in July 1995. The snout of North Hubley Glacier is in the background. (b) close- up showing offset anil up- glacie) and dip of the “fault”.The instrument indicating the dip is about 40 cm long.

Horizontal velocities: general characteristics

Figure (5) shows velocity vectors, averaged over the entire period of measurement 1993 95 (heavy arrows) and over a period of high melt, from 3 to 24 July 1993 (light arrows). In Figure (6), mean annual velocities from 1993 lo 1995 along the center line are shown together with those of 1970- 72. Annual velocity along the center line reaches a maximum of 16- 18 ma⁻¹ at S = 4000 m, just down- glacier of a pronounced narrowing and shallowing of the channel (Fig. 1). A second lfical maximum of about 14ma⁻¹ in occurs Over another bedrock threshold at s = 2300 m, while a shallow velocity minimum at s= 3000 m corresponds to thicker ice. The local velocity maximum at s = 2300 m seems to be associated with the “steady- state” equilibrium line approximated by the boundary of the ablation zone in 1970 (Stippled area in Figure (5)), when McCall Glacier had a mass balance close to zero (Reference Trabant and Benson,Trabant and Benson, 1986). Center- line velocities in the 1970s show the same pattern as in the 1990s but there was an additional local at maximum around s = 5200 m. In a later section, we calculate the change in velocity from the 1970s to the 1990s using changes in ice thickness and surface slope. Our calculations predict that the maximum at S = 5200 m would disappear in the 1990s and that there would be no up- glacier migration of the other two velocity maxima (open symbols in the 1990s velocity curve in Figure (6)).

Fig. 5. Velocity vectors: menu annual value 1993- 95 and summer value during the peak season. 3- 24 July 1993, Elevation contours are in m.a.s.l.: they refer to the 1958 topography. The stippled area is the approximate outline of the ablation zone in 1970. when the net balance was close to zero (Reference Trabant and Benson,Trabant and Benson, 1986).

Transverse profiles of mean annual velocity are presented in Figure (7). Solid and dotted lines represent the annual averages of 1993- 95 and 1970- 72. respectively. Both the actual velocity vectors, as well as their decomposition in components parallel and perpendicular to the central flow- line, are shown. Velocities were not measured at the glacier margins but wen- arbitrarily set to zero there. All velocity profiles, except T6, are parabolic in shape. This suggests that flow occurs primarily by ice deformation. The T6 velocity profile is somewhat more plug- shaped. hinting at significant sliding. Transverse velocities are small at T3,T6 and T11 At T20. where the combined flow of upper and middle cirques spreads into the confluence area, the transverse velocities show divergence. In contrast, there is pronounced convergence at T10. where the tee of the confluence is focused into the main trunk of the glacier.

Fig. 6. Mean annual velocity along the center line during 1970- 72 and 1993- 95. The heavy open symbols represent values that were calculated from Equation(11).

In Figure (8). we compare mean annual and seasonal velocities from 1993 to 1995. A detailed discussion of seasonal velocity fluctuations will follow in a later section; here, we make some general remarks. Annual velocities and those measured during winter and spring are all similar (solid symbols in Figure (8).) while velocities observed during the melt season open symbols are considerably higher. It is interesting that these elevated seasonal velocities lend to peak near s= 2300 m rather than at S = 4000 m. where the annual velocity has its maximum. Short- term velocities during periods of intense melt (large, open symbols) also show this tendency of an overall up- glacier shift of the region of maximum velocities. The highest measured velocity was 28 ma¹ at s = 2880 m, averaged over a period of 7 days during the summer of 1995. The summer period, during which velocities arc elevated, is too short to cause the annual average to be significantly higher than the winter average.

Longitudinal strain rate and climatic regime

The velocity field of McCall Glacier is characterized by small Longitudinal strain rates

these ate about an order of magnitude less than those typically found on temperate glaciers. Here, we show that these small strain rates reflect the Arctic climate regime of the glacier, as characterized by a small mass- balance gradient with elevation. .

The mass- balance rate at a longitudinal position s, is given by

where β is the tangent of the bed slope and s_EL is die position of the equilibrium line relative to the glacier head (s= 0). Here, we have neglected the small contribution to elevation difference that arises from changes in ice thickness along the glacier, In this approximation, t In steady- slate Length is L = 2_sEl, and the flux at the equilibrium line is

where W is the (constant width of the glacier. The mean longitudinal strain rate (absolute value) along the glacier is given by

being u_El surface velocity at the equilibrium line. I Using the relation

between surface and average velocity I Reference Paterson,Paterson, 1994), a parabolic cross- section of depth h_EL with a

and the above expression for q_El the mean longitudinal strain rate becomes

, This shows that the average longitudinal Strain rate is proportional to the mass- balance gradient (A similar result for the local value of longitudinal strain rate was obtained by Reference Nye,Nye (1957)). Continental glaciers, therefore, show smaller longitudinal strain rates than maritime glaciers with a similar geometry.

Fig. 7. Transverse profiles of mean annual velocity (ma⁻¹): shown for each profile are: longitudinal (heavy solid line with solid symbols (1993- 95), heavy dotted line with open symbols (1970- 72) ) and transverse component (heavy line (solid and dotted, respectively)) without symbols in the lower part of each plot, and velocity vectors for 1970- 72 (dotted) and 1993- 95 (solid). The vectors are sums of the two components.

A complication in Equation (1) is that the center- line depth _EL itself depends on L. β and

, and the flow- law parameter A We may take this into account by noting that, for slopes sleeper than about 5° β is approximately equal to the surface slope at the equilibrium line, α_El. Then, assuming only deformational flow, and, using Equation (1) we find

Lower ice temperatures of Arctic glaciers lead to smaller values of A and this will add to the effects of smaller balance gradients in decreasing the mean strain rate.

As a comparative example, we calculale longitudinal Strain rates using Equation (1) with n = 3 for two glaciers: McCall Glacier, with

, β=7.9° h_EL = 220 m and L = 7.8 km: and sub- Arctic Storglaciären, a polythermal glacier in northern Sweden which has a significantly larger mass- balance gradient a^-1 km^-1; evaluated from Reference Haeberli., Herren, and Hoelzle,Haeberli and others (1993)) and somewhat different geometry: β = 6.3°, h_EL, = 165 m.L = 3.8 km (evaluated from Reference Hooke., Calla,, Holmlund,, Nilsson, and Stroeven,Hook and other (1989)). For McCall Glacier, while for Storglaciären . These calculated values agree reasonably well with those observed (0.0036 and 0.0111 a¹, respectively, as determined from 2u _EL/L, and u_EL = 14 ma⁻¹ on McCall Glacior and u _El=21 m a⁻¹ on Storglaciären (Reference Hooke., Calla,, Holmlund,, Nilsson, and Stroeven,Hooke and others, 1989) The Larger strain rates on Storglaciären are a direct effect of a larger mass- balance gradient (i.e. a warmer Climate) there. The smaller bed slope and length of Storglaciärcn will decrease the mean strain rate some what relative to that on McCall Cilaeier but the climatic effect dominates.

Fig. 8. Intra- annual variations of center- line velocity. 1993- 95.

Emergence velocities

Emergence velocity, defined as the vertical rise or subsidence of the ice at a given horizontal position, can be calculated as

where u _H and v are horizontal and vertical components of the annual velocity of a marker pole and α₁ is the local surface slope in the flow direction. Measurements of emergence velocity can be used to evaluate effects of vertical advection on the ice- temperature distribution of a glacier or to determine the longiludinal ice flux q(s) by using the continutiy equation. Emergence velocities, calculated from Equation (3) for a number of poles along the center line, are shown as open circles in Figure (9). The optically surveyed velocity components U_H and v had errors of 0.02 and 0.05 ma⁻¹, respectively. Local Surface slopes over distances of 25- 40 m Were obtained from a kinematic GPS profile along the glacier that passed through the poles. The uncertainty in this slope dominates the overall error in the emergence velocity; maximum errors are 0.15 m a⁻¹.

Average emergence over the width of the glacier can also be obtained indirectly from the difference between annual elevation change and the mass- balance rate, each averaged over w:

The broken line in Figure (9) shows emergence velocities calculated from Equation (4). using mean elevation change and mass balance 1993- 95 at about 30 locations on the glacier surface. Corresponding errors are about 0.2 m a⁻¹ for both elevation change and mass balance, which leads to an overall error of 03 m a⁻¹ for this method.

Both methods suggest that the has a maximum of about 0.6 m a¹ at around s = 4000 m. Emergence at the center line is consistently higher than the width- averaged value. This is to be expected because the equilibrium line of McCall Glacier is not perpendicular to the center line (Fig. 5). This leads to (negative) accumulation- area contributions in the width average of emergence velocity up- glacier of S = 3000 m, while the center line itself is in the ablation area until about S= 500 m. In particular, the negative emergence velocity between TU and T20 calculated from Equation (4) is caused by including the entire accumulation bowl of the middle cirque in the width average. Furthermore, an up- glacier shift of the equilibrium line leads to emergence velocities that are temporarily higher at the center line than at the margins because the glacier needs to build a crown in the new ablation area (Reference Raymond,Raymond, 1971). For McCall Glacier, this effect might partly account for center- line emergence being higher than width- averaged emergence between cross- sections T9 and T6.

Fig. 9. Average emergence velocities 1993- 95: (i) on the center line, calculated from surveys of individual marker poles and local surface slopes in their vicinity, and (ii) width- averaged, calculated from width- averaged local mass balance and elevation change.

Analysis of surface velocity

It is generally assumed that there is little or no basal sliding beneath cold- based glaciers (e.g. Reference Echelmeyer, and Zhongxiang,Echelmeyer and Zhong- xiang, 1987). However, for a polythermal glacier whose base is at least locally at the melting point, basal sliding could be significant. In order to investigate this, we analyze the dcfor- mational speed, u _d, of McCall Glacier in terms of ice geometry.

Following Reference Nye,Nye (1965). the center- line surface speed of a glacier flowing without sliding in a channel of parabolic cross- section is given by

where p is the density of ice, g is the acceleration of gravity. n=3 the flow- law exponent, A is the flow- ice parameter. h is the local center- line ice thickness, α is the local surface slope and f is a shape factor that depends on the width to depth ratio of the channel. If ice temperature varies with depth. T = T(z), an effective depth- weighted flow parameter

that produces the correct surface value of u_d in Equation (5) can be used. Because we do not know how the temperature- depth profile varies along the channel, we use T(z) as determined at cross- section T6 by Reference Trabant,, Harrison,, Benson,, Weller, and Bowling,Trabant and others (1975). Using this temperature profile and A(T) as given by Reference Paterson,Paterson (1994, table 5.2), we obtain A_eff = 0.179 a⁻¹ bar⁻³, which corresponds to an effective ice temperature of −0.5° C In reality. A may vary along the glacier, especially as the thickness of the temperate layer at the bed varies, but ii should be between 0.145 (for an ice column at −1° C and 0.194 (for a 50 m temperate layer beneath ice at −5° C).

Surface speed along the glacier was calculated from Equation (5) using local ice thickness and surface slope as evaluated along the center line at 150 m intervals (Figs 2 and 3, respectively). The parameters are shown together with the channel- shape factor in Figure (10)a. The shape factor was interpolated as a function of s by assuming a parabola at each cross- section with half- width and depth taken from Figure (2). The fit between Equation (5) and the observed speed is poor (Fig. 10b). It is unlikely that the quality of the fit can be improved by reasonable variations in A and, as the deviation between the calculated and observed speed is highly variable, it is equally unlikely that the fit can be improved by including fluctuations in basal sliding at this stage.

Fig. 10. Modeled velocity of McCall Glacier compared to measurements, (solid symbols): (a) local input parameters; (b) velocity from Equation (5) using local parameters: (c) longitudinal averaging following Equation (7) with best fit coupling length l = 3.0h (solid line) and l = 1.8h (dashed line) and no sliding; (d) longitudinal averaging (l =3.0h) with a localized sliding anomaly included (u _B). Arrows in (c) and (d) indicate the datum location where measured and modeled velocities are equal.

The large variations in surface slope and ice thickness give rise to longitudinal stress gradients, which require modification of Equation (5). Following the analysis of Reference Kamb, and Echelmeyer,Kamb and Echelmeyer (1986), the deformational velocity of a glacier with small perturbations in slope and thickness is given by an exponentially weighted average of these perturbations along the glacier length L:

Here u _d ⁽⁰⁾ is the velocity defined for some average parallel- sided datum state and Δ specifies deviations from this datum state (e.g. Δln h ⁿ = ln hⁿ − ln hⁿ ₀ ) for ice thickness h(s) and the datum thickness h₀ ). The range of the stress coupling is governed by the characteristic length. l. of the exponential averaging kernel. For temperate glaciers with typical longitudinal strain rates of 0.01- 0.05 a¹, l is about 1.5- 2.5 times the ice thickness. The strain rates on McCall Glacier arc smaller (0.003- 0.005 a¹), being an indirect effect of the small mass- balance gradient of the glacier, as discussed above. The low strain rates lead to a longer averaging length for McCall Glacier of 3.0- 3.5 times the mean ice thickness.

In applying Equation (7) to McCall Glacier we recast it in a computationally more efficient form (Reference Kamb, and Echelmeyer,Kamb and Echelmeyer ,1986). equation (35). in which the datum state is defined by the average thickness, surface slope and shape factor at one representative location where the surface velocity u₀ is known. The velocity at any Other location along the glacier may then be determined by integration from this datum state. This procedure was applied m the local parameters given in Figure 10a for l= 3.0 ft results for a typical “temperate” coupling length, l = 1.8 h, are shown as a dashed line for comparison;. The fit to the observed speed (Fig. 10c) is much improved over that resulting from Equation (5), However, in the middle reaches of the glacier, between s = 4000 and 5000 m. the fit is not good, as the modeled speed is less than two- thirds the observed speed. Unlike other parts of the glacier, where decreasing depth is largely counterbalanced by increasing surface slope, here depth and slope decrease almost simultaneously from s = 3800 to 4800 m (Fig. 10a). This leads to a low predicted deformational speed in this region.

Use of the transverse- average surface slope, as defined before in the section surface slope, gives no significant improvement. The same is true if longitudinal variations of the flow- law parameter A are introduced. An up- glacier increase of A, which accounts for possible warming in the accumulation area caused by refreezing of meltwater. makes the fit worse.Using ,4 for T=0° C between . s = 1000 and 5000 m, and A at T = −2° C elsewhere (admitting that the observed layer of temperate basal ice could be a local phenomenon) . yields modeled velocities which are still about 25% to low. Finally, any observed motion on the thrust faults in this region is insufficient to explain the anomalous speed there.

The most likely explanation for the high velocities is a significant contribution of sliding at the temperate bed underlying the mid- glacier region. Iniroduction of such sliding a priori into the longitudinal averaging of Equation (7), following the procedures outlined by Reference Kamb, and Echelmeyer,Kamb and Echelmeyer (1986). produces excellent agreement between modeled and measured velocities along the entire glacier (Fig. 10d) .While the details of the sliding distribution are not uniquely defined, the main features of u _B are well constrained, tin must be concentrated at the center of the velocity anomaly. It must also account for the major part of ice motion there, not just the difference between the observed velocity and that calculated in the case of no sliding (Fig. 10c and d). Where there is a large amount of sliding, the surface profile of the glacier is different and the deformational velocity is decreased. This is discussed in the next section. As was mentioned before, the transverse- velocity profile at cross- section T6 (Fig. 7) shows more of a “plug” shape, which independently suggests significant sliding in this region. Thus, we infer that significant sliding occurs over this one region of the glacier, but little sliding need occur elsewhere.

Surface expression of a sliding anomaly

We next consider how localized basal sliding is reflected in the surface profile of a glacier. Simple arguments suggest that the glacier should be thinner when sliding speeds are high and this will cause a change in surface slope. To investigate this quantitatively, we determine h(s) and α(e) from the measured ice flux q(s), prescribed bed geometry z_b(s) and sliding velocity u_B(s). We then apply this relation to an idealized glacier with constant bed slope and flux gradient to illustrate the general topographic features associated with a sliding anomaly. This is followed by an application to the actual geometry and flux of McCall Glacier. The approach is complementary to the one in the previous section, in which measured h and α were given as inputs into Equations (7) to determine U_s .

The flux q through cross- section S(s) is not calculated from velocity measurements, which we only have at the ice surface. Rather, it is obtained from measured annual elevation change ∂h/ ∂t and mass balance b at a given time via mass continuity;

The righlhand side of Equation (8) equals the emergence velocity

averaged over the width of the glacier, w (Equation (4): Fig. 9). is insensitive to the particular iimc period dial is used to determine ∂h/ ∂t and because inter- annual variations of these data sets are synchronous and tend in cancel out. We use an average of 1994 and 1995 values. The flux, integrated from Equation (8), can be separated into contributions from cross- sectionally averaged tee deformation and sliding according to

For a parabolic cross- section we have

and for I he simple defbrmational flow given in Equation (5) we have with a second shape factor f ₁ from Reference Nye,Nye (1965). While f and f₁ depend weakly on die dimensions of the cross- section and the actual distribution of sliding, they are here approximated by their mean values for the cross- sections in Figure (2). For small surface slopes, , and from Equations (5), (8), and (9) we obtain

Equation 10 can be used to solve for the surface geometry. Z_s = z₁ + h, given ∂h/ ∂t and

This surface geometry will depend on the prescribed distribution of sliding velocity. In a different approach, Reference Balise, and Raymond,Balise and Raymond (1985) calculated how the sudden temporal onset of a basal- sliding anomaly banges the patterns of horizontal and emergence velocity for a linearly viscous parallel slab. Contrary to our method, they used glacier geometry as input to calculate emergence velocity. Their emergence velocity was no longer consistent with the given glacier geometry after initiation of sliding.

The ordinary differential Equation (10) was numerically integrated up- glacier from s = L, using a constant effective- flow parameter A_eff(Equation (6)) and n = 3. The glacier snout was treated as a rigid wedge with constant velocitv to avoid die singularity at s = L. where h goes to zero. Up- glacier integration is stable and the solutions h fulfill Equation (10) within the expected numerical error. Without sliding, the direction of Integration can be reversed at any s and the corresponding “down- glacier” solution agrees with the “up- glacier” solution. In theory this “reversibility” should hold in general as Equation (10) is order and has a unique solution for initial conditions h(S₀) =h₀ In practice, how ever, even for small, localized sliding, down- glacier integration was numerically unstable. It diverged from the up- glacier solution to produce unrealistically large down- glacier ice thicknesses that failed to decrease when approaching the desired location of the snout at s = L. This numerical instability is produced by the high sensitivity of the solutions up- glacier of any sliding anomaly as seen in Figure (11) and discussed below), and by the shape of the ter- minus region. The thickness profiles that we calculate from Equation (10) are uniquely determined by fixing the terminus at s = L; they are not an artifact of the up- glacier integration.

Fig. 11. Modeled profiles of ice thickness and surface slope on a linear bed. (a) Without sliding (0) and with narrow (1) and broad (2) Gaussian- shaped sliding anomalies; (b) with a rectangular- shaped sliding anomaly. For each case, the surface velocity (u _s) is also shown.

As a first example of the influence of localized sliding on ice thickness and surface topography, we consider a constant bed slope of 6.3°, equal to the mean bed slope of the modeled section of McCall Glacier, and a constant emergence velocity

. Changes in ice thickness and surface slope relative to the no- sliding case are denoted by Δh and Δα, respectively; as the bed z _b is fixed (Δz_b = 0), we have . The effect of Gaussian- shaped sliding anomalies of different widths is illustrated in Figure (11)a. There is a decrease in ice thickness within a limited region around the sliding anomaly relative to the no- sliding case (labeled as curve “0”). The maximum of |∆h| is located up- glacier of the maximum sliding speed and |∆h| is non- zero well up- glacier of the sliding anomaly but not down- glacier of it. Surface slope decreases down- glacier from the center of the anomalies but increases up- glacier. This up- glacier increase is significantly more pronounced for the wider anomaly. The rectangular sliding anomaly shown in Figure (11)b consists of two opposite step changes in u _B. These steps asymptotically cause constant up- glacier offsets in ice thickness (of opposite sign). In terms of their effect on deformational velocity, the observed changes in h and α restore the (prescribed) total mass flux within and around the region that transports mass by sliding. Smaller ice thicknesses up- glacier and outside the sliding anomaly require a higher overall ice velocity; this is achieved by an increased surface slope there. Within the sliding anomaly, deformational velocity must decrease and this is accompanied by decreasing ice thickness and surface slope.

Next we solved Equation (10) for the measured bed profile and emergence velocity of McCall Glacier (

in Figure (9)). The results for no- sliding anomaly are shown in Figure (12)a, together with measured ice thickness and surface slope (dashed lines). These solutions for h and α differ markedly from the measured quantities in the anomalous region s = 4000 m to s = 5000 m: h is too large there and α is too large in the lower part and too small in the upper part of the anomalous region. This resembles qualitatively the typical “fingerprint” of a sliding anomaly found above (Fig. 11). The surface expression of a trapezoidal sliding anomaly similar to the one proposed in the previous section is shown in Figure (12)b. There is excellent quantitative agreement between calculated and measured ice thicknesses and surface slopes.

Fig. 12. Modeled profiles (solid curves) of ice thickness and surface slope of McCall Glacier compared to measured values (dashed curves) for: (a) no basal sliding and (b) a trapezoidal- shaped sliding anomaly.

Longitudinal stresses and sliding were both equally important when flow was modeled directly from measured h and α in Equation (7) However, longitudinal stresses have a far less profound effect on the solutions for h and α from Equation (10) than does sliding. The reason for this “asymmetry” between Equations (7) and (10) is that flow depends sensitively on measured ice thickness and surface slope, q ∝ hⁿ⁺²αⁿ , while the change in ice thickness depends only weakly on flux,

. Solutions of Equation (10) are therefore relatively insensitive to small variations in q.

An important result of this section and the previous one is that there exists localized sliding at the base of this polythermal glacier over a 2 km long region and that this sliding accounts for up to 70% of the annual ice flux there. The measured velocity increase during the short summer season (Fig. 8) can account for at most 5–7% of the annual ice flux and thus the glacier must slide year- round at an almost uniform rate. Why does the glacier slide in this region but not in others? One reason could be that just up- glacier of the anomalous region the largest supraglacial stream on McCall Glacier disappears beneath the surface. From the analysis of Reference ShreveShreve (1972), the englacial channel should dip at about 11 times the surface slope, or around 70°. Surface water could therefore reach the bed near the beginning of the anomalous region and facilitate sliding there. However, there is no water flow from September to May. Perhaps there is significant water storage and correspondingly high water pressures throughout the winter.

Long- term changes of annual velocity between the 1970s and 1990s

Annual velocities in the lower part of McCall Glacier have changed markedly since the 1970s (Figure (6) and diamond symbols in Figure (13)). To determine whether these changes can be understood quantitatively from corresponding changes in the geometry of the glacier, we calculate the expected variations in deformational ice velocity (Equation (5)) that arise from small changes in surface slope, Δα, and ice thickness, Δh. To first order, this variation is given as

The local values of Δu/u, calculated from Equation (11) with n = 3 (heavy solid line) are compared with the observed velocity changes in Figure (13). Measured Δh/h and Δα/α, represented by solid and dotted lines, respectively, were both averaged over 300 m or about twice the ice thickness. The agreement over the entire lower ablation area is quite good. Between T9 and TlO observed Δu/u is close to zero while the calculated value is negative because the surface slope appears to have decreased. This apparent decrease is caused by the extrapolated 1970s slope being erroneously high due to known inaccuracies of the 1958 topographic map in the region.

Fig. 13. Comparison of measured (symbols) and modeled (curve) relative velocity change Δu/u since the 1970s.

While effects of localized sliding and longitudinal stress coupling are important to an understanding of the flow of McCall Glacier, changes in flow are apparently well represented by Equation (11), which neglects these effects. A tentative explanation is that changes in the spatial distributions of basal sliding and longitudinal stress since the 1970s were negligible. However, there is the alternative possibility that these changes were not negligible but merely obey equations similar to Equation (11) (Reference Echelmeyer,Echelmeyer, 1983).

Velocity profiles across T10, T11 and T20 (Fig.7) all show subtle changes from the 1970s to the 1990s. As the horizontal positions of poles in 1993 were chosen to be identical to those in 1972, these changes are likely to be real. The velocity at T11 actually increased by 5% despite the observed thinning of c. 5 m since 1972 at that location. Surface slopes from the 1970s are not available for this region, so we cannot determine whether it is a change in surface slope or is the seasonal sliding rate that offsets the effects of thickness change across T11.

Spatial pattern of seasonal velocity fluctuations on McCall Glacier

The higher velocities which we observed during the melt season (Fig. 8) are likely to be due to temporal increases in the rate of sliding as a result of meltwater production from early June to early August. To study the spatial pattern of enhanced ice motion, we calculate the differences Δu _s between seasonal velocities and average winter velocity for the years 1993–95 (Fig. 14a). Seasonal velocity is evaluated from the first to the last survey in a given year. The 1995 “summer season” (20 June- 11 August) spans almost the entire melt season, while the 1993 “summer season” (3–24 July) represents a warm period with higher than average ablation and the 1994 (11 June- 5 July) data span a cold period with very little melt.

Fig. 14. Seasonal velocity increases Δu _s = (u_s - u_w) for the 1993, 1994 and 1995 “summer” seasons: (b) per- cent increase in velocity relative to the winter velocity, Δu_s/u_w . These results are shown for the center- line markers.

In 1995, Δu _s had a maximum near cross- section T10, well up- glacier of the maximum in the winter velocity, u _w. The ratio Δu _s/u _w (Fig. 14b) shows a minimum close to the position of the maximum in u _w and it reaches maxima in the lower ablation area and above the equilibrium line. While the percentage increase in velocity is relatively large, the short duration of the “summer season” leads to only a small additional summer displacement. Thus, only 3–5% of the annual ice flux is due to the extra summer motion.

During the warm measurement period in 1993 the spatial patterns of Δu _s and Δu _s/u _w on the lower glacier were similar to those in 1995, but above T9 the short- term velocity was greater than in 1995. During the 1994 cold spell, the patterns of Δu _s and Δu _s/u _w were different from those in 1993 and 1995. A limited region between T3 and T6 showed a velocity increase, while the lowermost glacier and the confluence region had velocities that are lower than the winter average.

The pattern of Δu _s observed during the cold spell in June 1991 may be typical for the early stages of the ablation season. Speed- up starts between T3 and T6, where there is a temperate basal layer of finite thickness. As the melt season progresses, the region of strongest velocity increase is shifted up- glacier as indicated in the 1993 and 1995 patterns of Δu _s. We can only speculate on how this shift is accomplished. A possible mechanism could be an up- glacier migration of the zone of highest meltwater input during the progressing melt season. Up- glacier of s = 3500 m there are no obvious moulins but there are a few larger crevasses above s = 1500 m that might facilitate meltwater passage to the bed. Alternatively, Reference Andreasen,Andreason (1985) has proposed that a velocity increase in the upper reaches of polythermal glaciers might be caused by the backing- up of meltwater along the glacier sole due to weak basal drainage on the lower glacier. However, we feel that this mechanism is physically untenable for McCall Glacier.

Temporal variations in velocity and ablation for McCall Glacier and White Glacier

Short- term velocity measurements during the 1995 summer allow us to investigate the relation between temporal changes in ablation and velocity on McCall Glacier. Figure (15) shows velocity of the center- line poles at T11, T9, T6 and of a pole at the confluence of upper and middle cirque (UC/MC), along with representative ablation measured near T11. For all poles, the velocity shown prior to 20 June is an average from 5 July 1994, so it represents approximately the mean annual velocity. The ablation record shows that meltwater generation ceased glacier- wide during a cold spell from 25 June to 1 July.

Velocity at T11 and UC/MC, averaged over the 8 day period containing the cold spell, is greater than mean annual velocity. The velocity maintains this elevated value over the next 7 day period when a maximum in the ablation rate occurred and also during the following 36 day period. For the upper ablation area, this suggests that there is no short- term causal relation between ablation and velocity.

Velocity at cross- section T9, averaged over the period containing the cold spell, is slightly smaller than mean annual velocity, while velocity during the subsequent period of enhanced ablation shows an increase. There appears to be no time lag between ablation and velocity in this mid- ablation area, at least at the resolution of the measurements. Any reasonable interpolation of the velocity at this site indicates that the sensitivity of the velocity to ablation appears to decrease in the second half of the ablation season.

In order to bring out effects that are typical for polythermal glaciers, we compare our results with the record of velocity and ablation from White Glacier, Axel Heiberg Island, Canada (Reference Iken.Iken, 1974). White and McCall Glaciers share the same Arctic climate that limits significant ablation to early June through early August, with a prominent period of high melt in July. Ice- surface temperatures and ice depth are similar but White Glacier is longer (14.5 km), wider (1 km) and has a somewhat gentler slope (6°) than McCall Glacier. The measured annual center- line velocities on White Glacier are 30–35 m a⁻¹, and there does not appear to be any year- round sliding.

Fig. 15. Seasonal- velocity variations and ablation, summer 1995.

The high- resolution record of Reference Iken.Iken (1974) shows several interesting features in the temporal nature of the velocity; these are more clearly delineated than those from our coarse data set. Common features of the McCall and White Glacier records are (i) a substantial summer velocity increase on both glaciers that contributes only a small fraction to the annual ice motion, (ii) on the lower glacier, the sensitivity of the velocity to high ablation diminishes with the progress of the melt season, (iii) minimum velocities are reached at the onset of the melt season and velocity in autumn is somewhat higher than this minimum and (iv), especially on White Glacier, velocity may lag ablation by several days at some locations. However, on White Glacier the relative magnitude of seasonal sliding decreases up- glacier to the equilibrium line, while on McCall Glacier the highest contributions were found near the equilibrium line.

Speculations on the velocity fluctuations of polythermal glaciers

Temperate glaciers reach their maximum velocity in early summer, while minimum velocities coincide with the end of the melt season. This is probably a consequence of the seasonal evolution of the englacial drainage network (Reference Röthlisberger,Röthlisberger, 1972). In autumn, R- channel diameter is at its maximum and meltwater input is small; this leads to low water pressure and minimum sliding velocities. Over the course of the winter, R- channel closure causes a gradual increase in water pressure and consequent sliding velocity.

Polythermal glaciers (this paper; Reference Iken.Iken, 1974; Reference Andreasen,Andreason, 1985; Reference Hooke., Calla,, Holmlund,, Nilsson, and Stroeven,Hooke and others, 1989; Reference Willis.Willis, 1995) seem to behave somewhat differently: (i) highest velocities occur during highest melt, (ii) velocities in the late melt season are higher than those at the beginning of the melt season when they actually seem to be at a minimum and (iii) at some locations in the Upper ablation area there is evidence for a time lag of several days between ablation and velocity events. All three observations may be related to prolonged water storage in polythermal glaciers. This could also indicate that a mature R- channel network is absent during mid- to late- summer. The short Arctic ablation season, with surface melt rates comparable to those on temperate glaciers, may cause high basal meltwater input but insufficient time for R- channels to enlarge. The cold ice along the margins and at the terminus may lead to poor subglacial drainage in the lower ablation area, thus prolonging water Storage (although not totally inhibiting it). Some of this stored water is gradually released from storage during winter, as is evidenced by the growth of a large aufeis field in front of McCall Glacier in winter.