3.1. Ice deformation

Tiltmeters measure tilt from horizontal for two mutually perpendicular axes. The tiltmeter data are transformed from the tiltmeter frame to a map frame with three rotations (Fig. 4), enabling calculation of tilt from vertical, θ, and rotation angle, ψ (Reference BlakeBlake, 1992). It is impossible to resolve the azimuth, ϕ, from the tiltmeter data alone. Magnetometers were employed to determine this angle but they did not function properly (see section 2). Tiltmeter records are labeled by borehole name and height above the bed (in meters).

Fig. 4. Tiltmeter orientation in a map frame: are longitudinal, transverse and vertical up directions, respectively. The angles θ and ψ can be resolved by rotating the tiltmeter from the coordinate system to the coordinate system.

The tiltmeter records generally indicate steady deformation over periods of weeks to a year (Fig. 5). Superimposed on these general trends are large fluctuations and/or steps in tilt angle. Small gaps in the records are due to failure during data transmission; large gaps are the result of the tiltmeter turning off and restarting at a later date. This could be due to a loss of battery power, though we are unable to provide satisfactory explanations for many of the data gaps.

Fig. 5. Tilt angle as a function of time. The solid curves represent tiltmeter data; the dashed curves are model-derived synthetic tilt curves (see section 4.1). Dotted curves indicate time periods during which data were omitted during generation of the synthetic tilt curves. Tiltmeter records are labeled by borehole name and height above the bed (in meters).

N2-11 (Fig. 5a) was the only tiltmeter in the N2 borehole to operate for more than a few days. Its tilt rates were considerably higher than those observed by other tiltmeters. The negative initial slope of the tilt curve suggests that the tiltmeter was initially tilted up-glacier; vertical shearing probably caused the top of the tiltmeter to move towards vertical and the tiltmeter to achieve a minimum tilt angle as it passed through the vertical transverse plane.

Two tiltmeter records from the N1 borehole exceeded 100 days. N1-51 (Fig. 5b) experienced steady tilt rates throughout the year, with slightly higher rates in winter than in summer. N1-11 (Fig. 5c) exhibited very low tilt rates during summer, punctuated by steps in tilt of up to 1.5°.

Only limited tilt data are available from the CEN borehole. The short and sparse record from CEN-6.5 (Fig. 5d) indicates a decrease in tilt rate as summer progressed.

Tiltmeters in the S1 borehole exhibited highly variable behavior. S1-51 (Fig. 5e) experienced large fluctuations in tilt angle during early summer, including two large events on days 154 and 180. The tilt rate was constant between days 192 and 241. The S1-6 tilt record (Fig. 5f) spanned days 148–252. In general, its tilt angle steadily decreased during summer; superimposed on this trend is a positive departure between days 164 and 194. S1-6 was presumably tilted up-glacier. Two tiltmeter records from the S1 borehole were neglected from this analysis because their tilt angles often greatly exceeded design specifications.

Tilt data from the S2 borehole are also very limited: only S2-6 (Fig. 5g) yielded reliable data. It demonstrated steady tilt rates throughout summer, with a large drop in tilt on day 180. The tilt rates before and after this event were nearly identical.

Tiltmeters in the S3 borehole (Fig. 5h–k) provided the longest and most complete records. These tiltmeters were located farther from the bed than the tiltmeters in the other boreholes because they became stuck during deployment, as indicated by a drop in cable tension and a stabilization of the water-pressure readings at an elevation of 70.6 m above the bed. The uppermost tiltmeter, S3-170.6 (Fig. 5h), exhibited large tilt variations in summer. During winter (days 271–426), its tilt rate remained steady and maintained a slope similar to the mean slope of the summer data. In late winter its tilt rate rapidly increased. S3-120.6 (Fig. 5i) exhibited the steadiest deformation of all the tiltmeters. The tilt angle of S3-75.6 (Fig. 5j) oscillated during the first few weeks of summer and then slowly and steadily decreased. The lowest tiltmeter, S3-70.6 (Fig. 5k), demonstrated slow, steady changes in tilt throughout the study period, except for a large departure from the general trend between days 172 and 200. Tiltmeters S3-75.6 and S3-70.6 were, presumably, initially tilted up-glacier.

Finally, all of the tiltmeters exhibited small diurnal fluctuations in tilt during summer (Fig. 6) that cannot be seen at the scale of an entire tilt curve. These fluctuations are very similar to diurnal tilt fluctuations that have been observed by tiltmeters installed in till beneath other glaciers (Reference Iverson, Baker, Hooke, Hanson and JanssonIverson and others, 1999; Reference Kavanaugh and ClarkeKavanaugh and Clarke, 2006). In our measurements the maximum and minimum daily tilt typically occur at 0600 h and 1800 h, respectively, throughout the summer, though it should be noted that measurements were only obtained every 6 hours. The lack of drift in the timing of the maximum and minimum daily tilt indicates that the fluctuations cannot be attributed to the principal components of Earth tides. Furthermore, the fluctuations cease abruptly on or before day 270 (27 September 2002).

Fig. 6. Examples of diurnal fluctuations in tilt angle from three boreholes. All tilt curves showed in-phase diurnal variations during summer that stopped abruptly on or before day 270.

3.2. Surface velocity

The mean spring (days 120–136, 30 April to 16 May 2002) and summer (days 136–257, 16 May to 14 September 2002) surface velocities were approximately 1.2 and 1.5 times larger, respectively, than the mean winter (days 257–490, 14 September 2002 to 5 May 2003) surface velocities, across the drilling transect (Fig. 7). This is typical of velocities on Black Rapids Glacier (Reference Heinrichs, Mayo, Echelmeyer and HarrisonHeinrichs and others, 1996).

Fig. 7. Mean spring (30 April to 16 May 2002; squares), summer (16 May to 14 September 2002; circles) and winter (14 September 2002 to 5 May 2003; triangles) surface velocities across the drilling transect. The curves are model results (see section 4.1) for various relative contributions of the velocity and tiltmeter root-mean-square errors; the velocity data are not adequately reproduced when the inverse model depends only on tiltmeter data (not shown). There is only one curve for the spring velocity profile because no tiltmeter data were available for that time period, and so the inverse model depended solely on velocity data. The thick curves represent the solutions we used in our analysis: the velocity and tiltmeter errors were weighted by 0.25 and 0.75, respectively. The error in the velocity measurements is smaller than the line thickness.

Daily surface velocities near the N1 borehole during summer are presented in Figure 8. Velocities were between 0.15 and 0.2 m d⁻¹ (55 and 73 m a⁻¹) prior to the annual spring speed-up, which began on day 140. The maximum spring speed-up velocity of 0.55 m d⁻¹ (201 m a⁻¹) was reached on day 146. Velocity peaks were either small but long-lived (days 146, 153, 160 and 170) or large but shortlived (days 180, 196, 201 and 214). Similar velocity variations during summer have been observed previously on Black Rapids Glacier (Reference Nolan and EchelmeyerNolan and Echelmeyer, 1999a; Reference Truffer, Echelmeyer and HarrisonTruffer and others, 2001) and elsewhere (e.g. Reference WillisWillis, 1995). We do not know vertical velocities during this period because the GPS station was resting on the ice surface.

Fig. 8. Twelve-hour mean surface velocity (solid curve) compared to water level measured in the till (dotted curve) at the N1 borehole during summer. Velocities are given in m a⁻¹ for comparison with model results.

3.3. Water pressure

Although pressure transducers were installed in each borehole, records are only available from two boreholes. Data from this study were therefore supplemented with water-pressure measurements from a pressure transducer installed in the till near the N1 borehole (Truffer and Harrison, 2005).

The N1-till record (Fig. 9a) is lengthy, but it is also sparse due to wireless data transfer methods (Reference Harrison, Truffer, Echelmeyer, Pomraning, Abnett and RuhkickHarrison and others, 2004). Major peaks in water pressure were observed on days 146, 198 and 249, and a large drop in pressure occurred around day 175. Owing to the sparsity of the dataset, it is difficult to comment on magnitudes and rates of water-pressure fluctuations.

Fig. 9. Piezometric surface for (a) N1, (b) S1 and (c) S3 boreholes. The N1 pressure data were obtained from two pressure transducers installed in the till near the N1 borehole. The dashed lines represent the ice-overburden pressure at each borehole. The data in (a) were collected at irregular intervals due to wireless data transfer methods (Reference Harrison, Truffer, Echelmeyer, Pomraning, Abnett and RuhkickHarrison and others, 2004).

The S1 record (Fig. 9b) lasted from day 148 to day 170. Water pressure was generally high during this period, with major drops in pressure occurring on days 149, 156 and 163. Diurnal fluctuations in water pressure began on day 165.

The most complete record is from the S3 borehole (Fig. 9c). It indicates steady water pressure until day 160, followed by large diurnal fluctuations throughout the summer. Large drops in water pressure occurred on days 163 and 173; major peaks were observed on days 165, 170 and 180. Water pressure was generally low between days 190 and 200; in late summer it gradually increased until it reached a nearly constant value that was sustained throughout the winter.

4.1. Model description

Ice flow is modeled through a transverse cross-section (Fig. 2) 16 km from the glacier headwall (Fig. 1). The model is simplified by assuming that there is (1) no extension or compression along the flowline and (2) no vertical or transverse flow. The first assumption appears valid over seasonal timescales, during which longitudinal strain rates at the glacier surface are on the order of 10⁻³ a⁻¹ (Reference NolanNolan, 2003). This turns out to be more than one order of magnitude smaller than the shear strain rates estimated by our FE model at the location of the uppermost tiltmeter (S3–170.6) and should therefore have little effect on the model solutions, especially due to the non-linearity of the flow law (see Equation (2) below). The second assumption requires a flat transverse surface profile, a straight channel and little influence from tributary glaciers, all of which are reasonable approximations for this site. These assumptions reduce the ice-flow equations to a non-linear Poisson equation:

(1)

where η is the stress-dependent viscosity, ρ is ice density, α is the mean out-of-plane surface slope, g is acceleration due to gravity, u is the out-of-plane velocity, y is the transverse direction and z points upward and is perpendicular to the mean glacier slope. A power-law rheology (Glen’s flow law) is used to describe the non-linear viscosity of ice:

(2)

where A is the flow-law parameter, n is an empirical constant, is the second invariant of the strain-rate tensor and = 10⁻¹⁵ a⁻² is a finite viscosity parameter used to prevent infinite viscosity at low stresses. We set α = 1.8°, A = 0.1 a⁻¹ bar⁻³ (3.17 × 10⁻²⁴ s⁻¹ Pa⁻³) and n = 3. The surface slope, α, was estimated from GPS measurements. The flow-law parameter, A, was determined by adjusting its value until a flow model with no sliding along the bed produced the lowest measured velocities at this site during the past 30 years (Reference Heinrichs, Mayo, Echelmeyer and HarrisonHeinrichs and others, 1996).

The model domain is derived from a radio-echo sounding (RES) profile (Reference GadesGades, 1998) and borehole depths. The errors in the RES data are on the order of ±5 m, but the close agreement between our borehole depths, known to better than ±0.1 m, and the RES data suggests that the bed profile is known to better than ±5 m. The surface is assumed to be stress-free, which results in a Neumann boundary condition ( = 0, where is the outward normal vector), and we assume a Dirichlet condition (u = f(y)) for the bed.

We seek a basal velocity function, f(y), that minimizes the error between model results (velocities and model-derived synthetic tilt curves, described below) and measured surface velocities and tiltmeter data. In reality, basal velocity is a function of the stress field, bed roughness and other physical characteristics of the bed, which would require a mixed (Robin) boundary condition or, in the case of till deformation, a stress boundary condition. We adopt the Dirichlet boundary condition in order to derive the basal velocity distribution, and ultimately the basal stress distribution, without attempting to make a conclusive statement about the nature of the actual boundary condition. Note that, for the purposes of an inversion, it is irrelevant which basal boundary condition is assumed. Existence and uniqueness of the forward problem (i.e. solving Equations (1) and (2) with the stated boundary conditions) were proven by Reference Colinge and RappazColinge and Rappaz (1999); thus, the velocity field that best reproduces all the data is uniquely determined regardless of the boundary condition implemented in the model. We therefore choose to use a Dirichlet boundary condition along the bed because it provides for fast convergence and is easier to implement than mixed or stress boundary conditions.

Synthetic tilt curves are generated according to the procedures outlined in Reference Gudmundsson, Bauder, Lüthi, Fischer and FunkGudmundsson and others (1999). The modeled velocity gradients (∂u/∂y and ∂u/∂z), initial tilt angles, θ _o, and initial azimuths, ϕ ₀, are required to generate synthetic tilt curves. We use the initial tilt angle from a given tiltmeter record and the modeled velocity gradients at the tiltmeter’s location, and iteratively search for the value of ϕ _o that minimizes the root-mean-square (rms) error between the synthetic and measured tilt curves.

Determining the basal velocity distribution from measurements of ice deformation and surface velocity is an inverse problem. A basal velocity distribution can be found that reproduces the data exactly; potentially there are many such solutions. Due to errors in the data and uncertainties in the physical model, a solution that reproduces the data exactly would probably be unrealistic. In general inverse theory, a preferred solution is one that fits the data to within a given error and also minimizes some undesired quality of the unknown boundary condition (Reference ParkerParker, 1994). In the problem at hand, a reasonable solution might be one that forces the basal velocity distribution to be a smooth function with small first derivatives (Reference TrufferTruffer, 2004);the data are then fitted to within a given error by minimizing a norm that reflects the smoothness of the function.

Applying general inverse theory to a highly non-linear, 2-D FE model with an irregular boundary is non-trivial and computationally expensive. To simplify the problem, we seek a smooth basal velocity function that can be approximated by a fourth-order polynomial, with the additional assumption that the velocity equals zero at the margins. Although this choice of basal velocity function is somewhat arbitrary, it appears to be the simplest that can satisfactorily reproduce our measurements without restricting the location or magnitude of basal velocity maxima and minima. (We also tried higher-order polynomials and found the results to be very similar.) We are forced to specify the velocity at the margins because the actual velocity there is unknown, and leaving it unspecified results in very unrealistic solutions. However, the velocity specified at the margins appears to have very little effect on the model solution, especially near the glacier center line (see appendix A in Reference AmundsonAmundson, 2006). The basal shear stress distribution, which we ultimately desire, can be computed from the model solution.

Ideally the flow field would be determined as a continuous function of time. This is impossible with our velocity data because they are of insufficient temporal resolution to constrain the inverse model over short time periods. Furthermore, all our measurements were restricted to a transverse cross-section, thus forcing us to neglect longitudinal strain. This approximation is adequate for investigating strain rates over periods of months or longer in our study area (discussed above). However, over shorter time periods such as during the spring speed-up (Reference NolanNolan, 2003) and in summer (when velocity variations are large), longitudinal strain may be important. Our model also does not take into account viscoelastic effects that may be associated with short-lived motion events.

The data were therefore divided into spring (days 120–136, 30 April to 16 May 2002), summer (days 136–257, 16 May to 14 September 2002) and winter (days 257–490, 14 September 2002 to 5 May 2003), corresponding to trips to the field that occurred in spring and autumn 2002, and spring 2003. This division roughly agrees with changes from above-freezing to below-freezing air temperatures at nearby Gulkana Glacier (Fig. 3). In this manner we limit our investigations to seasonal changes in the mean stress field. Acknowledging these limitations, we will, however, qualitatively discuss short-term stress variations by considering tilt deviation from the synthetic tilt curves.

The coefficients of the seasonal basal velocity functions were determined simultaneously using a multidimensional unconstrained non-linear optimization (Nelder–Mead method, implemented with Matlab’s fminsearch) to minimize the error between data and model results. We coupled the summer and winter models by maintaining continuity in the tilt angle and azimuth of the synthetic tilt curves at the summer/winter transition; the spring model was run separately because it was based entirely on surface velocity data.

The modeling approach is summarized as follows:

1. 1. Guess the coefficients of the summer and winter basal velocity functions. Only three coefficients are required for each, since we are requiring the basal velocity functions to be fourth-order polynomials with zero velocity at the margins.

2. 2. Solve the spring, summer and winter FE models. Export strain rates and surface velocities.

3. 3. Calculate the percentage rms error between modeled surface velocities and observations (later referred to as the velocity error).

4. 4. For the summer and winter models, generate synthetic tilt curves for each tiltmeter record using the modeled velocity gradients and compute the percentage rms error between synthetic and measured tilt curves. The percentage rms error of each tiltmeter is weighted by the number of data points in the respective set to give the longer records more importance. Sum the weighted rms errors and divide by the total number of data points in all the records (later referred to as the tiltmeter error).

5. 5. Compute the total percentage error using some combination of the velocity and tiltmeter errors. Since it is not clear how to do this, we considered several possibilities (see section 4.2).

6. 6. Iteratively search for the coefficients of the basal velocity functions (summer and winter) that minimize the total percentage error using the Nelder–Mead method.

4.2. Model results

The model captured the seasonal tilt rates and surface velocities well (see Figs 5 and 7), although the agreement between synthetic and measured tilt curves could partly be attributed to the unknown initial azimuth of the synthetic tilt curves. Modeled velocity fields are shown in Figure 10a–c; the corresponding octahedral stress fields were determined from the modeled velocity gradients and are shown in Figure 10d–f. The results suggest that basal shear stresses in a zone about 500 m north of the deepest point in the channel were approximately 0.07 bar and 0.16 bar (7% and 16%) lower in spring and summer, respectively, than in winter (Fig. 11), and that basal shear stresses near the margins were correspondingly higher in summer than in winter (thus allowing the glacier to remain in stress equilibrium). This stress redistribution is qualitatively similar to that predicted from forward models assuming a Coulomb-friction rheology for the till (Reference Truffer, Echelmeyer and HarrisonTruffer and others, 2001). However, it should be noted that the ‘spring’ period was a short period that immediately preceded the spring speed-up, an event that is typically associated with longitudinal extension (e.g. Reference NolanNolan, 2003). Thus, the assumption of no longitudinal strain may be invalid during this period.

Fig. 10. Model results (shown in order of increasing velocity): (a–c) mean winter, spring and summer velocity fields (m a⁻¹); (d–f) mean winter, spring and summer octahedral stress fields (×100 kPa). Flow is into the plane.

Fig. 11. Stress difference plots. The winter octahedral stress field is subtracted from (a) the spring octahedral stress field and (b) the summer octahedral stress field. The shaded regions indicate areas where stresses were lower in spring and summer than in winter. The values given on the contours are in units of 100 kPa; note that the plots have different contour intervals.

The model results were largely independent of the relative weights assigned to the velocity and tiltmeter errors in the inverse model (see Fig. 7), although surface velocities had to be included in the analysis to ensure a good fit to all the data. Removing tiltmeter data from the analysis did not significantly affect the solutions. This provides tentative support for inverse models based solely on surface data, although in our case the agreement between the synthetic and measured tilt curves greatly increases our confidence in the model results. Hereafter we have weighted the velocity error by 0.25 and the tiltmeter error by 0.75, as this provided the best fit for both velocity and tiltmeter data.

Solution uniqueness cannot be proven for most nonlinear inverse problems (Reference ParkerParker, 1994); it is therefore difficult to determine whether our model solutions represent global or local error minima. However, the coefficients of the basal velocity functions are at least locally well constrained (Reference AmundsonAmundson, 2006).

5.1. Rapid motion events

Seasonal velocity variations at Black Rapids Glacier appear to be the result of basal stress redistribution between a well-defined region north of the center line and the margins. This well-defined region coincides with the area where Reference Nolan and EchelmeyerNolan and Echelmeyer (1999a, b) observed the till becoming seismically transparent during lake drainages and associated rapid motion events. The bed there seems to be highly susceptible to changes in water flux, maybe due to the several lakes along the northern margin of the glacier that might drain subglacially and significantly influence the subglacial drainage system.

The influence of the drainage events on the basal stresses can be investigated with our FE model by making a few assumptions:

1. 1. The velocity spike on day 180 (Fig. 8) can be attributed to a lake drainage event since it is similar in magnitude and duration to the velocity events that Reference Nolan and EchelmeyerNolan and Echelmeyer (1999a) and Reference Truffer, Echelmeyer and HarrisonTruffer and others (2001) related to lake drainages. These events occur every summer.

2. 2. The velocity prior to day 180 is very close to the mean velocity during the following winter, and therefore the basal velocity immediately preceding the motion event can be approximated by the previously determined winter basal velocity distribution.

3. 3. Water from lake drainages is routed through the region near the N1 borehole.

4. 4. The water influx from the drainage events is sufficiently large to locally decrease the normal stresses imparted to the bed by the glacier (Reference Nolan and EchelmeyerNolan and Echelmeyer, 1999b). This results in ice–bed decoupling, which can occur at pressures beneath the ice overburden pressure (e.g. Reference Iverson, Baker, Hooke, Hanson and JanssonIverson and others, 1999, Reference Iverson2003).

Thus, we numerically decouple the ice from the bed near the N1 borehole (by setting the local basal shear stresses equal to zero) and specify the basal velocity along the rest of the bed according to the previously determined winter basal velocity distribution. By adjusting the length and position of the decoupled zone until the modeled velocity at the N1 borehole agrees with the velocity observed on day 180, we can estimate the basal shear stress distribution during a lake drainage event.

Under these assumptions, the section of the boundary that must be decoupled in order to achieve the velocity on day 180 is very close to the length of the bed over which mean basal shear stresses are lower in spring and summer than in winter (Fig. 12). The same surface velocity can be obtained by shifting the region of the bed that is decoupled from the ice by a few hundred meters in either direction; increasing or decreasing its length by more than 100 m significantly reduces the agreement between the modeled velocity at the N1 borehole and our observations. Although this model was not well constrained (we used a surface velocity measurement at just one point) and did not take into account additional till deformation or basal sliding that could result from increasing the shear stresses adjacent to the decoupled zone (Reference Mair, Willis, Fischer, Hubbard, Nienow and HubbardMair and others, 2003), it does reinforce the contention that large increases in surface velocity (over a variety of timescales) are due to stress transfer from a well-defined zone near the N1 borehole toward the margins.

Fig. 12. Model-derived mean summer and winter basal shear stresses (dashed and dotted curves, respectively) and a possible basal shear stress distribution of a decoupling event (solid curve). The spring basal shear stress distribution is not shown because it is difficult to see the difference between the winter and spring stress distributions at this scale. The negative shear stresses and the large jump in stress observed at the zero basal shear stress/specified velocity transition are due to a mathematical singularity (Reference Hutter and OlunloyoHutter and Olunloyo, 1980).

Observations of till deformation (Reference Truffer and HarrisonTruffer and Harrison, 2006) agree with this model in that little till deformation appears to have occurred near the N1 borehole during the short-lived motion events. However, enhanced till deformation, which occurred almost exclusively during the spring speed-up, would also have produced stress redistribution toward the margins (Reference Truffer, Echelmeyer and HarrisonTruffer and others, 2001). The spring speed-up and the short-lived motion events that we attribute to lake drainages (Fig. 8, days 142–163, 180, 196, 201 and 214) significantly influenced the mean velocities during summer; the mean velocity at the N1 borehole during these periods was 192 m a⁻¹, while the mean velocity during the rest of the summer was 46 m a⁻¹, about 7 m a⁻¹ slower than the mean velocity during winter 2002/03. These rapid motion events, which we suggest are caused by rapid stress redistribution toward the margins, are thus responsible for the overall larger ice displacement in summer than in winter.

5.2. Diurnal tilt fluctuations

Tiltmeters react to a strain field, which cannot be determined uniquely from the tiltmeter measurements without additional data or assumptions. It is thus difficult to quantitatively interpret the diurnal fluctuations in tilt angle observed by all of the tiltmeters (Fig. 6), as we are unable to neglect any of the strain-rate components over short timescales (see section 4.1). By considering the effect of each strain-rate component on the tiltmeter readings we can, however, attempt to qualitatively explain the diurnal fluctuations in tilt. For example, diurnal fluctuations in shearing across the glacier would cause some tiltmeters to achieve maximum daily tilt at the same time that others are achieving minimum daily tilt, except in the very unlikely event that all tiltmeters were oriented to the same side of the vertical longitudinal plane. Thus, shearing across the glacier is an unlikely mechanism for the diurnal tilt fluctuations. Fluctuations in shearing along the flowline, for which we have no information, could partly account for the diurnal tilt fluctuations. However, the negative tilt rates during the diurnal tilt fluctuations (for tiltmeters oriented down-glacier) cannot be explained by diurnal fluctuations in shearing along the flowline, as that would imply up-glacier deformation. This suggests that longitudinal strain is also important.

The diurnal tilt fluctuations are consistent with observations and a mechanism proposed by Reference Sugiyama and GudmundssonSugiyama and Gudmundsson (2003). They obtained detailed observations of vertical strain in Unteraargletscher, from which they concluded that increased water influx to the bed during the day causes the glacier to accelerate; the acceleration is amplified up-glacier because the drainage efficiency is lower there, thus causing longitudinal compression. At night, when water influx decreases, the longitudinal flow speeds become more uniform and motion is dominated by shearing. This mechanism could cause tiltmeters to achieve maximum daily tilt in the morning and minimum daily tilt in the evening, similar to our observations (Fig. 6). However, the diurnal increases in tilt angle exhibited by tiltmeters oriented up-glacier (e.g. Fig. 6b) suggest that some extension occurred at night. This does not contradict the mechanism proposed by Reference Sugiyama and GudmundssonSugiyama and Gudmundsson (2003); compression could result in vertical thickening that is larger down-glacier, thus causing subsequent extension. Further supporting this mechanism on Black Rapids Glacier is the observation that the diurnal fluctuations in tilt angle abruptly cease on or before day 270, which probably coincides with the onset of winter conditions (Fig. 3) and decreased water influx.

5.3. Non-steady deformation

In our model we have assumed that ice deforms as a continuous medium. However, most tiltmeters experienced short-lived tilt events that are probably indicative of nonsteady deformation. Some of these events were localized and had little or no correlation with other tiltmeter, velocity or water-pressure records. For example, the large fluctuations in tilt angle during summer observed by tiltmeters N1-11, S1-51 and S3-170.6 (Fig. 5c, e and h) are poorly correlated with any of the instrument records. There were, however, at least two instances in which tilt events were well correlated between two instruments located in different boreholes but were not observed by more proximal tiltmeters: the large changes in tilt exhibited by S1–51 and S2-6 on day 180 (Fig. 5e and g) and the slope reversals between days 170 and 180 in the S1-6 and S3-70.6 tilt curves (Fig. 5f and k). Although the lack of correlation between most tilt events could be due to poor coupling between the tiltmeters and the ice, the strong synchronicity of the diurnal fluctuations in tilt (Fig. 6) and the strong correlation between the S1-6 and S3-70.6 tilt curves suggest that the tiltmeters were well coupled to the ice.

Other tiltmeter records exhibit similar behavior that cannot be explained by steady-state flow models. For example, tiltmeters installed near Jakobshavn Isbræ, Greenland, indicated overthrusting at depth (Reference Lüthi, Funk and IkenLüthi and others, 2003), and some (but not all) tiltmeters installed in a single borehole in Unteraargletscher, demonstrated quasi-repetitive tilt oscillations (Reference Gudmundsson, Bauder, Lüthi, Fischer and FunkGudmundsson and others, 1999). More knowledge of small-scale glacier mechanics and the coupling of tiltmeters to ice will be required if we are to be able to interpret all the fluctuations in tilt rate observed by a tiltmeter. Also, it remains unclear how much non-steady deformation, such as faulting, influences the net motion of a glacier.