Equations

We assume the "shallow-ice approximation" and a power-creep law for the behaviour of the ice which gives the depth-averaged deformation velocity, U_d, as a function of thickness, H, and basal shear stress, τ, to the power n. We will use here the analytical expression of Reference LliboutryLliboutry (1981):

(1)

(2)

(3)

B(T) is the function of how the ice deformation depends on ice temperature via an Arrhenius-type relation with and p = n - 1 + kGH. Here G is the bottom temperature gradient taken as 0.022°C m⁻¹ corresponding to a 50 mW m⁻² geothermic flux, T_b is the mean bottom temperature, T_m is the melting temperature at the bottom of the ice sheet (T_m = 273 - H/1503 expressed in R Q is the activation energy for creep taken as 60 kJ mol⁻¹, R is the universal gas constant (8.31 J mol⁻¹ K⁻¹), p is the density of ice, g is the gravitational constant and B₀ the flow-law coefficient which seems to depends on crystal structure, orientation of the c axes and impurity concentration.

Lliboutry obtained Equation (3) from an analytical integration of du/dz based on a relation between temperature and depth. The error due to this analytical integration with respect to a numerical integration is less than 1 %. The only major assumptions are that the rheological parameters B₀, n and Q are constant through the ice thickness.

For steady-state flow with neither convergence nor divergence, the balance velocity U_b is given by

(4)

where b is the accumulation rate.

Assuming equality of these two velocities we have, from previous equations, a relation between the thickness and the shear stress:

(5)

where f is the accumulation integral divided by B(T).

Most models are based on some form of these relations. Although the flow law is well supported by mechanical experiments, its use in models poses a few problems (Reference AlleyAlley, 1992). The different parameters n, B₀ and Q have a large set of empirical values which can vary with the state of the ice (temperature, viscosity, strain rate, types of creep, etc.). The scale of calculation of the shear stress is also a "parameter" which can have an influence on the validity of the law. In this paper we endeavour to use the precision of our data to estimate these different parameters.

Input values

Spatial scale for calculating the slope

The surface slope and its direction play a role in estimating the magnitude and direction of ice flow. Due to the presence of undulations at the 10 km scale, the length of smoothing of the topography is critical. Many authors (Reference PatersonPaterson, 1981; Reference Cooper, Mclntyre and RobinCooper and others 1982; Reference Young, Goodwin, Hazelton and ThwaitesYoung and others, 1989) suggest using a 50 km scale when calculating the shear stress to smooth the effects of the longitudinal stresses. Investigations on the influence of spatial smoothing with our data led to the same result. In this study we will use the surface slope calculated from the ERS-1 altimetric surface by linear regression at each point at a scale of 20 times the thickness. Using this smoothing, Figure 3b shows the variation of basal shear stress along the Terra Nova Bay-Dome C profile.

Fig. 3. Quantities involved in the theory from Terra Nova Bay (km 0) to Dome C (km 1160): (a) Basal temperature estimated from the thermodynamic equation (solid line) and the basal melting temperature correction (dashed line); ( b) shear stress calculated at a 20 times thickness scale from the ERS-1 geodetic topography; and (c) strain rate calculated from balance velocity and thickness.

Accumulation rate, basal temperature and shear strain rate

Accumulation rates are interpolated along the profile using maps from Reference Radok, Jenssen and MclnnesRadok and others (1987). The bottom temperature (Fig. 3a) is estimated from the thermodynamic equation under a steady-state assumption and considering dissipation, vertical diffusion, vertical and horizontal advection. More details can be found in Reference RitzRitz (1987), Reference Remy, Ritz and BrissetRemy and others (1996) or Reference Remy and MinsterRemy and Minster (1993). The quantity U/ H derived from balance velocities and the thickness is supposed to be proportional to the mean shear strain rate ε through the ice thickness (Fig. 3a). We will use a "temperature-dependent shear strain-rate" ε′ = (p + 2)U _b/H exp[k(T _b – T _m)] in this paper. With this notation the Glen flow law can be written:

(6)

For a flow law invariable along the profile

Assuming that the parameters involved in the flow law (n, B₀, Q, spatial scale of shear-stress calculation, etc.) are constant along the profile, linear regression between r and ε’ for the whole profile provides a statistical estimation of the parameters n and B₀, with n = 1.5 and B ₀ = 0.3 ba^−1.5 a"¹ and a correlation coefficient of r = 0.42. The influence of the n in the p parameter (see Equation (3)) is negligible for the regression. These n and B₀ values are consistent with those found in the literature.

We now compare the along-profile variations of the two terms of Equation (5) to examine the validity of our hypothesis and of the flow law with all the parameters (calculated or estimated). Comparison of these two quantities (Fig. 4) shows strong discrepancies. If we focus on the part of the profile which corresponds to the plateau (i.e. after km 200), we see that the righthand term decreases slowly up-flow, representing the decreasing value of the accumulation-rate integral. On the contrary, the lefthand term, H2Tn shows an up-flow trend of either increasing magnitude (n = 1) or constant trend (n = 3). Any feasible estimation of the accumulation rate, or of the bottom temperature, cannot explain this difference. To explain such discrepancies both a positive bottom temperature from 200 to >500km and a temperature decreasing from –10°C to –60°C from 700–1160 km would be needed. For the accumulation-rate integral similarly unrealistic values are needed, diminishing by a factor 10 from km 900 (near Dome C) to km 200 (edge of the plateau).

Fig. 4. Comparison of the left and right term of Equation (5) in 10 6 m² bar". Solid line corresponds to the righthand term and is equal to the integral of accumulation rate divided by B(T). The three others curves correspond to H²Tn for n=l,2 and 3, respectively.

This difference may be explained if the flow-law coefficient B₀ has a variation of a factor of 100 and the same global trend as the sum of accumulation rate (see Equation (5) and Fig. 4). It is difficult to asses the evolution of B₀ along the profile since its behaviour depends on various factors or complex processes (fabric, grain-size, impurities, anisotropy, orientation of the c axes, strain, etc.). However, the large-scale trend of B₀ is qualitatively in agreement with Reference Pimienta and DuvalPimienta and Duval (1987) who note that as the flow increases, a strongly oriented fabric develops, which in turn enhances the flow. Figure 5 shows the magnitude of B₀ required to satisfy Equation (5) (with n = 1.5).

Fig. 5. Magnitude ofB0 ( bar-n a–1) needed to satisfy Equation (5).

For a flow law which varies with the location

Figure 6a shows the shear stress, τ, versus the temperature-dependent shear strain rate ε’. Values of n and B₀ can be calculated from the slope of the points in this plot. Using a linear regression with a 100 km wide window, the parameters B₀ and n were calculated each kilometre from km 50 to km 1110 (Fig. 7). From these two figures we can see that the observed behaviour of the ice along the 1160 km long profile is more complicated than can be explained by a constant flow law. In spite of the apparent confusion of points in Figure 6a, some significant features can be observed. Most important is that the flow law seems to be valid over distinct segments. Five separate areas show a consistent value of n (Fig. 7), each of them corresponding to a maximum of correlation. The mean values of n and B₀ for each of these intervals with their associated correlation, r, are given below.

Fig. 6. (a) logT(bar) vs log ε’ (a⁻¹), (b) Bedrock geometry (solid line) and surface slop ε (dashed line). (c) Path of the profile on contour elevation map. Below 2300 m (bold contour) the contours are each 15 m, above 2300 m contour interval is 5 m.

Fig. 7. Linear regression between log τ and loge’ each km along the profile with a 100 km window. The solid line corresponds to the Glen flow exponent (i.e. regression coefficient) and the dashed line to the correlation coefficient, r, of the regression (right scale). Dotted lines indicate significant correlation at the 99% confidence level.

We consider that the first 350 km does not fit our hypothesis, because zone 1 is a convergence zone and at the beginning of the plateau the profile from 190 to 350 km is not aligned in the flow direction.

With the exception of zone 4 (discussed below), the n value seems to increase from the coast to Dome G. Reference Remy, Ritz and BrissetRemy and others (1996) found a value of n around 1 for temperatures less than –10°C and increasing n values for warmer temperatures. One can see from Figures 3a and 7 that the large-scale behaviour of n and T_b is similar from km 200 to Dome G.

Sliding

Reference Siegert and RidleySiegert and Ridley (1998), using an analysis of radio-echo sounding, show that the flat-surface region (133° E, 73–74° S) lying over the Adventure Subglacial Trench had a relatively bright return from the ice substrate, but no mirror-like reflection. They conclude that water is not trapped within a subglacial hollow and infer that the transport of water probably occurs through an unknown basal hydrological system. This is probably the channel that is crossed by our profile at km 850. This region corresponds to a maximum in thickness (>3700 m) and should have a positive basal temperature (Fig. 3a). A careful examination of the surface topography of this region (between D and E in Fig. 6c) shows a flatter zone. There is also a 10–20 km wide channel in the topography at the eastern border of the flatter zone which reaches the Adventure Subglacial Trench and could mark the position of the subglacial hydrological system. The occurrence of sliding produces a flatter surface and could explain the strong diminution of τ (between D and E in Fig. 6a) and the low value of n in zone 4.

Local sliding could also partly explain the anomalies observed between points C and D and between points Fand G where negative values of n are found (Fig. 7). In particular, these anomalies occur when there is a decrease of τ at constant ε’ which seems to correspond to basal sliding (see above between D and E). The presence of relatively high radar-echo reflectivity in both cases (in spite of negative temperature, Fig. 3a) suggests the presence of water at the bottom of the ice. Local phenomena such as enhanced geothermal flux or strain dissipation could generate the right conditions for the existence of water in the ice. This also applies for local anomalies from the flow law in the Dome C region (for instance in Figure 6a between km 1123 and 1135 where τ decreases between points A and B). This behaviour could be a consequence of the subglacial channels which occur beneath Dome C (Reference Remy, Shaeffer and LegresyRemy and Tabacco, 1999).

Anisotropy

Anisotropy is suggested by Figure 6a between points G →H and E → F (or between zones 2 and 3). At large scale, zones 2 and 3 have the same behaviour with a mean slope around 1.5 but we observe a clear shift between them. Although the ice is subject to shear stress of the same order of magnitude and the temperature is taken into account in ε’, these two zones show strong differences in the intensity of the shear strain rate. Changes in the ice-crystal-orientation fabric between these two regions can lead to such behaviour, and correspond to a change of B₀.

In our case, B₀ is enhanced by a factor of 4.4. Reference Jacka and BuddJacka and Budd (1989) found that in isotropic ice the development of anisotropy leads to an enhancement factor of 3 when isotropic ice is subject to compression and a factor 8 in shear. Reference Remy, Ritz and BrissetRemy and others (1996) found an enhancement factor of 3.

Longitudinal stress

A loop is observed between km 650 and km 800 which perturbs the estimation of n in zone 3 and explains the low value of r. This is due to the fact that in this region τ and ε’ vary in a similar way, but with a shift of 20 km. A look at this region on the surface map shows that the profile crosses a flat zone, but the estimation of basal temperature in this region is about –5°C, and the depression in the bedrock does not give a strong radar echo. A small trough followed by a bump of 2 m amplitude is observed at km 740 in the topography. This kind of signature is also observed in the transition zone between strong friction and weak friction in the Lake Vostok area and was quantitatively attributed to the effect of longitudinal stress (Reference Remy, Shaeffer and LegresyRemy and others, 1999). Therefore, this loop could be a residual effect of longitudinal-stress gradients which have not been entirely smoothed.