3.3 Basal Boundary Condition

The Section describes the formulation of the basal boundary condition. An interface for the upper boundary condition (atmospheric BC) is easily defined by the surface temperature and mass balance. Similarly, the basal boundary consists of mechanical and thermal boundary conditions. The complications arise because the thermal and mechanical boundary conditions depend on each other. The interface of the basal boundary can be described with the following fields (see also Fig.3.8):

  1. basal traction: this field specifies a parameter which is used to allow basal sliding.
  2. basal heat flux: heat flux entering the ice sheet from below.
  3. basal water depth: the presence of basal melt water affects the basal ice temperature

Additionally, the ice sheet model calculates a melt/freeze rate based on the temperature gradient and basal water depth. This is handled by GLIDE.


pict

Figure 3.8: Basal boundary condition.


3.3.1 Mechanical Boundary Conditions

If the ice is not frozen to the bed, basal décollement may occur. This can be parameterised by a traction factor, tb. Within the ice sheet model tb is used to either calculate basal sliding velocities, ub, in the case of zeroth order physics, i.e.

ub = tbτb
(3.64)

where τb is the basal shear stress. Alternatively, tb can be used as part of the stress–balance calculations when the model is used with higher order physics. In simple models tb may be uniform or prescribed as a spatial variable. More complex models may wish to make tb dependant on other variables, e.g. basal melt rate. Typically tb will depend on the presence of basal water.

The second mechanical boundary condition, basal melting/freeze–on , is handled within the ice sheet model. The details are described in Section 3.3.2.

3.3.2 Thermal Boundary Conditions

The thermal boundary condition at the ice base is more complicated than the mechanical BC. The ice is heated from below by the geothermal heat flux. Heat is generated by friction with the bed. Furthermore, the ice temperature is constrained to be smaller or equal to the pressure melting point of ice. The thermal boundary is set to the basal heat flux if there is no water present. If there is water, the thermal boundary condition is set to the pressure melting temperature1 .

Basal Melting and Freezing

At the ice base, z = h, we can define outgoing and incoming heat fluxes, Ho and Hi:

Ho = -kice   ||
∂T-||
∂zz=h+ (3.65a)
and
Hi = -krock   |
∂T-||
∂z |z=h- + ub τb + {
  ρice ˙B∕L when B˙< 0
  0       otherwise (3.65b)
where kice and krock are the thermal conductivities of ice and rock, ub τb is the heat generated by friction with the bed and L is the latent heat of fusion of pure water. The basal melt/freeze–on rate, can then be calculated from the difference between the incoming and outgoing heat fluxes:
B˙ = Ho---Hi
      ρiceL
(3.66)

Freeze–on occurs if is negative, basal melting occurs if is positive.

Geothermal Heat Flux

The heat flux accross the basal boundary depends on past temperature variations since temperature perturbations penetrate the bed rock if the ice is frozen to the ground (Ritz1987). The heat equation for the bed rock layer is given by the diffusion equation

                             ( 2     2      2 )
∂T- = --krock--∇2T = --krock--  ∂-T-+ ∂-T-+ ∂-T-  ,
 ∂t   ρrockcrock      ρrockcrock  ∂x2    ∂y2   ∂z2
(3.67)

where krock is the thermal conductivity, ρrock the density and crock the specific heat capacity of the bed rock layer.

Initial conditions for the temperature field T are found by applying the geothermal heat flux, G to an arbitrary surface temperature T0:

T(x,y,z) = T + -G--z.
           0   krock
(3.68)

This ensures that initially the geothermal heat flux experienced by the ice sheet is equal to the regional heat flux. The basal boundary condition of the bedrock layer is kept constant, i.e.

T(x,y,H    ) = T + -G--H   .
       rock    0   krock rock
(3.69)

Lateral boundary conditions are given by

   ||        ||          ||        ||
∂T-||   =  ∂T||     = ∂T-||   = ∂T-||    = 0.
 ∂x x=0   ∂x x=Lx   ∂y y=0   ∂y  y=Ly
(3.70)

At the upper boundary, the heat flux of the rock layer has to be matched with the heat flux in the basal ice layer when the ice is frozen to the bed, i.e.

       |             |
krock ∂T||     = kice ∂T-||   .
     ∂z|z=-0      ∂z |z=+0
(3.71)

Otherwise the temperature of the top bedrock layer is set to the surface temperature (if the cell has been occupied by ice, but there is no ice present) or the basal ice temperature (if there is ice). Equation (3.71) is automatically fulfilled if we set the top bedrock temperature to the basal ice temperature everywhere and then calculate the geothermal heat flux to be used as boundary condition for Equation (3.38).

3.3.3 Numerical Solution

The horizontal grid is described in Section 3.1.1. The vertical grid is irregular like the vertical grid of the ice sheet model. However, it is not scaled. Also for now, I have ignored topography or isostatic adjustment, i.e. the bedrock layer is assumed to be flat and constant.

The horizontal second derivative in Equation (3.67) becomes using finite–differences

   |
∂2T||                 Ti+1,j,k --2Ti,j,k-+-Ti--1,j,k
∂x2|      = Txx,i,j,k =          Δx
    xi,yi,zi
(3.72)

and similarly for 2T∕∂y2. The vertical second derivative 2T∕∂z2 is similar to Equation (3.43):

∂2T ||                       2Ti,j,k-1                 2Ti,j,k
-∂z2||     = Tzz,i,j,k = (zk --zk-1)(zk+1---zk-1) - (zk+1---zk)(zk --zk-1)
     xi,yi,zi
                                                    + -------2Ti,j,k+1-------
                                                      (zk+1 - zk)(zk+1 - zk-1)
(3.73)

Using the Crank-Nicholson scheme, Equation (3.67) becomes

Tt+1- Tt       { Tt+1  + T t      T t+1   + Tt       Tt+1  + T t    }
-i,j,k---i,j,k-= D   -xx,i,j,k---xx,i,j,k+ -yy,i,j,k---yy,i,j,k + -zz,i,j,k---zz,i,j,k  ,
    Δt                  2                2                 2
(3.74)

with D = krock(ρrockcrock). Equation (3.74) is solved by gathering all Tt+1 terms on the LHS and all other terms on the RHS. The index (i,j,k) is linearised using ι = i + (j - 1)N + (k - 1)NM. The resulting matrix system is solved using the same bi–conjugate gradient solver as for the ice thickness evolution.

3.3.4 Basal Hydrology

It is clear from the discussion above that the presence of basal water plays a crucial role in specifying both the mechanical and thermal boundary conditions. However, the treatment of basal water can vary greatly. Basal water is, therefore, left as an unspecified interface. GLIDE does provide a simple local water balance model which can be run in the absence of more complex models.

3.3.5 Putting It All Together

The basal boundary consists of the individual components described in the previous sections. All components are tightly linked with each other. Figure 3.9 illustrates how the modules are linked and in what order they are resolved.


pict

Figure 3.9: Flow diagram illustrating how the various modules communicate with each other by exchanging data fields.


The order of executions is then:

  1. Find the basal heat flux by either solving the equation describing the thermal evolution of the lithosphere, Eq. (3.67), or by using the geothermal heat flux directly. The upper boundary condition of (3.67) is the same as the lower boundary condition of the thermal evolution of the ice sheet.
  2. The lower boundary condition for the thermal evolution of the ice sheet is either given by the basal heat flux from Step 1; or if melt water is present the basal temperature is set to the pressure melting point of ice.
  3. Calculate the temperature distribution within the ice sheet given the boundary condition found during Step 2 and the atmospheric BC.
  4. Calculate a melt/freeze–on rate using Equation (3.66) given the outgoing heat flux calculated during Step 3, friction with the bed (calculated during the previous Step 7) and the incoming heat flux from Step 1. Freezing only occurs when there is basal water.
  5. Track basal water. This is a user supplied module which can take any complexity. Inputs will typically be the melt/freeze–on rate determined during Step 4.
  6. Calculate the basal traction parameter. Again, this is a user supplied module which typically will involve the presence of basal melt water (calculated during Step 5).
  7. Solve the mechanical ice equations given basal traction parameter from Step 6.

Clearly, this scheme has the problem that heat is lost if the basal heat flux is such that more water could be frozen than is available. This might be avoided by iterating the process. On the other hand if time steps are fairly small this might no matter to much.