Flow law, flow-law parameters and enhancement factor

The isotropic flow law for ice (Reference GlenGlen, 1955, Reference Glen1958) may be expressed in the relations between strain rates and stresses (τij) for components of shear and compression as

(1)

(2)

where the subscripts xz and zz denote horizontal shear and vertical compression, respectively, denotes the appropriate component of the deviatoric stress tensor, B is the flow-law parameter; n the flow-law exponent is taken equal to 3 in the present work; and the octahedral shear stress (τ_o) is calculated from the deviatoric stresses, and may be expressed as

(3)

for the case of ice flow corresponding to a confined vertical compression stress combined with a horizontal shear stress τ_xz. In ice sheets the shear stress τ_xz is frequently calculated by assuming a linear increase with the depth all the way to the bed as

(4)

in terms of the ice density ρ, the acceleration due to gravity g, the mean surface slope α (taken here as 0.007, averaged along the flowline over 5 times the ice thickness at the DSS site) and the ice depth z.

For the present study we generalize Equations (1) and (2) by the introduction of a single, variable enhancement factor E into the parameter B, which is taken here as

(5)

where the factor E may involve the effects of ice fabric, ice crystal size and impurities on the deformation of the ice. A_o(T) is the temperature-dependent parameter representing the minimum octahedral creep rate per unit octahedral shear stress for ice with an isotropic crystal fabric, expressed in terms of , a constant, the activation energy Q (40 kJ mol^–1; Reference Weertman, Whalley, Jones and GoldWeertman, 1973), the gas constant R(8.314 J mol^–1 K^–1) and the absolute temperature T. from Equations (1) and (2) it is clear that the single enhancement factor E represents an enhancement of the flow law relating octahedral strain rate to octahedral shear stress, relative to the isotropic case, and it can be compared with a similar quantity extracted from laboratory experiments involving combined shear and compression loads by Li and others (1996). A model treating separate enhancement functions for shear and compression was outlined by Reference Wang and WarnerWang and Warner (1999), but as one of the main purposes of the present paper is to extract an estimate of enhancement from analysis of measured shear strain rates, we consider here a single depth-dependent function E. This enhancement factor tends to the corresponding enhancement for the shear component of flow when that component is dominant (see Reference Wang and WarnerWang and Warner (1999) for details). Note that the assumption of a common E factor is essential to the simple elimination of the unknown compressive deviatoric stress in favour of the shear stress and the strain rates.

FromEquations (13), the flow relation can be rewritten as

(6)

The flow-law parameter B can be obtained from the above equation once the vertical strain rate is known, since the shear strain rate is available from the borehole inclination measurements (Fig.1b) and the shear stress can be calculated from Equation (4).

Here we use an approach to calculate vertical strain rate under steady-state conditions, based on the assumption that ice is incompressible and that the depth profile of horizontal velocity varies only slowly along the flowline:

(7)

where the coordinate x is taken along the flow direction, U(z) is the depth-dependent horizontal velocity integrated from the measured shear strain rate (Fig.1b), a is the surface accumulation rate (0.678ma^–1 at DSS, taking a 50 year average), D is the total ice thickness (1177m, averaged over 5 times the ice thickness, along the flowline at the DSS site), and is a function of depth which gives the shape of the horizontal velocity as a function of the depth z:

(8)

To demonstrate the extent of the enhanced flow, the enhancement factor E^obs is then calculated from Equation (5) using B obtained from Equation (6) and chosen to make the mean value of the enhancement factor (E^obs) equal to 1 for the top 400 m depth of the borehole.

Relationship of enhancement factor to fabric

The flow-law parameter B, the temperature-dependent parameter A_o(T) and the enhancement factor E^obs are shown in Figure 2a and b for comparison with various measures of the ice crystal fabric in Figure 2c–e. the fabric measurements (Li, 1995) involve 185 fabrics, typically at 5– 10m spacing, with over 100 measured c axes per fabric.

Fig. 2 The profiles of (a) flow-law parameter B, and temperature parameter A_o(T); (b) enhancement factor E^obs (solid line), and smoothed fit (dashed line); (c) fabric parameter ϕ(1/4); (d) fabric parameter ϕ_m, the modal value from the histogram plot of crystal number vs c-axis co-latitude angle (in 5˚ bins); (e) normalized ratios from fabric analysis. the ratios R_{2 0–3 0} and R_0–15 are, respectively, the percentage of the co-latitudes of all measured c axes falling within 20–30˚ and 0–15˚ of the vertical, normalized by the corresponding values for a randomly distributed fabric. Four crystallographic zones (I–IV) are indicated, with boundaries at 200, 450 and 1000 m.

Figure 2b shows that significant enhancement starts at around 500 m depth, where the shear and compression strain rates become comparable. It then increases rapidly (although the variability in the measurements produces scatter in the estimates of enhancement), and from 700 to 1000 m the values cluster around 6, with perhaps a slightly increasing trend, which is consistent with the strong single-central-maximum crystal fabrics (Fig. 2c–e). There is a substantial reduction in enhancement commencing below the depth of the maximum shear strain rate. This corresponds to the region where the fabrics compatible with simple shear (i.e. a single-maximum pattern) appear to become distorted, as signified by an increase in the fabric parameter ϕ(1/4) discussed above. This may reflect the influence of a more complicated pattern of stresses, associated with longitudinal stresses in the portion of the ice sheet near the bedrock, associated with high relief in the surrounding regional bedrock topography.

The variable ϕ(1/4) does not differentiate strongly between the small circle girdle and single-central-maximum crystal fabrics, which are regarded as compatible fabrics for compression and shear, respectively. Accordingly, as an indication of the fabric development towards a girdle (or single-maximum) pattern, the ratios R_20–30 and R_0–15 describing the relative concentrations of c axes in the co-latitude bands 20–30˚ and 0–15˚, have been calculated (personal communication from Li Jun,1999). They were calculated as averages over 50 m sections and are shown in Figure 2e. If the fabric develops towards a girdle (or single maximum), the value of R_20–30 (or R_0–15) will increase along with the fabric development. Based on this analysis of crystal-orientation fabrics, the evolution of crystal c-axis orientations in the DSS core has been divided into four zones (I–IV) as indicated in Figure 2e.The figure shows that the ice deformation produces crystal fabrics that become dominated by compression in the upper part of the borehole. In the first 200 m (zone I) both R_20–30 and R_0–15 increase above the isotropic value, due to the disappearance of crystals with c axes near the horizontal plane from the initially isotropic fabrics, although the effect on R_0–15 is small compared to the values attained for single maxima deeper in the core. In the second region, R_0–15 generally decreases while R_20–30 increases until 450 m depth. This zone could accordingly be interpreted as tending to group c axes around the girdle pattern more than about a central maximum. In zone III, R_20–30 holds its value until about 600m and generally takes lower values below that depth. There is a renewed increase in R_0–15 at around 450m which suggests that shear starts to become significant there, and there is a larger increase in R_0–15 around 650 m, associated with a proportionate decrease in R_20–30, indicating the transition to the fabrics compatible with a shear-dominated stress situation. This is consistent with the appearance of the fabrics at 511 and 666m shown in Figure 1e. the shear compatible fabrics come to dominate and this continues until around 1000 m. Below this point, in zone IV, the fabrics compatible with simple shear (i.e. single-maximum pattern) begin to become distorted as previously mentioned, and R_0–15 is somewhat decreased. This is also reflected in a corresponding increase in ϕ(1/4) (Fig. 2c) in this portion of the core (except within the spike in shear strain rate associated with the LGM) back to values typical of the region around 600 m depth. This indicates an opening of the crystal fabric, suggesting that the shear stress is no longer completely dominant, and that, as previously discussed, other (presumably longitudinal and transverse) stresses, associated with proximity to the bedrock and surrounding high points in bedrock topography, are coming into play.

Finally near the bedrock (well below the LGM), multi-maxima crystal fabrics are seen, such as the one from 1153.9 m displayed in Figure 1e. These are associated with larger crystal sizes, indicative of extensive recrystallization and possible stress relaxation.

A smoothed fit to the observed enhancement factor, E^obs (dashed line in Fig. 2b) is plotted in Figure 3 as a function of the (smoothed) measured fabric verticality represented by ϕ(1/4). the figure shows that the enhancement factor increases with the strengthening of fabric verticality, but also indicates that for the same fabric verticality, the enhancement factor below the maximum shear layer value is considerably lower than that from above the maximum value. Given that these crystal fabrics are equally central (by this measure), and compatible with shear stress, one would expect proportionately similar responses to applied shear stress. This suggests that the linear increase of shear stress with depth, which was assumed in the calculation of E^obs, overestimates the shear stresses below the depth of the maximum shear strain rate. Note that if we assume that other horizontal transverse and longitudinal stresses are starting to modify the fabrics below 1000 m, then there should also be additional contributions to the octahedral stress τ_o which would suggest that an even larger decrease in horizontal shear stress would be required. It is worth emphasizing that there is no physical requirement that the maximum shear stress τ_xz must occur at the base of the ice sheet. This is demonstrated by finite-element calculations of the stress fields in flowline models, such as the work of Reference Budd and Rowden-RichBudd and Rowden-Rich (1985) discussed in the review by Reference Budd and JackaBudd and Jacka (1989). the necessity to improve this aspect of application of the standard shallow-ice approximation regarding the linear depth dependence of the driving shear stress was also previously mentioned by Reference Wang and WarnerWang and Warner (1999).

The model

The previous discussion suggests that the relationship between enhancement and fabric verticality seen in Figure 3 should be regarded as more reliable above the maximum shear layer. Accordingly, the enhancement factors in the lower part of the borehole (below the maximum enhancement) were re-evaluated, based on the relation between E^obs and ϕ(1/4) above the maximum enhancement (indicated in Figure 3 by the circle symbols), using the ϕ(1/4) data from Figure 2c. These new enhancement factors are plotted as a dashed line in Figure 4b. This procedure assumes that the fabric verticality remains a suitable guide to enhancement for horizontal shear, and that the variation of the shear stress with depth near the bed is the uncertain quantity. for the present study we assume that the major discrepancy comes from the simple estimate of the shear stress. There would generally also be some contribution made by longitudinal stresses to the octahedral stress which would suggest that an even greater reduction in the imputed shear stress would be required, but this would be taking us out of the assumed stress configurations behind Equation (6) in any case.

Fig. 3 The smoothed fit to the enhancement factor (E^obs) plotted as a function of the smoothed fabric parameter (ϕ(1/4)). Points lying above the depth of the maximum shear strain rate are denoted (∘) and those from below the maximum shear strain rate are denoted (*).To avoid problems with smoothing, the data within the narrow spike of higher shear strain rates (see Fig. 1b), associated with the ice of the LGM, centred at 1113 m depth are not shown.

Fig. 4 Results from the model. Profiles of (a) shear strain rates from observation (thin line) and from the model (smooth line), and vertical compressive strain rate (dashed line); (b) enhancement factors E^obs and E(λ_S) (solid lines), and the modified E^obs in the basal zone (dashed line); (c) horizontal velocity (U) and vertical velocity (V_x); (d) modified shear stress (τ_xz), compression stress deviator and octahedral shear stress (τ_o); (e) modelled shear strain (ε_xz), compressive strain (ε_zz) and the age (dotted line) in 103years.

New ``corrected’’ or inferred shear stresses below the depth of maximum enhancement were accordingly determined by inserting the enhancements (E^obs) (obtained from the matching of enhancement and fabric verticality across the maximum shear layer as discussed above) and the observed shear strain rates (Fig. 1b) into Equation (6), leading to a marked shear stress reduction between the depth of the maximum shear strain rate and the depth of the LGM (solid line in Fig. 4d). This gives an indication of the magnitude of local stress variations that may be relevant close to the bedrock.

Observations of the vertical profiles of crystal-fabric properties are not generally available throughout ice sheets, and modelling of the evolution of such profiles in polar ice sheets is still in its infancy. Accordingly, we have used these modified shear stresses, together with a model enhancement function drawn from the analysis of laboratory ice-deformation experiments involving combined compression and shear stresses (Li and others, 1996; Warner and others, 1999; W. F. Budd and others, unpublished information) to model the expected shear strain-rate profile for DSS. This type of enhancement model (Reference Wang and WarnerWang and Warner, 1999) ignores the delays involved in the development of compatible crystal fabrics under changing stresses (the strain history), and relates enhancement to the nature of the applied stresses. Here we consider the simple relationship for shear flow enhancement

(9)

where the constant E_S represents the enhancement appropriate to a dominant shear stress, and the shear fraction (λ_S) (Li and others, 1996), which ranges between 0 for dominant compression and 1 for dominant shear, may be calculated (for the geometry associated with Equations (1) and (2)) as a function of the compressive and shear strain rates

(10)

The value of E_S was taken as 6 from consideration of a large number of laboratory experiments on ice with strong single-maximum fabrics, covering a range of different stresses and temperatures, which gave enhancements ranging from 4.5 to 12. the shear flow enhancement-factor Equation (9) is based on analyses of laboratory combined stress experiments, and while it represents well the marked modification of the isotropic Glen law it is not clear how far into the compression-dominated zone it should be extrapolated. Further-more, for the present application to Equation (6) we need to remember that our enhancement E in Equation (5) represents the enhancement of the octahedral stress–strain-rate relationship and should approach a factor describing the compressive component of deformation as λ_S goes to zero. As our main focus here is the shear strain-rate profile, for the present modelling study we have simply taken

(11)

to set the enhancement factor E to unity in the uppermost part of the ice sheet. This also gives reasonable values for E in the zone where enhanced compressive flow is expected. Using the ``corrected’’ shear stress, this prescription of enhancement works well in the present case, except that from about 45 mbelow the maximum shear layer the reduction in strain rate (due to the reduced shear stress) halts and the model begins to predict progressively larger strain rates than are observed (except in the LGM spike which clearly involves other influences on flow rate), leading to a prediction close to the bedrock of nearly 10 times the observed value. Previously Reference Wang and WarnerWang and Warner (1999) used a model of this general type successfully in the zone above the maximum shear layer. Here we see more clearly that the formulation of shear enhancement based on Equations (9) and (10) becomes inapplicable in these lower regions. Using these equations with the observed shear and inferred vertical compression strain rates gives a shear fraction (λ_S) close to unity (and consequently maximal enhancement) throughout the lowest 250m of the borehole. As discussed above, the observed enhancements and corresponding crystallographic information indicate that this is inappropriate, and laboratory experiments on large-crystal, multi-maximum basal ice also suggest that a shear enhancement factor near unity is appropriate near to the bedrock.

Accordingly, for the present study the model enhancement function is further modified, to reduce linearly to unity between the depth of maximum shear strain rate and a location about 50 m above the bottom of the borehole. This is intended to approximate the reduction of shear enhancement implied by the more open crystal fabrics (larger values of ϕ(1/4)) and the multi-maximum crystal fabrics seen within 50 mof the bedrock.

Given an initial estimate of vertical strain rate of a/D based on the observed accumulation and the total ice thickness, Equation (6) (combined with Equation (5), and the model enhancement factor described above involving Equations (8–11)) is iterated until the shear strain-rate profile, and hence the profiles of horizontal velocity, (z), vertical strain rate, λ_S and E(λ_S) converge.

Results from the model

The results from the model are summarized in Figure 4. the enhanced shear strain rates are shown in Figure 4a by the smooth solid line, and the calculated enhancement, E(λ_S), is shown in Figure 4b by a smooth line.

The shear strain-rate profile computed from the model using laboratory-determined rheology shows a good agreement with the observations from the borehole. the reduction of the shear strain rates near the bed is controlled by the counteracting combination of increased temperature, reduced shear stress and the imposed reduction in shear flow enhancement.

The base temperature at DSS is –6.9˚C, which indicates that no significant basal sliding can be expected over the bedrock, so the integration of the shear strain rate over the ice column yields the horizontal velocity at the surface. the computed surface velocity (Fig. 4c) of 2.4 ma^–1 is consistent with the value of 2.4ma^–1 obtained from the integration of the measured shear strain rate (Fig.1b) over the borehole depth.

The ratio of the depth-averaged velocity to the surface velocity is calculated from the model to be 0.73, which is close to the value of 0.74 determined from the borehole velocity measurements by Reference Morgan, van Ommen, Elcheikh and JunMorgan and others (1998), and their value of 0.76 calculated by assuming the local ice sheet is in mass balance (based on data from long-term surface elevation observations) and matching the surface accumulation rate (averaged over the last ~50 years) with the outflow around DSS determined by surface strain-rate data. the ratio of the depth-averaged velocity to the surface velocity represents the relationship between the column-integrated ice flux and measurements of the ice-sheet thickness and surface velocity. to determine the state of ice-sheet balance from the ice accumulation distributions and ice-sheet surface topography, this relationship must be considered (Reference Budd and WarnerBudd and Warner, 1996; Reference Wang and R.C.Wang and Warner, 1998). In areas where ice sliding can be neglected, this ratio strongly reflects the rheological properties of the ice and depends on the shear strain-rate profile through the ice column, which is influenced by the profiles of stress, temperature and ice crystal fabrics (Reference Wang and R.C.Wang and Warner, 1998) and may be strongly affected by the large-scale bedrock roughness. Comparing these values with the value of 0.8 for a cubic flow law with isothermal isotropic ice (see Reference PatersonPaterson, 1994, p. 251– 252), at DSS the value of 0.73 is significantly smaller. This probably involves the reduced shear stresses near the bed and the distortion of the fabrics, and is also complicated by the added questions of flow within a few ice thicknesses of the Dome Summit as well as the influence of quite rough bedrock topography upstream of the drilling site (see Reference Morgan, Wookey, Jun, van Ommen, Skinner and FitzpatrickMorgan and others, 1997, fig. 2).

The compressive and shear strains, and ice age are also calculated in the model and are shown in Figure 4e. Due to the low velocities and strain rates in the basal layers, the ice at the bottom of the DSS borehole has been computed to have undergone large strains (about 105% shear strain) and to be approaching 10⁶ years old. This model, however, assumes steady state, and for such long times the non-steady-state variations associated with the ice-age cycles would need to be taken into account.