ISOSTATIC BED ADJUSTMENT

The Antarctic ice sheet is in sorne places over 4500 m thick, with a rnean value of around 2000 rn. A considerable pressure is therefore exerted on the underlying bedrock, which bends downwards to restore the equilibrium of forces. The rnodel includes a time-dependent calculation of bedrock depression because the relaxation time for isostatic adjustment is of cornparable rnagnitude to the reaction time of the ice sheet to changing environmental conditions. This is particularly important for changes in the position of the grounding line, which are the result of a subtle interplay between variations in ice thickness and water depth, so that delayed isostatic displacernents rnay provide an additional cause for grounding-line rnigration.

In our approach, bed adjustrnents are calculated with a two-layer earth deforrnation model that consists of a viscous fluid asthenosphere enclosed by a uniform, thin and elastic lithospheric plate, which is in floating equilibrium with the underlying substratum of the asthenosphere. The elastic properties of the rigid lithosphere determine the deflection, which gives the ultirnate shape of the bedrock depression. The rate of bedrock adjustrnent, on the other hand, is governed by the rate of outflow in the asthenosphere and thus by its viscous properties.

Appropriate values for asthenospheric thickness and viscosity are not very well known, but data from observations of glacial rebound in Scandinavia and North America provide an estimate for an asthenospheric diffusion constant (Walcott, 1973). The time-dependent response of the underlying substratum may therefore be modelled by a diffusion equation for bedrock elevation h (e.g. Oerlernans and Van der Veen, 1984, chap.7), yielding:

where Da = 0.5 X 108 m*y-1 is the asthenospheric diffusivity, hn the undisturbed bed topography in the absence of loading and ^W' the deflection, taken positive downwards. Using equation (4.49) implies that the characteristic time scale for bedrock sinking depends On the size of the load.

With a typical length scale (L) of 1000 km this leads to a relaxation time for bedrock adjustment of the order T = L2/Da = 20000 years.

The rigidity or flexural stiffness of the crust causes the effect of a point load to be spread out beyond the boundary of the load itself. The lithosphere will not only bend downward under the load, but deflect upward in more remote areas, thereby creating a weak forebulge. These rheological properties of the lithosphere give rise to deviations from local isostatic equilibrium, that become especially rnarked near the ice sheet edge, where grounding line migration takes place (order of magnitude: 30 m). Thus, the rigidity of the lithosphere rnust be taken into account. The thin plate equation for flexural deformation may be written in its two-dimensional form as (Brotchie and Silvester,-! 969; Turcotte and Schubert.1982):

where D, is lithospheric flexural rigidity and pm mantle rock density [3300 kgm- 31. The right-hand side gives the applied load q minus the upward buoyancy force arising when the lithospere bends downwards into the asthenosphere;

the left-hand side represents the bending resistance opposing lithospheric flexure. D. is proportional to the elastic thickness of the plate raised to the third power. If the rigidity of the plate is low (implying that its thickness is small), then buoyancy forces dominate and compensation for the load is local. If Dr is large then the bending resistance term becomes dominant and the compensation is more regional in character. D,. is taken as 1025 Nm for the

entire Antarctic bedrock (cf. Drewry, 1983, sheet 6). This corresponds to a lithospheric thickness of 11 5 km.

Bending of a rigid plate is a linear process, so that the total isostatic displacement in one point on a numerical grid can be obtained by superimposition of contributions from neighbouring gridpoints. The solution for (4.50) reads (Brotchie and Silvester, 1969):

q Ax Ay L,

W . . (X) =

1J kei (X)

where wfij (X) is the deflection in an elernent ij caused by a central point load q at normalized distance X frorn the centre of loading, r the distance to the centre of loading [m], Lr the so-called 'radius of relative stiffness', which depends on the elastic properties of the lithosphere [= 132573 rn], Ax = Ay = 40 km the gridpoint spacing and kei (X) are Kelvin functions of Zero order. Tables of this function can be found e.g. in Abrarnowitz and Stegun (1 970).

Dividing both sides of (4.51) by the load q means that the right hand side now contains constants that can be conveniently tabulated in a matrix as a function of distance X. The deforrnation coefficients obtained define a bell-shaped function and are shown graphically in fig. 4.5. They represent the contribution of each gridpoint to the total deflection. In order to calculate lateral deflections caused by a central point load, a Square with sides equal to 800 km was used, comprising a 21 X 21 = 441 elernent rnatrix; all the deforrnation coefficients together add up to 1. Downward deflections can be Seen to extend over a distance of four tirnes the radius of relative stiffness (approx. 530 km) from the centre of loading at i = j = 11. Beyond this distance, a weak forebulge results.

The actual total downward displacernent in any point is then found by calculating the respective deflections for all 440 surrounding gridpoints, multiplied by the appropriate deformation coefficient. The resulting isostatic depression is then the sum of the contributions from all 441 elernents.

flg. 4.5: Pictorial representation of the crustal deflection caused by a central point load in (1 1 ,Tl). The central depression extends to a lateral distance of about 530 km. One side unit equals 40 km. In the model, deforrnation coefficients are rnultiplied by the central point load to obtain lateral isostatic displacements.

Glaciostatic depression bv present ice load lml

fig. 4.6: Isostatic uplift that would result frorn complete removal of the subsequent isostatic rebound.

current ice load and

Fig. 4.6 shows the crustal downwarping caused by the present Antarctic ice load as calculated with this rnodel. It also represents the uplift that would result if the current ice load were to be removed and full isostatic recovery were to take place. The effect turns out to be a very srnooth irnprint of the local ice loading distribution. Isostatic adjustrnents caused by changing sea level andlor ocean water rnerely displacing submarine ice are automatically incorporated in the calculations, because (4.50) takes into account the loading caused by both ice and water.

In a recent Paper, Stern and ten Brink (1989) have criticized the rnodel of lithospheric flexure used here, which assurnes a continuous, constant rigidity plate. They propose that both East and West Antarctica are separate plates with free edges at their cornmon boundary (the Transantarctic Mountains).

This results in a sornewhat different Pattern of induced flexure when subjected to an ice cap load. In particular, lithospheric displacernents in the Transantarctic Mountains are in their model smaller (they are up to 500 m in our approach). Stern and ten Brink agree on the Dr value for East Antarctica, but estirnate the flexural rigidity for the ROSS Ernbayment to be rnore than two orders of magnitudes less at 4 X 1022 Nm (corresponding to a lithospheric elastic thickness of only 19 km). Using this alternative model would rnake isostatic adjustrnents more local in West Antarctica.

In retrospect, it appears that also the diffusion approach for athenospheric viscous flow (eq. 4.49) may not be the most appropriate way to rnodel the transient upper mantle response (G. Deblonde, personal cornmunication). In particular, the irnplied response time scale rnay be too long and the resulting transient forebulge may not be entirely realistic. This is inferred frorn comparison of rnodelled postglacial rebound with observed relative sea-level data (e.g. Peltier, 1985, 1988). It is not clear whether these irnprovernents would influence the fundamental outcome of the ice sheet rnodel. however.