Treatment of the ice-shelf backpressure and buttressing in two horizontal dimensions

Abstract

The ice discharge from the grounded parts of marine ice sheets into the ocean is modulated by their floating extensions – ice shelves. The ice-shelf impact on the grounded ice is typically described as ‘backpressure’ or ‘buttressing’. Theoretical analyses of their effects have been restricted to one horizontal dimension. This study revisits the concepts of ‘backpressure’ introduced by Thomas (1977) and ‘buttressing’ numbers and ratios introduced by Gudmundsson (2013) and extends their theoretical analysis to two horizontal dimensions. Using the integral form of the momentum-balance formulation suitable for fast-flowing ice streams and ice shelves, our analysis provides a natural definition for the total backpressure force exerted by an ice shelf to the grounded ice upstream of its grounding line. The results of numerical analyses suggest that ice shelves whose second principal stress component is compressional over larger areas may provide more buttressing compared to ice shelves with smaller areas of compressional stresses or to ice shelves with both principal stresses being tensile.

1. Introduction

The dynamics of marine ice sheets’ grounding lines – locations where the grounded ice loses its contact with the underlying bedrock and starts to float forming ice shelves – control ice discharge into surrounding oceans and consequently, contributions of marine ice sheets to sea level. In turn, the grounding line dynamics depend on the geometric and dynamic conditions of the grounded portions of marine ice sheets and ice shelves. The geometric conditions are the presence or absence of lateral confinement and the variability of the bed topography under the grounded ice. The dynamic conditions are the stress regimes of the ice flow on the grounded and floating parts; these regimes are determined by the dominant components of the ice-flow momentum balance (Schoof, 2007a, 2007b; Haseloff and Sergienko, 2018; Sergienko and Wingham, 2019, 2022). The current conceptual understanding of the conditions at the grounding lines has been developed using one-dimensional flow-line models. Considering a laterally unconfined marine ice sheet resting on a flat bed (flat before the ice sheet was developed on top of it), Weertman (1974) suggested that such a configuration cannot attain a stable steady state if the bed slopes towards the interior of the ice sheet. This result, known as the ‘marine ice-sheet instability’ hypothesis, has been widely used to interpret the observed behavior of present-day ice sheets (e.g., Shepherd and others, 2018) and simulated behavior under future climate conditions (e.g., Cornford and others, 2015; Seroussi and others, 2017). The existing theoretical analyses of the grounding line behavior (e.g., Weertman, 1974; Schoof, 2007a, 2011, 2012; Tsai and others, 2015; Sergienko and Wingham, 2019, 2022, 2024; Sergienko, 2022b) have considered one horizontal dimension and, in the case of laterally confined configurations, parameterized the effects of lateral shear in the momentum balance of ice flow (e.g., Pegler, 2016; Schoof and others, 2017; Haseloff and Sergienko, 2018, 2022; Sergienko, 2022a).

Investigations of the effects of transverse variability on the conditions at the grounding line have been done using numerical models applied to idealized configurations (e.g., Goldberg, 2009, 2012a, 2012b; Gudmundsson and others, 2012; Gudmundsson, 2013) or realistic configurations (e.g., Seroussi and others, 2017; Reese and others, 2018; Sun and others, 2020). A few laboratory experiments and theoretical analyses built on the experimental results have been performed for laterally unconfined ice shelves (Pegler and Worster, 2012, 2013). Their results suggested that ice viscous deformation in the direction transverse to the main flow gives rise to hoop stresses that could potentially affect the stress regime at the grounding line. However, estimates for the unconfined parts of the Antarctic ice shelves and ice tongues suggest that the effects of hoop stresses are very small (Wearing and others, 2020).

Although about five decades ago, Thomas (1973, 1979) argued that the shear of the side walls of the ice shelves or the presence of ice rises can affect the stability of the grounding line, it is the results of fairly recent numerical studies (Gudmundsson and others, 2012; Gudmundsson, 2013) that demonstrated that Weertman's marine ice-sheet instability hypothesis does not hold if the marine ice sheet is laterally confined. Later theoretical studies in which the effects of lateral confinement have been parameterized confirmed this result by analyzing expressions of the ice flux through the grounding line (e.g., Schoof and others, 2017; Haseloff and Sergienko, 2018) and by linear stability analysis (Haseloff and Sergienko, 2022; Sergienko and Haseloff, 2023).

The process by which ice shelves impede ice discharge from the grounded part of marine ice sheets is termed ‘buttressing’. Its quantitative measure, also known as ‘backpressure’, was introduced by Thomas (1973), who defined it as the difference between the maximum depth-integrated driving stress experienced by unconfined ice shelves and the depth-integrated deviatoric stress. While concepts of the absence of buttressing and, as a consequence, zero back pressure, are straightforward for a laterally unconfined marine ice sheet, which is also uniform in the direction transverse to the ice flow and whose grounding line is a straight line, they are ambiguous in the presence of transverse variability and curved grounding lines. This is because its effect is non-local and arises as a result of interactions of the ice-shelf flow with obstacles – either the ice-shelf lateral boundaries or ice rises located far away from the grounding line – and is transmitted via the ice shelf deformation back to the grounding line.

To quantify backpressure, MacAyeal (1987) has introduced concepts of ‘form drag’ and ‘dynamic drag’, partitioning the total force at a given point of the grounding line into the ice deformation (the dynamic drag) and the hydrostatic (the form drag) components. Gudmundsson (2013) took a different approach to define the local effects of buttressing (i.e., at a given point of the grounding line). He introduced the normal and tangential buttressing numbers K _N,T and buttressing ratios Θ _N,T (K _N = 1 − Θ _N, K _T = Θ _T) that represent the ratio of the normal and transverse components of the force at the grounding line to the hydrostatic pressure. MacAyeal (1987), Gudmundsson (2013), and many subsequent studies aiming to quantify the effects of buttressing (e.g., Reese and others, 2018), used the results of numerical model simulations to compute the stress components at the grounding lines and evaluate the respective metrics.

Defined in terms of the components of stress at the grounding line, expressions for these metrics do not include any information about an ice shelf whose buttressing they are meant to quantify. The ice-shelf effects on these metrics are implicit: via its impacts on stress at the grounding line. This study aims to establish how the ice-shelf stress distribution and its boundary conditions affect buttressing, and makes their effects explicit in considerations of buttressing and backpressure. It revisits the concepts of backpressure introduced by Thomas (1977) and buttressing numbers introduced by Gudmundsson (2013) in the context of marine ice sheets that experience variability in the direction transverse to the dominant ice-flow direction. Using the integral form of the momentum balance typically used for fast flowing ice streams and ice shelves (the Shallow Shelf Approximation) (MacAyeal, 1989), we derive the expressions of the total forces provided by an ice shelf at the grounding line. These expressions can be naturally used as a definition of the total backpressure force provided by the ice shelf to its grounding line. It can be used as an integral metric characterizing the force balance of an ice shelf as a whole. Analysis of the point-wise backpressure shows that for two dimensional (i.e., non-uniform in the transverse direction) unconfined ice shelves it is non-zero, even though the total backpressure is zero. Such ice shelves do not provide buttressing to their grounding lines and the upstream ice flow as a whole, but the point-wise backpressure force may be non-zero. The results of numerical simulations show that spatial distributions of submarine melting have strong effects on the ice-shelf stress distribution, and as a result, on the grounding line and its buttressing. Analysis of the principal stress components suggests that ice shelves with larger spatial extent of the compressive second principal stress may provide more buttressing than those with less area experiencing compression or no compression at all. These results suggests that the second principal strain-rate component, which is proportional to the second principal stress, can be used as a proxy for the ice shelf buttressing and its evolution.

The manuscript is organized as follows: The model is described in section 2. The next section, section 3, provides a description of the total backpressure force. Derivations of the point-wise buttressing metrics are described in section 4. The results of numerical simulations are presented in section 5. Readers less interested in the mathematical aspects of the analysis can proceed to sections 5–7, which provide a physical interpretation of the results and their discussion.

2. Model description

Despite the complex geometry of the grounding lines of Antarctic ice shelves, the ice flow on the ice shelves exhibits a predominant direction – towards the calving front. As shown in Fig. 1, the streamlines on the ice shelves are nearly straight, even though the patterns on the grounded parts are very complex.

Figure 1. Ice flow of the Pine Island and Thwaites ice shelves. Gray lines represent streamlines and colors indicate ice speed (m yr⁻¹) (Rignot and Scheuchl, 2017). Red lines indicate the grounding lines.

Motivated by these observations and to simplify our analysis¹ , we choose a Cartesian coordinate system aligning the x-axis with the direction of dominant ice flow and the y-axis transverse to that direction (Fig. 2).

Figure 2. Model geometry: plane view L ^d–ice divide location, x ^g–grounding line location; L ^c–calving front location. Ice flows from left to right.

We use a vertically integrated momentum balance of ice flow typically used to describe ice-stream and ice-shelf flows (MacAyeal, 1989). In two horizontal dimensions, the momentum balance of ice flow on the grounded part G is given by:

(1a)$$\left[2\nu H( 2u_x + v_y) \right]_x + \left[\nu H\left(u_y + v_x\right)\right]_y -\tau^{bx} = \rho g HS_x,\; \quad \{ x,\; y\} \in G$$

(1b)$$\left[\nu H\left(u_y + v_x\right)\right]_x + \left[2\nu H( u_x + 2v_y) \right]_y-\tau^{by} = \rho g HS_y,\; \quad \{ x,\; y\} \in G$$

Here, subscripts indicate partial derivatives; H represents ice thickness, u and v are vertically averaged horizontal components of ice velocity v = {u,  v}, g is the acceleration due to gravity; S = B + H is the surface elevation, and B is the bed elevation; ν denotes a vertically averaged ice viscosity:

(2)$$\nu = {\bar B\over 2\left[u_x^2 + v_y^2 + \frac{1}{4}\left(u_y + v_x\right)^2 + u_xv_y\right]^{{(n-1)}/{2n}}},\; $$

with $\bar B$ as the constant ice-stiffness parameter ($\bar B =$ 1.68×10⁸ Pa s^1/3).

The basal shear τ ^b = {τ ^bx,  τ ^by} follows a power-law sliding law:

(3)$${\bf \tau}^b = -C_b\vert {\bf v}\vert ^{m-1}{\bf v},\; $$

where C _b = 7.6×10⁶ Pa m^−1/3 s^1/3 is the sliding coefficient, and m = 1/n = 1/3 is the sliding exponent.

The momentum balance of the floating ice shelf F is as follows:

(4a)$$\left[2\nu H( 2u_x + v_y) \right]_x + \left[\nu H\left(u_y + v_x\right)\right]_y = \rho g^\prime HH_x,\; \quad \{ x,\; y\} \in F$$

(4b)$$\left[\nu H\left(u_y + v_x\right)\right]_x + \left[2\nu H( u_x + 2v_y) \right]_y = \rho g^\prime HH_y,\; \quad \{ x,\; y\} \in F.$$

Here, g ^′ represents the reduced gravity:

(5)$$g^\prime = \delta g,\; $$

where δ denotes the buoyancy parameter:

(6)$$\delta = {\rho_w-\rho\over \rho_w},\; $$

and ρ and ρ _w are the densities of ice and sea water, respectively.

We define unit normal vectors to the grounding line $\vec {n}^g$ and to the calving front $\vec {n}^c$ as follows:

(7a)$$\vec{n}^g = \{ n^{g}_{x},\; n^{g}_{y}\} = {1\over \sqrt{1 + ( x^{g}_{y}) ^2}}\{ 1,\; -x^{g}_{y}\} $$

(7b)$$\vec{n}^c = \{ n^{c}_{x},\; n^{c}_{y}\} = {1\over \sqrt{1 + ( x^{c}_{y}) ^2}}\{ 1,\; -x^{c}_{y}\} ,\; $$

where, {x ^g,  y ^g} ∈ L ^g is the grounding line, {x ^c,  y ^c} ∈ L ^C is the calving front (Fig. 2), and $\displaystyle {x^{g, c}_{y} = {{\rm d}x^{g, c}}/{\, {\rm d}y}}$.

Boundary conditions at the upstream boundary, L ^D, can take different forms. We choose this boundary to represent an ice divide and use the following conditions:

(8a)$$u = v = 0,\; \quad \{ x,\; y\} \in L^D,\; $$

(8b)$$( H + B) _x = 0,\; \quad \{ x,\; y\} \in L^D.$$

At the calving front, the deviatoric stress in the ice shelf balances the hydrostatic pressure deficit caused by ice buoyancy:

(9a)$$\left[2\nu H( 2u_x + v_y) \right]n^{c}_{x} + \left[\nu H\left(u_y + v_x\right)\right]n^{c}_{y} = {\rho g^\prime\over 2}H^2n^{c}_{x},\; \quad \{ x,\; y\} \in L^c$$

(9b)$$\left[\nu H\left(u_y + v_x\right)\right]n^{c}_{x} + \left[2\nu H( u_x + 2v_y) \right]n^{c}_{y} = {\rho g^\prime\over 2}H^2n^{c}_{y},\; \quad \{ x,\; y\} \in L^c.$$

The boundary conditions at the lateral boundaries L ^N and L ^S will be specified in the following section.

The mass balance of the ice stream is

(10)$$H_t + \vec\nabla\cdot\vec Q = \dot a,\; \quad \{ x,\; y\} \in G,\; $$

where $\vec Q = \{ uH,\; \, vH\}$ is the ice flux, $\vec \nabla = \{ \partial _x,\; \, \partial _y\}$ is the divergence operator, and $\dot a$ is the net accumulation rate (positive for accumulation), mostly dominated by surface accumulation/ablation. The mass balance of the ice shelf is:

(11)$$H_t + \vec\nabla\cdot\vec Q = \dot a-\dot m,\; \quad \{ x,\; y\} \in F,\; $$

where, $\dot a- \dot m$ is the net ablation/accumulation rate. This could be dominated by ablation/accumulation at the ice-shelf surface or melting/refreezing at the ice-shelf base, or the two could balance each other.

At the grounding line, the ice thickness, velocity components, and normal and tangential stress components are continuous, and the flotation condition is

(12)$$H( x,\; y) = -{B( x,\; y) \over 1-\delta},\; \quad \{ x,\; y\} \in L^G$$

3. Total backpressure force

Backpressure is caused by the interactions of ice-shelf flow with obstacles – either lateral confinements or ice rises. To develop a conceptual understanding of backpressure in two horizontal dimensions, we focus on the effects of lateral boundaries, leaving considerations of the effects of ice rises for future studies.

To determine the total (or integral) backpressure that an ice shelf provides to the grounding line, we consider the integral form of the momentum balance (4). A vector/tensor form of the ice-shelf momentum balance (4) is given by:

(13)$$\nabla\cdot{\bf T} = \rho g^\prime H\vec\nabla H,\; $$

where $\nabla$ is the divergence operator in a given set of coordinates, and T is

(14)$${\bf T} = \left[\matrix{2\nu H( 2u_x + v_y) & \nu H\left(u_y + v_x\right)\cr \nu H\left(u_y + v_x\right)& 2\nu H( u_x + 2v_y) }\right],\; $$

which could be viewed as a ‘vertically integrated’ deviatoric stress-tensor (although it is not exactly that, as it takes into account the incompressibility equation, $\vec \nabla \cdot \vec v = 0$, and relies on the assumption that the vertical shear is negligible).

The right-hand side of the momentum balance (13) can be written as the gradient of the scalar field $\displaystyle {{H^2}/{2}}$

(15)$$\nabla\cdot{\bf T} = \rho g^\prime \vec\nabla{H^2\over 2},\; $$

Integrating both sides of (15) over the surface area of the ice shelf F, using the Gauss divergence theorem and the same justifications of its application to the Stokes or Navier-Stokes equations (e.g., Lamb, 1932), one obtains the integral form of the ice-shelf momentum balance (15):

(16)$$\oint_L\left({\bf T}\cdot\vec n-\rho g^\prime{H^2\over 2}\vec n\right)\, {\rm d}l = 0,\; $$

where $\vec n$ is an outward-pointing unit vector, ${\bf T}\cdot \vec n = T_{ij}n_j$ represents forces at the ice-shelf boundaries L, which include the calving front L ^C, the grounding line L ^G, and the lateral boundaries L ^N and L ^S (Fig. 2). Equation (16) represents a vertically integrated force balance of an ice shelf; it is satisfied for individual components of the force balance, such as normal and tangential components.

The boundary condition at the calving front (9), written in vector form, is

(17)$${\bf T}\vec n^c = \rho g^\prime{H^2\over 2}\vec n^c,\; \quad \{ x,\; y\} \in L^{C}.$$

Consequently, $\int _{L^{C}}( {\bf T}\cdot \vec n^c-\rho g^\prime {H^2}/{2}\vec n^c) \, {\rm d}l = 0$, and (16) become:

(18)$$\int_{L^{G}}\left({\bf T}\cdot\vec n^g-\rho g^\prime{H^2\over 2}\vec n^g\right)\, {\rm d}l = \int_{L^N}\left({\bf T}\cdot\vec n-\rho g^\prime{H^2\over 2}\vec n\right)\, {\rm d}l + \int_{L^{S}}\left({\bf T}\cdot\vec n-\rho g^\prime{H^2\over 2}\vec n\right)\, {\rm d}l.$$

Note the change in sign due to the direction of the normal vector $\vec n_g$ at the grounding line (it points in the same direction as the normal vector at the calving front). The quantity on the left-hand side is the backpressure integrated along the length of the grounding line, i.e., the force exerted by the ice shelf on the grounding line in addition to the force associated with the pressure deficit between ice and sea water. We denote this force as $\vec F^{BP}$:

(19)$$\vec F^{BP} = \int_{L^{G}}\left({\bf T}\cdot\vec n-\rho g^\prime{H^2\over 2}\vec n\right)\, {\rm d}l.$$

$\vec F^{BP}$ has two components corresponding to the coordinate system – either x − and y −components or normal and tangential components. The x − and y −components corresponding to the chosen geometry (Fig. 2) are

(20)$$\eqalign{& F^{BP}_x = \int_{L^{G}}\left[2\nu H( 2u_x + v_y) n_x + \nu H\left(u_y + v_x\right)n_y-\rho g^\prime{H^2\over 2} n_x\right]\, {\rm d}l\cr & F^{BP}_y = \int_{L^{G}}\left[\nu H\left(u_y + v_x\right)n_x + 2\nu H( u_x + 2v_y) n_y-\rho g^\prime{H^2\over 2} n_y\right]\, {\rm d}l.}$$

As is apparent from Eqn (18), $\vec F^{BP}$ depends on the conditions at the ice-shelf lateral boundaries and the length of these boundaries.

3.1 Laterally unconfined ice shelf

In this case, the boundary conditions at the lateral boundaries are the same as at the calving front (17), and (18)–(19) becomes

(21)$$\vec F^{BP} = \int_{L^{G}}\left({\bf T}\cdot\vec n-\rho g^\prime{H^2\over 2}\vec n\right)\, {\rm d}l = 0.$$

This indicates that if the ice-shelf lateral boundaries experience only the imbalance between hydrostatic pressures in ice and water due to the buoyancy of ice, then the ice shelf does not provide buttressing to the grounding line in the integral sense. However, this does not necessarily imply that $\displaystyle {{\bf T}\cdot \vec n = \rho g^\prime {H^2}/{2}\vec n}$ at each point along the grounding line, and locally the internal deformation may differ from the imbalance of the hydrostatic pressures in ice and water (this is discussed below in section 4.1). It is the total backpressure force of the unconfined ice shelf that is zero.

3.2 No flow at the lateral boundaries

If ice shelves are laterally confined and ice flow at their lateral boundaries is very slow (compared to the trunk of an ice shelf), it can be approximated by no-slip (or no-flow) conditions

(22)$$u = v = 0,\; \quad \{ x,\; y\} \in L^{N, S}.$$

For the chosen geometry (Fig. 2), this implies

(23)$$u_x = v_x = 0,\; $$

and

(24)$${\bf T}\cdot\vec n = \left[\matrix{\nu H u_y \cr 4\nu H v_y }\right]n_y,\; \quad \{ x,\; y\} \in L^{N, S}.$$

Physically, Eqn (24) represents friction between the ice shelf and its lateral boundaries. Consequently, the total backpressure force at the grounding line (18)–(19) is determined by the friction and the length of the lateral boundaries.

3.3 Shear at the lateral boundaries

As suggested by Thomas (1977), the friction at the ice-shelf lateral boundaries could be approximated by, for instance, a plastic yield stress of ice. If the magnitudes of lateral shear are known from direct observations or laboratory experiments, then instead of boundary conditions on velocities, boundary conditions on the stress could be prescribed:

(25)$$\vec t\cdot{\bf T}\vec n = -\vec\tau^w \quad \{ x,\; y\} \in L^{N, S},\; $$

where $\vec \tau ^w$ is a vertically integrated lateral shear, and $\vec t = \{ -n_y,\; \, n_x\}$ is a tangent unit vector such that $\vec t\cdot \vec n = 0$.

In this case, the total backpressure force is described by (18), where the components of ${\bf T}\vec n$ on the lateral boundaries L ^N and L ^S are determined by (25).

4. Local backpressure and buttressing numbers

The previous section has considered the total backpressure provided by the ice shelf to the grounding line and has demonstrated that in the absence of ice rises it can determined from the lateral boundary conditions only. This section focuses on the local buttressing effects.

As their measure, Gudmundsson (2013) introduced the buttressing numbers

(26a)$$K_N = 1-{N\over N_0},\; $$

(26b)$$K_T = {T\over N_0},\; $$

where

(27a)$$N = \vec n_g^\prime\cdot{\bf T}\vec n_g,\; $$

(27b)$$T = \vec t_g^\prime\cdot{\bf T}\vec n_g,\; $$

(27c)$$N_0 = {\rho g^\prime\over 2}H^2$$

$\vec n_g^\prime$ and $\vec t_g^\prime$ indicate transpose vectors. (Here, the definitions of N, T and N ₀ differ from those by Gudmundsson (2013) by a factor of H.) Using Eqn (14) the above expressions provide definitions of K _N,T and Θ _N,T (Θ _N = 1 − K _N, Θ _T = K _T) in terms of the stresses at the grounding line. However, as written, these definitions are oblivious to the ice shelves and depend on their properties and processes implicitly, i.e., via their effects on the grounding-line stresses. Since the physical meaning of buttressing numbers is to represent the effects of the ice shelves, it is expedient to express them via characteristics of the ice shelves. In order to do so, we largely follow an approach used in a one-dimensional analysis of laterally confined configurations of marine ice sheets (e.g., Pegler, 2016; Schoof and others, 2017; Haseloff and Sergienko, 2018, 2022; Sergienko and Haseloff, 2023).

4.1 Point-wise backpressure force

In order to determine the force balance at the grounding line, we integrate the ice-shelf momentum balance (4) from x ^g to x ^c and apply Leibniz's rule. The detailed derivations are described in Appendix A. Their result is the components of the force balance at the grounding line

(28a)$$2\nu H( 2u_x + v_y) n_{x}^{g} + \nu H\left(u_y + v_x\right)n_{y}^{g} = {\rho g^\prime\over 2}H^2n_{x}^{g} + {1\over \sqrt{1 + ( x_{y}^{g}) ^2}}\partial_y\int_{x^g}^{x^c}\nu H\left(u_y + v_x\right){\rm d}x,\; $$

(28b)$$\nu H\left(u_y + v_x\right)n_{x}^{g} + 2\nu H( u_x + 2v_y) n_{y}^{g} = {\rho g^\prime\over 2}H^2n_{y}^{g} + {1\over \sqrt{1 + ( x_{y}^{g}) ^2}}\partial_y\int_{x^g}^{x^c} \left[2\nu H( u_x + 2v_y) -{\rho g^\prime\over 2}H^2\right]{\rm d}x,\; $$

where $\displaystyle {x_{y}^{g} = {{\rm d}x^{g}( y) }/{\, {\rm d}y}}$ and $\displaystyle {\{ n^{g}_{x},\; \, n^{g}_{y}\} = {1}/{\sqrt {1 + ( x^{g}_{y}) ^2}}\{ 1,\; \, -x^{g}_{y}}\}$. On the left-hand side are components of the depth-integrated force due to internal deformation in the ice at the grounding line; on the right-hand side are components of the depth-integrated force provided by the ice shelf. The right-hand side components have two terms. The first of which are components of the buoyancy force, ρg ^′/2H ² and are the same if the ice shelf is absent. The second terms are components of the backpressure provided by the ice shelf at each point at the grounding line. These term are y − derivatives of the respective components of the depth-integrated ice-shelf deformation (shear (28a) and the deviation of the extension (or compression) from the ice buoyancy (28b)) integrated through the length of the ice shelf.

The above equations can be written as

(29a)$$f^{BP}_x = 2\nu H( 2u_x + v_y) n_{x}^{g} + \nu H\left(u_y + v_x\right)n_{y}^{g}-{\rho g^\prime\over 2}H^2n_{x}^{g} = {1\over \sqrt{1 + ( x_{y}^{g}) ^2}}\partial_y\int_{x^g}^{x^c}\nu H\left(u_y + v_x\right){\rm d}x,\; $$

(29b)$$f^{BP}_y = \nu H\left(u_y + v_x\right)n_{x}^{g} + 2\nu H( u_x + 2v_y) n_{y}^{g}-{\rho g^\prime\over 2}H^2n_{y}^{g} = {1\over \sqrt{1 + ( x_{y}^{g}) ^2}}\partial_y\int_{x^g}^{x^c} \left[2\nu H( u_x + 2v_y) -{\rho g^\prime\over 2}H^2\right]{\rm d}x,\; $$

where $\{ f^{BP}_x,\; \, f^{BP}_y\}$ are the components of the point-wise backpressure force. The relationship between the components of the point-wise and total backpressure force (19) is

(30a)$$F^{BP}_x = \int_{L^G}f^{BP}_x{\rm d}l$$

(30b)$$F^{BP}_y = \int_{L^G}f^{BP}_y{\rm d}l.$$

The right hand sides of (29) are determined by the y −derivatives. This implies that the point-wise backpressure is a two-dimensional (plane view) phenomenon and is determined by the transverse variability of the ice shelves; hence, the laterally uniform ice shelves provide no backpressure to their grounding lines. This also indicates that the point-wise backpressure of a laterally unconfined ice shelf with transverse variability is non-zero. Its components are determined by the transverse variability of the lateral shear (eqn. (29a)) and imbalance between the buoyancy force and the normal stress in the y-direction (eqn. (29b)) integrated through the length of the ice shelf. It also depends on the shape of the grounding line (i.e., on how it bends and curves), which in its turn depends on the variability of the bed topography in the direction transverse to the ice flow. The effects of the shape of the grounding line have been demonstrated numerically in idealized (Schoof, 2006, section 4.1) and realistic (e.g., Fürst and others, 2016; Gudmundsson and others, 2023) configurations. It should be empasized, however, unconfined ice shelves exert no total backpressure to their grounding lines, as indicated by Eqn (21).

4.2 Buttressing numbers and ratios

The grounding-line force balance (28) gives the following expressions for the buttressing numbers

(31a)$$\eqalign{& K_N = {1\over \displaystyle{\frac{\rho g^\prime}{2}H^2}\left(1 + ( x_{y}^{g}) ^2\right)} \cr& \left\{x_{y}^{g}\partial_y\int_{x^g}^{x^c} \left[2\nu H( u_x + 2v_y) -{\rho g^\prime\over 2}H^2\right]{\rm d}x -\partial_y\int_{x^g}^{x^c}\nu H\left(u_y + v_x\right){\rm d}x\right\},\; }$$

(31b)$$\eqalign{& K_T = {1\over \displaystyle{\frac{\rho g^\prime}{2}H^2}\left(1 + ( x_{y}^{g}) ^2\right)} \cr& \left\{\partial_y\int_{x^g}^{x^c} \left[2\nu H( u_x + 2v_y) -{\rho g^\prime\over 2}H^2\right]{\rm d}x + x_{y}^{g}\partial_y\int_{x^g}^{x^c}\nu H\left(u_y + v_x\right){\rm d}x\right\}}.$$

The corresponding buttressing ratios are Θ _N = 1 − K _N and Θ _T = K _T, respectively.

Expressions (31) show that in addition to the transverse variability through the ice shelf and the grounding-line shape that control the point-wise backpressure components, the buttressing characteristics depend on the ice thickness at the grounding line, and hence the bed topography.

5. Impact of the lateral boundary conditions and submarine melting on backpressure and buttressing

To get a quantitative sense of the effects of lateral boundary conditions and submarine melting on the backpressure of a steady-state configuration, we consider an idealized marine ice sheet flowing over bed topography that varies along and across the direction of ice flow

(32)$$B( x,\; y) = B_0 + B_1\cos{{\pi x\over L_x}} + B_2\cos{{12\pi x\over L_x}}\cos{{6\pi y\over L_y}}$$

All model parameters are listed in Table 1. Figure 3 illustrates the shape of the ice sheet with no slip at the lateral boundaries (Fig. 3a), bed topography and the grounding line positions for a spatially variable melt rate (cyan line) and a spatially uniform melt rate (magenta line). The spatially variable melt rate is

(33)$$\displaystyle{\dot m( x,\; y) = \dot m_0\left[1-\left({x-x_g\over L_x-x_g}\right)^{{1}/{3}}\right]\left[1 + \left({y\over L_y}\right)^2\right] + \dot a},\; $$

where $\dot m_0$ is a constant with different values for different boundary conditions. The spatially uniform melt rate is the area averaged value of (33)

(34)$$\displaystyle{< \dot m( x,\; y) > = {1\over F}\iint_{F}\dot m( x,\; y) \, {\rm d}x{\rm d}y},\; $$

where F is the ice-shelf area. The functional form of melt rate, Eqn (33) is chosen purely for its simplicity, however it mimics the observed and simulated increase of melt rates along the northern boundary due to the effects of sub-ice-shelf cavity circulation (Goldberg and others, 2012a; Adusumilli and others, 2020).

Table 1. Model parameters

($^\star$1 m yr⁻¹ for the unconfined ice shelf.).

Figure 3. (a) Steady-state shape of a marine ice sheet with no slip at the lateral boundaries. (b) Steady-state grounding-line positions obtained with a spatially variable melt rate (33) (cyan line) and a spatially uniform melt rate (magenta line).

The undulated bed topography B (Eqn (32)) results in meandering grounding lines (cyan and magenta lines in Fig. 3b).

We consider three kinds of boundary conditions at the lateral boundaries – no slip (22), lateral shear (25), and a laterally unconfined ice shelf

(35)$${\bf T\cdot}\vec n^{S, N} = \rho g^\prime{H^2\over 2}\vec n^{S, N},\; \quad \{ x,\; y\} \in L^{S, N}.$$

For the no slip and lateral shear we assume the same conditions on the grounded and floating parts; for the laterally unconfined ice shelf we use slip conditions (no shear) at the lateral boundaries of the grounded part.

For each kind of the lateral boundary condition and melt rate we obtain a steady-state configuration as a solution of an optimization problem. To do so, we use the finite-element solver Comsol^TM(COMSOL, 2024) and optimize the grounding line position in such a way that the momentum (1)-(4) and steady-state forms of the mass (10)-(11) balances together with the boundary conditions at the divide (8), calving front (9) and the grounding line (flotation condition (12)) are simultaneously satisfied. For this procedure we use an optimization solver based on the Sparse Nonlinear OPTimizer (SNOPT) algorithm (Gill and others, 2005). The mesh resolution is 5 km away from the grounding line and 500 m in 10 km zone of the grounding line (5 km upstream and downstream).

For each steady-state configuration, we analyze the effective stress τ _eff (the second invariant of the three-dimensional stress tensor) and the principal stress components τ _I and τ _II, both their orientation and magnitude; the buttressing ratios Θ _N and Θ _T; the point-wise backpressure force (Eqn (29)); and the total backpressure at the grounding line $\vec F^{BP}$ (Eqn (19)). In all simulations, we assume that the calving front is fixed.

5.1 No slip

In the case of lateral confinement with no-slip conditions at the lateral boundaries, the ice flow has a characteristic pattern of slow flow near the lateral boundaries and faster flow in the trunk of the grounded and floating portions (Fig. 4a). The presence of undulations on the bed (gray contour lines in Fig. 4a) upstream of the grounding line (the white line in Fig. 4a) and also the spatial variability of the melt rate (33) in the transverse direction cause slight deviations of ice flow from being parallel to its lateral boundaries (the black vectors in Fig. 4a).

Figure 4. Ice flow and stress characteristics for no-slip lateral conditions and spatially variable melt rates $\dot m( x,\; \, y)$ ($\dot a = 1$ m yr⁻¹). (a) ice speed (m yr⁻¹) (color) contour lines are bed elevation; (b) effective stress (kPa) (color), white and black vectors are principal stress components (white–extensional, black–compressional); (c) first principal stress τ _I (kPa) (horizontal color bar) and normal buttressing ratio Θ _N (vertical color bar); (d) second principal stress τ _II (kPa) (horizontal color bar) and tangential buttressing ration Θ _T (vertical color bar); white lines a contours of τ _II =0.

The boundary layers, or shear margins, ~10 km wide are formed on the grounded and floating parts near the lateral boundaries due to the no-slip condition. In the shear margins, the effective stress is of the order of ~80 kPa (Fig. 4b). The principal stress components (white (extensional) and black (compressional) vectors in Fig. 4b) are aligned at ~45$^\circ$ with respect to the direction of ice flow. Both principal stress components are of the order 100–120 kPa (Figs. 4c,d). The first principal stress is always tensile (Fig. 4c) and the second is predominantly compressional (Fig. 4d; the white contour line indicates τ _II = 0) Away from the shear margins, the magnitudes of the effective stress as well as the principal stress components are substantially lower (~20 kPa) (Figs. 4b,d). The presence of the bed undulations results in a slight compression when ice flows around them (Fig. 4d).

At the grounding line (the white line in Fig. 4a), the effective stress is of the order of 80 kPa (Fig. 4b) and is primarily determined by the first principal stress component, which is extensional there (Fig. 4c). The curve of the grounding line is primarily caused by the bed undulations and also by the spatial variability of melt rates (Eqn (33)). As a result of meander of the grounding line both normal Θ _N and tangential Θ _T buttressing ratios (and buttressing numbers K _N and K _T) are non-zero (gray colors in Figs. 4c, d). The magnitude of Θ _N is lager than the magnitude of Θ _T (~0.6 − 1 vs ~0.2).

Comparison of the results of simulations with the spatially variable melt rate (Eqn (33)) to those with the spatially uniform melt rate (Eqn (34)) allows to assess the influence of the melt rate spatial variability on the marine ice-sheet state – its geometry (the ice-thickness distribution and the grounding line position), flow and stress regimes. In the case of the spatially variable melt rate, the grounding line is slightly upstream of the grounding line in the case of spatially uniform melt rate (Figs. 5e,d). Because of the melt-rate variability in the y-direction, the grounding line is not symmetric with respect to the center-line, and its upstream displacement from the grounding line with the spatially uniform melt rate progressively increases from ~1.5 km at the southern boundary L ^S to ~5 km at the northern boundary L ^N. This displacement results in a faster ice flow immediately upstream of the grounding line by ~30 − 40 m yr⁻¹ and also over the whole grounded part by ~5 − 10 m yr⁻¹ (Fig. 5a). The spatial patterns of the speed difference are more complicated on the ice shelf: the flow is faster in the immediate vicinity of the grounding line because of its overall upstream position and also in the shear margins up to ~50 km from the calving front, and it is slower in the rest of the ice shelf.

Figure 5. The effects of spatial variability of melt rates for the case of no-slip lateral conditions. Panels (a)-(d) show differences between configurations obtained with spatially variable $\dot m( x,\; \, y)$ (Eqn (33) and spatially uniform $< \dot m( x,\; \, y) \gt$ melt rates (Eqn (34)); (a) speed (m yr⁻¹); (b) ice thickness (m); (c) first principal stress τ _I (kPa); (d) second principal stress τ _II (kPa); (e) Normal buttressing ratio Θ _N; (f) Tangential buttressing ratio Θ _T. In panels (a) and (b) the white and black lines are the grounding line. In the panel (d) the white lines are contour lines of τ _II =0 (solid with $\dot m( x,\; \, y)$ and dashed with $\langle \dot m( x,\; \, y) \rangle$). In panels (e) and (f) the left, upstream, lines are the grounding lines with $\dot m( x,\; \, y)$ and the right, downstream lines are the grounding lines with $< \dot m( x,\; \, y) \gt$; gray lines are contour lines of bed elevation.

The large-scale patterns in the ice-thickness differences are similar to those of the speed differences. Overall, the ice is slightly thinner (~10 m) on the grounded part (Fig. 5b). It is significantly thinner (more than 100 m) in the immediate vicinity of the grounding line and particularly closer to the northern boundary where the melt rate is the largest (eqn. (33)). On the ice shelf, the ice thickness is smaller almost everywhere except from the vicinity of the calving front where the ice thickness becomes larger compared to that with the spatially uniform melt rate (eqn. (34)) with magnitudes up to 100 m in the shear margins (Fig. 5b). The differences in the ice-shelf thickness are largest where the melt rates are largest.

The magnitudes of the differences in the respective principal stress components obtained in two simulations are largest in the ice-shelf shear margins. In the case of spatially uniform melting the shear-margin spatial extent is smaller compared to that in the case of the spatially variable melting (Figs. 5c, d). The displacement of the grounding line due to spatially variable melting upstream of its position in the case of spatially uniform melting changes the magnitude of the first principal component by ~20 kPa. The slightly different locations of the grounding lines and stress regimes around them result in slightly different magnitudes of the buttressing ratios and numbers, however, their spatial patterns and magnitudes are quite similar for the two spatial distributions of melt rates (Figs. 5e, f).

The spatial patterns of the point-wise backpressure components reflect variations of the bed topography at the grounding line in the transverse direction (Figs. 6a, b). The magnitudes of the backpressure force are the largest near the lateral boundaries, where the lateral shear is the largest. In the case of spatially variable melt rates (dark solid lines), the magnitudes of force components are slightly smaller compared to those produced by the spatially uniform melt rates (dark dashed lines). The differences increase as the impact of the transverse variability of the melt rates increases towards the northern lateral boundary (towards larger values of y, on the left in Figs. 6a, b). Because the buttressing ratios are normalized by the ice thickness at the grounding line (hence the bed elevation), their patterns are less reflective of the bed topography. The magnitudes of the buttressing numbers have larger deviation from unity, in the case of Θ _N, and zero, in the case of Θ _T, towards the lateral boundaries (Figs. 6a, b, light blue and green lines, right vertical axes).

Figure 6. Point-wise backpressure force and buttressing ratios as a function of y for various lateral boundary conditions.(a)–(b) no-slip; (c)–(d) lateral shear; (e)–(f) unconfined ice shelf. The left column shows f _x and Θ _N; the right column shows f _y and Θ _T. The left axes are for f _x,y, the right axes are for Θ _T. Solid lines correspond to the case of the spatially variable melt rates; dashed lines correspond to the spatially uniform melt rates. Note the reverse direction of the horizontal axes, y.

The scalar characteristics, such as the magnitudes of the backpressure force and its components are summarized in Table 2. For the case of the spatially variable melt rate (Eqn (33)), the total backpressure force components computed with Eqn (20) are $F^{BP}_x$=4.33×10¹²N and $F^{BP}_y = 2.48\times 10^{12}$N. The difference between these values and those computed with Eqn (18), i.e., as a sum of integrals along the lateral boundaries L ^N and L ^S is less then 0.1%, and is due to the numerical errors associated with computing the stress components and integrals numerically. For the case of spatially uniform melt rate (Eqn (34)) these values are $F^{BP}_x = 5.4\times 10^{12}$N and $F^{BP}_y = 1.72\times 10^{12}$N. The difference between computations with expressions (20) and (18) is similar – less then 0.1%. The magnitude of the total backpressure, $\vert \vec F^{BP}\vert$, in the case of the spatially variable melt rate is 5×10¹²N, which is smaller than that in the case of the spatially uniform melt rate, 5.7×10¹²N. These values can be compared to the force provided by the basal shear upstream of the grounding line. This force in a two-kilometer zone is 4.69 × 10¹³N (4.66 × 10¹³N in the case of the spatially uniform melt rates) – almost an order of magnitude larger than the total backpressure force.

Table 2. Scalar metrics of the ice-shelf buttressing

Values in parentheses correspond to spatially uniform melt rates.

5.2 Lateral shear

When shear is prescribed at the lateral boundaries, we assume that in the boundary conditions (25) $\vec \tau ^w = -C_wH{\bf v}$, where C _w = 10¹⁰ Pa  m⁻¹s. This boundary condition is thought to mimic the effects of ice softening in the shear margins that develops with time due to fracturing and crevassing – processes that are not represented in the used model. With such a formulation and the chosen parameters, the lateral shear is of the order of 15–20 kPa on the grounded part and 50–60  kPa on the ice shelf. As a result, the ice flow is only about 35–40% slower at the lateral boundaries than the fastest flow in the trunk of the ice stream/ice shelf (Fig. 7a). Away from the lateral boundaries, the ice flow is fairly similar in the cases of no-slip (Fig. 4a). The direction of ice flow is affected by the presence of undulations and spatial variability of the melt rates.

Figure 7. Ice flow and stress characteristics for prescribed shear stress at the lateral boundaries and spatially variable melt rates $\dot m( x,\; \, y)$ ($\dot {a} = 0.5$ m  yr⁻¹). Panels are the same as in Fig. 4.

Apart from the vicinity of the grounding line, where the magnitudes of the effective stress are similar for the two cases of the lateral boundary conditions, the effective stress is substantially lower in the case of the lateral shear boundary conditions (Fig. 7b). Distinct shear zones in which principal stress components change their sign (Figs. 7c, d) are still present, but they are narrower and the magnitudes of the principal stress components are lower than those in the case of no slip at the lateral boundaries (Figs. 4c, d). At the grounding line, the effective stress is of the order of 70  kPa; it is dominated by the first principal stress component (Fig. 7c). The curvature of the grounding line that is formed due to spatial variability of the bed topography, the presence of lateral boundaries and also due to the spatially variable melt rate results in both the normal and tangential buttressing ratios and numbers (the gray colorbars in Figs. 7c, d). Their magnitudes do not substantially differ from those in the case of the no-slip lateral boundary conditions.

In the case of spatially variable melt rate (Eqn (33)), the grounding line slightly diverts to the left (white line in Fig. 7a). Compared to that obtained with the spatially uniform melt rate (Eqn (34)), the ice flow is slightly faster at the northern boundary of the ice shelf (Fig. 8a), and slightly slower through the rest of the ice shelf. The spatial patterns in the differences in the ice thickness are such that that the ice is thinner almost everywhere on the ice shelf with larger thinning in its northern part, and slightly thicker on the southern part near the calving front (Fig. 8b). The magnitudes of the differences of the ice speed and ice thickness are smaller on the ice shelf and are similar on the grounded part to those in the case of no slip at the lateral boundaries (Figs. 5a, b). Differences in the principal stress components (Figs. 8c, d) indicate narrower ice-shelf shear zones in the case of spatially uniform melt rate. Although the magnitudes of the buttressing ratios and numbers are similar for the both kinds of the lateral boundary conditions (Figs. 5e, f and 8e, f), in the case of the lateral shear and spatially variable melt rate, the grounding line position relative to its position in the case of the spatially uniform melt rate is farther upstream by ~2 km in the northern part of the domain, compared to its relative upstream position in the case of no slip at the lateral boundaries (Figs. 5e, f).

Figure 8. The effects of spatial variability of melt rates for the case of the prescribed shear at the lateral boundaries. Panels are the same as in Fig. 5.

The point-wise backpressure components (Figs. 6c, d) are similar to those for the no-slip lateral conditions (Figs. 6a, b). However, in the case of the lateral shear, the magnitudes of the force components, as well as the difference of these magnitudes near the lateral boundaries and away from them, diminish. The effect of the spatially variable melt rates is similar to that of the case of no slip at the lateral boundaries.

The components of the total backpressure force are $F^{BP}_x = 1.7\times 10^{12}$ N, $F^{BP}_y = 0.25\times 10^{12}$ N with the force magnitude of 1.72 × 10¹² N, for the case of the spatially variable melt rate and $F^{BP}_x = 2.6\times 10^{12}$ N, $F^{BP}_y = 0.1\times 10^{12}$ N with the force magnitude is 2.6 × 10¹² N, for the case of the spatially uniform melt rate. The magnitude of the total backpressure is smaller in the case of the spatially variable melt rate than in the case of the spatially uniform melt rate. Both values are much smaller (by a factor of 2 to 3) compared to those in the case of no slip at the lateral boundaries. The force provided by the basal shear in the two-kilometer zone upstream of the grounding line is 5 × 10¹³N, which is slightly larger than the magnitude of this force in the case of no slip at the lateral boundaries.

5.3 Unconfined ice shelf

In the case of a laterally unconfined ice shelf, with conditions Eqn (17) prescribed at the ice-shelf lateral boundaries and calving front, the ice flow is almost uniform downstream of the grounding line (Fig. 9a). The only slight variations in it are caused by the undulated bed topography upstream of it, and the spatially variable melt rate (black vectors in Fig. 9a). The spatial variability of the ice-shelf flow is significantly less compared to the other cases of the lateral boundary conditions (Figs. 4a and 7a). The effective stress is of the order 10–20 kPa through both the grounded and floating parts, except the grounding line and its immediate vicinity, where it is of the order of 70 kPa (Fig. 9b). The first principal stress is extensional and oriented along the ice flow, its magnitude is larger than the magnitude of the second principal stress (black and white vectors in Fig. 9b). There are spatial variations in the first and second principal stresses on the grounded part and downstream of the grounding line, and the second principal stress is both extensional and compressional in these regions (Figs. 9c, d). This variability in the principal stresses is caused by the undulated bed topography and its effects on ice flow upstream of the grounding line.

Figure 9. Ice flow and stress characteristics for an unconfined ice shelf and spatially variable melt rates $\dot m( x,\; \, y)$. Panels are the same as in Fig. 4.

In contrast to the no-slip and shear at the lateral boundaries, for the unconfined ice shelf, the effects of the spatially variable melt rates have no impact on the ice sheet upstream of the grounding line and are confined to the ice shelf only (Fig. 10). With the spatially variable melt rate, the ice flow is slightly faster immediately downstream of the grounding line and slower for the most part of the ice shelf (Fig. 10a). The ice-shelf is thinner except from a zone near the calving front, which is larger near the southern boundary, where the ice is thicker compared to that with the spatially uniform melt rates (Fig. 10b). The spatial patterns of the principal stress components are somewhat similar to those of the ice thickness – the magnitudes of the principal stresses are lower by ~10 kPa through the ice shelf, and slightly larger near the calving front near the southern boundary (Figs. 10c, d). The buttressing ratios and numbers are very similar for the spatially variable and spatially uniform melt rates (Figs. 10e, f). The same is true for the point-wise backpressure forces (Figs. 6e, f). (The small differences are due to numerical artifacts.) Additionally, the force components and buttressing numbers near the lateral boundaries have similar magnitudes to those away from the boundaries. The spatial patters in the force components and the buttressing numbers reflect topographic variability at the grounding line.

Figure 10. The effects of spatial variability of melt rates for an unconfined ice shelf. Panels are the same as in Fig. 5.

The components of the total backpressure force computed with Eqns (20) are $F^{BP}_x = 1.9\times 10^{4}$N and $F^{BP}_y = 0.6\times 10^{4}$N, the integrals on the right-hand side of (18) are zero. The nonzero values obtained with Eqns. (20) are due to numerical errors associated with the numerical nature of integration of these expressions. The force resulting from basal shear in the 2  km zone upstream of the grounding line is 5.99 × 10¹³ N.

6. Discussion

In our analysis we have revisited the concepts of backpressure introduced by Thomas (1977) and buttressing numbers and ratios introduced by Gudmundsson (2013).

6.1 The total and point-wise backpressure force

Starting with the Shallow Stream/Shelf Approximation (SSA) of the momentum balance appropriate for ice-stream and ice-shelf flows (MacAyeal, 1989) and focusing on the effects of the conditions at the lateral boundaries of ice shelves, we have written the ice-shelf momentum balance in an integral form (Eqn (16)), which represents a force balance of the whole ice shelf. This form gives a natural definition of the total backpressure force - a force exerted by the ice shelf on ice at the grounding line (Eqns (19)–(20)). According to the ice-shelf force balance (Eqn (16)), it depends on the conditions at the lateral boundaries and the length of these boundaries.

The integral form of the momentum balance provides an explanation for a widely accepted fact that in the absence of pinning points or ice rises a laterally unconfined ice shelf, as a whole, does not provide any buttressing to the grounded ice upstream of the grounding line. This is not necessarily the case in the point-wise sense, and the point-wise components of the backpressure force (Eqns (29)) are non-zero along the grounding line (6e, f). Their spatial variability is determined by the bed topography variations along the grounding line. Applications of different distributions of the melt rate (Eqns (33) and (34)) to the laterally unconfined ice shelf have no impact on the grounding line and the upstream ice flow; it affects only the ice shelf – its flow, ice-thickness and stress distribution (Fig. 10).

Contrary to the laterally unconfined ice shelves, the lateral confinement with no slip or prescribed shear at the lateral boundaries gives rise to buttressing in both the point-wise and total sense. The magnitudes of the pointwise components increase towards the lateral boundaries and are significantly larger at the boundaries in the case of no slip (6a, d). The total backpressure force is of the order of 10¹²N (sections 5.1–5.2). In comparison, the force exerted by basal shear in a two-kilometer zone upstream of the grounding line is an order of magnitude larger (Table 2). Its magnitude is determined by the magnitude of the sliding coefficient C _b (Eqn (3)). With the chosen value (Table 1) that has been used in many theoretical and numerical studies (Schoof, 2007a, 2007b; Pattyn and others, 2012), the stress-balance of the ice flow upstream of the grounding line is dominated by the basal shear and the driving stress. However, this is not the only possible stress regime, and in a regime of low basal and driving stress (e.g., Sergienko and Wingham, 2019), the magnitude of the total backpressure may be of the same order or exceed the magnitude of the basal shear force. In such circumstances, buttressing may have the dominant effect on the grounding line dynamics.

The boundary conditions with the prescribed shear aim to mimic the effects of ice softening due to fracturing and crevassing in the shear margins of the ice shelves. The lateral shear with the magnitudes of 15–20 kPa leads to a more than twofold reduction in the total backpressure force compared to the case of no-slip at the lateral boundaries, which assumes no changes in the ice stiffness associated with damage of ice in the shear zones.

The expressions for the components of the backpresseure (29), which lead to (31), illustrate that the stress components at the grounding line cannot be approximated by the expression of the ice flux for a laterally uniform ice stream with an unconfined ice shelf derived for a one-dimensional geometry by Schoof (2007a, 2007b) as done in several large-scale ice-sheet models (Ritz and others, 2015; DeConto and Pollard, 2016; Pattyn, 2017; Quiquet and others, 2018; DeConto and others, 2021; Coulon and others, 2023). This is because equating the normal stress component to the vertically integrated pressure deficit at the grounding line, as required by the Schoof (2007a, 2007b) expression, implies that the point-wise backpressure force is zero (29), and the ice shelf has no effect on the stress at the grounding line in the point-wise sense.

6.2 Buttressing numbers and ratios

The buttressing numbers and ratios introduced by Gudmundsson (2013) are defined in terms of the stress components at the grounding line. However, these expressions do not provide any information how stress and its variability through the ice shelf affect stress at the grounding line. The derived expressions for buttressing numbers (31) demonstrate that in addition to the shape of the grounding line, the buttressing numbers and ratios are determined by the transverse variability of the imbalance between normal stress in the across-flow direction and lateral shear integrated over the ice-shelf length. This demonstrates once again that the local backpressure force is a two-dimensional effect (plane view) and without variability in the transverse direction the local backpressure is zero.

A strong dependence of expressions (31) on the ice-shelf transverse variability suggests that the accurate knowledge of transverse variability in the ice-shelf properties (e.g., the ice-shelf stiffness) and processes (e.g., submarine melting) is necessary to accurately capture their impacts on the grounding lines in ice-sheet models. Similar to the results of numerical sensitivity studies (e.g., Feldmann and others, 2022) Eqns (31) imply that to adequately account for the effects of submarine melting on buttressing, numerical models should accurately represent the spatial distribution of melt rates, which are determined by the interactions of ice shelves with ocean circulation in sub-ice-shelf cavities. This requires the use of coupled ice-sheet/ocean models (e.g., Goldberg and Holland, 2022), or parameterizations that could accurately mimic their behavior and account for the dependence of melt rates on the ocean circulation in the cavity, in which the ocean pressure gradients in the direction transverse to the ice-shelf flow and Coriolis force play equally important roles as those of the ocean pressure gradients in the direction along the ice-shelf flow (e.g., Goldberg and others, 2012a; Sergienko, 2013; Goldberg and Holland, 2022).

6.3 Ice-shelf stress distribution

Analysis of the principal stress components obtained in numerical simulations of the laterally confined ice shelves shows that the first principal stress component (defined as the largest eigenvalue) is tensile for all boundary conditions (panel (c) in Figs. 4, 7, 9). The second principal stress can be compressive as well as tensile (panel (d) in Figs. 4, 7, 9; the white contour represents τ _II = 0). The spatial extent of compressive stress depends on the lateral boundary conditions. As Table 2 shows, the fraction of the ice-shelf area in which the second principal stress is compressive ($R^{\tau _{II}< 0}$) is the largest in the case of no slip (67%) and the smallest in the case of the unconfined ice shelf (8%). In the latter case, stress is compressive in the immediate vicinity of the grounding line, and most likely is due to the impact of bed topography on ice flow immediately upstream of the grounding line (Fig. 9(d)).

The largest impact of the spatial variability of the melt rates on the transition of the second principal stress from compressive to tensile is observed for the lateral shear boundary conditions (Fig. 8(d) and Table 2). This is in contrast to the case of no slip, in which the spatial pattern of the compressive stress is slightly different for the two melt-rate distributions (the solid and dashed white contour lines in Fig. 5(d)), however, the area fraction with compressive stress is similar, 67–68% (Table 2). This integral metric, along with others considered in this study (the total backpressure force at the grounding line and the force provided by the basal shear upstream of it) are useful indicators of the ice stress regimes. They could be used to diagnose its temporal evolution in numerical models as well as in observational analyses.

Fürst and others (2016) have used the direction of the second principal stress component to establish the ‘passive shelf ice’, or the ice-shelf ‘safety band’. Their choice was inspired by the ‘compressive arch’ – the compressive principal strain rate at the calving front used by Doake and others (1998) as a criterion for the calving-front stability. The signs of the principal stresses are the same as the signs of the principal strain rates. In contrast to the stresses that require accurate knowledge of the ice viscosity or the ice stiffness, the principal strain-rate components can be estimated from remote sensing observations of the ice-shelf surface velocity. As Fig. 11 illustrates, on the Pine Island Ice Shelf, the second principal strain-rate component $\dot e_{II}$ is predominantly compressive with very large magnitudes at its shear margins (Fig. 11b). In contrast, on the Thwaites Eastern Ice Shelf, $\dot e_{II}$ is predominantly tensile. In the immediate vicinity of the grounding line it is compressional (Fig. 11d). It appears to be caused by the effects of the bed topography on the ice flow upstream of it, similar to the compressive pattern of the second main stress of the unconfined ice shelf (Fig. 9(d)).

Figure 11. Principal strain-rate components of the Pine Island Glacier and Thwaites ice shelves (Rignot and Scheuchl, 2017). Magenta lines are contour lines of $\dot e_{I, II} = 0$. Green lines indicate the grounding lines.

The results of numerical analysis by Fürst and others (2016) indicate that except for a small area near the calving front, almost the entire Pine Island Glacier Ice Shelf provides buttressing, and the removal of large parts of the ice shelf leads to rapid retreat of the grounding line. As Fig. 11b illustrates, its second principal strain rate is predominantly compressive. In contrast to the Pine Island Glacier, the recent disintegration of the Thwaites Eastern Ice Shelf (e.g., Benn and others, 2022) has not caused substantial changes in the dynamics of its grounding line. Before the ice-shelf disintegration its principal strain rate was tensile (Fig. 11c). This leads to a hypothesis that the compressive second principal stress or strain-rate components could be used as a proxy of the amount of buttressing provided by an ice shelf to ice upstream of its grounding line. A physical justification for this hypothesis is similar to the idea of ‘compressive arch’ at the calving front proposed by Doake and others (1998). The negative second principal strain-rate component (hence the negative second principal stress) on the ice shelf indicates that ice is under compression, and its flow is inhibited compared to the case when the ice-shelf flow is purely extensional (both principal stresses are positive). The larger horizontal extent of the compressive stresses may indicate the larger backpressure force provided by the ice shelf to its grounding line. More detailed analysis of this hypothesis is needed; numerical investigations of the relationship between the extent of an ice shelf experiencing compressional stresses and the backpressure at the grounding line will be the subject of future studies.

7. Conclusions

We have revisited the concepts of backpressure introduced by Thomas (1977) and buttressing numbers introduced by Gudmundsson (2013) for marine ice sheets without pinning points or ice rises on their ice shelves. Our results show that backpressure and point-wise buttressing are two-dimensional effects that arise due to transverse variability of the grounded and floating parts of the marine ice sheets. The integral form of the ice-stream and the ice-shelf momentum balance (SSA) provides an innate definition of the total backpressure force at the grounding line. For laterally confined ice shelves, it depends on the stress at the lateral boundaries and their length. For laterally unconfined ice shelves it is zero. However, the point-wise backpressure force for such ice shelves can be non-zero.

The results of numerical analysis show that buttressing of confined ice shelves is highly sensitive to the spatial distributions of submarine melting. They also show that ice shelves with more buttressing tend to have larger areas with a compressive second principle stress. This suggests that the spatial extent of the compressive second principle strain rate can be used as a proxy for buttressing, and changes in this spatial extent may be indicative of the temporal variability of the ice-shelf buttressing.