Introduction

The speed of glaciers responds to many physical processes at different space and time-scales. Here, variations in speed of a tide-water glacier are examined with respect to tidal period variations in sea-level at the terminus, and diurnal period variations in melt-water input at the surface. The spatial domain considered is the lower 10 km reach of Columbia Glacier, a region of high speed and high basal water pressure.

The purpose of this work is to interpret the calculations of speed in the terminal region that were derived from field experiments in 1984, 1985, and 1986. The flow of Columbia Glacier over weekly to seasonal time-scales is inherently unsteady because of important changes in the basal water system and calving flux. These changes are difficult to understand quantitatively. However, the glacier is quasi-steady for time-scales of a day or shorter and the data set for glacier speed and a variety of physical variables is relatively rich for these time-scales. Thus, an analysis of these motions can provide a window into the dynamics of the lower glacier.

Study Area

The field site is Columbia Glacier, a large tide-water glacier in south-central Alaska (Fig. 1). The glacier has an area of approximately 1100 km², and a length of about 65 km along the main branch from the Chugach Mountains to the edge of Prince William Sound. Reference PostPost (1975) determined that an irreversible retreat would occur if the glacier lost its footing on the terminal moraine. In anticipation of this retreat, a major research effort was begun in 1974 (Reference MeierMeier and others, 1978) and there have been numerous field experiments since then. The retreat began in 1977 and accelerated in 1984. By early 1986, the terminus was 2 km behind the terminal moraine and considerable glacier thinning had occurred, as shown in Figure 2 (Reference Krimmel and VaughnKrimmel and Vaughn, 1987). The water depth at the terminus is approximately 250 m and is increasing as the glacier retreats. In a recent paper, Reference Meier and PostMeier and Post (1987) summarized the state of knowledge of tide-water glaciers with Columbia Glacier as an example.

Fig. 1.

Photomosaic of Columbia Glacier for 4 August 1985, showing the location of the terminus in 1984, 1985, and 1986. The electronic distance-meter (EDM) station is located on Heather Island at H, the meteorological station is 100 m north-west of H, and the tide gage is on the east side of Heather Island. The reflectors for 1984 are at T84; for 1985 at T85A, T85B, and T85C; and for 1986 at T86. The numbers from 52 km to 66 km are distances along the center line from the head of the glacier. (Adapted from Reference Krimmel and VaughnKrimmel and Vaughn, 1977.)

Fig. 2.

Longitudinal profile of the terminus region of Columbia Glacier with a 5: 1 vertical exaggeration. (Adapted from Reference Krimmel and VaughnKrimmel and Vaughn, 1987.)

Data and Analysis

As part of a larger field program, distance was measured to reflectors on the surface of Columbia Glacier during the late summer of 1984 and 1985, and March 1986. The program in 1984 has been described in detail by Reference VaughnVaughn and others (1985). The 1985 program was similar, whereas the 1986 study only included meteorological and electronic distance-meter (EDM) data (see Reference Krimmel and VaughnKrimmel and Vaughn, 1987; Reference Walters and DunlapWalters and Dunlap, 1987). The EDM data sets vary from 2 weeks to 1 month in length and contain numerous gaps because of inclement weather. The reflector speed is calculated using cubic splines to interpolate the data sets to a standard time reference and to approximate the derivative of distance with respect to time. The speed is subjected to a low-pass digital filter with a cut-off of about 1.5 h, slightly longer than the sub-sampling interval of 1 h. Detailed information on the data-processing methods, and the analysis of the 1984 and 1985 data sets is in Reference Walters and DunlapWalters and Dunlap (1987). The focus here will be on the tidal period variations in speed.

The analysis of the 1984 field data (Reference Walters and DunlapWalters and Dunlap, 1987) showed that speed averaged 10 m/d at a single reflector and contained three noteworthy features: (1) a low-frequency variation with a time-scale of several days and amplitude of up to 6 m/d related to precipitation events; (2) diurnal and semi-diurnal variations with amplitude of about 1 m/d related to ocean tides at the terminus; and (3) diurnal variations with amplitude of less than 0.5 m/d related to surface melt-water input (Fig. 3). The tide height and the glacier speed are inversely related; that is, the speed is highest when the tide is lowest. A pressurization of the basal water system due to the variation of sea-level would have the opposite effect. If there is, in fact, a combination of processes, the retarding effect of a rise in sea-level is dominant.

Fig. 3.

Data for 1984. a. Ice speed showing data gaps; b. Solar insolation at the meteorological station; c. Tide stage at Heather Island; d. Precipitation at the meteorological station. Time is in Julian days from 1 January 1984. The tick marks are at 00.00 h of the day labeled.

In 1985, distance was measured to three reflectors spaced longitudinally up-glacier. The average speed varied from about 10 m/d at the lower reflector to about 7 m/d at the upper reflector (spacing of 3.3 km). This speed difference is explained by thinning of the glacier toward the terminus. The speed contained a low-frequency variation from the effect of the changes in precipitation from successive storms, a diurnal variation from the effect of surface melt, and diurnal and semi-diurnal variations from the effect of ocean tides (Fig. 4). The tidal signal was damped with distance from the terminus with an e-folding length of 2 km (where the amplitude is reduced by e ^-1). The signal related to melt-water input has an amplitude of about 5% of the mean speed and seems to be relatively constant in space, although the data are too sparse at the upper reflector for a quantitative analysis. Field experiments in 1987 (unpublished information of M.F. Meier and B. Kamb) indicate that there is about a 5% variation in the mean speed associated with the melt-water signal at a location near that of the upper reflector in 1985. In 1986, the speed at a single reflector averaged about 12 m/d with a conspicuous tidal period signal (Reference Krimmel and VaughnKrimmel and Vaughn, 1987).

Fig. 4.

Data for 1985. a. Ice speed at the upper reflector (T85C) showing data gaps; b. Ice speed at the middle reflector (T85B); c. Ice speed at the lower reflector (T85A); d. Tide stage at Valdez, Alaska. Time is in Julian days from 1 January 1985. The tick marks are at 00.00 h of the day labeled.

Results indicate that diurnal variations in melt-water input and tidal period variations in sea-level at the terminus are the dominant forcing processes in this frequency range. The melt-water signal has a continuous spectrum with a broad peak centered at 1 cycle/d. The tidal signal has a line spectrum with the strongest line at 1.93 cycles/d (M₂ tidal constituent). Because of the broad peak for melt water, most of the tidal lines are severely contaminated, particularly the diurnal lines. Thus, only the M₂ line is used in the calculations. Harmonic analysis was applied to the ice speed and the sea-level data for the three field experiments in order to extract the amplitude and phase of the M₂ constituent. The results indicate that the relative magnitude of the speed variations decreases with distance from the terminus, the magnitude at the lower reflector has some year-to-year variation, and the phase remains nearly constant.

An exponential least-squares fit was made to the 1985 data to derive an e-folding length of 2.0 km (Fig. 5). The amplitude at the farthest reflector was below the noise level in the data as estimated during the harmonic analysis. Hence the length scale depends upon the fit with the amplitudes at the other two reflectors, subject to the constraint that the amplitude is smaller than the noise level at the third reflector. The results for the length scale were not sensitive to the details of the fitting procedure.

Fig. 5.

Results for the harmonic analysis of ice speed and tide. Distance is measured from the width-averaged location of the terminus. The amplitude for ice speeds is at the M₂ tidal frequency and has been normalized by the tidal stage to form an admittance (m/d of ice speed per meter of tide stage). The triangle is the 1984 data point, circles are the 1985 data points, and the square is the 1986 data point. The solid line is the exponential fit to the 1985 data. The station locations are shown in Figure 1.

Next, the magnitude of the surface displacements is calculated from the spatial distribution of the speed variations (Fig. 5) and the bottom topography of the glacier (Fig. 2). These estimates are necessary for the development of the flow equation. Expanding the thickness and the longitudinal speed in terms of a mean plus a perturbation, and introducing these quantities into the continuity equation, one derives

where h ₁ and u _l are the perturbations in thickness and depth-averaged speed, respectively; h ₀ is the mean thickness of the glacier at a distance x from the terminus, and t is time. A second term in the gradient, the perturbation in thickness times the depth-averaged mean speed, h ₁ u ₀, is negligibly small. Using characteristic values for Columbia Glacier, h ₁ has a vertical displacement of about 0.003 m and a horizontal displacement of about 0.06 m at the terminus. The displacements decay up-glacier with a length scale approximately the same as that of u ₁, about 2 km. Thus, the speed perturbation results in small displacements in the geometry of the glacier. This leads to a perturbation expansion of the flow equation where perturbations in surface elevation and slope can be neglected.

Simple Model

In the first article in a recent series of papers, Reference Kamb and EchelmeyerKamb and Echelmeyer (1986) examined the influence of longitudinal stress coupling on flow. Based on previous work, they developed a longitudinal flow equation and expanded the dependent variables in terms of mean values and small deviations. The flow equation for the deviations is approximately linear.

With some modification, this work can be applied to the speed varlattons at Columbia Glacier. Certain assumptions are made for the purpose of simplifying the analysis and bringing out the essence of the behavior of the perturbations. These assumptions generally follow those of Kamb and Echelmeyer. The glacier is assumed to be adequately described by plane-strain flow of a wide glacier down an inclined plane. The local ice thickness h(x) and the surface slope α(x), which is positive when the surface slopes downward in the direction of flow, are functions of the distance x which is positive down-glacier. The slope α is small, such that sin α = α. The resulting vertically integrated flow equation (stress-equilibrium equation) is

where τ ′ _xx is the vertically averaged stress deviator, ρ is ice density, g is gravitational acceleration, and τ _B is bottom stress (Reference Kamb and EchelmeyerKamb and Echelmeyer, 1986, equation (4)).

Because the flow at Columbia Glacier is primarily by basal sliding, the relation between depth-averaged speed and basal stress is approximated as (Reference BindschadlerBindschadler, 1983)

where p _e = p _i – p _w is the effective basal pressure, p _i is the ice-overburden pressure, p _w is the basal water pressure, and c ₁ is an empirical constant. Note that p _i is a constant in time for diurnal periods. Using the same relation for the stress deviator as that used by Kamb and Echelmeyer (Reference Kamb and Echelmeyer1986, equation (5)), and substituting Equation (3) into Equation (2), the flow equation becomes

where η is a depth-averaged viscosity. With the assumptions presented here, and the fact that Columbia Glacier is dominated by sliding, the gradient in u is relatively constant with depth so that η is approximately equal to the effective viscosity.

At this point, the development departs from that of Kamb and Echelmeyer because of the differences in the flow problems under consideration. Kamb and Echelmeyer examined the effects of longitudinal coupling on perturbations in slope and ice thickness where the perturbations are smoothed by the stress. Here, the perturbations in h and α are negligible, but there is a perturbation in u introduced through the normal stress as a boundary condition at the terminus. As a result, u and p _e are expanded in terms of a mean and a perturbation, but h and α are not expanded.

where u ₀ is the mean speed, u _l is the perturbation in speed, p ₀ is the mean effective pressure, p ₁ is the perturbation in basal water pressure, ∊ is a small parameter, and 0( ) denotes order. Expanding Equation (3) for τ _B to order ∊

Expanding Equation (4), collecting terms of order E, and neglecting perturbations in the viscosity,

where τ ₀ = ( u ₀ p ₀ /c ₁)^1/m . For this simple case, let h and η be constant. Then Equation (7) can be written as

where the length scale L is given by (9)

The primary uncertainty at this point is the value for η . Taking m = 3, using a guess of 1 bar a for η (Reference Kamb and EchelmeyerKamb and Echelmeyer, 1986, fig. 3), and characteristic values τ ₀ = 1.5 bar, u ₀ = 10 m/d, and h = 250 m for Columbia Glacier, then L = 2.7 km, which is in reasonable agreement with an observationally based estimate of 2.0 km.

Because Equation (8) is linear, the two solutions arising from a perturbation in normal stress at the terminus and from a perturbation in basal water pressure superimpose, and are examined separately. The solution to Equation (8) for the case where p ₁ = 0 is straightforward:

where x ₀ is the location of the terminus and u ₁(x ₀) is determined from the normal stress balance at the terminus, which can be written as

where the second term on the right-hand side is the normal stress. The perturbation equation derived from Equation (11) is

where p _T is the pressure variation arising from the tides. Using Equation (10) in Equation (12) results in u ₁(x ₀) = p _T L/ η . Using the same parameter values as used in the estimate of L, the speed perturbation at the terminus is approximately 0.8 m/d, which is in good agreement with the value derived from the data analysis. In effect, the perturbation introduced at the terminus by changes in the normal stress is damped over a length scale L by the longitudinal stress. A similar damping in speed was observed by Reference Kamb and EchelmeyerKamb and Echelmeyer (1986) for a speed perturbation below an ice fall on Blue Glacier.

The seasonal and inter-annual changes in tidally induced glacier-speed variations arise through a variety of mechanisms. The tidal range changes over time as the various tidal constituents reinforce and cancel each other. However, this can be corrected for by choosing a sufficiently long time period for analysis (>28 d) or by normalizing by the tidal amplitude to derive an admittance which is plotted in Figure 5. A much more dominant cause of changes in speed variations is changes in the glacier mean speed and geometry. As can be seen in Equation (9), the length scale L depends on the speed, the ice thickness, the viscosity, and the basal stress. All of these change to some degree over the seasonal cycle, and inter-annually as the glacier retreats (Fig. 2). Moreover, the coupling between the normal pressure variations caused by the tide and the normal stress in the glacier also depends on the geometry and viscosity of the glacier. The amplitudes plotted in Figure 5 show the effect of these variations for the three data sets. The amplitude for 1984 (when the terminus shape was complicated) is relatively small, whereas the amplitude for 1986 (when the glacier had retreated) is relatively large. A quantitative analysis of these changes cannot be made with this simple model; rather, it must wait for a more generalized model that takes into account the spatial variation in the geometry and flow .

Whereas the tidal part of the speed variations is localized in time period and spatial scale, the melt-water-related speed variations encompass large time and space scales. Seasonal variations in water input to the glacier cause large seasonal changes in speed and geometry, particularly at the terminus. The seasonal variation in water input presumably modifies the basal water system and changes the relation between water input, basal pressure, and hence basal stress. There are insufficient data for Columbia Glacier even to begin analysis of these seasonal variations. The approach here is to consider those periods when the diurnal variation in melt-water input represents a perturbation in a larger base flow. Under these conditions, the perturbations in basal water pressure are also small when compared with the average basal water pressure. This is particularly true for the lower reach of Columbia Glacier, since it lies below sea-level and is therefore pressurized to at least this level.

To examine the effects of perturbations in basal pressure, Equation (7) is solved for the case where tidal pressure variations at the terminus are neglected (p ₁ = 0) and there are finite perturbations in basal water pressure (p ₁ ≠ 0). This equation is equivalent to equation (13) in Reference Kamb and EchelmeyerKamb and Echelmeyer (1986) with F(x) replaced by F(x) = p ₁/ p ₀. The solution is given by their equation (15), namely

If the pressure perturbation is a function of x, it is smoothed over a length scale L up-and down-glacier from the location x.

The analysis of observational data from 1985 indicates that the diurnal variations in ice speed that are correlated with melt-water input are a relatively constant percentage of the mean speed for locations in the lower reach of Columbia Glacier, with perhaps a slight increase toward the terminus. The amplitude of the variation amounts to about 5% of the mean speed. If the pressure variations are also spatially constant, the solution to Equation (13) reduces to

From speed calculations, the predicted variations in basal water pressure from Equation (14) would be about 5% of the effective pressure. However, the observed percentage was much larger than this, of the order of 100%. The amplitude of the diurnal basal water-pressure variations observed in the 1987 field experiment is about 10 m and the observed p _e is very small, also of the order of 10 m. This value for p _e is of similar magnitude to its error because of the uncertainty in ice density. These results imply that there is a large part of the basal stress that is not correlated with effective pressure as given by Equation (3). Other sources of stress that would not be affected by basal water pressure are form drag associated with irregular bed topography and lateral stress. In addition, Equation (3) could be modified by introducing a reference pressure p _r for effective pressure, where p _r = p _e + a constant that depends upon glacier geometry and flow. However, there does not appear to be a plausible physical reason to make this change. In essence, the basal water pressure is so high that sliding friction is not a major component of the stress at the bed.

This issue will hopefully be resolved from two directions: a detailed analysis of the speed and pressure observations from the 1987 field experiment, and numerical model experiments for the mean flow and the perturbations including full geometric variations.

Concluding Remarks

This simple model contains the conceptual essence of the effects of normal pressure variations at the terminus of a tide-water glacier. Much of the behavior of the speed perturbations observed in the lower reach of Columbia Glacier can be explained with a tidally induced normal pressure perturbation that is coupled into a speed perturbation by the longitudinal stress. The speed perturbation is damped over a coupling length of about 2 km. An analysis of the diurnal signals indicates that there is a large part of the bed stress that is uncorrelated with basal water pressure. Possible sources for the additional stress are form drag from irregular bed topography and lateral stress. The full quantitative assessment must await the numerical solution of the flow equation for the mean flow, and that for the perturbations, including all variations in geometry.