# Controlling basal strength¶

When using option `-stress_balance ssa+sia`

, the SIA+SSA hybrid stress balance, a
model for basal resistance is required. This model for basal resistance is based, at least
conceptually, on the hypothesis that the ice sheet is underlain by a layer of till
[83]. The user can control the parts of this model:

- the so-called sliding law, typically a power law, which relates the ice base (sliding) velocity to the basal shear stress, and which has a coefficient which is or has the units of a yield stress,
- the model relating the effective pressure on the till layer to the yield stress of that layer, and
- the model for relating the amount of water stored in the till to the effective pressure on the till.

This subsection explains the relevant options.

The primary example of `-stress_balance ssa+sia`

usage is in section Getting started: a Greenland ice sheet example of
this Manual, but the option is also used in sections MISMIP,
MISMIP3d, and Example: A regional model of the Jakobshavn outlet glacier in Greenland.

In PISM the key coefficient in the sliding is always denoted as yield stress \(\tau_c\),
which is `tauc`

in PISM output files. This parameter represents the strength of the
aggregate material at the base of an ice sheet, a poorly-observed mixture of pressurized
liquid water, ice, granular till, and bedrock bumps. The yield stress concept also extends
to the power law form, and thus most standard sliding laws can be chosen by user options
(below). One reason that the yield stress is a useful parameter is that it can be
compared, when looking at PISM output files, to the driving stress (`taud_mag`

in
PISM output files). Specifically, where `tauc`

\(<\) `taud_mag`

you are likely to
see sliding if option `-stress_balance ssa+sia`

is used.

A historical note on modeling basal sliding is in order. Sliding can be added directly to
a SIA stress balance model by making the sliding velocity a local function of the basal
value of the driving stress. Such an SIA sliding mechanism appears in ISMIP-HEINO
[84] and in EISMINT II experiment H [19], among
other places. This kind of sliding is *not* recommended, as it does not make sense to
regard the driving stress as the local generator of flow if the bed is not holding all of
that stress [22], [41]. Within PISM, for historical reasons,
there is an implementation of SIA-based sliding only for verification test E; see section
Verification. PISM does *not* support this SIA-based sliding mode in other contexts.

## Choosing the sliding law¶

In PISM the sliding law can be chosen to be a purely-plastic (Coulomb) model, namely,

Equation (8) says that the (vector) basal shear stress \(\boldsymbol{\tau}_b\) is at most the yield stress \(\tau_c\), and that only when the shear stress reaches the yield value can there be sliding. The sliding law can, however, also be chosen to be the power law

where \(u_{\text{threshold}}\) is a parameter with units of velocity (see below). Condition (8) is studied in [38] and [85] in particular, while power laws for sliding are common across the glaciological literature (e.g. see [71], [49]). Notice that the coefficient \(\tau_c\) in (9) has units of stress, regardless of the power \(q\).

In both of the above equations (8) and (9) we call
\(\tau_c\) the *yield stress*. It corresponds to the variable `tauc`

in PISM output files.
We call the power law (9) a “pseudo-plastic” law with power \(q\) and
threshold velocity \(u_{\text{threshold}}\). At the threshold velocity the basal shear
stress \(\boldsymbol{\tau}_b\) has exact magnitude \(\tau_c\). In equation
(9), \(q\) is the power controlled by `-pseudo_plastic_q`

, and the
threshold velocity \(u_{\text{threshold}}\) is controlled by `-pseudo_plastic_uthreshold`

.
The plastic model (8) is the \(q=0\) case of (9).

See Table 16 for options controlling the choice of sliding law. The
purely plastic case is the default; just use `-stress_balance ssa+sia`

to turn it on.
(Or use `-stress_balance ssa`

if a model with no vertical shear is desired.)

Warning

Options `-pseudo_plastic_q`

and `-pseudo_plastic_uthreshold`

have no effect if
`-pseudo_plastic`

is not set.

Option | Description |
---|---|

`-pseudo_plastic` |
Enables the pseudo-plastic power law model. If this is not set the sliding law is
purely-plastic, so `pseudo_plastic_q` and `pseudo_plastic_uthreshold` are
inactive. |

`-plastic_reg` (m/a) |
Set the value of \(\epsilon\) regularization of the plastic law, in the formula
\(\boldsymbol{\tau}_b = - \tau_c \mathbf{u}/\sqrt{|\mathbf{u}|^2 + \epsilon^2}\). The
default is \(0.01\) m/a. This parameter is inactive if `-pseudo_plastic` is set. |

`-pseudo_plastic_q` |
Set the exponent \(q\) in (9). The default is \(0.25\). |

`-pseudo_plastic_uthreshold` (m/a) |
Set \(u_{\text{threshold}}\) in (9). The default is \(100\) m/a. |

Equation (9) is a very flexible power law form. For example, the linear case is \(q=1\), in which case if \(\beta=\tau_c/u_{\text{threshold}}\) then the law is of the form

(The “\(\beta\)” coefficient is also called \(\beta^2\) in some sources (see [37],
for example).) If you want such a linear sliding law, and you have a value
\(\beta=\) `beta`

in \(\text{Pa}\,\text{s}\,\text{m}^{-1}\), then use this option
combination:

```
-pseudo_plastic \
-pseudo_plastic_q 1.0 \
-pseudo_plastic_uthreshold 3.1556926e7 \
-yield_stress constant -tauc beta
```

This sets \(u_{\text{threshold}}\) to 1 \(\text{m}\,\text{s}^{-1}\) but using units \(\text{m}\,\text{a}^{-1}\).

More generally, it is common in the literature to see power-law sliding relations in the form

where \(C\) is a constant, as for example in sections MISMIP and MISMIP3d. In that case, use this option combination:

```
-pseudo_plastic \
-pseudo_plastic_q m \
-pseudo_plastic_uthreshold 3.1556926e7 \
-yield_stress constant \
-tauc C
```

## Determining the yield stress¶

Other than setting it to a constant, which only applies in some special cases, the above discussion does not determine the yield stress \(\tau_c\). As shown in Table 17, there are two schemes for determining \(\tau_c\) in a spatially-variable manner:

`-yield_stress mohr_coulomb`

(the default) determines the yields stress by models of till material property (the till friction angle) and of the effective pressure on the saturated till, or`-yield_stress constant`

allows the yield stress to be supplied as time-independent data, read from the input file.

In normal modelling cases, variations in yield stress are part of the explanation of the
locations of ice streams [38]. The default model ```
-yield_stress
mohr_coulomb
```

determines these variations in time and space. The value of \(\tau_c\) is
determined in part by a subglacial hydrology model, including the modeled till-pore water
amount `tillwat`

(section Subglacial hydrology), which then determines the effective
pressure \(N_{till}\) (see below). The value of \(\tau_c\) is also determined in part by a
material property field \(\phi=\) `tillphi`

, the “till friction angle”. These
quantities are related by the Mohr-Coulomb criterion [71]:

Here \(c_0\) is called the “till cohesion”, whose default value in PISM is zero (see
[38], formula (2.4)) but which can be set by option `-till_cohesion`

.

Option combination `-yield_stress constant -tauc X`

can be used to fix the yield stress
to have value \(\tau_c = X\) at all grounded locations and all times if desired. This is
unlikely to be a good modelling choice for real ice sheets.

Option | Description |
---|---|

`-yield_stress mohr_coulomb` |
The default. Use equation (10) to determine \(\tau_c\). Only
effective if `-stress_balance ssa` or `-stress_balance ssa+sia` is also set. |

`-till_cohesion` |
Set the value of the till cohesion (\(c_{0}\)) in the plastic till model. The value is a pressure, given in Pa. |

`-tauc_slippery_grounding_lines` |
If set, reduces the basal yield stress at grounded-below-sea-level grid points one cell away from floating ice or ocean. Specifically, it replaces the normally-computed \(\tau_c\) from the Mohr-Coulomb relation, which uses the effective pressure from the modeled amount of water in the till, by the minimum value of \(\tau_c\) from Mohr-Coulomb, i.e. using the effective pressure corresponding to the maximum amount of till-stored water. Does not alter the reported amount of till water, nor does this mechanism affect water conservation. |

`-plastic_phi` (degrees) |
Use a constant till friction angle. The default is \(30^{\circ}\). |

`-topg_to_phi` (list of 4 numbers) |
Compute \(\phi\) using equation (11). |

`-yield_stress constant` |
Keep the current values of the till yield stress \(\tau_c\). That is, do not update
them by the default model using the stored basal melt water. Only effective if
`-stress_balance ssa` or `-stress_balance ssa+sia` is also set. |

`-tauc` |
Directly set the till yield stress \(\tau_c\), in units Pa, at all grounded locations
and all times. Only effective if used with `-yield_stress constant` , because
otherwise \(\tau_c\) is updated dynamically. |

We find that an effective, though heuristic, way to determine \(\phi=\) `tillphi`

in
(10) is to make it a function of bed elevation
[59], [3],
[30]. This heuristic is motivated by hypothesis that basal material
with a marine history should be weak [4]. PISM has a mechanism
setting \(\phi =\) `tillphi`

to be a *piecewise-linear* function of bed elevation. The
option is

```
-topg_to_phi phimin,phimax,bmin,bmax
```

Thus the user supplies 4 parameters: \(\phimin\), \(\phimax\), \(\bmin\), \(\bmax\), where \(b\) stands for the bed elevation. To explain these, we define \(M = (\phimax - \phimin) / (\bmax - \bmin)\). Then

It is worth noting that an earth deformation model (see section Earth deformation models) changes \(b(x,y)=\mathrm{topg}\) used in (11), so that a sequence of runs such as

```
pismr -i foo.nc -bed_def lc -stress_balance ssa+sia -topg_to_phi 10,30,-50,0 ... -o bar.nc
pismr -i bar.nc -bed_def lc -stress_balance ssa+sia -topg_to_phi 10,30,-50,0 ... -o baz.nc
```

will use *different* `tillphi`

fields in the first and second runs. PISM will print a
warning during initialization of the second run:

```
* Initializing the default basal yield stress model...
option -topg_to_phi seen; creating tillphi map from bed elev ...
PISM WARNING: -topg_to_phi computation will override the 'tillphi' field
present in the input file 'bar.nc'!
```

Omitting the `-topg_to_phi`

option in the second run would make PISM continue with the
same `tillphi`

field which was set in the first run.

## Determining the effective pressure¶

When using the default option `-yield_stress mohr_coulomb`

, the effective pressure on
the till \(N_{till}\) is determined by the modeled amount of water in the till. Lower
effective pressure means that more of the weight of the ice is carried by the pressurized
water in the till and thus the ice can slide more easily. That is, equation
(10) sets the value of \(\tau_c\) proportionately to \(N_{till}\). The amount
of water in the till is, however, a nontrivial output of the hydrology (section
Subglacial hydrology) and conservation-of-energy (section Modeling conservation of energy) submodels in
PISM.

Following [86], based on laboratory experiments with till extracted
from an ice stream in Antarctica, [87] propose the following
parameterization which is used in PISM. It is based on the ratio \(s=W_{till}/W_{till}^{max}\)
where \(W_{till}=\) `tillwat`

is the effective thickness of water in the till and
\(W_{till}^{max}=\) `hydrology.tillwat_max`

is the maximum amount of water in the
till (see section Subglacial hydrology):

Here \(P_o\) is the ice overburden pressure, which is determined entirely by the ice
thickness and density, and the remaining parameters are set by options in
Table 18. While there is experimental support for the default
values of \(C_c\), \(e_0\), and \(N_0\), the value of \(\delta=\)
`basal_yield_stress.mohr_coulomb.till_effective_fraction_overburden`

should be
regarded as uncertain, important, and subject to parameter studies to assess its effect.

Option | Description |
---|---|

`-till_reference_void_ratio` |
\(= e_0\) in (12), dimensionless, with default value 0.69 [86] |

`-till_compressibility_coefficient` |
\(= C_c\) in (12), dimensionless, with default value 0.12 [86] |

`-till_effective_fraction_overburden` |
\(= \delta\) in (12), dimensionless, with default value 0.02 [87] |

`-till_reference_effective_pressure` |
\(= N_0\) in (12), in Pa, with default value 1000.0 [86] |

Previous | Up | Next |