Hostname: page-component-848d4c4894-2pzkn Total loading time: 0 Render date: 2024-05-06T17:01:10.171Z Has data issue: false hasContentIssue false

On the nonlinear viscosity of the orthotropic bulk rheology

Published online by Cambridge University Press:  18 May 2022

Nicholas M. Rathmann*
Affiliation:
Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark
David A. Lilien
Affiliation:
Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark Centre for Earth Observation Science, University of Manitoba, Winnipeg, MB, Canada
*
Author for correspondence: Nicholas M. Rathmann, E-mail: rathmann@nbi.ku.dk
Rights & Permissions [Opens in a new window]

Abstract

We compare different ways the bulk flow nonlinearity of glacier ice can be captured in an orthotropic rheology. Specifically, we compare the unapproximated orthotropic rheology, derived from plastic potential theory, to existing approximations that assume either the nonlinear viscosity or fluidity is identical to that of Glen's isotropic flow law. We find, overall, a reasonable agreement between the three orthotropic rheologies, and with existing Dye 3 ice-core deformation tests, although assuming Glen's viscosity provides the best approximation to the unapproximated rheology. Our results therefore suggest that previous studies based on either approximation to the orthotropic rheology are on relatively safe ground in the sense that both approximations generally agree with the unapproximated rheology and experimental data. Finally, we provide the forward and inverse analytical forms of all three rheologies for use in future numerical ice-flow modelling.

Type
Letter
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (https://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution, and reproduction in any medium, provided the original work is properly cited.
Copyright
Copyright © The Author(s), 2022. Published by Cambridge University Press

Introduction

The crystal orientation fabric of glacier ice co-evolves with deformation and can locally enhance the rate of deformation by orders of magnitude (Shoji and Langway, Reference Shoji and Langway1985; Pimienta et al., Reference Pimienta, Duval and Lipenkov1987). How exactly fabric anisotropy evolves and affects the large-scale flow of ice masses has received a lot of attention in the literature, covering a broad range of topics. Questions examined include: How is the flow and age–depth relationship at ice divides affected by fabric anisotropy (Durand and others, Reference Durand2007; Martín and others, Reference Martín, Gudmundsson, Pritchard and Gagliardini2009; Pettit and others, Reference Pettit2011; Martín and Gudmundsson, Reference Martín and Gudmundsson2012)? How does fabric anisotropy affect the flow of grounded ice and ice shelves (Ma and others, Reference Ma2010)? How does fabric anisotropy affect the dynamics of mountain glaciers and ice streams compared to temperature variations (Hruby and others, Reference Hruby2020)? Do near-surface variations in fabric contain climatic information (Kennedy and others, Reference Kennedy, Pettit and Di Prinzio2013; Kennedy and Pettit, Reference Kennedy and Pettit2015)? Can fabric be used as a proxy with memory of past flow conditions (Thorsteinsson and others, Reference Thorsteinsson, Waddington and Fletcher2003; Wilson and Peternell, Reference Wilson and Peternell2011) for discovering e.g. palaeo ice streams (Lilien and others, Reference Lilien, Rathmann, Hvidberg and Dahl-Jensen2021; Llorens and others, Reference Llorens2021)? How might fabric affect inferred basal sliding and hence mass fluxes (Rathmann and Lilien, Reference Rathmann and Lilien2021)?

Addressing such questions requires being able to accurately model the large-scale flow of anisotropic ice. This demands, in turn, a bulk anisotropic rheology that can capture the relevant viscous effects of fabric anisotropy for a given application. Many anisotropic rheologies have been proposed (see e.g. Montagnat and others, Reference Montagnat2014) and often fall into one of two categories:

  1. (1) The bulk rheology is taken to be the grain-averaged, anisotropic monocrystal rheology, for some appropriately defined average (e.g. Martín and others, Reference Martín, Gudmundsson, Pritchard and Gagliardini2009; Rathmann and others, Reference Rathmann, Hvidberg, Grinsted, Lilien and Dahl-Jensen2021, and references therein).

  2. (2) The grain c-axis orientation distribution function (ODF) is assumed to possess certain symmetries, thereby allowing a bulk rheology to be derived from e.g. plastic potential theory (e.g. Gillet-Chaulet and others, Reference Gillet-Chaulet, Gagliardini, Meyssonnier, Montagnat and Castelnau2005; Rathmann and Lilien, Reference Rathmann and Lilien2021, and references therein).

Both approaches have pros and cons; choosing one over the other is a trade off. The first approach places no constraints on what ODFs are permitted, whereas the second approach is, strictly speaking, only valid for a subset of possible ODF patterns (symmetries). Although this makes the first approach attractive, it depends on high-order fabric structure tensors (ODF moments) for flow exponents n > 1 – e.g. structure tensors through order eight for the popular flow exponent of n = 3 – which makes deriving the inverse rheology (needed for numerical modelling) too involved even for simple cases like assuming a homogeneous stress field over the polycrystal scale (Rathmann and others, Reference Rathmann, Hvidberg, Grinsted, Lilien and Dahl-Jensen2021). While a linear-viscous monocrystal rheology (n = 1) or a nonlinear monocrystal fluidity that is independent of crystal orientation (Pettit and others, Reference Pettit, Thorsteinsson, Jacobson and Waddington2007) does not necessarily suffer from this problem, neglecting the observed directionally-dependent nonlinear fluidity of monocrystals (Duval and others, Reference Duval, Ashby and Anderman1983) can lead to bulk strain-rate enhancements being underestimated by at least an order of magnitude (Rathmann and others, Reference Rathmann, Hvidberg, Grinsted, Lilien and Dahl-Jensen2021).

The second approach is widely adopted in ice-flow modelling but requires specifying the directional viscosities (due to fabric) in the directions of the ODF symmetry axes (elaborated on below). The problem is therefore twofold in the sense that the directional viscosities must be calculated using a separate model that depends on the orientation fabric (a so-called ‘micro-macro’ (μ-M) model, that may or may not function as an average over the grains). For example, the bulk anisotropic ‘GOLF’ rheology in the model Elmer/Ice (Gillet-Chaulet and others, Reference Gillet-Chaulet, Gagliardini, Meyssonnier, Montagnat and Castelnau2005) uses tabulated directional viscosities calculated from a separate viscoplastic self-consistent model.

The GOLF rheology in Elmer/Ice is currently the only operational anisotropic rheology of the second category above. The GOLF rheology assumes that the fabric is orthotropic; that is, the ODF has three planes of reflection symmetry with normals mi (i = 1,  2,  3) (Fig. 1a). The initial version of the rheology (Gillet-Chaulet and others, Reference Gillet-Chaulet, Gagliardini, Meyssonnier, Montagnat and Castelnau2005, Reference Gillet-Chaulet, Gagliardini, Meyssonnier, Zwinger and Ruokolainen2006) considered a linear bulk viscosity, but has (e.g. Ma and others, Reference Ma2010; Licciulli, Reference Licciulli2018; Lilien and others, Reference Lilien, Rathmann, Hvidberg and Dahl-Jensen2021) been extended following Martín and others (Reference Martín, Gudmundsson, Pritchard and Gagliardini2009) (henceforth M09 approximation) and Pettit and others (Reference Pettit, Thorsteinsson, Jacobson and Waddington2007) (henceforth P07 approximation) to include a nonlinear viscosity or fluidity, respectively, taken to depend only on the second isotropic invariant of the strain-rate or stress tensor, respectively, by analogy to Glen's isotropic flow law.

Fig. 1. Panel a: An orthotropic c-axis ODF and the three axes of reflection symmetry, mi (principal directions). Panel b: Directional enhancement factors introduced by an orthotropic crystal orientation fabric.

The validity of assuming Glen's viscosity or fluidity is, however, untested, and does not conform with plastic potential theory where the fluidity (viscosity) should depend nonlinearly on the six orthotropic invariants of the stress (strain-rate) tensor. That is, when fabric is not isotropic but rather orthotropic, the second invariant in Glen's flow law decomposes into a weighted sum over six invariants, where the weights depend on the directional enhancement factors along the orthotropic fabric symmetry directions (elaborated on below). Indeed, this was pointed out by Ma and others (Reference Ma2010) who noted that no theoretical or experimental results are currently available to discriminate between the approximations.

In this letter, we provide such a comparison between the M09 approximation (assuming Glen's viscosity), the P07 approximation (assuming Glen's fluidity) and the unapproximated orthotropic rheology derived from plastic potential theory. In doing so, we find that the three orthotropic rheologies largely agree and are generally consistent with Dye 3 ice-core deformation tests. In addition, we give the forward and inverse analytical forms of the three orthotropic rheologies useful for future ice-flow modelling.

In the following, we start by introducing the orthotropic rheology, including the two approximations, followed by a validation and comparison of the three rheologies.

Notation

Throughout, inner, double inner and outer products, are denoted by ${\bf A}\cdot {\bf B} = \sum _{j}A_{ij} B_{jk}$, ${\bf A}\cdot \cdot \, {\bf B} = \sum _{i}\sum _{j} A_{ij} B_{ji}$ and ab = a ib j, respectively, where lower and upper case bold denote vectors and second-order tensors, respectively.

Orthotropic bulk rheology

The orthotropic bulk rheology can be constructed from plastic potential theory by demanding that the rheology be invariant under reflections with respect to mi. Objectivity implies that the forward rheology must depend on the six stress-tensor invariants of the symmetry (reflection) transformations, i.e. τ · · mimi and (τ · τ) · · mimi, or equivalently (Naumenko and Altenbach, Reference Naumenko and Altenbach2007)

$$\eqalign{I_{1}( {\boldsymbol \tau}) & = {\boldsymbol \tau}\cdot\cdot \, {{\bf m}_{2}{\bf m}_{2} - {\bf m}_{3}{\bf m}_{3}\over 2},\; \quad I_{4}( {\boldsymbol \tau}) = {\boldsymbol \tau}\cdot\cdot \, {{\bf m}_{2}{\bf m}_{3} + {\bf m}_{3}{\bf m}_{2}\over 2},\; \cr I_{2}( {\boldsymbol \tau}) & = {\boldsymbol \tau}\cdot\cdot \, {{\bf m}_{3}{\bf m}_{3} - {\bf m}_{1}{\bf m}_{1}\over 2},\; \quad I_{5}( {\boldsymbol \tau}) = {\boldsymbol \tau}\cdot\cdot \, {{\bf m}_{3}{\bf m}_{1} + {\bf m}_{1}{\bf m}_{3}\over 2} ,\; \cr I_{3}( {\boldsymbol \tau}) & = {\boldsymbol \tau}\cdot\cdot \, {{\bf m}_{1}{\bf m}_{1} - {\bf m}_{2}{\bf m}_{2}\over 2},\; \quad I_{6}( {\boldsymbol \tau}) = {\boldsymbol \tau}\cdot\cdot \, {{\bf m}_{1}{\bf m}_{2} + {\bf m}_{2}{\bf m}_{1}\over 2}.}$$

Written compactly, the forward rheology is (see Naumenko and Altenbach (Reference Naumenko and Altenbach2007) for a derivation)

(1)$${\dot{\bf \epsilon}} = \eta^{-1}\sum_{i = 1}^{3} \left[\lambda_{i} I_{i}{{\bf m}_{\,j_i}{\bf m}_{\,j_i} - {\bf m}_{k_i}{\bf m}_{k_i}\over 2} + \lambda_{i + 3} I_{i + 3} {{\bf m}_{\,j_i}{\bf m}_{k_i} + {\bf m}_{k_i}{\bf m}_{\,j_i}\over 2} \right],\; $$

where the nonlinear fluidity is

(2)$$\eta^{-1} = A \Bigg(\sum_{i = 1}^{3} \Big[\lambda_{i}I_{i}^{2} + \lambda_{i + 3}I_{i + 3}^{2} \Big]\Bigg)^{( n-1) /2}.$$

Here, A is the flow-rate factor, λ i are six free material parameters, I i = I i(τ) is assumed implicit, and the index tuples are defined as

(3)$$( j_1,\; j_2,\; j_3) = ( 2,\; 3,\; 1) ,\; $$
(4)$$( k_1,\; k_2,\; k_3) = ( 3,\; 1,\; 2) .$$

Notice that a factor of 1/2(n−1)/2 has been absorbed into A compared to the conventional definition of A in the literature.

The rheology (1)–(2) may be posed in a form that is more relevant to glaciology by expressing λi in terms of directional strain-rate enhancement factors (caused by the crystal orientation fabric), defined relative to Glen's isotropic law ${\dot {\bf \epsilon }}^{{\rm Glen}} = A( {\boldsymbol \tau }\cdot \cdot \, {\boldsymbol \tau }) ^{( n-1) /2} {\boldsymbol \tau }$. Specifically, we consider the shear and longitudinal strain-rate enhancements w.r.t. the three symmetry axes m1, m2, and m3 (Fig. 1b):

(5)$$E_{{i}{\,j}} = {{\dot{\bf \epsilon}}( \hat{{\boldsymbol \tau}}( {\bf m}_{i},\; {\bf m}_{\,j}) ) \cdot\cdot \, {\bf m}_{i}{\bf m}_{\,j}\over {\dot{\bf \epsilon}}^{{\rm Glen}}( \hat{{\boldsymbol \tau}}( {\bf m}_{i},\; {\bf m}_{\,j}) ) \cdot\cdot \, {\bf m}_{i}{\bf m}_{\,j}}\hbox{ for }i,\; j = 1,\; 2,\; 3 ,\; $$

where $\hat {{\boldsymbol \tau }}( {\bf m}_{i},\; \, {\bf m}_{j})$ is an idealized compression (i = j) or shear (i ≠ j) stress-tensor function aligned with the fabric symmetry directions:

$$\hat{{\boldsymbol \tau}}( {\bf m}_{i},\; {\bf m}_{\,j}) = \tau_0 \left\{\matrix{ {\bf I}/3-{\bf m}_{i}{\bf m}_{i} \hfill & \hbox{if }i = j \hfill \cr {\bf m}_{i}{\bf m}_{\,j} + {\bf m}_{\,j}{\bf m}_{i} \hfill & \hbox{if } i \neq j \hfill }\right. .$$

Calculating the six independent constraints from Eqn (5), gives

$$\lambda_{i} = {4\over 3} \Big(E^{{{2/( n + 1) }}}_{{\,j_i}{\,j_i}} + E^{{{2/( n + 1) }}}_{{k_i}{k_i}} - E^{{2/( n + 1) }}_{{i}{i}}\Big),\; \quad \lambda_{i + 3} = 2E^{{2/( n + 1) }}_{{\,j_i}{k_i}} ,\;$$

where E ij remain to be specified using a separate model that depends on the local orientation fabric (elaborated on below).

Inverse rheology

Posing Eqns (1)–(2) in its inverse form, ${\boldsymbol \tau} ({\dot {\epsilon}})$, amounts to solving a nonlinear matrix equation. Unlike other anisotropic rheologies, such as the transversely isotropic rheology (e.g. Rathmann and others, Reference Rathmann, Hvidberg, Grinsted, Lilien and Dahl-Jensen2021), the orthotropic rheology does not tensorially depend on τ but on its projections τ · · mimj. Inverting the rheology is therefore possible by calculating the invariants $I_{l}( { \dot {\bf \epsilon }})$ (l = 1,  2,  3,  4,  5,  6) using Eqns (1)–(2) and solving for I l(τ), combined with the assumption that the inverse law should, too, depend tensorially only on ${\bf m}_{j_i}{\bf m}_{j_i}-{\bf m}_{k_i}{\bf m}_{k_i}$ and ${\bf m}_{j_i}{\bf m}_{k_i} + {\bf m}_{k_i}{\bf m}_{j_i}$. It follows from long but arithmetically straight-forward calculations that

(6)$${\boldsymbol \tau} = \eta \sum_{i = 1}^{3} \Bigg[{\lambda_{i}\over \gamma} (I_{\,j_i}-I_{k_i}) {{\bf I} - 3{\bf m}_{i}{\bf m}_{i}\over 2} + {4\over \lambda_{i + 3}}I_{i + 3}{{\bf m}_{\,j_i}{\bf m}_{k_i} + {\bf m}_{k_i}{\bf m}_{\,j_i}\over 2} \Bigg],\; $$

where $I_{i} = I_{i}( { \dot {\bf \epsilon }})$ is assumed implicit, and the nonlinear viscosity is

(7)$$\eta = A^{-1/n} \, \Bigg(\sum_{i = 1}^{3} \Bigg[{\lambda_{i}\over \gamma} ( I_{\,j_i}-I_{k_i}) ^{2} + {4\over \lambda_{i + 3}}I_{i + 3}^{2} \Bigg]\Bigg)^{( 1-n) /2n} .$$

For convenience, the shorthand γ is defined as

(8)$$\gamma = \sum_{i = 1}^{3}\Big[2E^{{{2/( n + 1) }}}_{{\,j_i}{\,j_i}}E^{{{2/( n + 1) }}}_{{k_i}{k_i}} - E^{{4/( n + 1) }}_{{i}{i}} \Big].$$

Indeed, we have numerically verified that the inverse satisfies ${ \dot {\bf \epsilon }}( {\boldsymbol \tau }( { \dot {\bf \epsilon }}_{{\rm 0}}) ) = { \dot {\bf \epsilon }}_{{\rm 0}}$ for various ${ \dot {\bf \epsilon }}_{{\rm 0}}$. Note that the inverse rheology was simplified by applying the identity m1m1 + m2m2 + m3m3 = I, and that $I_{j_i}-I_{k_i} = -3/2 { \dot {\bf \epsilon }}\cdot \cdot \, {\bf m}_{i}{\bf m}_{i}$.

Pettit and others (Reference Pettit, Thorsteinsson, Jacobson and Waddington2007) approximation

Pettit and others (Reference Pettit, Thorsteinsson, Jacobson and Waddington2007) proposed capturing the bulk flow nonlinearity in a linear-viscous anisotropic rheology by replacing the fluidity with that of Glen's law:

(9)$$\eta^{-1} = A( {\boldsymbol \tau}\cdot\cdot \, {\boldsymbol \tau}) ^{( n-1) /2} .$$

If Eqn (2) is replaced by (9), the bulk flow enhancements (5) imply that

(10)$$\lambda_{i} = {4\over 3} \Big(E_{{\,j_i}{\,j_i}} + E_{{k_i}{k_i}} - E_{{i}{i}}\Big),\; \quad \lambda_{i + 3} = 2E_{{\,j_i}{k_i}} .$$

The inverse rheology can be determined following the same procedure outlined above, yielding Eqn (6) but with the viscosity

(11)$$\eqalign{\eta & = A^{-1/n} \, \Bigg(\sum_{i = 1}^{3} \Bigg[{3\over 2}{\lambda_{i}^{2}\over \gamma^{2}} ( I_{\,j_i}-I_{k_i}) ^{2} + {8\over \lambda_{i + 3}^{2}}I_{i + 3}^{2} \cr & \quad+ {3\over 2}{\lambda_{\,j_i}\lambda_{k_i}\over \gamma^{2}} ( I_{i}-I_{\,j_i}) ( I_{i}-I_{k_i}) \Bigg]\Bigg)^{( 1-n) /2n} ,\; }$$

where

(12)$$\gamma = \sum_{i = 1}^{3} \Big[2E_{{\,j_i}{\,j_i}}E_{{k_i}{k_i}} - E^{2}_{{i}{i}} \Big].$$

Martín and others (Reference Martín, Gudmundsson, Pritchard and Gagliardini2009) approximation

Martín and others (Reference Martín, Gudmundsson, Pritchard and Gagliardini2009) proposed capturing the bulk flow nonlinearity in a linear-viscous anisotropic rheology by instead replacing the viscosity with that of Glen's law:

(13)$$\eta = A^{-1/n} ( {\dot{\bf \epsilon}}\cdot\cdot \, {\dot{\bf \epsilon}}) ^{( 1-n) /2n}.$$

In order to derive the forward and inverse orthotropic rheology consistent with the viscosity (13), consider replacing the fluidity of the forward rheology (2) with

(14)$$\eta^{-1} = A \, \Bigg(\sum_{i = 1}^{3} \Big[\lambda_{i}'I_{i}^{2} + \lambda_{i + 3}' I_{i + 3}^{2} \Big]\Bigg)^{( n-1) /2} ,\; $$

where λi′ ≠ λi. The corresponding inverse rheology, ${\boldsymbol \tau }( { \dot {\bf \epsilon }})$, is then given by Eqn (6) with Glen's viscosity (13) if

(15)$$\lambda_{i}' = {1\over 2} \Bigg(\lambda_{i}^{2} + \lambda_{i}{\lambda_{\,j_i} + \lambda_{k_i}\over 2} - {\lambda_{\,j_i}\lambda_{k_i}\over 2} \Bigg),\quad\lambda_{i + 3}' = {1\over 2} \lambda_{i + 3}^{2},\; $$

and

(16)$$\gamma = {9\over 16} ( \lambda_{1}\lambda_{2} + \lambda_{1}\lambda_{3} + \lambda_{2}\lambda_{3}) .$$

Relating λi to the enhancements E ij using Eqn (5), yields

(17)$$\lambda_{i + 3} = 2 E^{1/n}_{{\,j_i}{k_i}},\; $$

and

(18)$$E_{{i}{i}} = \Bigg({3\over 16} \Bigg)^{( n + 1) /2} 2\, \Big(\lambda_{\,j_i} + \lambda_{k_i}\Big)\Big(\lambda_{\,j_i}^{2} + \lambda_{k_i}^{2} + \lambda_{\,j_i}\lambda_{k_i}\Big)^{( n-1) /2} ,\; $$

where λi(E 11,  E 22,  E 33) can be inverted for numerically.

Comparison of rheologies

It is not immediately clear how the three orthotropic rheologies compare. Although the rheologies are constructed to reproduce the same directional enhancement factors when the principal stress and fabric directions are aligned by virtue of Eqn (5), the rheologies are not guaranteed to agree when the stress and fabric directions depart from these idealized ‘calibration states’ – except for n = 1 where the rheologies are identical by construction. To quantify the potential discrepancies, we consider two experiments: (1) numerically reproducing existing deformation tests made on Dye 3 ice-core samples, and (2) calculating strain-rates for synthetic ODFs that are modelled separately.

Dye 3 deformation tests

Shoji and Langway (Reference Shoji and Langway1985) investigated how strain rates depend on the misalignment between the axis of uniaxial compression (stress direction) and the preferred direction of a strong single maximum fabric in Dye 3 ice-core samples. Specifically, they considered the longitudinal strain-rate parallel to the stress direction, divided by that predicted by Glen's law for an isotropic fabric. In terms of the notation above, they calculated the strain-rate ratio

(19)$${\dot{\epsilon}_{vv}\over \dot{\epsilon}_{vv}^{{\rm Glen}}} = {{\dot{\bf \epsilon}}( \hat{{\boldsymbol \tau}}( {\bf v},\; {\bf v}) ) \cdot\cdot \, {\bf v}{\bf v}\over { \dot{\bf\epsilon}}^{{\rm Glen}}( \hat{{\boldsymbol \tau}}( {\bf v},\; {\bf v}) ) \cdot\cdot \, {\bf v}{\bf v}} ,\; $$

where $\dot {\epsilon }_{vv}$ is the experimentally determined quantity, v(θ) = [sin(θ), 0,  cos(θ)] is the compressive axis, and θ is the angle between the compressive axis and the single-maximum direction. Only when θ = 0° are the three rheologies guaranteed (constructed) to agree (by virtue of Eqn (5)), and when θ = 90° by virtue of symmetry in the rheology. Considering intermediate angles allows, therefore, to make the differences between the rheologies clear, and to validate the rheologies against Dye 3 deformation tests.

Following Shoji and Langway (Reference Shoji and Langway1985), we take n = 3 and note that both A and the stress tensor magnitude cancel by the division in (19) (the same applies, in principle, to other isotropic contributions, such as impurities or strain-hardening, insofar as their effect can be captured in A). We furthermore assume that the orientation fabric is effectively unidirectional (all c-axes aligned with ${\hat {\boldsymbol z}}$) which is approximately true (Herron and others, Reference Herron, Langway and Brugger1985). The six bulk enhancement factors E ij (the remaining free parameters of the rheologies) are calculated using a grain-averaged, transversely isotropic monocrystal rheology that depends on the ODF, calibrated to reproduce directional enhancement factors observed from simple-shear deformation experiments (see Rathmann and Lilien (Reference Rathmann and Lilien2021) for details). Notice that the viscoplastic self-consistent approach adopted by Elmer/Ice (Gillet-Chaulet and others, Reference Gillet-Chaulet, Gagliardini, Meyssonnier, Montagnat and Castelnau2005) could equally well have been used to calculate E ij, and that E ij are constants independent of θ as far as the orthotropic rheologies are concerned.

Figure 2 shows the predicted relative strain-rates (lines) compared to those measured by Shoji and Langway (Reference Shoji and Langway1985) (markers). We find that the M09 approximation provides the best approximation to the unapproximated orthotropic rheology, but that all three rheologies compare relatively well with observations. However, the symmetry around θ = 45°, predicted by the orthotropic rheologies, is not supported by the deformation tests; that is, the deformation tests suggest that compression perpendicular to the single-maximum direction is slightly softer than compression parallel to it.

Fig. 2. Strain rates, relative to Glen's law, predicted by the three rheologies (black and grey lines) for a vertical unidirectional fabric (strong single maximum) as a function of the angle θ made between the stress direction (axis of compression, v) and the preferred fabric direction. Markers indicate measurements reported in Shoji and Langway (Reference Shoji and Langway1985).

The fluidity always cancels in the division (19) for the P07 approximation, implying that the strain-rates predicted by the P07 approximation, relative to Glen's flow law (i.e. Fig. 2 and results below), are identical for all n. Hence, in the linear limit (n = 1) all three rheologies produce relative strain-rates identical to those shown here for the P07 approximation (recall the three rheologies are constructed to conform in the linear limit).

We mention in passing that the community appears lately to have become more receptive to evidence that alternative flow exponents, n ≃ 4, might be relevant in some circumstances (Bons and others, Reference Bons2018; Qi and Goldsby, Reference Qi and Goldsby2021; Millstein and others, Reference Millstein, Minchew and Pegler2022). The grey lines in Figure 2 show the corresponding behaviour for n = 4. Although the data provided by Shoji and Langway (Reference Shoji and Langway1985) is calculated assuming n = 3, and therefore cannot be used to discriminate between different n, we note that the difference between n = 3 and n = 4 is relatively small compared to n = 1 (dotted black line).

Synthetic fabrics

We also quantified the potential discrepancies between the three rheologies by calculating the strain rates predicted for a synthetic ice parcel with an evolving fabric. Specifically, we used our spectral fabric model (Rathmann and Lilien, Reference Rathmann and Lilien2021; Rathmann and others, Reference Rathmann, Hvidberg, Grinsted, Lilien and Dahl-Jensen2021) to generate synthetic fabric (ODF) histories of an ice parcel subject to confined vertical compression or vertical simple shear, prescribed in terms of a constant strain-rate tensor. The two experiments consider exclusively the effect of lattice rotation (strain-induced rotation of c-axes) for ~3/4 of the total modelled parcel strains, whereas discontinuous dynamic recrystallization (DDRX) is assumed to dominate throughout the remaining deformation. We refer the reader to Rathmann and Lilien (Reference Rathmann and Lilien2021) for details on the model, but note that the (normalized) stress tensor – which determines the orientation of nucleated grains during DDRX – is assumed to be coaxial with the strain-rate tensor, and in this way no ice-flow modelling is involved. Examples of the simulated ODFs for different parcel strains, ε, are shown in Figures 3b–e and 4b–e for the two deformation experiments (red dots denote the fabric principal directions).

Fig. 3. Synthetic confined compression experiment. Panel a: Strain rate predicted by the three orthotropic rheologies, relative to Glen's isotropic flow law, as a function of an evolving ODF (i.e. modelled parcel strain, ε) for a constant stress tensor. Blue (lattice rotation) and red (DDRX) regions indicate the process dominating fabric evolution as a function of modelled parcel strain, ε. Panel b–e: ODFs at selected parcel strains. Red dots (mi) indicate the principal directions of the second-order structure tensor.

Fig. 4. Synthetic simple shear experiment. Figure text same as in Figure 3.

Given the synthetic ODFs at each degree of parcel strain, we calculate E ij in the same way as done above for the Dye 3 experiments.

Figures 3a and 4a show the longitudinal ($\dot {\epsilon }_{zz}$) and shear ($\dot {\epsilon }_{xz}$) strain rates predicted by the three rheologies, relative to Glen's law, for a constant (shared) stress-tensor that is coaxial with the strain-rate tensor used to model the parcel deformation. Note again that both the rate factor, A, and the stress-tensor magnitude cancel by considering strain rates relative to (divided by) Glen's law. We find that the largest difference between the M09 approximation and the unapproximated rheology is 5% (across both experiments), whereas for the P07 approximation the largest difference is 30%.

We point out that orthotropy might be a poor approximation for fabrics (ODFs) strongly affected by DDRX, such as seen in Figure 3e where a tetragonal symmetry is found. In this case, the double maximum pattern (two maxima in each hemisphere) introduce two additional symmetry axes along which the directional viscosities cannot be specified in an orthotropic rheology.

Discussion and conclusion

Overall, we find that both M09 and P07 approximations provide a reasonable nonlinear extension to the linear orthotropic rheology when compared to the unapproximated rheology derived from plastic potential theory. We suggest, therefore, that previous studies relying on an orthotropic bulk rheology with either the M09 or P07 approximation (see Introduction) are on relatively safe grounds even though not rigorously justified. Of course, carefully constructed deformation tests are needed to determine which version of the orthotropic rheology best agrees with the behaviour of real polycrystalline ice. For this purpose, deformation tests like those conducted on Dye 3 ice-core samples (Shoji and Langway, Reference Shoji and Langway1985) might be useful (i.e. varying the misalignment between the compressive stress direction and the preferred fabric direction). However, existing experimental results (e.g. Shoji and Langway, Reference Shoji and Langway1985, considered here) seem too ambiguous for discriminating between the rheologies (Fig. 2). In this light, the fact that the rheologies are found to produce relatively similar responses is a comforting result.

We emphasize, however, that the P07 approximation to the unapproximated orthotropic rheology is less accurate than the M09 approximation. This is particularly true for ice under compression with a stress direction that is misaligned with the preferred fabric direction (Fig. 2), possibly relevant in warm, highly stressed areas of ice sheets where DDRX is active (Figs. 3 and 4) (Duval and Castelnau, Reference Duval and Castelnau1995). As numerical ice-flow models develop support for DDRX in the future, our results provide an important caveat: the choice of a bulk flow nonlinearity in anisotropic rheologies (of the second kind described in the introduction) might affect modelled ice velocities for DDRX-induced fabrics. The two approximations might additionally lead to larger discrepancies in coupled ice-flow modelling as reported by Martín and Gudmundsson (Reference Martín and Gudmundsson2012); if small differences in viscosity lead to different deformation and fabric, minor differences could, potentially, be amplified to become significant.

We end by noting that both the forward and inverse analytical forms of all three rheologies were newly provided here in coordinate-independent form (i.e. not the fabric eigen basis) for easy use in future numerical ice-flow modelling.

Code availability

The model code (rheologies, fabric evolution and enhancement-factor calculations) is available at https://github.com/nicholasmr/specfab.

Acknowledgements

We wish to thank Editors Alan Rempel and Ralf Greve, and two anonymous reviewers for providing comments that greatly improved this paper. Notably, this includes more clearly demonstrating the differences between the three orthotropic rheologies by the numerical experiments presented here. The research leading to these results has received funding from the Villum Investigator grant nr. 16572 as part of the IceFlow project.

References

Bons, PD, and 6 others (2018) Greenland ice sheet: higher nonlinearity of ice flow significantly reduces estimated basal motion. Geophysical Research Letters 45(13), 65426548. doi: 10.1029/2018GL078356CrossRefGoogle Scholar
Durand, G, and 8 others (2007) Change in ice rheology during climate variations – implications for ice flow modelling and dating of the EPICA Dome C core. Climate of the Past 3(1), 155167. doi: 10.5194/cp-3-155-2007CrossRefGoogle Scholar
Duval, P, Ashby, MF and Anderman, I (1983) Rate-controlling processes in the creep of polycrystalline ice. The Journal of Physical Chemistry 87(21), 40664074. doi: 10.1021/j100244a014CrossRefGoogle Scholar
Duval, P and Castelnau, O (1995) Dynamic recrystallization of ice in polar ice sheets. Le Journal de Physique IV 5(C3), C3197. doi: 10.1051/jp4:1995317Google Scholar
Gillet-Chaulet, F, Gagliardini, O, Meyssonnier, J, Montagnat, M and Castelnau, O (2005) A user-friendly anisotropic flow law for ice-sheet modeling. Journal of Glaciology 51(172), 314. doi: 10.3189/172756505781829584CrossRefGoogle Scholar
Gillet-Chaulet, F, Gagliardini, O, Meyssonnier, J, Zwinger, T and Ruokolainen, J (2006) Flow-induced anisotropy in polar ice and related ice-sheet flow modelling. Journal of Non-Newtonian Fluid Mechanics 134(1), 3343.CrossRefGoogle Scholar
Herron, SL, Langway, CC Jr and Brugger, KA (1985) Ultrasonic Velocities and Crystalline Anisotropy in the Ice Core from Dye 3, Greenland, 23–31. American Geophysical Union (AGU). doi: 10.1029/GM033p0023.CrossRefGoogle Scholar
Hruby, K, and 5 others (2020) The impact of temperature and crystal orientation fabric on the dynamics of mountain glaciers and ice streams. Journal of Glaciology 66(259), 755765. doi: 10.1017/jog.2020.44CrossRefGoogle Scholar
Kennedy, JH and Pettit, EC (2015) The response of fabric variations to simple shear and migration recrystallization. Journal of Glaciology 61(227), 537550. doi: 10.3189/2015JoG14J156CrossRefGoogle Scholar
Kennedy, JH, Pettit, EC and Di Prinzio, CL (2013) The evolution of crystal fabric in ice sheets and its link to climate history. Journal of Glaciology 59(214), 357373. doi: 10.3189/2013JoG12J159CrossRefGoogle Scholar
Licciulli, C (2018) Full Stokes ice-flow modeling of the high-Alpine glacier saddle Colle Gnifetti, Monte Rosa: flow field characterization for an improved interpretation of the ice-core records. Ph.D. thesis, University of Heidelberg, Heidelberg, Germany. doi: 10.11588/heidok.00023981.CrossRefGoogle Scholar
Lilien, DA, Rathmann, NM, Hvidberg, CS and Dahl-Jensen, D (2021) Modeling ice-crystal fabric as a proxy for ice-stream stability. Journal of Geophysical Research: Earth Surface 126(9), e2021JF006306. doi: 10.1029/2021JF006306.Google Scholar
Llorens, MG, and 9 others (2021) Can changes in ice-sheet flow be inferred from crystallographic preferred orientations?. The Cryosphere Discussions 2021, 124. doi: 10.5194/tc-2021-224Google Scholar
Ma, Y, and 5 others (2010) Enhancement factors for grounded ice and ice shelves inferred from an anisotropic ice-flow model. Journal of Glaciology 56(199), 805812. doi: 10.3189/002214310794457209CrossRefGoogle Scholar
Martín, C and Gudmundsson, GH (2012) Effects of nonlinear rheology, temperature and anisotropy on the relationship between age and depth at ice divides. The Cryosphere 6(5), 12211229. doi: 10.5194/tc-6-1221-2012CrossRefGoogle Scholar
Martín, C, Gudmundsson, GH, Pritchard, HD and Gagliardini, O (2009) On the effects of anisotropic rheology on ice flow, internal structure, and the age–depth relationship at ice divides. Journal of Geophysical Research: Earth Surface 114(F4), F04001. doi: 10.1029/2008JF001204CrossRefGoogle Scholar
Millstein, JD, Minchew, BM and Pegler, SS (2022) Ice viscosity is more sensitive to stress than commonly assumed. Communications Earth & Environment 3(1), 17. doi: 10.1038/s43247-022-00385-xCrossRefGoogle Scholar
Montagnat, M, and 11 others (2014) Multiscale modeling of ice deformation behavior. Journal of Structural Geology 61, 78108.CrossRefGoogle Scholar
Naumenko, K and Altenbach, H (2007) Modeling of creep for structural analysis. Berlin, Heidelberg: Springer. doi: 10.1007/978-3-540-70839-1.CrossRefGoogle Scholar
Pettit, EC, and 6 others (2011) The crossover stress, anisotropy and the ice flow law at Siple Dome, West Antarctica. Journal of Glaciology 57(201), 3952. doi: 10.3189/002214311795306619CrossRefGoogle Scholar
Pettit, EC, Thorsteinsson, T, Jacobson, HP and Waddington, ED (2007) The role of crystal fabric in flow near an ice divide. Journal of Glaciology 53(181), 277288. doi: 10.3189/172756507782202766CrossRefGoogle Scholar
Pimienta, P, Duval, P and Lipenkov, VY (1987) Mechanical behaviour of anisotropic polar ice. In The Physical Basis of Ice Sheet Modelling, IAHS Publication No. 170, 57–66, IAHS Press, Wallingford, UK.Google Scholar
Qi, C and Goldsby, DL (2021) An experimental investigation of the effect of grain size on ‘dislocation creep’ of ice. Journal of Geophysical Research: Solid Earth 126(9), e2021JB021824. doi: 10.1029/2021JB021824. e2021JB021824 2021JB021824.Google Scholar
Rathmann, NM, Hvidberg, CS, Grinsted, A, Lilien, DA and Dahl-Jensen, D (2021) Effect of an orientation-dependent non-linear grain fluidity on bulk directional enhancement factors. Journal of Glaciology 67(263), 569575. doi: 10.1017/jog.2020.117CrossRefGoogle Scholar
Rathmann, NM and Lilien, DA (2021) Inferred basal friction and mass flux affected by crystal-orientation fabrics. Journal of Glaciology 1, 117. doi: 10.1017/jog.2021.88Google Scholar
Shoji, H and Langway, CC (1985) Mechanical Properties of Fresh Ice Core from Dye 3, Greenland, 39–48. American Geophysical Union (AGU). doi: 10.1029/GM033p0039.CrossRefGoogle Scholar
Thorsteinsson, T, Waddington, ED and Fletcher, RC (2003) Spatial and temporal scales of anisotropic effects in ice-sheet flow. Annals of Glaciology 37, 4048. doi: 10.3189/172756403781815429CrossRefGoogle Scholar
Wilson, CJ and Peternell, M (2011) Evaluating ice fabrics using fabric analyser techniques in Sørsdal Glacier, East Antarctica. Journal of Glaciology 57(205), 881894. doi: 10.3189/002214311798043744CrossRefGoogle Scholar
Figure 0

Fig. 1. Panel a: An orthotropic c-axis ODF and the three axes of reflection symmetry, mi (principal directions). Panel b: Directional enhancement factors introduced by an orthotropic crystal orientation fabric.

Figure 1

Fig. 2. Strain rates, relative to Glen's law, predicted by the three rheologies (black and grey lines) for a vertical unidirectional fabric (strong single maximum) as a function of the angle θ made between the stress direction (axis of compression, v) and the preferred fabric direction. Markers indicate measurements reported in Shoji and Langway (1985).

Figure 2

Fig. 3. Synthetic confined compression experiment. Panel a: Strain rate predicted by the three orthotropic rheologies, relative to Glen's isotropic flow law, as a function of an evolving ODF (i.e. modelled parcel strain, ε) for a constant stress tensor. Blue (lattice rotation) and red (DDRX) regions indicate the process dominating fabric evolution as a function of modelled parcel strain, ε. Panel b–e: ODFs at selected parcel strains. Red dots (mi) indicate the principal directions of the second-order structure tensor.

Figure 3

Fig. 4. Synthetic simple shear experiment. Figure text same as in Figure 3.