Hostname: page-component-7c8c6479df-ws8qp Total loading time: 0 Render date: 2024-03-27T03:18:43.466Z Has data issue: false hasContentIssue false

Two orthotropic models for strain-induced anisotropy of polar ice

Published online by Cambridge University Press:  20 January 2017

R. Staroszczyk
Affiliation:
School of Mathematics, University of East Anglia,Norwich NR47TJ, England
O. Gagliardini
Affiliation:
Laboratoire de Glaciologie et Geophysique de l'Environement du CNRS, BP96, 38402 Saint-Martin-d’Heres Cedex,France
Rights & Permissions [Opens in a new window]

Abstract

As polar ice descends from the free surface to depth in a large ice sheet, it undergoes deformations which give rise to the formation and subsequent evolution of a fabric and associated anisotropy. In this paper two orthotropic models of such strain-induced anisotropy are considered. Model A is based on analysis of the microscopic behaviour of an individual ice crystal with transversely isotropic response and assumed uniform stress in a polycrystal. The macroscopic response of the ice aggregate is then derived by applying the concept of an orientation distribution function, and the resulting viscous law relates the strain rate to the stress and three structure tensors. In model B it is assumed that the macroscopic response of ice is determined by the fabric induced entirely by macroscopic deformations, and all microprocesses taking place at the grain level are ignored. A constitutive relation is derived from a general frame-indifferent law for orthotropic materials, and expresses the stress in terms of the strain rate, strain and throe structure tensors. The two models are applied to determine the viscous response of ice to continued uniaxial compression and simple shearing in order to compare the predictions of both theories.

Type
Research Article
Copyright
Copyright © International Glaciological Society 1999

1. Introduction

Polar ice cores drilled at different sites in Antarctica and Greenland (Reference Gow and WilliamsonGOW and Williamson, 1976: Reference Russell-Head and BuddRussell-Head and Budd, 1979; Reference Herron and LangwayHerron and Langway, 1982; Lpenkov and others, 1989; reveal strong fabrics, with significant alignment off axes of individual crystals along some preferential directions, induced by strains to which the ice is subject as it descends from the free surface to depth in an ice sheet. The fabric created during the ice deformation gives rise to macroscopic anisotropy of the medium which considerably affects polar ice-sheet dynamic behaviour over long time-scales. This has been confirmed by numerical simulations carried out by Reference MangeneyMangeney and others (1996), who have applied a transversely isotropic flow law. with a rotational .symmetry axis assumed to be vertical everywhere, to the ice-sheet flow problem with a non-evolving ice fabric derived from the measurements made along the GRIP ice core (Thorstcinsson and others, 1997). It has been found that for a given, fixed in time, free-surface elevation the assumption of the ice anisotropy results in much faster ice (low compared to the isotropic case, and an estimated accumulation rate required to maintain a steady-state flow exceeds by more than 1.5 times the rate for isotropic ice. Similar conclusions have been drawn by Reference Mangeney, Califano and HutterMangeney and others 1997), who considered a more realistic problem in which the free-surface elevation was calculated lor a fixed accumulation rate. Their numerical results showed that the ice anisotropy, as well as leading to a globally faster flow, significantly increases shear stresses near the bedrock, makes the free surface flatter in the ice-divide region, and smooths out the effects of the bed topography. This clearly indicates that the strain-induced anisotropy must be included in any large-scale numerical model of polar ice sheets if realistic results are to be obtained. However, there is still no generally accepted constitutive law describing the anisotropic behaviour of polar ice. This is partly because of the complexity of physical phenomena taking place on the microscale of single crystals during the fabric formation, but also because of difficulties associated with t he mathematical modelling of the strain-induced, evolving in time, anisotropy.

Following Reference Lliboutry and DuvalLliboutry and Duval (1985) and Aliev 1992 . we distinguish three main regions along the depth of a typical polar ice sheet in accordance with the microscopic processes which dominate the ice-fabric evolution.

  • (1) In the upper part of a large ice sheet, extending from approximately 100m under the free surface to about one-third of its thickness, shear stresses are negligibly small compared to normal stresses. Deformation is mainly due to dislocation glide on basal planes, and caxes of grains rotate towards compressional axes. The grain-sizc increases linearly with the age of ice, which is here nearly proportional lo its depth, and no new grains are created.

  • (2) Between approximately one-third and two-thirds of the ice cap thickness, the shear stresses gradually increase, though they are still smaller than the normal stresses. With increasing deformation and strain energy stored in grains, the fabric continues to strengthen and the process of polygomsation(also called rotation mrvstallisationj occurs. In this process new grains, with orientations similar to old grains which are not consumed by new nuclei, are produced. The average grain-size changes little with depth throughout this region,

  • (3) In the region directly over the bedrock, the shear stresses dominate the normal stresses, which combined with high temperatures initiates the process of migration(or dynamic) recrystallisation.New, strain-free, grains with their caxes at high angles (~45°) to compressional axes are created at the expense of old grains which are consumed by new nuclei. Abrupt changes in the average grain-size are usually observed in this region, but generally the grains are much larger than in region 2.

In order to construct a macroscopic constitutive law for anisotropic polycrystalline ice, a basic and physically motivated approach is to derive an average response of ice aggregate from the properties of individual grains and assumptions on crystal interactions, To date, most of the existing anisotropic models include only the mechanism of the grain c-axis rotation, which dominates the fabric formation and evolution in the upper part of an ice sheet, although some attempts to incorporate the ecrystallisation process, prevailing in deeper regions of ice sheets, have already been made. Reference AzuraaAzurma (1994) and Reference Azuma and Goto-AzumaAzuma and Goto-Azuma (1996) assume that individual crystals deform only by basal glide, the glide direction is determined by that of the maximum macroscopic shear stress in the polycrystal, and the crystal (microscopic) and polycrystal (macroscopic) stresses, which may be different in this model, arc related by a geometric tensor associated with the c-axis and glide directions. Reference LliboutryLliboutry (1993) assumes that the microscopic stress acting on each individual grain is equal to the bulk macroscopic stress applied to the polycrystal, and formulates a flow law for a transversely isotropic aggregate. Reference Castelnau and DuvalCastelnau and Duval (1994) and Van derVeen and Whillans (1994) extended Lliboutry’s homogeneous stress model (also called the static model) to any type of evolving anisotropy by considering a polycrystal consisting of a finite number of grains. Reference Van der Veen and WhillansVan der Veen and Whillans (1994) are the first to include the recrystallisation mechanism in their model, and consider two alternative criteria to determine the onset of recrystallisation in terms of threshold accumulated strains of crystals. Another, more general, approach is the visco-plastic self-consistcnt (VPSC) model developed by Reference Castelnau, Duval, Lebensohn and CanovaCastelnau and others (1996). In this formulation, the single crystal is treated as an embedded idealised geometric inclusion in an infinite medium with properties of an assumed form supposed to represent the macroscopic behaviour. In contrast to the models discussed above, in which individual crystals can glide only on basal planes, in the VPSG model crystal slips on basal, prismatic and pyramidal planes arc considered, and stresses and strain rates are assumed to depend on the crystal lographic orientation.

In all the above theories, which can be called discrete grain models, usually about 200—400 grains are needed to properly describe the ice fabric at a given point. Since in typical large-scale numerical models for ice-sheet flows the number of mesh nodes can exceed 105, it is obvious that discrete rain models are not suitable for simulations run on currently available computers. Therefore, in order to significantly reduce the number of variables involved in the description of ice fabric, another approach, in which the polycrystalline aggregate is treated as a continuum, has been adopted. In this approach a so-called orientation distribution function (ODF), defining continuous weightings to the grain C-axis orientation, has been applied. Although the eonccpt of the ODF is already well established in material science (Reference BungeBunge, 1982), it was first introduced to the field of theoretical glaciology only recently by Reference LliboutryLliboutry (1993), Reference Meyssonnier and PhilipMeyssonnier and Philip (1996) and Reference Svendsen and HutterSvendsen and Hotter (1996). Reference Meyssonnier and PhilipMeyssonnier and Philip 1996 formulated a transversely isotropic flow law for a polycrystalline ice using the ODF concept and the VPSC model by adopting some simplifying assumptions on the single grain behaviour, namely, that the crystal is transversely isotropic and offers small resistance to shearing parallel to its basal planes, and the response is governed by a linearly viscous law. Reference Svendsen and HutterSvendsen and Hutter (1996) employed the ODF approach to derive analytically a frame-indifferent viscous law which incorporates the fabric through a single structure tensor defined by an axis of assumed transverse isotropy. This theory has been considerably extended by Gödert and Hutter(1998, in press). However, the complicated calculations required to follow the evolving properties of individual ice elements will add considerably to numerical treatments ol large ice-sheet flows. A transversely isotropic flow law that avoids the use of an orientation distribution function has been proposed by Reference Van der Veen and WhillansVan der Veen and Whillans (1990). They modify Reference JohnsonJohnson’s (1977) law for a transversely’ isotropic viscoelastic solid by replacing material measures of stress and strain rate by spatial measures. However, they include the vertical (gravity) direction in the material structure, so this model cannot be treated as a valid constitutive relation for the response to general loading.

Dependence of the fabric on ice deformation implies that transverse isotropy of the medium can occur only if the flow in the ice sheet induces uniform stretches in some plane, say a horizontal plane in a gravity-dominated flow under a central dome, in which case a transversely isotropic fabric with the vertical as its rotational symmetry axis develops. Elsewhere in the ice sheet, however, such particular symmetry of the iceflow does not take place in general, and consquently the icedeformation, and hence the ice fabric, does not have rotational symmetry about any axis. Therefore, in order to more realistically describe the strain-induced ice fabric, a more general form of anisotropy is needed. Such increased generality is offered by orthotropic models. Although, strictly, the orthotropy does not occur for any arbitrary loading because of the fabric evolution which destroys any (already developed) orthotropic symmetries in the material, the assumption of orthotropy is strongly supported by observational data showing that the fabric in ice sheets does, in fact, usually reveal this type of material symmetry (this is in part due to the dynamic recrystallisation, not considered here, which helps to sustain the earlier created orthotropv despite the changes in the strain configuration).

In the paper, two orthotropic viscous models arc presented. The first model, formulated by Reference Gagliardini and MeyssonnierGagliardini and Meyssonnier (in press), stems from analysis of the behaviour of a single ice grain in a polycrystal under the assumption of stress homogeneity, following Reference LliboutryLliboutry (1993) and Reference Van der Veen and WhillansVan der Veen and Whillans (1994). Assuming transverse isotropy for each individual crystal and a linear response to deviatoric stresses, the macroscopic strain rate of the aggregate is derived by applying the ODF concept. Compared to Reference Meyssonnier and PhilipMeyssonnier and Philip (1996), Reference Svendsen and HutterSvendsen and Hutter (1996) and Reference Godert and HotterGodert and Hutter (1998), additional parameters arc introduced into the ODF. Owing to the extra parameters, the derivation of a macroscopically orthotropic (instead of transversely isotropic) viscous law relating the strain rate to the deviatoric stress and three structure tensors is possible.

The second orthotropic model considered here has been formulated by Reference Morland, Staroszczyk and StaroszczykMorland and Staroszczyk (1998) and further extended by Staroszczyk and Morland (in press). In this approach it is assumed that the macroscopic mechanical respouse of ice can be described in terms of the fabric induced purely by macroscopic deformation, and all microscopic processes occurring at the grain level are ignored. The assumption that the induced anisotropy depends only on the current strains and not on the deformation history is a considerable simplification since, in general, the fabric evolution is a path-dependent process. It is believed, however,that this approximation is the simplest approach to an evolving anisotropic viscous law that could be tractable in large-scale sheet dynamics, since it requires that only current deformation gradients are calculated in addition to the velocity and pressure fields. The constitutive relation is derived from the general frame-indifferent orthotropic representation (Reference Boehler and BoehlerBoehler, 1987), and expresses the dcviatoric stress in terms of the current strain rate, deformation and three structure tensors. The relation is separable in the isotropic dependence on strain rate and fabric dependence on deformation, and in its simplified form has only one independent fabric function that fully describes the orlhotropic viscous response of ice. Although in this approach local interactions between individual crystals are excluded from the analysis, this method allows good qualitative agreement with observations, and flexibility to correlate with detailed experimental results.

The two orthotropic models are used to determine the viscous response of icc to simple stress and strain configurations, corresponding to those occurring in the uniaxial compression and simple shear tests carried out in a laboratory. The predictions of both theories are compared, and the influence of some model parameters on calculated responses is investigated. Additionally, the results for simple shear given by the micro macroscopic model are compared with the results obtained from a discrete grain model in which no assumptions are made about material symmetries) in order to verify the validity of the assumed orthotropic behaviour of polycrystalline ice.

2. Micro-Macroscopic Model

This model, formulated by Reference Gagliardini and MeyssonnierGagliardini and Meyssonnier (in press), incorporates the basic micromechanism taking place on the grain level during the ice deformation, namely, the rotation of crystal caxes towards the axes of compression and away from the axes of extension. The macroscopic viscous outflow law which expresses the strain rates in terms of the deviatoric stresses is derived from the behaviour of individual grains by applying the homogenisation method based on the ODF approach,

In the following, three Cartesian reference frames are used to describe the behaviour of ice on the microscale of a single grain and the macroscale of a polycrystal:

{R} with axes xiis a fixed global reference frame;

{Ro}is a privileged frame of an orthotropic polycrystal, whose axes (xi)o coincide with the orthogonal privileged directions in the material;

{Rg} with axes (xi)g is a local frame associated with an individual grain, whose (x3)g axis coincides with the c axis of this grain.

Microscopic quantities associated with an individual grain are indicated by a tilde, and superscripts “g”and “o” are used to denote non-scalar quantities expressed, respectively, in the focal {Rg} and the privileged orthotropic {Ro} frames. Where no superscript is applied, respective symbols refer to macroscopic quantities expressed in the global coordinate system {R}. Since the individual ice grain is assumed here to be transversely isotropic, with its c axis being the rotational symmetry axis, the grain position relative to the global reference frame {R} can be uniquely described by means of two angles: the co-latitude for zenith angle) θ and the longitude Φ (see Fig. 1) These two angles determine the rotation matrix R. (cf. Reference Meyssonnier and PhilipMeyssonnier and Philip. 1996) which connects components of vectors and tensors in {R} and {Rg}.

Fig. 1. Global and local reference frames, with angles θ and Φ difining the c-axis orientation if a grain.

Following Reference LliboutryLliboutry (1993) and Reference Van der Veen and WhillansVan der Veen and Whillans (1994), the hypothesis of tfic stress homogeneity in ice aggregate has been adopted

(1)

stating that the microscopic stress in each grain, irrespective of the grain orientation, is equal to the macroscopic (bulk) stress σ applied to the polycrystal. The deviatoric stress σ´ is defined by

(2)

where pis the mean pressure, I is the unit tensor, and tr σdenotes the trace of σ. Since the ice is assumed incompressi ble, pis a workless constraint not given by a constitutive law, but determined by the momentum balance and boundary conditions.

Further, it is assumed that the transversely isotropic crystal deforms mainly by simple shear parallel Lo its basal plane,and its response to stress is linearly viscous. Adopting the simple relation (Reference Meyssonnier and PhilipMeyssonnier and Philip, 1996) between the microscopic strain rate D~ and the microscopic deviatoric stress which is equal to the macroscopic stress σ´ in view of Equation (1), the transversely isotropic How law can be expressed in the form

(3)

where M~=c ® c is the structure tensor defined in the global reference frame {R} by the unit vector c = (sinθcosΦ, sinθsinΦ, cosθ) associated with the grain caxis, and ψis the fluidity (reciprocal viscosity) for shearing parallel to the crystal basal plane. The parameter β is the ratio of the shear viscosity in a plane parallel to the caxis to the shear viscosity in a plane of isotropy (normal to the c axis) and can be regarded as a measure of the grain anisotropy. When β=0, the grain can deform only by basal glide, as assumed in the Reference LliboutryLliboutry (1993) and Reference Van der Veen and WhillansVan der Veen and Whillans (1994) models, while β=1 means that the grain is isotropic. In what follows it is also assumed that each grain in a polycrystal occupies the same volume and the number of grains does not change during the deformation, i.e. the grain growth, polygonisation and dynamic recrystallisation phenomena are not accounted for in this model.

The macroscopic (bulk) strain rateD of the polycrystal, in response to the macroscopic deviatoric stress σ´ is defined as the average of the strain rates of its constituent grains. Various homogenisation techniques can be employed to calculate the average strain rates in the ice aggregate. In the discrete grain models (Reference Azuma and Goto-AzumaAzuma, 1994; Reference Van der Veen and WhillansVan der Veen and Whillans, 1994; Reference Castelnau, Duval, Lebensohn and CanovaCastelnau and others, 1996) with a finite number of grains, the components of D are simply arithmetic means of the corresponding components of D. In our continuum model with an infinite number of grains, we make use of the ODFconcept, which describes the ice fabric in terms of the relative density of grains whose caxes have the orientation (θ,Φ)in the global reference frame {R}. In the ODF approach, the weighted average of a quantity A(θ,Φ) is defined by

(4)

where f(θ,Φ) is the proportion of grains with orientation (θ,Φ) in the element “area” θdθdΦ, and by definition

(5)

Note that in the case of isotropy we have f(θ,Φ) = 1, and for transverse isotropy, with x3 being the rotational symmetry axis, the ODFdoes not depend on Φ i.e. f(θ,Φ)== f(θ).

Hence, by the weighted average definition (Equation (4)), the components of the macroscopic strain rates in the global reference frame {R} are given by

(6)

Once the ODF has been determined, Equations (3) and (6) yield the relation between the macroscopic strain rate D and the deviatoric stress σ´In order to describe the evolution of the ODF, consider the motion of the single grain, i.e. the rotation of its caxis in the global coordinate system {R}. Following Reference Meyssonnier and PhilipMeyssonnier and Philip (1996), the rotation matrix R that defines t he orientation of the grain in {R} is governed by the relation

(7)

where RT is the transpose of R, the superposed dot denotes the time derivative, Wg and W~ are spins (rates of rotation) of the grain in the local and global reference frames, respectively, and cg = (0, 0,1). Since during the grain deformation the original parallel glide planes remain parallel to each other, the velocity component along the crystal caxis, when expressed in the rotating frame {Rg} attached to the grain, is a function ofXj only. This leads to the kinematic relations

(8)

which are a direct consequence of the adopted system of reference frames. In order to close the system of governing equations, three additional conditions arc required, and these arc adopted by assuming that the microscopic and macroscopic spins are equal, i.e.

(9)

The above equation is closely related to the Reference TaylorTaylor (1938) assumption which postulates the equality of micro- and macroscopic velocity gradients. The combination of constraints (Equations (1) and (9)) adopted here is exactly the same as that used by Reference Godert and HotterGodert and Hutter (1998). With Equations (8) and (9), Equation (7) provides two relations which describe the change in (lie grain orientation by

(10)
(11)

Since recrystallisation is not taken into consideration, and hence the total number of grains in the polycrystal is conserved, it follows from the continuity equation (5) that the ODF satisfies the relation

(12)

Equations (10-12) completely describe the evolution of ice fabric for any type of macroscopic anisotropy that is based on the assumed behaviour of ice on the microscale of individual crystals. In the case of orthotropic anisotropy, the medium possesses three planes of reflexional symmetry, which in the reference frame {R°} are the planes ((x1)o,(x2)o),((x1)o),(x3)o)and ((x2)o,(x3)o) these material symmetries must be accounted for in the ODF. Following Reference Meyssonnier and PhilipMeyssonnier and Philip (1996), and using analytical results obtained by Reference Gagliardini and MeyssonnierGagliardini and Meyssonnier (in press), we adopt the following form of the ODF;

(13)

where k1, k2,k3 and Φ° are parameters, the latter being the angle of rotation of the orthotropic frame { Ro} with respect to the global reference frame -{r}. Since the grain-conservation relation (Equation (5)) implies that k1, k2, k3 = 1,only three parameters in Equation (13) are independent. Now, by substituting Equation (3) into Equation (6), and using Equation (13) in Equation (4), we obtain the orthotropic, linearly viscous law relating the macroscopic strain-rates to the macroscopic deviatoric stresses by

(14)

where Mr o = er o ® er o arc three structure tensors defined by the orthotropic axes unit vectors er o, and Jr = tr(Mr oσ´)are three invariants. The six response coefficients αi(i=1,...,6) are functions of the grain rheological parameters ψ and β

(15)

and N30,N32,N50,N52,N54 are five moments defined by

(16)

In the case of ice isotropy, i.e. when f(θ,Φ) = 1 and kl = k2 = k3 = 1, the above moments are

(17)

and the six material coefficients defined by Equation (15) become

(18)

With the above parameters, and in view of the identity M1 o+M2 o+M3 o= I , the orthotropie constitutive law (Equation (14)) reduces to the linear Glen’s flow law for isotropic ice

(19)

where ψ0 is the fluidit y of the mac rosco pica lly isotropic p olycrysta l, which is related to the individual crys tal fluidity ψby

(20)

It follows from the la tter formula th at in the case of isotropic crystals (β = 1) the macroscopic fluidity ψ0 equals the grain fluidity Φ while in the case of the most ani so tropic crystals (β= 0), when the grain deform ation occ urs only by basal glide, ψ0= 0.4ψ.

3. Continuum model

In this continuum approach, proposed by Morl and Staroszczyk (1998) and further developed by Staroszczyk and Morland and (in press), the phenomena taking place on the microsca le of individual grains, as opposed to the model considered in section 2, are ignored , and it is assumed that the macroscopic response of ice depends only on the fabric induced solely by the mac roscopic deformations. Although this is a considerable simplification, it is beli eved th at such an approach to an eyoh-ing an iso tropic constitutive law is wel l suited to the large-scale ice-sheet modelling,since it req uires only that the deformation g radi ent field has to be determined during the ice flow in orde r to describe the fa bric development. The chosen (arm of the orth otropic viscous law is the relation between the frame-indifferent deviatoric Cauchy stress σ´, c urrent strain rate D , Cauchy-Green stra in tensor B, and three structure tensors Mr(r = 1,2,3) defined by the outer products of the current principal stretch axes unit vectors er(r= 1,2,3).

Let OXi(i = 1, 2,3) be spati al recta ng ular Cartesian coordinates with OXi(i = 1, 2,3) particle reference coordinates,and Vi the velocity components, then the deformation gradient F and strain rate D have the components

(21)

The strain B , unit vectors er(i=1,2,3)and squares of the principa l stretches br(r= 1, 2, 3) are defined by

(22)

where λr(r = 1, 2,3) are principal stretches a long the principal axes er, and the three structure tensors Mr are given by

(23)

Ice is assumed incompressible, so

(24)

In order to derive the macroscopic flow law for anisotropic ice, we follow Boehlcr’s (1987) theory for frame-indifferent relations between tensors and vectors, which ensures that material properties are independent of the observer.In the case or a symme tric tensor (here σ´) being a function of two oth er symmetric tensors (here D and B ), the general orthotropic representation is

(25)

where the 12 response coeffcients Φ(i = 1, ...., 12) are functions of the 19 invariants

(26)

subject to the constraints that the deviatoric Cauchy tress has zero trace, and the material is incompressible, so that only II coefficients Φi are independent, and only 18 invariants Ij are non-trivial, since I19 = 1. An alternative constitutive law for the strain rate D in terms or the deviatoric stress σ´ and the strain B , corresponding to the usual glaciological approach for the isotropic fluid model, can be formulated simila rly by interchanging σ´ and D in Equation (25) and (26).

We require that Equation (25) reduces to an isotropic visco us-fluid law

(27)

where Φ1, Φ2 depend on two invariants tr D2 and det D,when there is no fabric; that is, in the initial undeformed state F = I when the principal stretches are equal, necessarily λl = λ2 = λ3= 1 by the incompressibility co ndition (Equation (24)), or subsequently when F = I or when F corresponds to a rigid rotation of the ice clement. The conventional glaciological model is Φ2= 0 and Φ1 depends only on tr . The above prescription asserts that there is fabric - some alignment of initi ally randomly distributed c axes - only when there are differential principal stretches,i.e. when B≠I.

Following Staroszczyk and Morland (in press), we simplify the general law (Equation (25)) by ignoring tensor dependence on D2, B2 and DB, and the terms involving Mr without D or B, i.e. we restrict attention to a reduced model with a linear tensor dependence of σ´ on D and B. Hence, we adopt the relation

(28)

where the Φr+3 and Φ12 terms have been modi fi ed to recover zero trace, noting that the included scalar tr (MrD ) = Ir,and the scalar tr (DB) = I13+ I14 +I15. Equation (28) still represents a non-linearly viscous law since, in general, the response functions Φr+3 and Φ12 depend on D. Wc further express these fun ctions in sep a rable forms which factor out im'a riants depending only on the deformation B and retain a common dependence on inva riants involving the strain rate D:

(29)
(30)

Where

(31)

The response functions (Eq ua tions (29) and (30)) have to satisfy the isotropic fluid law (Equation (27)) when B=I and hence I21= I22= 3; thus

(32)

where the functions h and g are normalised by

(33)

We restrict attention to a simple model with the fabric response function h depending only on the principal stretches λr(through Ir+3= brr 2), and the function g depending only on the invariant measure of total deformation I21 = tr B; thus

(34)

where μ0= 1/Φ0 is the (consta nt) isotropic fluid viscosity (when bl = b2 = b3 = 1), and

(35)

Now, with Equations (34) and (32), the constitutive law (Equation (28)) expressed in terms of the two resp onse functions h and 9 and the viscosity μ0 takes the linear form

(36)

Staroszczyk and Morland (in press) derived equalities and inequ a lities which have to be sati sfi ed by instal1laneous viscositiesties ui)j (i, j = 1, 2,3; i =1= j ), depending on the relative magnitudes of the principal stretches. With the ordering b1≥b2≥b3 there are six distinct sets of relative values of (br)(r= 1, 2, 3), and for each of them corresponding relations order μ1213 and μ23 in the coordinate frame of the principal stretch axes. By using the viscosity relation corresponding to the pla ne flow, i. e. when λ2 = λ2 = 1, hence b3 = l/b1 and K = b1 + 1 + b1 -1, it is possible to express g(K ) in terms of h(br), namely,

(37)

where

(38)

The limit of Equation (37) as b1→1, K → 3, combined with the normalisation (34)4 gives

(39)

which is a restriction on h(br) at b=1. Rewriting the function g(K) as

(40)

where G(K) is bounded, it can be shown that the limit of Equation (37) as b1→α, K~b1, yields the relation

(41)

That is, in view of Equation (37), only one fabric response function, h(br), remains free for prescription, subject to Equation (39). In order to further restrict this function, the model defined by Equation (36) is used to predictthe viscous response of ice at large deformations in the axial compression and simple shear conditions. For such configurations, Reference Budd and JackaBudd and Jacka (1989) present experimental data obtained for a steady-state flow (when the microprocesses of grain rotation, polygonisation and dynamic recrystallisation balance one another) and determine the limit ratios off abric-induced viscosity to (maximum) isotropic viscosity. These empirical ratios, expressed in terms of so-called enhancement fctors for compression and shear, provide two further relations connecting h(α), h(0) and G(α)(Staroszczyk and MOrland, in press ). For the uniaxial compression the viscous law (Equation (36)) yields

(42)

where A is the reciprocal of the axial enhancement factor, and for the simple shear in the plane flow we obtain

(43)

where S is the reciprocal of the shear enhancement factor. Equations (4143) provide the following values of h(0), h(α) and G(α):

(44)

Despite several simplifying assumptions adopted to derive the model (36) with the single fabric function h(br) from the general orthotropic law (Equation (25)), this approach still retains considerable flexibility to correlate with observed data. In fact, only three specific restrictions are imposed on the function h: Equation (39) to yield a valid response in the isotropic state, and Equations (44) 1 and (44)2 to match the enhancement factors for compression and simple shearing at large strains. However, it is also required that the function h(br) yield a response which satisfies the viscosity relations derived in Straroszcyk and Morland (in press), otherwise physically invalid responses can be obtained.

4. Model Comparisons

The two orthotropic models described in sections 2 and 3 arc now applied to explore the viscous response of ice to simple deformation histories, corresponding to those taking place in typical laboratory tests. For the sake of brevity, throughout this section we will refer to the micro-macroscopic model as model A, and to the continuum model as model B.

First we consider an axial compression in the x2 direction, defined by the axial stretch λ<1, with equal lateral stretches λ1 and λ3 along the x1 and x3 coordinate axes, for which the deformation field is described by

(45)

and the associated velocities, in view of Equations (45) 4 and (45)5 are

(46)

With Equations (45) and (46), the Cauchy-Green strain tensor B and the strain-rate tensor D are given by

(47)

and the three structure tensors, due to the coincidence of the principal stretch axes er with the coordinate axes xr, are defined by

(48)

The deviatoric stresses (Equations (2)) are given by the diagonal tensor

(49)

The response of ice predicted by model A, defined by Equation(14), is illustrated in Figure 2, which shows the evolution of the normalised axial viscosity with increasing principal stretch λ1 for different values of the grain-anisotropy parameter β. The curves corresponding to β = 0.10, 0.15, 0.20, 0.25 are labelled Al, A2, A3, A4, respectively, and the same labelling applies in subsequent plots illustrating the results given by this model. It is seen from the figure that the micro-macroscopic model predicts very slight softening of ice (decrease in viscosity) during the first stage of uniaxial loading, for the stretches 1≤λ11.2 (0.7λ1). The softening stage is followed by a phase of considerable hardening of ice, particularly pronounced for the stretches in the range 1.3 ^A|^ 1.7 (0.35 ^ Aa ^ 0.6) and small values of the parameter β (i.e.for very anisotropic ice grains). As the deformation continues, all the grain c axes rotate towards the compression axis (here the x2 axis), and at large strains the macroscopic viscosity of a polycrystal approaches the axial viscosity of a single crystal. Since the latter viscosity is (1/β) larger than the viscosity for shearing on crystal basal planes, it follows from Equation (20) that the limit macroscopic viscosity for compression is given in terms of β by the relation

Fig. 2. Evolution qftlze ratio σ´ 22/(2μ0D22) with illmasing stretch λl increasing compression for different values of β in model A.

(50)

The viscous response of ice yielded by model B is determined by the function h(br). There is a variety of possible functions, as long as they satisfy some very general conditions (Staroszczyk and Morland, in press), and this gives the model ample flexibility to correlate with observations. For illustrations, we apply simple monotornic increasing functions

(51)
(52)

where h0 and hα=h(α) Equation (51) is a free parameter, and α in Equations (51) and (52) is determined by Equation (39). The limit values h(0) and h(α) are related through Equations (44) to A and S,, the reciprocalsof the enhancement factors. In Staroszczyk and Morland (in press) both A and S were chosen less than unity (both enhancement factors greater than unity), which corresponds to the values measured by Reference Budd and JackaBudd and Jacka (1989) for warm ice near melting. Here, in order to compare the predictions of the two orthotropic models, we adopt A > 1 (an enhancement factor for compression less than unity), meaning the increase in viscosity with increasing deformation, which according to Reference Pimienta, Duval and LipenkovPimienta and others (1987) is the case for the response of cold ice subjected to stress levels typically occurring in polar ice sheets. The value of A can be as high as 10 for a single-maximum fabric, as has been found experimentally by Reference Pimienta, Duval and LipenkovPimienta and others (1987), although recently Reference MangeneyMangeney and others (1996) calculated the value of about 3 for the ice near the bottom of the GRIP ice core, deduced from the data provided by Reference Thorsteinsson, Kipfstuhl and MillerThorsteinsson and others (1997). We carry out the simulations for two values of this factor, namely, A = 2.2, which is smaller, and A = 4.6, which is larger than the value obtained for the GRIP ice, though both orthotropic models are flexible in this respect and allow much larger values of A to be implemented. Since Expression (50), defining the asymptotic value of viscosity, is simply the reciprocal of the enhancement factor, A, it relates the parameter A in model B to the parameter β in model A. Hence, we find that the selected values A = 2.2 and A = 4.6 correspond to β=0.25 and β=0.10, respectively. Associated with the chosen values of A are the factors S=0.55 and S0.46= being the limit values of the normalised viscosity in simple shear, which correspond, respectively, to β=0.25 and β=0.10 in model A, as will be shown shortly in Figure 6. In view of Equations (44), we have the following connections:

(53)
(54)

The response functions hdetermined by the limit values in Equations (54) are illustrated in Figure 3, and very similar plots can be obtained for Equations (53). Henceforth, she curves labelled Bl, B2, B3 correspond to Equation (51) with m = 1,1.5,2, respectively, and the curves labelled B4 correspond to Equation (52).

Fig. 3. Fig. 3. Adopted forms of the fabric response function h(br) in model B.

With Equations (477#x2013;49), the constitutive law (Equation (36)) leads to the following relation describing the behaviour of ice in uniaxial compression:

(55)

where K = 2b1 + b −2 1 . The evolution of the axial viscosity d22/(2Φ0D22) with increasing stretch At for the adopted forms of the fabric response function h(br)is illustrated in Figure 4 for the parameters in Equation (53), and in Figure 5 for the parameters in Equations (54). We note that the inlluence of the function h(b r) on the predicted response of ice is most significant during the first phase of deformation. We see that the functions in Equation (51) with m = 1.5 and m = 2, curves B2 and B3, respectively, yield responses which agree quite well with the responses given by model A, curves A1 and A4, particularly for the set of parameters in Equations (54).

Fig. 4. Evolution of the ratio σ´ 22/(2μ0D22) with increasing stretch Al in uniaxial compressionJor dif.ferentJabric response functions h( br) in model B. The results Jor A = 4.6 and S = 0.46 are compared with the prediction of model AJor μ= 0.10.

Fig. 5. >Evolution of the ratio σ´ 22/(2μ0D22) with increasing stretch Al in uniaxial compressionJor differentfabric response Junctions h( br) in model B. The results Jor A = 2.2 and S = 0.55 are compared with the prediction of model A for β = 0.25.

Next we consider a simple shear at constant strain rate D12 = 1/2 γ which follows an initial plane compression and stretch that has been frozen at constant λ2 = λ1 -1 蠄 1 by the removal of the stress and strain rate. The deformation field is now described by

(56)

the corresponding velocities are

(57)

and the strain and strain-rate tensors arc defined by

(58)

The principal stretch squares bi (i=1,2,3), the eigenvalues of B, are given by

(59)

and the associated principal vectors er are determined by Equation (22) 2. In terms of these vectors, the three structure tensors are expressed by

(60)

The response of ice to simple shearing given by model A, defined by Equation (14), is illustrated in Figure 6, which shows the evolution of the normalised viscosity σ´ 12/(μ0 γ) With increasing shear k started from the isotropic state (λ12= 1) for different values of the grain-anisotropy parameter β We note that the model predicts initial hardening of ice (for strains K≲1), with a more significant increase in viscosity taking place for smaller values of β(i.e. for more anisotropic ice grains). As shearing continues, the initial hardening is followed by the softening phase, with a monotonic decrease in viscosity until an asymptotic value is reached at large strains, This limit value is given by

(61)

and follows from Equation (20), since at very large shear deformations all crystals are aligned for easy glide on basal planes I their caxes are approximately parallel to one another; and hence in the limit the macroscopic viscosity of a polycrystal approaches the viscosity of an individual grain.

In Figure 6 we also show the results predicted by a discrete grain model, with 800 grains, for β=0.10(curve Dl), and β=0.255 (curve D4). The discrete model is based on the same assumptions as model A, except that no restrictions on material symmetries are imposed, so il has a more general character than the orthotropic model. Therefore, comparisons between the predictions of both models should give some indication whether the assumption of orthotropy of ice is justified. The results obtained seem to support the validity of this assumption, since the maximum relative discrepancies between curves A1 and Di, and A4 and D4, do not exceed 15% for k ~ 1, and for larger shear strains (k≲5) die agreement between the results given by both models is very good.

The behaviour of ice predicted by model B, defined by

Fig. 6. Evolution of the ratioσ´ 12/(μ0 γ) with increasing strain k in simple shear started from an isotropic state (λ2= 1) for different values of βin model A. Also shown are the results given by the discrete grain model for β= 0.10 (curve Dl) and β=0.25(curve D4).

Equation (36) with Equations (5860), is described

(62)

where now k = b1+ b1 −1 + 1. Figures 7 and 8 show, for different fabric response functions (Equations (51) and (52)), the evolution of the viscosity ratio σ´ 12/(μ0γ) with shear _ started from the isotropic stare (λ1 = λ22 = 1). The same sets of parameters as in the uniaxial compression simulations are used, i.e. those in Equation (53) for Figure 7 and those in Equation (54) for Figure 8, and again the results yielded by model B are compared with the predictions of the micro-macroscopic model A (curves A1 and A4). It cats be observed that now, for simple shear, the discrepancies between the responses given by the two models are more considerable than in the case of uniaxial compression, especially for the parameters in Equation (53) (Fig. 7), corresponding to more anisotropic grains (smaller value of β) in model A, or “stiffer” ice (larger value of the viscosity factor A) in model B. For the parameters in Equation (54) (Fig. 8), a good qualitative agreement between the two models is still readied, although model A predicts more significant increase in ice viscosity during the initial phase of shearing. However, the discrepancies between the two models are still limited to relatively small (compared to those occurring in ice-sheet flows) shear strains, so it is anticipated that in realistic ice-sheet flow simu lations the two constitutive theories should yield similar predictions for most of an ice sheet, apart from a region near the ice divide, where (relatively) small shear deformations occur.

Fig. 7. Fig. 6. Evolution of the ratio σ´ 12/(μ0 γ) with increasing strain k in simple shear started from an isotropic state (λ2 = 1) for different response functions h(br) in model B. The resultsfor A = 4.6 and S = 0.46 are compared with the prediction qf model A jor β= 0.10.

Fig. 8. Evoluation of the ratio σ´ 12/(μ0 γ) with increasing strain r;, in simple shear started from an isotra/lic state (λ2 = 1) jor different fabric response functions h(br) in model B. The resultsfor A = 2.2 and S = 0.55 are compared with the prediction qf model A jor β= 0.25.

5. Conclusions

We have presented two orthotropic models for viscous response of ice based on two fundamentally different approaches. The first, a micro-macroscopic model, has been derived from the behaviour of individual ice grains and incorporates the mechanism of crystal rotation in response to shear stresses. The second, a macroscopic model, was constructed by applying the general theory of frame-indifferent(objective) relations for materials with orthotropic anisotropy. The results of numerical simulations for continued uniaxial compression and simple shearing have shown that for appropriately chosen model parameters a good correlation between the responses given by the two theories can be achieved. Both models predict significant hardening of ice in compression, and initial hardening followed by softening in simple shear, which seems to agree with the observed behaviour of cold ice at low stress levels. The theoretical results, however, should be verified against detailed experimental data; unfortunately, very few relevant data on cold ice are available, since practically all laboratory tests performed to date have been carried out on warm ice near melting and at relatively high deviatoric stresses. It is also anticipated that the numerical simulations of realistic ice-sheet flow’s and the comparison of results obtained in this way with in situ data should throw some light on the validity of the assumptions made in the paper: first of all, whether the assumption of orthotropic behaviour of ice is justified, or a more general form of anisotropic constitutive law is required.

Further work should concentrate on including in the models other micromechanisms, such as rotation and migration (dynamic) recrystallisation. This will, however, add considerably to the complexity of theoretical formulations. This is particularly true of the micro-macroscopic model, and, to the authors’ knowledge, very few attempts to incorporate recrystallisation have yet been made. In terms of extending its generality, the continuum model is more flexible, since only the modification of response functions is required, although in order to account for a more complex mechanical behaviour more tensor generators may be needed in the constitutive law. Finally, we would like to emphasise that in order lo construct a model that properly describes the behaviour of ice in polar ice sheets, more experimental work on cold ice needs to be done to provide reliable input data for the theory.

Acknowledgements

The authors would like to thank two anonymous referees for thorough reviews that have helped to improve the paper. Part, of this work was done while the first author was visiting the Laboratoire de Glaciologie, Saint-Martin-d’Heres, France, supported by EISMINT Exchange Grant No. 9709 funded by the European Science Foundation.

References

References

Alley, R. B. 1992. Flow-law hypotheses for ice-sheet modeling. J. Glaciol., 38 (129), 245256 Google Scholar
Azuraa, N 1994. A flow law for anisotropic ice and its application to icesheets. Earth Planet. Set, Lett. 128 3-4 601614 Google Scholar
Azuma, N. Goto-Azuma, K 1996. An anisotropic flow law For ice-sheet ice and its implications. Ann. Glaciol. 23, 202208 Google Scholar
Boehler, J. P. 1987. Representations for isotropic and anisotropic non-polynomial tensor functions In Boehler, J. P., ed Application of tensor functions in solid mechanics., 3153 (Course 292.) Berlin, etc.. Springer-Verlag. International Centre for the Mechanical.SciencesCrossRefGoogle Scholar
Budd, W.F. Jacka, T. H. 1989. A review of ice rheology for ice sheet modelling Cold Reg. Sci. Tecknol., 16,(2), 107144 Google Scholar
Bunge, H. J. 1982. Texture analysis in material science., London, Butterworths Google Scholar
Castelnau, O. Duval, P. 1994. Simulations of anisotropy and fabric development in polar ices. Ann. Glaciol., 20, 277282 Google Scholar
Castelnau, O. Duval, P. Lebensohn, R. Canova, G. R. 1996. Viscoplastic modeling of texture development in polyerystalline ice with a self-consistent approach: comparison with bound estimates. J. Geophys. Res.., 101, (B6), 13,85113,868 CrossRefGoogle Scholar
Gagliardini, O. Meyssonnier, J. In press. Analytical derivations lor the behaviour and fabric evolution of a linear orthotropic ice polycrystal. J.Geophys. Res.,Google Scholar
Godert, G. Hotter, K. 1998. Induced anisotropy in large ice shields: the ory and its homogenization. Continuum Mech. Thermodyn., 10, (5) 293318 Google Scholar
Godert, G. Hutter, In press.A kinematic model lor induced anisotropy in large ice shields. Proc R Soc London Google Scholar
Gow, A.J. Williamson, T. 1976. Rhcological implications of the internal structure and crystal fabrics of the West Antarctic ice sheet as revealed by deep core drilling at Byrd Station. Geol.Soc.Am.Bull., 87 (12) 16651677 Google Scholar
Herron, S. L. Langway, G. C. Jr 1982.A comparison of ice fabrics and textures at Gamp Century, Greenland and Byrd Station, Antarctica. Ann. Glaciol., 3, 118124 Google Scholar
Johnson, A. F. 1977. Creep characterization of transversely-isotropic metallic material. J.Mech. Phys. Solids., 25, 1117—;126 Google Scholar
Lipcnkov, V.Ya. Barkov, N. I. Duval, P. Pimienta, P. 1989. Crystalline texture of the 2083 m ice core at Vostok Station, Antarctica. J. Glaciol.., 35 (121) 392398 Google Scholar
Lliboutry, L. 1993. Anisotropic, transversely isotropic nonlinear viscosity of rock ice and rheological parameters inferred from homogenization. Int.J Plasticity., 9, (5), 619632 Google Scholar
Lliboutry, L Duval, P. 1985.Various isotropic and anisotropic ices found in glaciers and polar ice caps and their corresponding rheologies. Am. Geophys.., 3 (2) 207224 Google Scholar
Mangeney, A.F. Califano and O. Gastelnau. 1996. Isothermal flow of an anisotropic ice sheet in the vicinity of an ice divide. J.Geophys. Res., 101 (B12),28,18928,204 Google Scholar
Mangeney, A., Califano, F. Hutter, K. 1997. A numerical sutdy of anisotropic, low Reynolds number, free surface flow for ice sheet modeling. J.Geophys. Res.., 102, (B12) 22,74922,764 Google Scholar
Meyssonnier, J. Philip, A 1996. A model for the tangent viscous behaviour of anisotropic polar ice. Ann. Glacial 23,253261 CrossRefGoogle Scholar
Morland, L Staroszczyk, R Staroszczyk, 1998 .Viscous response of polar ice with evolving fabric. Continuum Mech.Themodyn., 10, (3), 135152 Google Scholar
Pimienta, P. Duval, P. Lipenkov, V.Ya 1987. Mechanical behavior of anisotropic polar ice. International Association of Hydrological Sciences Publication., 170 5766 (Symposium at Vancouver 1987 — The Physical Basis of Ice Sheet Modelling)Google Scholar
Russell-Head, D.S. Budd, W. F. 1979. Ice-sheet flow properties derived from bore-hole shear measurements combined with ice-core studies. J. GlacioL., 24, (90), 117130 Google Scholar
Staroszczyk, R. Morland, L. In press. Orthotropic viscous response of polar ice. J.Eng..Math.,Google Scholar
Svendsen, B. Hutter, K. 1996.A continuum approach for modelling induced anisotropy in glaciers and ice sheets. Ann. Glacial., 23, 262269 Google Scholar
Taylor, G.I 1938 Plastic strain in metals J.Inst. Met., 62 307324 Google Scholar
Thorsteinsson, Th. Kipfstuhl, J. Miller, 1997. Textures and fabrics in the GRIP ice core. J. Geophys. Res., 102 (Cl2), 26,58326,600 Google Scholar
Van der Veen, C.J. Whillans, I. M. 1990. Flow laws for glacier ice: comparison of numerical predictions and field measurements. J.Glaciol., 36 (124), 324339 Google Scholar
Van der Veen, C.J Whillans, I.M. 1994. Development of fabric in ice. Cold. Reg. Sci.Technol., 22,(2), 171195 Google Scholar
Figure 0

Fig. 1. Global and local reference frames, with angles θ and Φ difining the c-axis orientation if a grain.

Figure 1

Fig. 2. Evolution qftlze ratio σ´22/(2μ0D22) with illmasing stretch λl increasing compression for different values of β in model A.

Figure 2

Fig. 3. Fig. 3. Adopted forms of the fabric response function h(br) in model B.

Figure 3

Fig. 4. Evolution of the ratio σ´22/(2μ0D22) with increasingstretch Al in uniaxial compressionJor dif.ferentJabric response functions h( br) in model B. The results Jor A = 4.6 andS = 0.46 are compared with the prediction of model AJor μ= 0.10.

Figure 4

Fig. 5. >Evolution of the ratio σ´22/(2μ0D22) with increasing stretch Al in uniaxial compressionJor differentfabric response Junctions h( br) in model B. The results Jor A = 2.2 and S = 0.55 are compared with the prediction of model A for β = 0.25.

Figure 5

Fig. 6. Evolution of the ratioσ´12/(μ0 γ) with increasing strain k in simple shear started from an isotropic state (λ2= 1) for different values of βin model A. Also shown are the results given by the discrete grain model for β= 0.10 (curve Dl) and β=0.25(curve D4).

Figure 6

Fig. 7. Fig. 6. Evolution of the ratio σ´12/(μ0 γ) with increasing strain k in simple shear started from an isotropic state (λ2 = 1) for different response functions h(br) in model B. The resultsfor A = 4.6 and S = 0.46 are compared with the prediction qf model A jor β= 0.10.

Figure 7

Fig. 8. Evoluation of the ratio σ´12/(μ0 γ) with increasing strain r;, in simple shear started from an isotra/lic state (λ2 = 1) jor different fabric response functions h(br) in model B. The resultsfor A = 2.2 and S = 0.55 are compared with the prediction qf model A jor β= 0.25.