## Abstract

The region of sea ice near the edge of the sea ice pack is known as the marginal ice zone (MIZ), and its dynamics are complicated by ocean wave interaction with the ice cover, strong gradients in the atmosphere and ocean and variations in sea ice rheology. This paper focuses on the role of sea ice rheology in determining the dynamics of the MIZ. Here, sea ice is treated as a granular material with a composite rheology describing collisional ice floe interaction and plastic interaction. The collisional component of sea ice rheology depends upon the granular temperature, a measure of the kinetic energy of flow fluctuations. A simplified model of the MIZ is introduced consisting of the along and across momentum balance of the sea ice and the balance equation of fluctuation kinetic energy. The steady solution of these equations is found to leading order using elementary methods. This reveals a concentrated region of rapid ice flow parallel to the ice edge, which is in accordance with field observations, and previously called the *ice jet*. Previous explanations of the ice jet relied upon the existence of ocean currents beneath the ice cover. We show that an ice jet results as a natural consequence of the granular nature of sea ice.

## 1. Introduction

The marginal ice zone (MIZ) at the edge of the sea ice pack is approximately 5–100 km wide and marks the region in which incident wave energy plays an important role in sea ice dynamics. In the MIZ, sea ice floes have a smaller area than in the central pack, being typically 50–500 m wide, 0.1–5 m deep, and with an approximately convex polygonal shape, see figure 1 (taken from Squire 1998). The MIZ is often the boundary between temperate and polar climate systems, where strong horizontal and vertical gradients in the atmosphere and ocean may be present. Understanding MIZ dynamics is particularly relevant for climate modelling and ice forecasting for navigation. There are many complex interactions involved in the dynamics of the MIZ, including melting and freezing, and it is useful to focus on particular processes to improve fundamental understanding.

This paper addresses the role of sea ice rheology in the MIZ. In many climate models and regional studies, sea ice interaction forces are described using the viscous-plastic rheology introduced by Hibler (1979): sea ice is treated as an isotropic plastic material with a yield surface that can harden or weaken according to local ice thickness and concentration, with an associated, normal flow rule, and a viscous sub-yield response. This rheology was originally motivated by observations of sea ice deformation in the central pack, which is characterized by closely packed, interlocking floes, during the Arctic ice dynamics joint experiment, University of Washington, 1970–1978 (Pritchard 1980). As one approaches the ice edge from the central pack, however, the interlocking ice floes give way to a more discrete, dilute system, which was studied during the marginal ice zone experiment (MIZEX) campaigns of the 1980s (e.g. Squire 1998). In the MIZ, a consideration of collisions between individual floes can be used to describe the sea ice rheology as a viscous material, e.g. Shen *et al*. (1987). The plastic and viscous rheological behaviours can be thought of as the extreme behaviours of a granular material: sea ice motions in the central pack correspond to slow deformation of closely packed granules with granule interaction being mainly of plastic form (frictional rubbing at floe edges and overriding) and motions in the MIZ correspond to more rapid, collisional interactions rather similar to a two-dimensional gas, except that the collisions are dissipative. The discriminating feature between these two rheological behaviours, aside from ice concentration and floe size, is the degree of random motion of the floes that can lead to floe collisions. In the central pack, a given floe is typically in constant contact with its near neighbours. In deformation without shear, for example, after some initial rearrangement of the floes due to the peculiarity of the initial configuration of floe–floe contacts, the random component of floe velocity will be negligible compared with the mean velocity. (In reality, local variations of air and ocean drags would cause some small random motion of the floes.) With shear, significant random floe motion is restricted to relatively narrow shear bands with nearly constant (plug) flow elsewhere; see Kwok (2001) for a description of sea ice motion in the central pack inferred from satellite imagery. By contrast, in the MIZ, the motion of floes is relatively unconstrained by floe contacts, and the random component of floe velocity is comparatively large compared with the mean floe velocity. It is the random component of floe velocity that causes floes to collide in the MIZ. The magnitude of the fluctuation kinetic energy, i.e. the kinetic energy associated with the random motion of the floes, can be used to define the so-called *granular temperature*, e.g. Haff (1983). The granular temperature satisfies a balance equation and is a dependent variable upon which the constitutive behaviour of a granular material can depend. For a somewhat fuller description of granular flows, particularly appropriate to this paper, the reader is referred to Haff (1983) and Savage (1998).

The aim of this paper is to apply granular flow theory to the MIZ of sea ice using a composite sea ice rheology consisting of a sum of plastic and collisional (viscous) components. The relative importance of each rheological component depends upon the granular temperature and ice concentration. This approach builds upon two significant previous studies. Firstly, the role of a plastic sea ice rheology in the MIZ was investigated by Leppäranta & Hibler (1985). They considered a steady-state, one-dimensional MIZ and discovered plug flow in the sea ice. Secondly, Lu *et al*. (1989) presented calculations using the same geometry in which the collisional rheology of Shen *et al*. (1987; slightly modified to account for two different floe sizes) was used up to 100 km from the MIZ edge with a constant fluctuation velocity to simulate wave interaction with the ice cover (i.e. constant granular temperature) and the viscous-plastic rheology of Hibler (1979) was used beyond this. The main advances made in this paper are that both rheological behaviours are applied throughout the MIZ and their domain of significance calculated rather than imposed *a priori*; the collisional stresses are related to granular temperature, which is a dependent variable; and a useful approximate, analytical solution is obtained using elementary methods which promotes understanding and allows easy calculation of relevant quantities.

In §2, the momentum balance of sea ice, granular heat balance equation, and proposed plastic and collisional rheologies are described. In §3, these equations are applied to an idealized, steady MIZ, some necessary constitutive postulates are proposed for ocean drag, granular diffusion, dissipation and generation, and boundary conditions are introduced. In §4, the equations are scaled with typical magnitudes of variation and appropriate sea ice parameters used to recast the governing equations into non-dimensional form. In §5, perturbation theory and elementary methods are used to obtain the leading order inner and outer solutions (appropriate to the dynamic behaviour near and far from the ice edge, respectively) and a composite solution is constructed. The composite solution is used in §6 to calculate various features of interest, especially at the MIZ edge. In particular, the *ice jet* is introduced, discussed, and compared with observations. The sensitivity of the ice jet to model parameters is explored. In §7, the main conclusions are summarized.

## 2. Momentum balance of sea ice with composite rheology

We take the vertically integrated, two-dimensional momentum balance of the sea-ice–ocean mixture layer to be given by(2.1)as derived in some detail by Gray & Morland (1993). The mixture layer is defined to contain the sea ice floes and the ocean between the sea ice floes with the upper boundary being given by the continuous extension of the floe upper surfaces and the lower boundary given by the continuous extension of the lower floe surfaces, see figure 2. Strictly, the mixed layer also contains the air between sea level and the continuous extension of the floe upper surfaces, but this can be safely ignored in the following development because air density is approximately 1000 times smaller than the density of the ocean or sea ice. A similar form of the momentum balance was derived by Hibler (1979), which identifies the sea ice component only and is identical to that above when one ignores sea-ice–ocean thrusts and replaces mixture mass *m* per unit area with ice mass per unit area. The sea-ice–ocean mixture mass per unit area , where *ρ*, *ρ*_{o} are ice and ocean density, respectively, *h* is ice depth, *h*_{o} is the ocean mixture depth, *A* is the areal ice concentration (0<*A*<1), and *ρ*_{o}*h*_{o}=*ρh* follows by hydrostasy. The mean ice velocity is * u*≡(

*u*

_{1},

*u*

_{2})=(

*u*,

*v*), which is also taken to be the mean velocity of the ocean in the mixture layer (Gray & Morland 1993),

*f*

_{C}is the Coriolis parameter,

*is a unit vector normal to the ice surface directed into the atmosphere, and are drag forces applied to the sea ice layer (only) by the atmosphere and ocean, respectively (the tractions on the upper and lower surfaces of the ocean component of the mixture cancel; Gray & Morland 1993), is an ocean tilt force (which will be identically zero for the situation considered later), and*

**k***is the ice interaction stress tensor. The stress tensor is defined extra to the vertically integrated hydrostatic water pressure.*

**σ**Here, we neglect thermodynamic processes leading to melting and freezing of sea ice or formation of new ice. However, the ice thickness may increase through mechanical redistribution of ice during pressure ridging. During convergence, ice stresses initially increase as the ice concentration increases (owing to greater floe contact). As the ice concentration approaches unity, however, larger ice stresses are supported through pressure ridging. The ice floes first break into blocks through failure in flexion, then these ice blocks are forced up and down against gravity and friction to form the sail and keel of a pressure ridge. Ice may become heavily ridged, with increased mean ice thickness. Here, we adopt the simple ridging scheme used by Hibler (1979). For ice concentrations less than unity, ice thickness is constant and stresses increase solely through increased ice concentration; once the ice concentration reaches unity, further stress increase is maintained through pressure ridging and the ice thickness increases. More detailed discussion of mechanical redistribution schemes may be found in Rothrock (1975), Hibler (1985), Shinohara (1990) and Gray & Killworth (1996).

In addition to the momentum equations, in order to specify the ice interaction force, we require a description of the ‘jostling’ of the floes, since this will determine the contribution of collisions to the internal ice stress tensor. Within a horizontal representative area element used for averaging in the momentum balance, the velocity of a floe is given by its mean plus a rapidly fluctuating departure velocity * u*′. By definition, over a representative area element, the average of this departure velocity is zero, 〈

*′〉=*

**u****0**. The total kinetic energy per unit area of the sea-ice–ocean mixture is given by(2.2)where we have used the approximation that the ice and ocean departure velocities in the mixture are equal in magnitude (any refinement of this would require a detailed model of the ocean flow in the mixture layer, together with a model of mass exchange with the upper mixed layer of the ocean). By analogy with the jostling motion of molecules in a gas, the departure velocity has been associated with the granular temperature (e.g. Haff 1983); here we define the granular temperature as(2.3)

Following Savage (1998), although he did not specifically consider sea ice, the rate of change of kinetic energy per unit area of the sea-ice–ocean mixture is given by(2.4)where * q* is the fluctuation kinetic energy flux vector, ∇.(

*.*

**σ***) is the rate at which internal stresses do work,*

**u***Γ*is the rate at which fluctuation kinetic energy is dissipated through frictional rubbing between floe edges and through imperfectly elastic collisions between floes, and the last term is the rate at which body forces do work, where is the body force. (The air and ocean drags are body forces here because we consider the vertically integrated momentum balance.) The rate of work by body forces is(2.5)with zero energy contribution from the Coriolis force and by definition ′=

**0**, . We have assumed that 〈

*A*′〉=0, since the fluctuation flux of floes leaving a given region will typically be balanced by those entering a given region as 〈

*′〉=*

**u****0**(this may not hold true if there is a gradient of floe sizes across the region; we do not consider this possibility in the following). We thus have(2.6)where the atmospheric and oceanic drags, and , are defined to be the mean average over the representative area element, and is a fluctuation kinetic energy source that arises through the positive correlation of transient turbulent air (and ocean) drag with the departure velocity within the representative area element. The turbulent energy source, which does not result from Reynolds averaging, is likely to be largest near the edge of the MIZ, where turbulent drags are produced during the encounter of air and ocean flows with the bluff sea ice layer, but could also be particularly important where the atmospheric boundary layer is unstably stratified. Subtracting from equation (2.6) the inner product of the momentum equation with

**, and making use of the tensor relation , we obtain an equation determining the evolution of the granular temperature:(2.7)**

*u*For the sake of completeness and concreteness, we present the plastic and collisional rheologies in full, although it will be seen that we only require shear stress and isotropic stress with zero dilation. Both of these rheologies are of reduced Reiner–Rivlin form. We only require the plastic portion of Hibler's (1979) viscous-plastic rheology, which is(2.8)where(2.9)(2.10)and(2.11)

The tensor quantity is the strain rate, *η*^{F} and *ζ*^{F} are shear and areal coefficients of viscosity, *P*^{F} is an ice strength term and *δ*_{ij} is the Kronecker delta function. The viscosities and ice strength are nonlinear functions of strain rate and scalars and are so defined that the stress lies on an elliptical yield curve in invariant stress space with eccentricity *e*, with position on the yield surface being determined by the strain rate and associated normal flow rule. The yield surface hardens (weakens) as the ice concentration and thickness increase (decrease). The parameter *P*^{*} is a constant, and the weighting function *g*(*A*) is closely related to the function introduced by Gray & Morland (1993) that describes the mean proportion of ice–ice contact length on floe boundaries. We shall use Hibler's expression *g*(*A*)=e^{−C(1−A)}, where *C* is a tuneable parameter typically set to *C*=20, and note that 0<*g*<1.

The collisional contribution to stress we take to be given by the analytical rheology of Shen *et al*. (1987), which is(2.12)where(2.13)and(2.14)where *L*_{f} is the floe diameter, *e*′ is the coefficient of restitution, and we have replaced the fluctuation velocity magnitude with . Shen *et al*. (1987) derived this rheology for the collision of identical cylindrical floes and the parameter *A*_{max} corresponds to the maximum ice concentration possible with close packing. Lu *et al*. (1989) considered a multiple size distribution, but the calculations showed little sensitivity to this refinement. Here, we further assume the same form of the rheology with irregularly shaped floes (where *L*_{f} is now floe span) and set *A*_{max}=1. In deriving this rheology, Shen *et al*. (1987) made the assumption that . In the calculations that follow, this requirement is easily met wherever the collisional rheology makes a non-negligible contribution to stress.

We consider a composite rheology with collisional and plastic stress contributions that allow us to gauge the relative importance of the contributions to stress for the boundary-layer flow described below. We make the constitutive postulate that the composite stress is given by the simple, linear sum of collisional and plastic components,(2.15)so that, for example, *η*=*η*^{C}+*η*^{F}. Such a postulate has been used successfully in earlier studies of granular flow, e.g. Johnson & Jackson (1987), but has so far not been justified on the basis of experiments or analysis.

## 3. Boundary-layer flow

We now apply the above evolution equations and rheologies to the geometry depicted in figure 3, which is a simplification of the MIZ. We use an orthogonal, Cartesian coordinate system (*x*, *y*); open ocean occupies the half-space *x*<0 and sea ice occupies *x*>0. In the tensor descriptions, the positive *x*-direction is labelled ‘1’ and the positive *y*-direction is labelled ‘2’ and we have {*u*_{1}, *u*_{2}}={*u*, *v*}. We consider the steady situation that arises when on-ice air drag forces compress the ice cover until they are opposed by ice forces and a stationary ice edge is obtained. (With an off-ice air drag force, the edge of the MIZ would diverge.) We assume no variation in the *y*-direction, and, since *u*=0, we have and . For simplicity, and in common with other treatments of the MIZ, e.g. Leppäranta & Hibler (1985), we additionally neglect ocean currents that imply that the sea surface tilt force is zero. Although this assumption is made for mathematical convenience, scaling estimates and measurements also show it to be a good approximation, e.g. Steele *et al*. (1997). With these assumptions, the momentum balance is(3.1)and(3.2)and the granular temperature equation is(3.3)

The energy flux *q*_{1} transports kinetic energy through collisions and the random motion of floes and, in common with treatments of heat and mass transport, this transport is parametrized as diffusion,(3.4)where following Haff (1983) the diffusion coefficient is given by(3.5)where *s* is the distance separating floes, and *λ* is an (1) dimensionless number. A feature of this type of nonlinear diffusion, which is additionally used in modelling porous media flows and population dispersion, is that the diffusion problem exhibits compact support, i.e. the solution is identically zero outside of some specified domain, e.g. Grindrod (1996). We shall see this feature in the solutions derived later. Fluctuation kinetic energy is dissipated through frictional rubbing between floes and through imperfectly elastic collisions between floes. Following Savage (1998), the dissipation from frictional rubbing is given by(3.6)where *Γ*′_{r}∼0.25*P*^{*}/*L*_{f} and *gP*^{*} is a measure of plastic strength. Dissipation from imperfectly elastic collisions is given by Haff (1983) to be(3.7)(*r* is an (1) dimensionless number). The total dissipation may be expressed as the sum of these separate contributions,(3.8)

To use these parametrizations for diffusion coefficient and dissipation, we require an expression for the distance between floes. On geometrical grounds, we choose(3.9)

Note that this expression satisfies *s*→0 as *A*→1 and *s*→∞ as *A*→0. The source of fluctuation kinetic energy due to correlation of the turbulent air and ocean drags with ice floe departure velocity is typically highly variable. Here, we adopt the following approximation(3.10)where *F*′ is a constant whose value characterizes the strength of the coupling of the turbulent air/ocean drag fluctuations with the floe fluctuation velocity.

The ocean and atmosphere drag forces are routinely parametrized using quadratic boundary‐layer laws. Throughout this paper, we treat the mean atmospheric drag as a forcing parameter. The mean ocean drag is(3.11)where *C*_{o} is dimensionless drag coefficient, is geostrophic ocean velocity, and *θ* is the turning angle. Since we neglect ocean currents, .

Substituting the expressions for the stress tensor, energy flux, dissipation, turbulent energy source and ocean drag into the evolution equations (3.1)–(3.3), we obtain(3.12)(3.13)and(3.14)

We now discuss boundary conditions. At the open ocean interface we require continuity of normal traction, and, since the ice stress tensor is defined *extra* to the hydrostatic stress, we write this as and . The non-hydrostatic compressive wave stress arises from ocean waves transferring a fraction of their momentum to the ice cover as they enter the ice cover. The plastic hydrostatic pressure, *P*^{F}, gives rise to an unrealistic divergent stress discontinuity at an open ice edge that will cause the ice cover to diverge in the absence of any applied forces. This violates the thermodynamic requirement that non-negative work is done by the ice stresses and is now well recognized and treated in numerical simulations, e.g. Ip (1993). Since the ice edge is steady, we subtract a constant isotropic stress from the plastic stress so that and thus at the edge since we consider the case . This does not affect the sea ice momentum balance, which depends only on the divergence of stress. Continuity of normal traction at the interface yields(3.15)and(3.16)

Formally, the first boundary condition provides a constraint on *T*(0)=*T*^{*} and *A*(0), which enters through *γ*, and the second provides a constraint on (∂*v*/∂*x*)(0), *T*(0) and *A*(0). Conservation of granular heat across the MIZ edge provides a further constraint upon (*k*∂*T*/∂*x*)(0), which completes the specification of the boundary conditions. However, both and are difficult to determine from observations or scaling analysis, and the manner in which wave energy is transferred into granular heat at the interface is both poorly understood and unmeasured. Observations of departure velocities (Martin & Becker 1987), ice velocity profiles using buoy data (Johannessen *et al*. 1983) and ice concentrations enable us to prescribe *T*(0), (∂*v*/∂*x*)(0) and *A*(0) directly, and this is sufficient to constrain the problem. We thus impose(3.17)

The behaviour far from the ice edge is determined by the momentum and granular energy balances and is calculated in the next two sections.

## 4. Scaling analysis

We begin by substituting the expressions for the viscosities and ice strength into the evolution equations (3.12)–(3.14) and non-dimensionalizing by using a length-scale *L*^{*}, a granular temperature *T*^{*} a velocity *v*^{*} and ice thickness *h*^{*}. We obtain(4.1)(4.2)and(4.3)where the variables {*v*(*x*), *T*(*x*), *A*(*x*), *h*(*x*)} are dimensionless and (1),(4.4)are (1) in the range of interest where *g*^{*}=*g*(*A*^{*}), *s*^{*}=(*A*^{*}^{−1/2}−1) and *A*^{*}=0.8, , and the non-dimensional constants − are (1) and given by(4.5)(4.6)(4.7)(4.8)where *ϵ* is much less than unity.

We chose the length-scale *L*^{*}=100 km as this is the typical length-scale ove which the observed mean sea ice velocity is observed to vary by its own magnitude, and chose typical sea ice values *v*^{*}=0.1 m s^{−1} and *h*^{*}=1 m (for example, Gray & Morland 1993). The granular temperature *T*^{*} is chosen to be the value at the edge of the MIZ, and provides an upper bound; we choose a value of 10^{−2} m^{2} s^{−2} which corresponds to equal mean and variance of ice velocity at the MIZ edge. This is a rather high value but later calculations show that the mean value near the MIZ edge is in fair agreement with observations made during the Greenland MIZEX campaign (Martin & Becker 1987). The chosen length-scale also corresponds to the typical grid-spacing used in basin-scale simulations of the sea ice cover. Further than one grid spacing away from the MIZ edge, we expect the plastic rheology to dominate and thus the effect of collisional rheology on sea ice dynamics constitutes a ‘sub-grid’ phenomena.

Appropriate parameter values for the MIZ are given in table 1 (for example, Hibler 1979; Haff 1983; Shen *et al*. 1987). From Savage (1998), The value of *F*′ is difficult to determine, but since far from the MIZ edge there will be a balance between generation and dissipation, we set *F*′≈*g*^{*}*Γ*′_{r}. Using these values and the typical magnitudes for the dependent variables, we can calculate the values of the non-dimensional parameters −, also shown in table 1, where we have taken *ϵ*=10^{−2}.

## 5. Asymptotic solution

To identify the regions in which the various terms in equations (4.1)–(4.3) play an important role in determining the solution, we introduce the non-dimensional coordinate stretching *x*=*ϵ*^{α}*η*, where we consider *α*≥0. Inspection of the resulting equations reveals the existence of two distinguished balances of terms that correspond to coordinate scalings of *α*=0, the outer, and *α*=1, the inner. Following the usual procedure, we determine the solutions in the outer and inner regions separately, and then apply matching conditions that determine the unknown constants of integration (which, more generally, may be functions). We then construct the simplest matched asymptotic solution, which is the union of the inner and outer solutions minus their interception.

### (a) Outer solution

We first consider the outer solution. We set *α*=0 and pose the asymptotic solution . To leading order, from (4.1) to (4.3), we have(5.1)(5.2)(5.3)where and . The first two equations describe the across and along MIZ momentum balances where collisional stresses are unimportant, and the last equation states that turbulent energy gained by the floes through correlation of fluctuation velocity with air/ocean drag fluctuations (third term) is dissipated through frictional rubbing and collisions (first and second terms). These equations may be solved in inverse, integral form analytically. We eliminate from (5.1) using (5.2) to reveal two possible quadratics in ,(5.4)with solutions(5.5)where(5.6)and the only physical solution is typically given by taking both positive signs in the right-hand side of (5.5) (otherwise *V*<0, which is inconsistent). Substituting this solution for velocity as a function of ice concentration into (5.2), separating the (ordinary) derivative and integrating yields the solution for ice concentration and ice thickness(5.7)where(5.8)may be inverted numerically, and *Ω* is a constant of integration that we later determine by matching with the inner solution to be *Ω*=*A*_{edge} so that(5.9)is the width of the region in which pressure ridging does not occur. From (5.3), the granular temperature is found algebraically to be zero or(5.10)and we construct the composite solution(5.11)since granular temperature cannot be negative and may be non-zero if there is a balance between the turbulent source and dissipation. Note that as .

### (b) Inner solution

We now proceed to the inner solution. Setting *α*=1, we pose asymptotic solutions of the form given above, with the superscript ‘O’ replaced by ‘I’. To leading order, we have(5.12)(5.13)(5.14)where , and . In this last equation, the first term is the contribution to granular heat from work done by the mean plastic stress, the second term describes the diffusion of granular heat, the third and fourth terms describe dissipation through frictional rubbing and collisions, respectively, and the fifth term is the turbulent source. Equations (5.12) and (5.13) immediately imply(5.15)and(5.16)where *β* is a constant of integration to be discussed later. Note that at leading order the plastic stress does not enter the momentum balance, but does enter the granular heat balance. Application of the ice concentration boundary condition at the ice edge (3.17) gives(5.17)

Substituting (5.16) into (5.14) and rearranging, we obtain(5.18)where(5.19)are positive constants. Equation (5.18) is nonlinear, but we are able to determine an analytical, inverse solution as follows. First, we repose the problem in terms of an inverse function *Ψ*,(5.20)where, since we expect to be strictly monotonically decreasing, is uniquely defined. After manipulations, (5.18) transforms into(5.21)where . Noting that and using the substitution , it may be shown that(5.22)solves equation (5.21), where *ω*_{0} is an arbitrary constant of integration. The indefinite integral of (5.22) is(5.23)where both negative signs are taken in (5.22) (see §6),(5.24)the constants {*X*_{1}, *X*_{2}, *X*_{3}} are the roots of the polynomial(5.25)*E*_{1}(*n*,*ϕ*) is the elliptic integral of the first kind, and *E*_{3}(*n*,*ϕ*,*m*) is the incomplete elliptic integral of the third kind (see Abramowitz & Stegun 1965). The implicit solution for is then , where the constant *ω*_{1} is determined by the requirement that(5.26)so that the implicit solution is(5.27)

Turning our attention to the velocity solution , we note that(5.28)so that we may re-write (5.16) as(5.29)

The right-hand side may be integrated to yield the indefinite integral(5.30)with the velocity solution given by(5.31)where *ω*_{2} is the constant of integration. It is straightforward to determine from using (5.27).

### (c) Matched solution

Matching of the inner and outer solutions requires that(5.32)where *η* is fixed in the limits on the left-hand side and *x* is fixed in the limits on the right-hand side (for example, Hinch 1994). These matching conditions are now used to determine the constants of integration in the inner solution and *Ω* in the outer solution.

Substitution of the inner and outer ice concentration solutions, (5.8) and (5.17) into (5.32), determines(5.33)so that the outer solution is now fully specified.

Since by assertion is strictly monotonically decreasing in , substituting for the outer temperature solution (5.11) in (5.32) yields(5.34)which implies(5.35)

This condition applied to equation (5.22) determines the constant of integration(5.36)and justifies taking the negative signs in (5.22). Note that specification of *ω*_{0} determines *Ψ*′(1) and hence , which, combined with *Ψ*(1)=0 or , allows either (5.21) or (5.18) to be solved numerically by backward or forward integration (respectively).

Substitution of the outer velocity solution (5.5) into (5.32) yields(5.37)where we have used monotonicity of , and the above matching conditions on temperature and ice concentration. Application of this condition to (5.31) determines the integration constant(5.38)

The composite approximations, uniformly valid in the inner and outer regions, may be constructed in the usual way (e.g. Hinch 1994) to be(5.39)(5.40)and ice concentration and thickness are given by (5.7) with the outer superscript removed.

## 6. Calculations and discussion

During the Norwegian remote sensing experiment north of Svalbard in 1979, Johannessen *et al*. (1983) made various observations in the MIZ. Their buoy measurements, and subsequent observations by many other researchers, revealed the existence of an ice jet under conditions of down ice winds with a maximum ice velocity at the edge of the MIZ. For direct comparison with observations, we plot the composite asymptotic solution developed here, (5.7), (5.39) and (5.40), in dimensional form. Where ice jets are commonly observed, *A*_{edge}≈0.7 and the ice velocity decreases by 0.1 m s^{−1} over 1–5 km, which implies from equation (5.16) that(6.1)

In figure 4, we show plots of the granular temperature, ice velocity and ice concentration for *β*=−0.25 and the sea ice parameter values introduced in §4. The plots on the left highlight the behaviour in the vicinity of the ice edge; the plots on the right show the outer behaviour without pressure ridging. The most striking feature is the sharp discontinuity in the plots of granular temperature and ice velocity. This is a direct result of the nonlinear diffusivity of granular heat and well describes the narrowness of the mobile layer of granules (ice floes). (A diffusivity which is not a function of temperature will give the sort of smooth profile more familiar in diffusion problems.) The discontinuity marks the extent of the ice jet and the transition to the outer solution. Within the ice jet, a typical value of granular temperature is *T*≈0.005 m^{2} s^{−2} in accordance with the observations of Martin & Becker (1987) made during MIZEX. Further from the ice edge, the ice velocity is slowly increasing, reaching a constant value where the ice concentration reaches unity. Beyond this point (not shown) pressure ridging causes the ice thickness to increase linearly in response to the constant atmospheric drag.

Comparison of the velocity plot with figure 5 (taken from Johannessen *et al*. 1983) showing the velocities of a series of buoys 50 m, 5 km, 50 km and 250 km from the edge of the MIZ reveals a good qualitative and reasonable quantitative match of our analytical solution with observations. This is especially true for figure 5*b*, where the wind direction ensures along-ice and on-ice air drag forces. The rapid drop in buoy velocity between 50 m and 5 km defines the extent of the ice jet as being 5 km or less in width, which supports our analytical solution. From equation (5.27), we define the width of the ice jet to be(6.2)and we define the intensity of the jet as(6.3)

These two equations provide explicit, closed-form relationships between ice jet width and intensity as functions of the controlling parameters. A relatively simple expression for the ice jet width can be derived with the assumption that so that the integral of (5.22) does not involve elliptical integrals. This expression is(6.4)and is found to be an excellent approximation in the parameter range investigated (accurate to within approximately 5%).

In figure 6*a–c*, we plot the ice jet width, intensity, and product of width and intensity against the parameter *β*, which is proportional to ice velocity shear at the edge of the MIZ. The greater the shear at the edge of the MIZ (*β* is more negative), the narrower the ice jet, and this dependence is monotonic. The velocity profiles are of similar shape for all realistic values of *β*, and the product *WI* is approximately twice the extra mass flux in the ice jet. This extra mass flux is seen to monotonically increase as the magnitude of ice shear at the ice edge increases.

Johannessen *et al*. (1983) used the free-drift relation between ice velocity and geostrophic wind with no ice stresses presented in Thorndike & Colony (1982) to calculate the effective drag coefficients and turning angles in the atmospheric drag stress **τ**_{a} required to explain the buoy observations. Since the effective drag coefficient was larger in the ice jet than far from the MIZ edge, they suggested that the ice jet was the result of a stronger coupling between the wind and ice drift at the ice edge. However, observational estimates of drag coefficients revealed that the increase of air drag required to account for the ice jet was implausibly large. Johannessen *et al*. suggested that another possible explanation of the ice jet would be the existence of an ocean current shear beneath the ice cover associated with wind-driven upwelling and concluded that such a current, combined with increased wind-ice drift coupling, could explain the existence of an ice jet. The analysis presented here shows that the ice jet naturally results from the granular nature of the ice cover: at the edge of the MIZ, ice floe collisions are important and shear stress is proportional to velocity shear so that the ice jet results; further from the MIZ edge, ice floe interaction becomes mainly plastic, so that stress is independent of strain rate magnitude and the velocity tends to a constant (plug flow).

In figure 6*d–f*, we set *β*=−0.25 and investigate the dependence of ice jet width on the controlling parameters , and . The magnitude of the parameter reflects the relative strength of diffusive change of fluctuation kinetic energy to frictional dissipation; as expected, as diffusion dominates, the ice jet width increases and, in fact, exactly. The parameter reflects the relative strength of collisional dissipation to frictional dissipation; for example, for constant frictional dissipation, the ice jet width sharply increases as decreases below unity. The parameter reflects the relative strength of fluctuation kinetic energy generation through local correlations of turbulent air and wind flows with ice motion to frictional dissipation; as expected, as dissipation exceeds generation, the narrower the width of the ice jet.

The inner solution describing the ice jet depends upon the air drag only through the requirement of matching with the outer solution (provided the air drag does not reach a size of (*ϵ*^{−1}), which would be unrealistic). However, the direction and intensity of the jet also strongly depends upon the ocean shear stress at the edge of the MIZ, (this determines the value of *β*), and so indirectly depends on the air drag driving the ocean shear at the edge of the MIZ. We have chosen typical on-ice and along-ice edge air drags. If the air drag becomes off-ice, the ice cover would diverge so that, within a day or so, the ice concentration becomes small and ice floe interaction forces become negligible in the momentum balance. In this case, our model would predict no ice jet. The inner solution is also to leading order independent of Coriolis forces since these are rather small in the ice jet.

In figure 7*a*,*b*, we plot the plastic and collisional contributions to compressive sea ice stress (solid line) and (dashed line), and contributions to shear stress (solid line) and (dashed line), for the same parameters as figure 4. Except for a region of the order of five floe widths, the plastic interaction between the floes dominates the compressive stress. Conversely, and only within the MIZ ice jet, the collisional contributions to shear stress exceed the plastic contribution. Moving away from the MIZ edge within the jet, the collisional shear stress contribution is approximately constant, then decreases quickly to a negligible value as the effect of the decrease in granular temperature dominates. Beyond the ice jet (further than about 5 km), the collisional contribution to stress is negligible, and the plastic compressive and shear ice stress continue to increase monotonically as ice concentration increases. Note that the plastic shear stress changes sign at the edge of the ice jet; this is because the velocity gradient changes sign here. The ice velocity slowly increases beyond the ice jet in response to the along-ice air drag. These stress calculations suggest that it is reasonable to apply zero normal traction at the sea ice edge in a model of sea ice that ignores collisional rheology, which is consistent with current practice in numerical modelling.

The results of the analysis presented here differ somewhat from the numerical calculations of Lu *et al*. (1989) who used a collisional rheology, based on that determined by Shen *et al*. (1987), up to 100 km from the MIZ edge, and the viscous-plastic rheology of Hibler (1979) beyond this. In the simulations most closely corresponding to the situation considered here, Lu *et al*. (1989) found a significant velocity shear in the whole of the region in which the collisional rheology was used, which one would interpret as an unrealistic ice jet width of 100 km. This almost certainly results from their assumption of a constant granular temperature in this domain. Also, our calculations show that the plastic ice stress dominates throughout the MIZ except in the ice jet where the collisional shear stress dominates, necessarily in contrast to Lu *et al*. (1989).

## 7. Concluding remarks

This paper has applied granular flow theory to an idealized marginal ice zone of sea ice. The model consists of the along and across momentum balance of the sea-ice–ocean mixture, a balance equation of fluctuation kinetic energy, and various constitutive postulates including a composite collisional-plastic sea ice rheology. The rather complicated governing equations were scaled and an approximate leading‐order solution obtained analytically. This solution revealed an ice jet, an observed region of strong ice flow at, and parallel to, the MIZ edge. Previous explanations of the ice jet relied on strong ocean currents beneath the sea ice; the calculations presented here reveal an ice jet is a natural result of the granular nature of sea ice. Sensitivity to model parameters was explored.

Several refinements of the approach presented here are feasible, e.g. modifications of the mixture relation for the composite sea ice rheology, calculation of the higher‐order correction terms to the asymptotic solution, or even a more explicit calculation of the ocean wave interaction with the sea-ice–ocean mixture layer, although this would be hard, and previous attempts have been inconclusive (e.g. Squire 1998). An interesting refinement might be to consider an anisotropic granular temperature. This might be relevant because the main source of random motion of the floes at the MIZ edge is interaction with the wave field and the radiation of wave power might be expected to be stronger perpendicular than parallel to the MIZ edge due to the damping effect of the sea ice field. This would require an extra balance equation for fluctuation energy, but should not affect the main conclusions.

## Acknowledgments

I thank Dr Wilchinsky, Professor Hutter and two anonymous referees for useful comments on an earlier draft.

## Footnotes

One contribution of 11 to a Theme ‘Geophysical granular and particle-laden flows’.

- © 2005 The Royal Society