Skip to main content Accessibility help


  • Access


      • Send article to Kindle

        To send this article to your Kindle, first ensure is added to your Approved Personal Document E-mail List under your Personal Document Settings on the Manage Your Content and Devices page of your Amazon account. Then enter the ‘name’ part of your Kindle email address below. Find out more about sending to your Kindle. Find out more about sending to your Kindle.

        Note you can select to send to either the or variations. ‘’ emails are free but can only be sent to your device when it is connected to wi-fi. ‘’ emails can be delivered even when you are not connected to wi-fi, but note that service fees apply.

        Find out more about the Kindle Personal Document Service.

        A computer scheme for rapid calculations of balance-flux distributions
        Available formats

        Send article to Dropbox

        To send this article to your Dropbox account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you use this feature, you will be asked to authorise Cambridge Core to connect with your <service> account. Find out more about sending content to Dropbox.

        A computer scheme for rapid calculations of balance-flux distributions
        Available formats

        Send article to Google Drive

        To send this article to your Google Drive account, please select one or more formats and confirm that you agree to abide by our usage policies. If this is the first time you use this feature, you will be asked to authorise Cambridge Core to connect with your <service> account. Find out more about sending content to Google Drive.

        A computer scheme for rapid calculations of balance-flux distributions
        Available formats
Export citation


A simple computer scheme developed by Budd and Smith (1985) and modified by D. Jenssen has been further developed to provide a rapid computation of steady-state balance fluxes over arbitrary ice masses, given the surface elevations and net accumulation distribution. The scheme provides a powerful diagnostic tool to examine the flux and state of balance over whole ice masses or limited regions to interpret field observations for dynamics or the state of balance.

In many cases the uncertainty in the state of balance may be much less than the uncertainty in the deformation and sliding properties of the ice and so the flux and velocities derived from balance could provide a useful guide for the dynamics where direct observations are sparse.

The scheme assumes that, on a horizontal scale of many ice thicknesses, the ice-flow direction is approximately down the steepest surface slope. The continuity equation is used to compute steady-state implied downslope fluxes at each grid point from integrations of the net accumulation over the area from the summits to the edges. The algorithm ensures the exact integral balance of the surface net flux over the area with flow through boundaries.

Applications are demonstrated for the whole of Antarctica and for regional areas. Comparisons are made between fluxes computed from observed ice thicknesses and velocities and those computed from balance. The observed ice thicknesses can also be used to compute surface velocities from assumed column-to-surface velocity ratios. The combined fluxes from observations and balance can be used to compute rates of change of elevation with time.

1. Background

The concepts of “balance flux” and “balance velocity” go back a long way (e.g. Budd and others, 1971) and refer to the conditions for steady state in a glacier or ice sheet such that the outward-flow distribution exactly matches the net accumulation, and consequently there is no ice-thickness change with time at any point. To be more precise, if over a horizontal (x, y) domain, A(x, y) is the net surface accumulation (or ablation rate) in ice-equivalent volume per unit area per unit time and M(x, y) is the net basal melting (or freezing rate), then at any point, with ice thickness Z (again in ice-equivalent units) and depth-averaged horizontal velocity the steady-state balance condition implies that


In this context, the divergence operator acts over the x, y domain and may be represented as


where (i, j) are the unit vectors in the directions x and y.

An appropriate version of Gauss’s divergence theorem may be used to express Equation (1) in integral form. For any closed area S in the (x, y) domain where dl is the elemental line length along the boundary C and n is the outward-directed unit vector, normal to the boundary, the total gain (or loss) through the upper and lower surfaces is balanced by the net outward flow through the boundary, i.e.


The quantity Φ = Z represents a vector horizontal-flux density for unit width (in km2 a−1) in the direction of flow. The quantity


is a scalar flux integral representing the volume flux through the section C i of the boundary. The term Φ(x, y) will be referred to as the balance-flux distribution. The corresponding balance-velocity distribution is given by


For many applications, there is interest in the variation of the velocity and flux along a flowline. For this case, if A(ξ) is the net accumulation (in absence of significant basal melt) at distance ξ along a flowline and η(ξ) is the width between a pair of neighbouring flowlines, one each side of the central flowline ξ, then the flux integral, from the start of the flowline to ξ, is given by


This technique was used by Budd and others (1971) to compute balance-flux and velocity distributions over the Antarctic ice sheet from a number of selected flowlines and an assumed accumulation distribution. The technique used by Budd and others (1971) was largely manually based with the time-consuming digitization of individual flowlines. Since then a variety of applications have been developed, as discussed below, using different schemes for automatic processing of digital gridded data sets to derive implied distributions of steady-state balance fluxes and velocities from prescribed distributions of net accumulation, surface elevation and ice thickness. So far, however, the details of an efficient numerical scheme for the procedure have not been presented. The purpose of this paper is to outline the basis and application of the current numerical scheme being used, which is rapid to compute, has an accurate integral-closure constraint and can be readily applied to a wide range of arbitrary ice-mass configurations specified by gridded data sets. The present computer code is freely available for general use from the authors.

2. Uses of Balance-Flux Calculations

Uncertainties in the flow properties of ice sheets stem from many factors, including uncertainties in the temperature distribution and in the nature and intensity of the ice-crystal anisotropy. Budd and Jacka (1989) showed enhancement of flow rates for anisotropic ice in tertiary flow under compression of about a factor of 3 relative to the minimum under compression of isotropic ice. For simple shear, the corresponding enhancement is about an order of magnitude. Strong anisotropy associated with horizontal shear appears to develop in the lowest third and fourth quartiles of depth but this may be greatly disturbed by flow over and around rough basal topography. Even greater uncertainties exist for the rate of basal sliding for ice.

As a consequence, the state of balance may be known in many cases to a greater degree of accuracy than the ice-flow or sliding properties. Also, since the balance-flux rates scale linearly with uniform fractional changes in the net surface balance over the domain, the balance-flux calculations can be used to assimilate a wide range of observations to construct complete distributions of derived features throughout the ice mass. This principle was used by Budd and others (1971) to compute a comprehensive range of derived features for the Antarctic ice sheet from the data available at that time. Similar results for Greenland have been presented by Budd and others (1982) and Radok and others (1982).

For glaciers, Budd and Allison (1975) showed that balance concepts could also provide the basis for the derivation of an extensive range of characteristics. Smith and Budd (1981) demonstrated that dynamic modelling of the response of glaciers to climate change could also benefit from rapid computations of balance fluxes and steady-state profiles. All of these studies used flowline-type techniques. An automatic numerical scheme for gridded data was devised and used for assessing the Antarctic mass budget by Budd and Smith (1985). Comparisons with observed velocities showed that the majority of points were indicative of a positive balance for the present accumulation-rate data. This same automatic gridded technique was used as an inverse technique to derive surface balances for the northern ice cap of Mars by Budd and others (1986). The use of Balance velocities also provided a useful technique, in the absence of observed velocities, to derive large-scale empirical relations for ice sliding of ice streams as given by McInnes and Budd (1984).

The power of the technique has been further illustrated by the application to basal-meltwater flux by Budd and Jenssen (1987). The basal-whaler flow rates are so high, compared with the ice flow, that the balance technique can be used for basal water in conjunction with time-dependent ice-flow modelling. An updated version of derived physical characteristics for the Antarctic on a 20 km grid was carried out by Radok and others (1986, 1987). Finally, the balance fluxes have been used to help tune the dynamic flow computations along with the limited observed surface velocities, e.g. Budd and Jenssen (1989) and Mavrakis (1993).

The balance-flux calculations now provide a powerful diagnostic tool for assimilating observations by adjustment of the fluxes at the observed points and then via continuity over the whole domain by best-fit statistical techniques. Since the basic procedure underlying the scheme has not been published, this paper presents the essential outline of the scheme and the basis of the numerical code.

3. Basic assumptions underlying the Balance-Flux Calculations

For the balance flux in its simplest form, if the distribution of surface elevation and net surface mass balance are given then, for negligible basal-melt flux, the horizontal-flux distribution can be derived from Equation (1) With certain assumptions:

  • i. The direction of ice flow is assumed to be orthogonal to the elevation contours. In practice, this assumption has been found to apply reasonably well for horizontal scales which are large (say 20 times) relative to the ice thickness, cf. e.g. Budd (1968), Budd and Young (1979) and Hamley and others (1985). For smaller scales, irregular bedrock can cause deviations in direction comparable with the corresponding fluctuations in slope, and basal, longitudinal and transverse stresses. Nevertheless, this assumption is also often used for large-scale dynamics modelling and it can be further tested by the comparisons with observations.

  • ii. The flux rector in the direction of flaw can be resolved into rectangular (x, y) components which can be suitably represented by finite-difference schemes to the resolution of interest. This entails that the grid-point density chosen is sufficient to represent reasonably the surface topography in enough detail to capture the essential flow characteristics. On the other hand, the small fluctuations in topography mentioned in (i) may need to be smoothed to represent properly the large-scale flow.

4. Data Sets Required

In the simplest form of application of balance fluxes, it is assumed that a surface-elevation data set is available and can be digitized into a regular rectangular (x, y) grid with grid spacing Δx, Δy, The spacing needs to be chosen to give sufficient detail, depending on the smoothness of the surface. No interior surface hollows (or exactly flat regions) are allowed (in the absence of net surface ablation). A corresponding digital data set for surface net accumulation is required over the same grid. These two data sets are sufficient to compute the balance-flux distribution implied by the condition of flow orthogonal to the elevation contours.

If data also exist for ice thickness or base elevation, then the balance velocities can also be computed from Equation (5). However, velocities can be calculated only at those locations where ice-thickness data exist, or are specified.

For more complex applications, if only some parts of these data sets are available, then the scheme could be used with inversion techniques or statistical assimilation methods to derive balanced, compatible, complete data sets.

5. Outline of Basic Algorithm

Only a brief outline of the basic principles of the scheme is given here with details available with the code from the authors. Many different numerical schemes can be used to represent the basic continuity Equations (1) and (4) with sufficient accuracy for the integral fluxes. In general, it has been found that for smooth data sets the results tend to converge for increasing resolution of the different schemes. The computer code has options for a number of these schemes. For simplicity, we consider here first a centred difference scheme for any point x i,j where the nearest neighbours are x i−1,j , x i+1,j x i,j−1 and x i,j+1.

The surface elevations at these corresponding points are denoted E i,j , E i−1,j , E i+1,j , E i,j−1 and E i,j+1.

The surface-slope vector components at x i,j are taken as


The magnitude of the slope is


and the direction of the slope is at an angle θ to the x axis given by


This also represents the direction of flow and the orthogonal to the elevation contour through xi,j .

We take A i,j as the net accumulation at x i,j , Φ i,j = ( Z) i,j as the flux-density vector, and (n δl) i,j as any line-element length vector intercepting the flow.

The resolved components of flux in the x, y grid directions may then be taken as


With the flow directions determined from the topography, a direct discretization of Equation (1) at each grid point would lead to a system of linear equations for the magnitudes of the vector-flux density at each grid point. For centred differences, this takes the form


This system could be solved by sparse matrix methods but this might be problematic as fluxes differ spatially by several orders of magnitude.

The key to the present scheme for determining the balance fluxes is to change the focus from the magnitude of the vector-flux density at each grid point Φ i,j to the total scalar flux of ice out of the (i, j) grid cell and rewrite the discretized continuity Equation (3) as


where the total scalar influx of ice is simply the sum of out-flowing contributions from neighbouring upstream grid points, and to apportion the outflow from the (i, j) cell into downstream neighbouring cells using the slope directions θ i,j .

For example:


where these X and Y tractions of the scalar outflux, , , contribute to the total influxes of the downslope neighbouring points in the x- and y-grid directions, respectively.

By sorting and treating the grid points in order of descending elevation, so that the total influx is always known from the contributions of previously treated upslope points, the balance-flux integral out of each (i, j) cell can be straighforwardly calculated. This scheme allows the treatment of an arbitrary number of summits, ridges and basins in the elevation data set. For points of equal elevation, the relative order results directly from the order in the data array.

Using the scalar outflow and the slope directions, the vector flux-density field can be produced. We find that for some of the more Involved scalar outflux-apportioning schemes which we have used, simply dividing the outflow magnitude by the length of the grid-cell edge provides a good measure of the magnitude of the balance-flux density field.

This partitioning of the flux integral into rectangular fractions on the x, y grid is an approximation dependent on the finite-difference scheme used but it allows an exact accounting of the flux over the domain. The computation of the balance fluxes has also been found not to be very sensitive to the numerical form of the partitioning, provided the data fields are smooth and the resolution is adequate. As an example, an alternative normalized form of the flux-integral components and has been tried and has been found to give approximately equivalent results, i.e. with flux-integral fractions apportioned as


Similarly, more complex forms of the numerical grid or finite-difference schemes have been found to be not as important as the conditions of smooth fields and high resolution. For example, the staggered-grid formulation has the x, y components of flow towards and away from each point governed by the mid-interval slopes:


with corresponding relations for the mid-point flux integrals across the corresponding cell faces, , , , and . The continuity equation for balance then allows the downstream flux to be computed from the upstream fluxes and the net accumulation rate through Equation (12) as before. For a staggered grid, the outgoing scalar-flux integral is partitioned between the cell faces in proportion to these positive (downslope) mid-interval slopes, e.g.


where N is the sum of all the positive (downward) mid-interval slopes from the point (i, j).

The staggered-grid scheme above can be modified to include the diagonal downslope neighbouring grid points by extending the role of the mid-interval slopes to include all elevation différences to the neighbouring points and scaling the diagonal directions by . A corresponding modified form of Equation (16) is used to partition the outflows. This scheme has been found to relied curved topography more accurately when using coarser grids.

The routine provides that the net accumulation is distributed as flux progressively downslope to the outward boundary. For any closed domain, the integral over the region of accumulation (S) necessarily equates to the outward flux summed around the boundary (C) to the precision of round-off error, i.e.


This corresponds to the discrete form of Equation (3) neglecting basal melt. A more detailed analysis of the approximations and validity of the various numerical schemes will be presented elsewhere. Here, the emphasis is placed on an illustration of the application. It suffices to say that the discrepancies between the numerical schemes result primarily in small differences between neighbouring cells but very little difference in the fluxes across the outflow sections compared with the differences resulting from the alternative data sets as shown below.

6. Application to the Whole Antarctic Continent

Previous applications of this type of analysis have been carried out for the whole of Antarctica by Budd and Smith (1985) and with subdivisions for the major regional drainage basins by Radok and others (1986). The method is useful for assimilating new data and evaluating the state of balance and the dynamics. As an example, Figure 1 shows the distribution of balance fluxes computed from the accumulation distribution of Budd and Smith (1985) and their corresponding surface-elevation data sel digitized on a 20 km grid from Drewry’s (1983) folio map. This type of output was used by Budd and Jenssen (1989) to assess the state of balance by comparison with observations of fluxes obtained from products of surface velocities (V s) and thickness (Z) at points over the continent. There are now many more points for which such an evaluation can be made but this may be more appropriate for separate regional studies as shown below. The results can also be used to assess the ice-dynamics formulations as shown by Budd and Jenssen (1989). The balance and observations intercomparisons indicated that most of the continent is within about ±50% of balance and large areas within ±20%. Because of the uncertainties in basal temperatures, the ice-crystal anisotropy and the flow and sliding properties of ice, it is difficult to compute the velocities or fluxes (a priori) to this degree of precision, The balance fluxes, normalized to the observations, can be used to evaluate the dynamics formulations in relation to the ice-sheet stresses and other features. For this paper, we wish to emphasize that with new improved data sets for surface elevation and accumulation rate, the balance fluxes provide an even more accurate and powerful means of assessing both the state of balance and the ice dynamics. This is done by a re-analysis of the Wilkes Land observed fluxes previously examined by Hamley and others (1985) but using the new surface elevations derived from satellite-radar altimetry and additional surface-traverse observations from Medhurst (unpublished) and Hazelton (unpublished).

Fig. 1. The pattern of ice-balance flux density over the Antarctic continent is shown as calculated from the balance-flux scheme using the SPRI-based elevations on a 20 km grid and the present-accumulation data set of Budd and Smith (1985). The graduated shading is shown using a log scale in units of km2 a−1. Over the ice shelves, the downslope flow is less applicable.

7. Wilkes land comparison of balance fluxes with observations

Hamley and others (1985) used the SPRI-based elevation data set to assess the state of balance of the Wilkes Land region by comparison of balance fluxes with observed values of surface velocity (V s) and ice thickness for V s Z. They concluded that the region was close to balance with the accumulation data set used, and using an implied column-average velocity to surface-velocity ratio of 0.89. Nevertheless, their traverse profile results showed a number of large discrepancies between computed and observed fluxes. Since then, additional traverse data have become available for the line between Casey and Mirny, and relatively high-resolution accurate surface elevations with practically complete coverage over the basin have become available from satellite-radar altimetry (Bamber, 1994). Bamber supplied these data analysed on the same 20 km grid as that used by Budd and Smith (1985) for the SLRI-based data. The satellite data provide a far more comprehensive basis for flow directions and computed fluxes. To illustrate this, a re-analysis of the interior and coastal (near 2000 m elevation) Wilkes Land traverse lines was carried out using the same accumulation data set and the same balance-flux algorithm but with the two different elevation data sets. The traverse lines are shown in Figure 2 with the points marked where the observed velocities used to compute fluxes were available. Comparisons between the two different balance-flux calculations and the observations along the traverse lines are shown in figure 3a and b.

Fig. 2. An enlargement of part of Figure 1 over the Wilkes Land region is shown with the traverse routes inland of Casey: (1) near the 2000 m contour from about 93° E to 128° E; and (2) from Mirny to Pioneerskaya and Dome C. The locations where both surface velocities and ice thickness were measured are shown by the circles. The balance fluxes are shown in units of 10−2 km2 a−1.

Fig. 3. The balance fluxes (Z) computed along the traverse routes as derived from the SPRI-based elevations (dotted lines) and the ERS-1 elevations (full lines) compared with the surface-velocity and ice-thickness products (VsZ) from the observations (circles). The satellite elevations give rise to much closer spatial matching with the pattern from the observations. (a) Dome C to Mirny; (b) near the 2000 m elevation line from east to west. Units are in km2 a−1.

A summary of the mean flux calculations and correlations with the observations are given in Table 1. It is clear that the fluxes computed from the satellite-derived elevations represent much better the observed character of the flow. The strong flux in the major ice streams penetrates a long way inland and is noticeable to elevations over 3000 m. In both cases, the observed fluxes for the coastal line appear close to balance, depending on the ratio of column average to surface velocity. However, the spatial pattern of the satellite-elevation fluxes is much closer to the observed with much higher correlation coefficients between calculated and observed fluxes.

Table I. Comparison of fluxes (VsZ) from observations (Obs) with those computed from balance () using satellite (ERS) elevations or the digitized (SPRI) map, shown for the inland line (In) and coastal line (Co)

For the more inland region, a positive balance of about 35% of the surface-accumulation rate is implied using an average to surface-velocity ratio of 0.87. This value is chosen as a broad-scale interior average from the range given by Budd and others (1971) of 0.85–0.92, going from the interior to the coast in this region. This would imply a rate of rise of the region’s interior to this line around 1–2 cm a−1. More detailed studies of the balance and dynamics will be treated in a separate paper. Here, it suffices to say that the improved surface-elevation data set will lead directly to improved estimates for the velocity and flux distributions over Antarctica.

8. Future Prospects and Conclusions

So far, the satellite-elevation data set only reaches to about 82° S. Nevertheless, many of the large accumulation basins of East Antarctica such as Wilkes Land and the Lambert and Shirase Glacier basins are completely covered. Eventually, complete satellite-elevation data may be obtained but in the interim much can still be done with a composite data set consisting of the satellite data to 82° S and the interpolated ground-based elevations further south. Additional surface velocities measured at strategic locations can then go a long way towards assessing the current distribution of the state of balance over the ice sheet or the expected rate of change of surface elevation (∂E/∂t). This then puts the emphasis on obtaining improved data sets for surface net accumulation and for the column-average to surface-velocity ratio (λ = /V s ), which is now contributing significantly to the remaining uncertainty. Eventually, from satellite altimetry (radar and laser), it should become possible to derive the rate of elevation change with time. Then, it will also become possible to derive the actual fluxes, from the existing routine, by distributing the quantity (A– ∂E/∂t) downslope over the surface.


The authors wish to thank Dr J. L. Bamber of the Mullard Space Laboratory, University College London, for providing his data set for surface elevation from the ERS-1 satellite altimetry. The Wilkes Land traverse data were supplied from reports of the ANARE glaciologists: T. Hamley, T. Medhurst and W. Hazelton, and from N. Young, who supervised the traverse programmes, and to whom the authors wish to acknowledge their thanks.


Bamber, J. L. 1994. A digital elevation model of the Antarctic ice sheet derived from ERS-1 altimeter data and comparison with terrestrial measurements. Ann. Glaciol., 20, 4854.
Budd, W. 1968. The longitudinal velocity profile of large ice masses. International Association of Scientific Hydrology Publication 79 (General Assembly of Bern 1967Snow and Ice), 5877.
Budd, W. F. and Allison, I. F.. 1975. An empirical scheme for estimating the dynamics of unmeasured glaciers. International Association of Hydrological Sciences Publication 104 (Symposium at Moscow 1971 — Snow and Ice in Mountainous Areas), 246256.
Budd, W. F. and Jacka, T. H.. 1989. A review of ice theology for ice sheet modelling. Cold Reg. Sci. Technol., 16(2), 107144.
Budd, W. F. and Jenssen, D.. 1987. Numerical modelling of the large-scale basal water flux under the West Antarctic ice sheet. In Van der Veen, C. J. and Oerlemans, J., eds. Dynamics of the West Antarctic ice sheet. Dordrecht, etc., Kluwer Academic Publishers, 293–320. (Glaciology and Quaternary Geology 4.)
Budd, W. F. and Jenssen, D.. 1989. The dynamics of the Antarctic ice sheet. Ann. Glaciol., 12, 1622.
Budd, W. F. and Smith, I. N.. 1985. The state of balance of the Antarctic ice sheet. Glaciers, ice sheets, and sea level: effect of a CO2-induced climatic change, Report of a Workshop held in Seattle, Washington, September 13–15, 1984. Washington, DC, US Department of Energy. Office of Energy Research, 172–177. (Attachment 9, Report DOE/ER/6O235–1.)
Budd, W. F. and Young, N. W.. 1979. Results from the I. A. G. P. flow-line study inland of Casey, Wilkes Land, Antarctica. J. Glaciol., 24(90), 89101.
Budd, W. F., Jenssen, D. and Radok, U.. 1971. Derived physical characteristics of the Antarctic ice sheet. ANARE Interim Rep. 124.
Budd, W. F., Jacka, T. H., Jenssen, D., Radok, U. and Young, N. W.. 1982. Derived physical characteristics of the Greenland ice sheet. Mark I. Parkville, Victoria, University of Melbourne. Meteorology Department. (Publication 23.)
Budd, W. F., Jenssen, D., Leach, J. H. I., Smith, I. N. and Radok, U.. 1986. The north polar ice cap of Mars as a steady-state system. Polar-forschung, 56(1–2), 4363.
Drewry, D. J., ed. 1983. Antarctica: glaciological and geophysical folia. Cambridge, Scott Polar Research Institute.
Hamley, T. C., Smith, I. N. and Young, N. W., 1985. Mass-balance and ice-flow-law parameters for East Antarctica. J. Glaciol., 31(109), 334339.
Hazelton, W. Unpublished. Data report. Vol.1: 1986 Casey Glaciology Programme. Hobart, Antarctic Division. Glaciology Section, (Internal report, 1987.)
McInnes, B. J., and Budd, W. F.. 1984. A cross-sectional model for West Antarctica. Ann. Glaciol., 5, 9599.
Mavrakis, E. 1993. Time dependent, three-dimensional modelling of the dynamics and thermodynamics of large ice masses. (M. Se. thesis, University of Melbourne.)
Medhurst, T. Unpublished. Data report: 1983 Casey Glaciology Programme. Hobart, Antarctic Division. Glaciology Section. Internal report, 1984.)
Radok, U., Barry, R. G., Jenssen, D., Keen, R. A., Kiladis, G. N. and McInnes, B.. 1982. Climatic and physical characteristics of the Greenland ice sheet. Part I and II. Boulder, CO, University of Colorado. Cooperative Institute for Research in Environmental Sciences.
Radok, U., Brown, T. J., Jenssen, D., Smith, I. N. and Budd, W. F.. 1986. On the singing potential of polar ice streams. Port IV, Antarctic ice accumulation basins and their main discharge regions. Boulder, CO, University of Colorado. Cooperative Institute for Research in Environmental Sciences: Melbourne, University of Melbourne. Meteorology Department. (US Department of Energy Report DE/ER/60197–5.)
Salamatin, A. N., Lipenkov, V. Y. and Blinov, K. V.. 1994. Vostnk (Antarctica) climate record time- scale deduced from the analysis of a borehole-temperature profile. Ann. Glaciol., 20, 207214.
Radok, U., Jenssen, D. and McInnes, B.. 1987. On the surging potential of polar ice streams. Antarctic surges — a clear and present danger? Boulder. CO, University of Colorado. Cooperative Institute for Research in Environmental Sciences; Melbourne, University of Melbourne. Meteorology Department. (US Department of Energy Report DOF./ER/60197-H1.)
Smith, I. N. and Budd, W. F.. 1981. The derivation of past climate changes from observed changes of glaciers. International Association of Hydrological Sciences Publication 131 (Symposium at Canberra 1979 — Sea Level, Ice and Climatic Change), 3152.