On the nonlinear viscosity of the orthotropic bulk rheology

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.

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, 1985; Pimienta et al., 1987). 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, 2007; Martín and others, 2009; Pettit and others, 2011; Martín and Gudmundsson, 2012)? How does fabric anisotropy affect the flow of grounded ice and ice shelves (Ma and others, 2010)? How does fabric anisotropy affect the dynamics of mountain glaciers and ice streams compared to temperature variations (Hruby and others, 2020)? Do near-surface variations in fabric contain climatic information (Kennedy and others, 2013; Kennedy and Pettit, 2015)? Can fabric be used as a proxy with memory of past flow conditions (Thorsteinsson and others, 2003; Wilson and Peternell, 2011) for discovering e.g. palaeo ice streams (Lilien and others, 2021; Llorens and others, 2021)? How might fabric affect inferred basal sliding and hence mass fluxes (Rathmann and Lilien, 2021)?

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, 2014) 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, 2009; Rathmann and others, 2021, 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, 2005; Rathmann and Lilien, 2021, 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, 2021). While a linear-viscous monocrystal rheology (n = 1) or a nonlinear monocrystal fluidity that is independent of crystal orientation (Pettit and others, 2007) does not necessarily suffer from this problem, neglecting the observed directionally-dependent nonlinear fluidity of monocrystals (Duval and others, 1983) can lead to bulk strain-rate enhancements being underestimated by at least an order of magnitude (Rathmann and others, 2021).

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, 2005) 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 m_i (i = 1,  2,  3) (Fig. 1a). The initial version of the rheology (Gillet-Chaulet and others, 2005, 2006) considered a linear bulk viscosity, but has (e.g. Ma and others, 2010; Licciulli, 2018; Lilien and others, 2021) been extended following Martín and others (2009) (henceforth M09 approximation) and Pettit and others (2007) (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, m_i (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 (2010) 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 m_i. Objectivity implies that the forward rheology must depend on the six stress-tensor invariants of the symmetry (reflection) transformations, i.e. τ · · m_im_i and (τ · τ) · · m_im_i, or equivalently (Naumenko and Altenbach, 2007)

$$\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 (2007) 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 m₁, m₂, and m₃ (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, 2021), the orthotropic rheology does not tensorially depend on τ but on its projections τ · · m_im_j. 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 m₁m₁ + m₂m₂ + m₃m₃ = 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 (2007) approximation

Pettit and others (2007) 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 (2009) approximation

Martín and others (2009) 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 ₁₁,  E ₂₂,  E ₃₃) 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 (1985) 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 (1985), 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, 1985). 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 (2021) for details). Notice that the viscoplastic self-consistent approach adopted by Elmer/Ice (Gillet-Chaulet and others, 2005) 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 (1985) (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 (1985).

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, 2018; Qi and Goldsby, 2021; Millstein and others, 2022). The grey lines in Figure 2 show the corresponding behaviour for n = 4. Although the data provided by Shoji and Langway (1985) 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, 2021; Rathmann and others, 2021) 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 (2021) 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 (m_i) 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, 1985) 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, 1985, 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, 1995). 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 (2012); 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.