Notations and main hypothesis
The problem to be solved consists of the gravity-driven flow of isothermal, incompressible and non-linear viscous ice. The geometry is restricted to a two-dimensional plane flow perpendicular to the y direction. Ice flows along the x direction, and the z axis is the vertical upward-pointing axis. Notation is detailed in Figure 1. The ice sheet flows over a rigid bedrock, z = b(x), assuming a non-linear friction law for the grounded part, and extends further as an ice shelf over the ocean. In the case of the ice shelf, the ice slides perfectly over the sea. The last grounded point defines the grounding line and is denoted x
G hereafter. The left-hand side of the domain is assumed to be symmetric and the shelf ends at the right-hand side of the domain.
Fig. 1. Notation for the problem to be solved.
The constitutive law for the ice behavior is a Norton–Hoff type law (Glen’s flow law in glaciology):
where S is the deviatoric stress tensor, Dij
)/2 are the components of the strain-rate tensor and u is the velocity. The effective viscosity, η, can be expressed as
where the strain-rate invariant, γ
e, is defined as
In the following applications, n is set to 3 and ice is assumed isothermal, such that the fluidity parameter, B, remains constant (see values in Table 1).
Table 1. Values of the parameters used in this study, which correspond to steps 5 and 6 of the Marine Ice-Sheet Intercomparison Project (MISMIP) benchmark, but are expressed differently, as here the fluidity parameter is B = 2A. However, numerically the constitutive relations are rigorously the same.
State equations and boundary conditions
The ice flow is computed by solving the Stokes problem, expressed by the mass-conservation equation in the case of incompressibility
and the momentum conservation equations
where σ = S − pI is the Cauchy stress with p the isotropic pressure; ρ
i is the ice density and g the gravity vector.
In the present problem, the lateral boundaries of the domain are a dome and a calving front. The dome is assumed to be axisymmetric for the flow problem, which implies that ux(0,z) = 0. The calving front is an artificial cutting of the shelf, which can be seen as the point where icebergs are calved. The exact position of the calving front is not relevant to our problem because it does not influence the upstream flow. This surface undergoes the sea-ice pressure, p
w(z,t), which evolves vertically as
w is the sea-water density and l
w the sea level.
In addition to the lateral boundaries, the ice body is bounded by two free surfaces, namely the stress-free upper surface, z = z
s(x, t), and the bottom surface, z = z
b(x, t), at the interface between the bed or sea and the ice. The evolution of the two free surfaces is determined by solving a local transport equation. Note also that the length of the sea/ice interface, starting from the grounding line, x
G(t), is not known in advance and is part of the solution.
The upper surface is a stress-free surface, which implies that n · (σ · n)|
s = p
atm ≈ 0, where n is the unit normal vector of the surface pointing outward and the subscript ‘s’ denotes the value taken at the ice/air interface. The equation governing the upper stress-free surface evolution reads:
s) denotes the upper surface velocity in the horizontal (i = x) and vertical (i = z) directions and a(x, t) is the accumulation/ablation function given as a vertical flux at the upper surface. In what follows, the accumulation is supposed constant both in space and in time (see Table 1).
The bottom sea stress-free surface obeys the following equation:
Table 1. Values of the parameters used in this study, which correspond to steps 5 and 6 of the Marine Ice-Sheet Intercomparison Project (MISMIP) benchmark, but are expressed differently, as here the fluidity parameter is B = 2A. However, numerically the constitutive relations are rigorously the same
where accretion of sea water by refreezing or melt of bottom ice is neglected. Note that this equation is still valid for the points on the bottom surface which are in contact with the bedrock. Assuming a rigid, impenetrable bedrock, z = b(x), the following topological conditions must be fulfilled by z
s and z
As a consequence, the unilateral link between the ice and the bedrock can be treated as a contact problem: the ice cannot penetrate the bedrock but is allowed to move away from it (Lestringant, 1994). Resolutions of the contact problem have been inspired by previous studies on subglacial cavities, namely the works of Schoof (2005) and Gagliardini and others (2007). At a given point, x, ice is assumed to be in ‘true’ contact with the bedrock (and corresponding boundary conditions are applied; see below) if the ice touches the bed and the stress exerted by the ice is larger than the sea-water pressure. Conversely, the ice is assumed to be in contact with the sea if the bottom surface is above the bed or if the ice touches the bed and the sea-water pressure remains larger than the normal stress, σ
nn. In other words:
(1) the ice/bedrock boundary condition applies if
(2) the ice/sea boundary condition applies if
where the subscript |
b denotes the value taken at the bottom surface. For the flow problem, two different boundary conditions have to be applied, depending on whether the ice is in contact with the bedrock or with the sea. For the ice/bedrock boundary condition (i.e. condition (1) above), a non-linear friction law is applied:
b is the basal shear stress, u
b = u · t is the sliding velocity at the base and t is the tangent vector to the surface z
b. The parameters C and m entering the friction law are given in Table 1. For the ice/sea boundary condition (i.e. condition (2) above), the ice is sliding perfectly over the sea, i.e. t · (σ · n)b = 0, and the normal stress is equal to the buoyancy sea-water pressure (Equation (6)).
The equations presented above have been implemented in the finite-element code Elmer (http://www.csc.fi/elmer/). For gravity-driven ice flow we solve the set of the Stokes equations neglecting inertia terms. More details on the numerics are given by Durand and others (in press). Note that a similar approach, i.e. full Stokes finite-element modeling of marine ice sheets, was initiated by Lestringant (1994) and used more recently by Nowicki and Wingham (2008).