The ice sheet model includes simple approximations for calculating isostatic adjustment. These approximations depend on how the lithosphere and the mantle are treated. For each subsystem there are two models. The lithosphere can be described as a
while the mantle is treated as a
At each isostasy time-step, the load of ice and water is calculated, as an equivalent mantle-depth (L). If the basal elevation is above sea-level, then the load is simply due to the ice:
where H is the ice thickness, with ρi and ρm being the densities of the ice and mantle respectively. In the case where the bedrock is below sea-level, the load is calculated is that due to a change in sea-level rise and/or the presence of non-floating ice. When the ice is floating (ρiH < ρo(z0 - h)), the load is only due to sea-level changes
whereas when the ice is grounded, it displaces the water, and adds an additional load:
here, ρo is the density of sea water, z0 is the change in sea-level relative to a reference level and h is the bedrock elevation relative to the same reference level. The value of h will be negative for submerged bedrock, hence the plus sign in (3.77).
This is model is selected by setting lithosphere = 1 in the configuration file. By simulatuing the deformation of the lithosphere, the deformation seen by the aesthenosphere beneath is calculated. In the absence of this model, the deformation is that due to Archimedes’ Principle, as though the load were floating on the aesthenosphere.
The elastic lithosphere model is based on work by Lambeck and Nakiboglu (1980), and its implementation is fully described in Hagdorn (2003). The lithosphere model only affects the geometry of the deformation — the timescale for isostatic adjustment is controlled by the aesthenosphere model.
The load due to a single (rectangular) grid point is approximated as being applied to a disc of the same area. The deformation due to a disc of ice of radius A and thickness H is given by these expressions. For r < A:
and for r ≥ A:
where Ber(x), Bei(x), Ker(x) and Kei(x) are Kelvin functions of zero order, Lr = (D∕ρmg))1∕4 is the radius of relative stiffness, and D is the flexural rigidity. The constants Ci and Di are given by
Here, the prime indicates the first spatial derivative of the Kelvin functions.
If a fluid mantle is selected, it adjusts instantly to changes in lithospheric loading. However, a relaxing mantle is also available.