Surface motion

To investigate the potential influence of subglacial hydrology on glacier dynamics, surface motion was recorded across an area centred ∼1500 m from the glacier snout (Fig. 1) through which the eastern PDA is predicted to pass (Reference SharpSharp and others, 1993). Glacier surface motion was determined by repeatedly surveying to prisms mounted on stakes drilled into the glacier surface. In early June 1994 a network of 20 motion stakes was established. Rows of stakes were staggered to create a network of surface strain triangles (Fig. 1). Surveys were carried out using a Geodimeter 410 total station system with an accuracy in the distance measurement of ±(2 mm + 3 ppm). This instrument error combined with human measurement error led to horizontal and vertical accuracy in stake positions of 9 and 4 mm, respectively, during the study period over the distances surveyed. Further details of survey set-up and errors can be found elsewhere (Reference Mair, Nienow, Willis and SharpMair and others, 2001).

Surface strain rates

Surface strain rates were calculated from the deformation of strain triangles formed by the staggered rows of motion stakes (Fig. 1). For each strain triangle, the position of the centre of the strain triangle, the two perpendicular principal strains (ε _1,2) and the angle of orientation of the first principal strain relative to north (θ) were calculated. Strain rate was assumed to be homogeneous within the strain triangles. Reference RamsayRamsay (1967) describes in full the method by which the magnitude and direction of the principal strains were calculated. Longitudinal strain (ε_xx ), transverse strain (ε_yy ) and shear strain (ε_xy ) were calculated where x is the flow direction, which is north, and y is cross-flow direction increasing to the east. Vertical strain (ε_zz ) was calculated, assuming incompressibility of ice (the vertical dimension, z, is measured in m a.s.l). Strain rates were obtained by dividing the strains by the time period between surveys.

Limitations of vertical motion analysis

The greatest limitation in using Equation (1) with field data lies in deciding upon an appropriate spatial scale over which these measurements may be representative and in the uncertainties in the measured quantities. In most studies of ice dynamics, bed and surface slopes and surface velocities are averaged over spatial scales of about 10 ice thicknesses. However, the surface vertical velocity pattern during the spring event shows dramatically anomalous surface uplift concentrated over areas with a spatial scale of about one or two ice thicknesses (Fig. 4b). The coincidence of the event with a period of high meltwater input and discharge suggests that some of this uplift may be induced by hydraulic processes at the glacier bed. Previous investigations of the subglacial hydrology of the glacier have shown dramatic spatial heterogeneity in the glacier drainage system over spatial scales of about one ice thickness (Reference Hubbard, Sharp, Willis, Nielsen and SmartHubbard and others, 1995). The dilemma is that horizontal averaging of surface velocities and strain rates over large dimensions, of the order of 10 ice thicknesses, may increase our confidence in interpreting the basal flow regime from surface measurements, but that such a scale of horizontal averaging is too large to be directly associated with expected patterns of subglacial hydrological forcing. Conversely, interpreting the basal flow regime from surface measurements over the scale of expected patterns of subglacial hydrological forcing, i.e. about one ice thickness, means that measurements and interpretation are at a horizontal scale over which there is an expectation for strain rates to vary with depth. Therefore, to determine whether there is a relationship between subglacial hydrology and surface uplift, we adopt a compromise. Patterns of surface velocities within a study area with side dimensions of about five ice thicknesses are analyzed to identify areas with horizontal dimensions of one ice thickness in which basal flow anomalies are most likely to occur. Analysis over this intermediate spatial scale is theoretically justified by Reference Balise and RaymondBalise and Raymond’s (1985) model which suggested that over spatial scales of between one and five ice thicknesses, a basal vertical flow anomaly may induce a surface response with the same spatial pattern but with reduced amplitude. We then examine the extent to which strain rates are expected to vary with depth within areas of anomalous basal flow, and assess the likelihood that cavity opening at the bed contributes to the anomaly.

The spatial scale over which the bed slope, β, is averaged was chosen to be less than the spatial resolution of the stake network in case any local variation in surface velocity between stakes may have been due to a localized bed-slope anomaly. However, the down-glacier bed slope in this region of the valley is low and fairly constant over length scales an order of magnitude greater than this. Small-scale averaging of β does not change its value significantly and is therefore not considered a major limitation. Interpolated patterns of surface horizontal and vertical velocity are highly influenced by the spatial resolution of the stake network. The higher the spatial resolution, the more confident we can be of interpolated grid values. Since patterns of horizontal velocity (Fig. 3) and of surface uplift (Fig. 4) generally extend over areas larger than the distance between adjacent stakes (∼100–150 m) the spatial resolution of the stake network appears to be good enough to constrain the interpolated grid values effectively. Although strain rates are calculated on the assumption that strain is homogeneous across strain triangles with dimensions of ∼100–150 m, the grids are constructed by interpolating between the values of homogeneous strain rate located at the centre of each triangle. This interpolation therefore treats values which are really areal averages as point data. This has the effect of focusing strain-rate values around specific points even though each value is calculated from the relative displacement of three points distributed over a much larger area. However, this problem will be reduced by the relatively high density of strain-rate measurements within the stake network. The surface expression of the basal velocity anomaly may be spread over a larger area than the basal anomaly due to the attenuation of vertical velocities between the bed and the surface over short length scales, i.e. of the order of one ice thickness (Reference Balise and RaymondBalise and Raymond, 1985). Thus basal cavity opening may be concentrated in only a part of the areas identified by the following analysis. However, this effect may be partly cancelled out by the effect of interpolating strain rates between the centre points of strain triangles as discussed above. Given these limitations, the following quantitative analysis can at best qualitatively identify areas where subglacial cavity enlargement was most likely to be within a spatial scale of about one ice thickness.

Using strain rates to constrain identification of most likely areas of subglacial cavity opening

Again, following Reference Hooke, Calla, Holmlund, Nilsson and StroevenHooke and others (1989), if and u _b are both 0, . If and u _b = u _s, . Thus, in the absence of cavity opening, and assuming that 0 ≤ u _b ≤ u _s (i.e. excluding the possibility of extrusive flow), w _s and w _s − u _s tan β represent limiting values of . Furthermore, if were independent of depth, would equal . Thus when does not lie between the limiting values for , it is likely that is varying with depth.

For ease of interpretation we can categorize the relationships between w _s, u _s tan β and as follows. If lies between w _s and w _s − u _s tan β then it is possible that all surface vertical velocity could be due to vertical strain and basal sliding without any cavity opening or enlargement. We call this category 1 (Cat.1). If is higher than w _s − u _s tan β then w _s could only be due to vertical strain and basal sliding without cavity opening if strain rates decrease (i.e. more vertical compression) with depth. We call this category 2 (Cat.2). If is lower than w _s then surface velocity could only be due to vertical strain, basal sliding and no cavity opening if strain rates increase (i.e. more vertical extension) with depth. We call this category 3 (Cat.3). For each of the time periods, before the event (JD169–173), during widespread surface uplift (JD173–176) and during widespread surface lowering (JD176–180), the study area can be divided according to these categories (Fig. 6).

Fig. 6. Spatial categorization of the relationships between w_s, u_s tan β and . (a) JD169–173, (b) JD173–176, (c) JD176–180. Definition of each category is explained in the text. Differences (mc⁻¹) between and w_s − u_s tan β are overlain where Cat.2 is most prevalent (a,c). Differences (m d⁻¹) between and w_s are overlain where Cat.3 is most prevalent (b). Location of supraglacial stream-fed moulins (i.e.“active”moulins) marked.

If cavity opening occurred, it is most likely to have done so during the initial stages of the spring event (JD 173–176) within the Cat.3 area. At this time, was most strongly negative whilst w _s was most strongly positive over much of the study area (cf. Figs 4b and 5b(ii)). For this time period, the differences between and w _s are overlain on the Cat.3 region (Fig. 6b). The more negative these differences are, the greater the increase in strain rates with depth must be if no cavity opening occurs. Reference Hooke, Calla, Holmlund, Nilsson and StroevenHooke and others (1989) suggested that where their results showed a value of that was ∼0.003 m d⁻¹ less than w _s then this may be interpreted as indicating that was positive. Ice thickness at the study location on Storglaciären, Sweden, was similar to the maximum ice thickness in our study area. At Haut Glacier d’Arolla, differences between and w _s are >0.003 m d⁻¹ across almost the entire Cat.3 area during the period JD173–176. The biggest difference measured is >0.03 m d⁻¹, an order of magnitude greater than at Storglaciären.

Identification of Cat.3 areas will be affected by errors in measurements to and w _s. Estimates of sources of error are presented for the period JD173–176, i.e. the time of the highest surface uplift (Table 1: note that the values of the errors differ from Reference Mair, Nienow, Willis and SharpMair and others (2001), since this earlier paper discussed errors from the whole melt season whilst in this paper we quote errors from a specific short time period during the spring event; this whole analysis of vertical motion was carried out partly because vertical velocity errors were lower than average during the spring event). The total measurement error quoted in Table 1 assumes that all the sources of measurement error are cumulative. Since this is unlikely, the value of 0.016 m d⁻¹ is probably a high estimate for the measurement error. If the measurement error in Table 1 acted only to reduce the differences between and w _s uniformly across the Cat.3 area shown in Figure 6b, which is highly unlikely, then would still be less than w _s over a large area in the southeast of the study area (reducing the maximum difference to 0.016 m d⁻¹) and over a significant area further down-glacier and towards the glacier centre (reducing the difference here to 0.008 m d⁻¹).

Table 1. Estimation of errors in terms affecting the identification of areas of cavity opening or enlargement, , JD173–176

Since vertical strain-rate gradients are likely to express themselves as similar to horizontal gradients, horizontal gradients in can be used to help constrain vertical variations in . During the period JD173–176, in the southeast area of the stake network, and Z = 130 m, giving (Fig. 5b(ii)). Horizontal gradients in were greatest at this time, reaching a maximum of ∼0.00035 a⁻¹ m⁻¹. If increased linearly with depth at this rate then . If increased linearly to reach zero at the bed (with a depth gradient of ∼0.00049 a⁻¹ m⁻¹), then . In addition to the errors above, such increases in with depth could reduce the differences within Cat.3 areas in Figure 6b by a further 0.008 or 0.0112 m d⁻¹, respectively. Even an unlikely combination of maximum error and these high depth gradients would leave differences of >0.003 m d⁻¹ in this southeast region. To completely account for the difference between and w _s in this area, must increase with a depth gradient of at least 0.00166 a⁻¹ m⁻¹, i.e. almost five times greater than the highest measured horizontal strain-rate gradient. Furthermore, these drastic measurement-error and strain-rate corrections are required simply to bring the value for up to its lower limiting value which assumes no basal sliding. Since surface velocities during this time period are 400–500% higher than annual average surface motion (Fig. 3b), this assumption is very unlikely. It seems more likely that cavity opening or enlargement did occur within a broad axis of the study area, centred along ∼6650E, during the period JD173–176. By similar analysis, it is likely that similarly high positive differences between and w _s − u _s tan β, evident within Cat.2 regions during the period of surface lowering (JD176–180) (i.e. ∼0.02–0.03 m d⁻¹, Fig. 6c), are at least in part caused by cavity closure. This analysis strongly suggests that cavity opening or enlargement is widespread during the spring event and that the surface uplift cannot simply be explained as a consequence of measurement errors or by depth variation of vertical strain rates. The possible spatial patterns of cavity enlargement (and closure) are now discussed in more detail.

Relationship between spatial pattern of cavity opening and subglacial drainage-system configuration

The conclusion of the above analysis is further strengthened by the fact that, during the period of surface uplift (JD173–176), the differences between and w _s are highest in those areas located immediately down-glacier of clusters of active moulins (cf. Figs 1 and 6b). Very high subglacial water pressures, capable of promoting widespread bed separation, would be expected in these areas where increasing surface meltwater discharges are conveyed rapidly to the subglacial drainage system. More detailed analysis of the potential spatial patterns of cavity opening and closure identifies further evidence for a relationship with subglacial drainage-system configuration.

Areas where cavity closure is most likely during surface lowering (i.e. areas of largest differences in Fig. 6c) are located adjacent to areas where cavity opening is most likely during surface uplift (i.e. areas of largest differences in Fig. 6b). The fact that the areas of greatest opening and greatest closure are not coincident indicates that by the end of the spring event, on JD180, areas of the bed may have existed where cavities remained open, i.e. where net cavity enlargement had occurred. Areas of most likely net cavity enlargement can be identified by integrating the magnitude of the differences between and w _s − u _s tan β over the period JD176–180, then subtracting this from the time-integrated differences between and w _s for the period JD173–176. Net cavity enlargement may have occurred within areas where the residual at the end of the spring event is positive (Fig. 7). There is a strong association between the spatial pattern of positive values and the location of the eastern PDA predicted from the map of the up-glacier contributing area on the subglacial hydraulic potential surface. The locations of areas of most likely net cavity enlargement also lie close to the VPAs that have been identified from borehole water-level measurements (Reference Hubbard, Sharp, Willis, Nielsen and SmartHubbard and others, 1995; Gordon, 1996). The similarities between these completely independently derived datasets add further support to the idea that subglacial cavities were enlarged by high subglacial water pressures during the spring event.

Fig. 7. Identification of areas of most likely net cavity enlargement. The magnitude of the differences between and w_s − u_s tan β during surface lowering, JD176–180, are integrated through time and subtracted from the time-integrated differences between and w _s during surface uplift, JD173–176. Areas where the residual (measured in m) is positive will be those most likely to have experienced net cavity enlargement by the end of the spring event. Location of the eastern PDA is marked with dashed lines.

Relationship between subglacial equipotentials, net cavity enlargement and evolution of channelized drainage?

It appears that during the spring event supraglacial water accessed the bed of the glacier via moulins and was forced by hydraulic gradients to flow across those areas of the bed predicted by the mapping of subglacial equipotentials. Most of the water flowed along the eastern PDA. Large increases in water input to a distributed system may have led to very high subglacial water pressures. These in turn may have promoted bed separation and the enlargement of basal cavities as the glacier sole was displaced vertically upwards and rapidly down-glacier by hydraulic jacking. Net cavity enlargement would have increased the hydraulic efficiency of subglacial drainage along these axes after the spring event. Water would therefore subsequently have drained at lower pressures along these axes than elsewhere, causing an increase in the hydraulic gradient forcing water towards these axes from adjacent areas of the glacier bed. A greater proportion of subglacial discharge may therefore have drained along the axes created by the cavity enlargement. Thus the vertical and horizontal movement of the base of the glacier during the spring event, causing axes of enlarged basal cavities, may have created a positive feedback whereby greater discharges were subsequently able to drain along these axes.

Later in the melt season, removal of the surface snow-pack results in increases in the magnitude and amplitude of diurnal discharge cycles entering the moulins (Reference Nienow, Sharp and WillisNienow and others, 1998). Higher discharges will increase the size of the cavities further and enlarge drainage pathways between cavities by frictional melting (Reference KambKamb, 1987). It is likely that the drainage system will develop so that water flows through a network of small channels located within the broad axis of net cavity enlargement. Larger channels within the axis will capture water from smaller ones, causing drainage to eventually be concentrated in a main channel concentrated in a relatively small area of the bed (Reference RöthlisbergerRöthlisberger, 1972). This mechanism for generating a zone of concentrated flow may solve the problem of having to assume an initial channel when trying to model channel evolution over a melt season.

Many previous field and theoretical studies indicate that cavity enlargement during periods of high subglacial water pressure and enhanced basal motion plays a crucial role in the subsequent development of more efficient, low-pressure subglacial drainage. Reference KambKamb (1987) developed a theoretical model to show that if water-pressure perturbations within linked-cavity drainage systems are high enough and sustained for long enough, linked cavities may evolve rapidly by frictional melting into large channels. This theory has support from field evidence. Rapid increases in subglacial discharge during early-melt-season rainstorms have been suggested to play a key role in initiating the opening of subglacial channels (Reference Hock and HookeHock and Hooke, 1993; Reference Gordon, Sharp, Hubbard, Smart, Ketterling and WillisGordon and others, 1998; Reference Nienow, Sharp and WillisNienow and others, 1998). Such early-melt-season high subglacial discharge events have been frequently associated with short-term glacier speed-up events (Reference Röthlisberger, Lang, Gurnell and ClarkRöthlisberger and Lang, 1987). Reference Raymond, Benedict, Harrison, Echelmeyer and SturmRaymond and others (1995) found that at both a surge-type glacier (Black Rapids Glacier, Alaska, U.S.A.) and a non-surge-type glacier (Fels Glacier, Alaska) spring events signalled a rapid evolution of the subglacial drainage system. Changes in mean water and sediment discharges indicated the enlargement of existing, and the opening up of new, drainage pathways (Reference Raymond, Benedict, Harrison, Echelmeyer and SturmRaymond and others, 1995). Our results add strong support to these suggestions. We demonstrate that it may be possible to actually determine the locations of cavity growth and that these locations are coincident with areas in which channelized drainage is shown to exist later in the melt season.

Reference KambKamb (1987) considered that the enlargement of highly pressurized orifices connecting subglacial cavities was dominated by the mechanism of frictional heating rather than by bed separation. However, as cavities grow larger during the spring, they may start to link up, with the result that orifices no longer provide the main means of connection. Cavity enlargement, due to the vertical and horizontal movement of the glacier sole, may enable the drainage system to accommodate the large discharges required for the eventual development of subglacial channels without requiring such high water pressures as envisaged by Kamb’s model. Hence separation, rather than melting, is the key to rapid drainage development.