3.1.1. Equations

The melting-refreezing process is ignored. The ice is assumed to be a continuous, homogeneous, isotropic, incompressible, isothermal medium. Disregarding the effects of gravity and inertia, the equations of equilibrium and mass conservation are:

where σ_ij are the stress tensor components and u is the velocity vector.

3.1.2. Constitutive law for ice

The assumed constitutive law is Norton-Hoffs power law which relates the strain-rates έ _ij to the deviatoric stresses σ’_ij:

where τ is the effective shear stress defined by: This constitutive law leads to the relation known as Glen’s flow law (Reference GlenGlen, 1955): where is the effective shear strain-rate given by

The steady-creep law, equation (3), was verified by Reference DuvalDuval (1976) in torsion-compression tests on natural ices (n was found equal to 3 for τ in the range 0.1–0.4 MPa). It was used in our computations, although when passing over a bump a material particle does not undergo steady-creep conditions (Reference PatersonPaterson, 1981). The isothermal ice hypothesis led to the use of a constant for B.

3.1.3. Coordinate systems, flow geometry, notation

The numerical solution is achieved using a Cartesian reference frame with the x-axis joining the summits S-S’ of the bumps and the z-axis along the rotation axis of the machine (see Figs 1 and 5). Plane flow is assumed in the x—y plane. For easier reading, the results are shown in polar coordinates (ρ,θ) as defined in Figure 5.

Fig. 5.

Mesh used for computation, coordinates, and notation (shadowed area: ice-bump contact)

The geometry of the flow domain D is simplified by replacing the actual notched outer boundary with a circle of equation ρ = R_e (see Fig. 5). As steady flow is assumed, the inner boundary formed by the two cavity roofs and ice-bump interfaces is a stream line and a trajectory with an equation ρ = ρ₀(θ). In the following, the point at which an ice particle moving along this trajectory loses contact with a bump is named the “separation point”, the point at which it re-contacts the next bump is named the “rejoining point”.

The following notation and conventions are used:

u Velocity vector.

u_x, u_y Velocity components in rectangular coordinates.

u _ρ , u _θ Velocity components in polar coordinates.

n Unit normal to the boundary of D, pointing outward from the ice.

t Unity tangent obtained by anti-clockwise rotation of n by π/2.

Σ Stress vector acting on the boundary of D.

Σ_n, Σ_t Normal and tangential components of Σ, along n and t (positive when tensile)

p Isotropic pressure (positive when tensile).

Ω Angular (driving) velocity of the apparatus.

3.1.4. Boundary conditions

Taking into account the symmetry with respect to the z-axis, the solution, expressed in polar coordinates, is periodical of period π. Consequently, the modelled region is halved (see Fig. 5).

The prescribed boundary conditions are as follows:

• (a) periodicity π expressed by setting u(x,0) = −u(−x,0) and Σ(х,0) = Σ(−x,0) (on ST and S’ T’ boundaries in Figure 5);

• (b) constant angular velocity Ω, resulting in tangential

  velocity U of constant magnitude on the outer ring TT’ : u_ρ = 0; u _θ = U = R _eΩ

• (c) free-surface condition Σ = 0 on the cavity roof S’ C;

• (d) frictionless ice-bump contact expressed by: un = 0 and Σ_t = 0 on SC. It must be verified a posteriori that the ice-bump interface is subject only to compression (Σ_n > 0).

3.1.5. Comments

Condition (c) above corresponds to a zero air pressure P_c inside the cavity. The solution for p_c = H = constant is simply obtained by adding H to the isotropic pressure (the velocity and deviatoric stress fields are unchanged). As ice incompressibility is assumed, the volume of the cavity, fixed by the initial shape given to the domain D, does not depend on p_c. The problem is thus different from the one treated by Brepson (Reference Brepson1979, Reference Brepson and Tryde1980): in his experiments the upper plate was free to move vertically and its upper side was submitted to water pressure p_w. During steady plane flow, the shape of the cavity adapted itself so that the action of the vertical stresses exerted by ice balanced the action of p_w.

3.1.6. Choice of units

The computation is made using dimensionless variables. The length unit is denoted by L. The magnitude U of the tangential driving velocity is taken as the velocity unit. Consequently, the strain-rates are expressed in units of U/L. For the sake of convenience, the Glen’s flow law coefficient B is a fixed a priori. The stress unit is therefore σ = (U/BL)^1/n .

As the reduced equations of the problem are identical to Equations (1) and (2), the notation used for the actual and reduced variables is the same.

3.2.1. Variational formulation

The problem is put into a form suitable for numerical solution by using the Reference BirdBird (1960) Reference JohnsonJohnson (1960) variational principle for steady flow of a non-Newtonian fluid with dissipation potential Φ. It has been proved (Reference Brepson and TrydeMeyssonnier, 1983) that the solution for a given cavity shape is obtained by making the functional

stationary over the set of kinematically admissible velocity fields, i.e. those which satisfy here:

• (a) u(x, 0) = −u(−x, 0);

• (b) u = U on the outer ring (ρ = R _e, 0 ≼ θ ≼ π);

  where uⁱ _j, pⁱ are the values of u_x, u_y, and p at the Ne and Me nodes i and are polynomial functions of x and y.

  The stationarity of J is expressed by solving the set of simultaneous equations with the unknowns uⁱ _j, pⁱ

• (c) un = 0 on the ice-bump contact.

3.2.2. Finite-element method solution

This formulation leads naturally to the choice of velocity and pressure as primitive variables (Reference Nickell, Tanner and CaswellNickell and others, 1974; Reference ThompsonThompson, 1975). The solution to δJ = 0 is carried out by discretizing domain D into a finite number of geometrically simple elements. In each element e the velocity U^e and pressure p^e are interpolated as:

where N is the total number of “velocity nodes” and M is the total number of “pressure nodes”, while prescribing boundary conditions on the velocities.

In this work Lagrangian triangular six-node elements, with a quadratic interpolation of the velocity and a linear interpolation of the pressure, were chosen (Reference ThompsonThompson, 1975; Reference and and CliffeJackson and Cliffe, 1981). Continuity of the velocity across triangle sides is ensured, but not that of the pressure so as to force complete incompressibility at each point within each element. So each triangle has six “velocity nodes” (three vertices + three mid-side nodes) and three internal nodes corresponding to the “pressure at triangle vertices” unknowns. The mesh was made up of 142 elements and 351 “velocity” nodes, of which 105 were triangle vertices (see Fig. 5). Its regular structure was broken so as to refine the element size near the rejoining point. Curved triangles were used to describe the inner and outer boundaries of D (strict incompressibility in these elements is then no longer ensured). It is worth pointing out that this mesh discretizes the space occupied by ice and not the ice itself, corresponding to an Eulerian description of the flow.

3.2.3. Assembly and solution of the finite-element system

Equations (4) are assembled using the contributions from each element to integral J. The boundary conditions on u are specified after assembly. The driving velocity U is prescribed using the Payne—Irons method. The sliding and periodicity conditions are prescribed by making the appropriate combinations of rows and columns of the system, and eliminating the excess unknowns. The resulting non-linear finite-element system is:

where Y is the nodal velocity and pressure-unknown vector, the matrix K depends on the velocities through the viscosity and B ₀ is the “load” vector resulting from the prescribed velocity condition. It is solved by successive updating, for which the convergence, for n ≥ 1 in Glen’s flow law, was proved by Friâa (1979). At step i the viscosity at each integrating point is calculated using the strain-rates obtained at step i − 1 as The resulting linear system, [K(u _{i − 1})]Y_i = B ₀, is solved by the Gauss algorithm.

In this work, the iterations were stopped at step i when ||Y_i − Y _{i − 1}|| / ||Y_i < 10⁻³ and ||B ₀ − [K(u_j )]Y_i ||² < 10^–5 (σL ²)², where ||V||² is the squared length of vector V.

3.2.4. Free-surface fitting

At each step the solution verifies only the natural condition Σ = 0 on the cavity roof S’ C (see Fig. 5): the extra requirement for steady flow is that S’ C must be a stream line. This is carried out as follows:

• (a) the initial shape of the cavity is fixed a priori as the arc of circle tangent to one bump at its summit S’ (ρ₀, π) and intersecting the other bump at point C (ρ₀, π/4), close to the observed position;

• (b) after calculating the solution to system (5), the trajectory of the ice particles passing through S’ is obtained by integrating the equation dy/dx = u_y/u_x;

• (c) this trajectory is then taken as the new free surface and the corresponding velocity field is re-computed;

• (d) the different computations are re-iterated until the condition un = 0 is satisfied on the free surface.

At each stage one must verify that the normal stresses at the ice-bump interface nodes are all compressive, and that the new computed free surface does not overlap the bump profile. This represents the main difficulty of the process, due to the changes in nature of the boundary conditions occurring near the ends of the cavity.

From the start, it appeared that the ice tended to contact with the bump down-stream, looking like the swelling of jets in extrusion die outlets (Reference Nickell, Tanner and CaswellNickell and others, 1974). The iterations were stopped when the ratio |un/uu| was lower than 10⁻² at all the nodes on the free surface. Their displacements, normal to the cavity roof, were then below 10⁻⁴ L.

4.1.1. Velocity field

Due to the large extent of the cavity, the velocity field is little disturbed in the zone situated above it. The u _θ profiles are close to those of a rigid body in the sector π/2 > θ > π(see Fig. 7). The disturbance becomes important in the vicinity of the rejoining point.

Fig. 7.

Profiles of the u_θ tangential velocity component along radii of polar angles θ = 0, θ = π/4, θ = π/2, computed with n = 3 (dotted lines: rigid motion u_θ = ρΩ profiles).

4.1.2. Strain-rates

Contours of derived from the computed velocity field are shown in Figure 8. This figure can be read as a map of τ contours by applying Glen’s flow law. The strain-rate gradient is very high in the region near the rejoining point where γ passes from 2 to 17.6 (U/L)its maximum value, over a distance equal to 0.015L, giving a gradient of about 10³(U/L ²)•

Fig. 8.

Effective shear strain-rate contours (Glen’s flow law, n = 3; dimensionless values).

4.1.3. Stresses

The nodal isotropic pressures are all compressive, lying between −0.2σ and −7.7σ. The minimum compression is reached on the cavity roof, just after separation. The maximum compression is reached at the rejoining point. The istropic pressure remains high in the whole region of convergence of the flow (see Fig. 9), being very high at the rejoining point and stretching over the ice-bump interface. Thus, the isotropic pressure contributes largely to the normal stress acting on the bump (see Fig. 6). A similar result was found by Reference Lliboutry and RitzLliboutry and Ritz (1978) for ice flow past a sphere. Figure 10 shows the principal stresses computed at the centroids of the triangles. All of these stresses are compressive. The influence of the free surface is significant and remains perceptible in the core of the ice above the cavity. On the other hand, the extent of the cavity-bump transition zones within the ice, characterized by the rotation of the principal stress directions, seems very limited.

Fig. 9.

Isotropic pressure contours (Glen’s flow law, n = 3; dimensionless values).

Fig. 10.

Principal stresses at the element centroids (n = 3).

4.2.1. Units and value of B in Glen’s flow law (3)

The length unit such as the dimensions of the finite-element mesh correspond to that of the apparatus is L = 1 m. The angular velocity O was measured very accurately thanks to the photographs and checked to remain constant for the duration of the experiment. It corresponded to a driving velocity U = 154.1 ± 0.5 m a⁻¹. The resulting strain-rate unit is U/L = 154 a⁻¹. The stress unit is derived from the comparison between the measured and computed values of the torque exerted by the ice on the torsion rod. The “theoretical” torque calculated from the Σ_n distribution shown in Figure 6, with a sample thickness of 0.13L is C _th = 48.5 × 10⁻⁴ (σL³). The “experimental” torque is taken as the average of the observed values (see Fig. 4), T being the period during which quasi-steady state was achieved. According to Reference Mellor and ColeMellor and Cole (1981), a better value C _e, accounting for the non-linearity of Glen’s flow law, would be such that

In the present case, the fluctuations of the torque, compared to C _e or C _exp, are small enough so that the two estimates are very close. The average value, taking into account a friction torque of 600 N m, is C _exp = 4150 Nm The C _exp/C _th ratio gives the stress unit: σ = 0.85 Mpa (8.5 bar).

The resulting viscosity unit is σL/U = 1.74 ×10¹¹ Pa s (= 0.055 bar a), and the corresponding value of coefficient B in Glen’s flow law (3) is: B = U/(La³) = 7.92 × 10⁻²⁴ Pa_-3 s^-1 (= 0.25 bar_-3 a^-1). This value will be discussed later.

4.2.2. Shape of the cavity

Diagrams of the cavity profile have been made using the photogaphs taken at time t₀ then after one and four revolutions of the machine (times t₁ and t₄ in Figure 4). Figure 11 shows the variations of ρ (measured to the millimetre) versus θ for the points on the cavity roof. A slight increase of the cavity volume can be noted between the extreme measurements. The computed cavity corresponds well to observations, especially for the one obtained after one revolution (t = t₁ ).

Fig. 11.

Radius ρ(θ) of the cavity roof: (A) cavity observed at t = t₀ (beginning of observations), (B) cavity observed at t = t₁ (after one revolution), (C) cavity observed at t = t₄ (after four revolutions). Solid line: cavity computed with n = 3; dotted line: polar radius of the bump profile.

4.2.3. Deformation of the marker network

The trajectories are solutions of the differential system:

Each trajectory is defined by its starting point M ₀ on the radius θ = π. The solution to Equations (6) is found starting from M ₀ at time t = 0 = t₀ and using the Runge-Kutta method. At time t, the current point M belongs to element e. The position of M at t + Δt is calculated using the nodal velocities solution to finite-element system (5) and the interpolation functions specific to element e. The time step Δt is taken small enough so that, if M crossses the boundary of e during Δt, its new position can be kept. This computation was performed for 19 stream lines, equally spaced at 5 mm on the radius θ = п. Figure 12 shows the deformation of the material line θ = 0, initially radial, after one and four revolutions.

Fig. 12.

Computed deformation (n = 3) of the radial line θ = 0: (0) initial position; (1) after one revolution; (4) after four revolutions (this line is not part of the network which is already deformed at time t₀).

During the experiment, the ice trajectories are observed by following the displacement of each marker.

The comparison of the calculated and observed deformations of the network was made after one and four revolutions. The initial position of the network (instant t₀ in Figure 4) is shown in Figure 13. In order to reduce the computational burden, each individual marker trajectory was linearly interpolated between the two nearest among the 19 already computed trajectories. The results are shown in Figures 14 and 15. After one revolution the computed and observed networks are practically the same. After four revolutions the results remain quite good: the mean deviation in the angular position θ of the markers is 1°.

Fig. 13.

initial position of the marker network (time t₀ in Figure 4).

4.2.4. Time description of the trajectories

More detailed checking of the flow can be carried out by studying how a marker moves along its trajectory with time. The study is restricted to the description of half a trajectory (0 > θ > π). About 50 photographs are available to measure the angular position θ(t) of a marker. The computed values are 12 times more numerous. The results for a marker on the innermost stream line are shown in Figure 16. As the differences between the observed and calculated θ are very slight after half a revolution, instead of θ(t) this figure shows the angular difference θ(t) to a point which would turn at the constant speed of rotation Ω of the machine, i.e. θ(t) = Ωt − θ(t). The comparison is made by taking the rejoining point on a bump as the starting point of the trajectory. It can be noted that when passing over the cavity, the differences θ are slightly greater for the marker than for the theoretical material point, which indicates a smaller increase in the observed

Fig. 16.

Angular difference θ(t) = Ωt – θ(t) for a marker moving along the ice-bump contact and the cavity roof, o: observed; solid line: computed (n = 3) values.

velocity than that calculated with n = 3 in Glen’s flow law.

5.1.1. Steady-state conditions

The cold-room temperature is very difficult to maintain at a constant value over long periods. These temperature variations give rise to heat exchanges with the ice and are largely responsible for the torque variations (see Fig. 4).

On the other hand, the steady-flow assumption supposes that the shape of the cavities does not change. The measured increase in volume of a cavity between zero and four revolutions was about 100cm³, corresponding to a 1% decrease in the volume of ice, half of this variation taking place during the first revolution (volume of the cavity after four revolutions: 680 cm³). This decrease in the volume of ice was found to be mainly caused by melting at the bump contact, as heat was provided by the outer environment, further densification of the ice after time t = t ₀ being undetectable (Reference MeyssonnierMeyssonnier, 1983). Two relatively important melting phases correspond to the decreases in torque recorded at times t = t ₀ and t = t ₂ in Figure 4.

Finally, it must be noted that the polar angle of the rejoining point varied little from the beginning to the end of the observations, passing from 44° to 43.5°.

5.1.2. Boundary conditions

The boundary conditions used for the computation and for which experimental validity must be verified are:

• (a) frictionless ice-bump contact: in order to ensure the existence of a lubricating water film at the ice-bump interface, the melting temperature T _m must be below the ambient temperature. According to Reference LliboutryLliboutry (1976): T_m = T₀ − βp_m − AN/w, where T₀ is the melting temperature at normal atmospheric pressure, p_m is the pressure at the ice-

  water interface, β ≈ 0.1°C MPa⁻¹ (Reference HarrisonHarrison, 1972), N is the number of saline equivalents per kg of ice, w is the water content, and A = 1.86°C kg mol⁻¹. With p _m = 6.5 Mpa, the mean normal stress corresponding to the measured torque, the drop in T_m is 0.65 C To reach the cold-room temperature N/w should be near 0.2 mol kg⁻¹. This is quite a probable value as w was probably very low. During the experiment no disappearance of air bubbles was noticed which would have indicated melting at grain boundaries (Reference Barnes, Tabor and WalkerBarnes and others, 1971) and which was observed by Reference BrepsonBrepson (1979), who obtained a layer of blue ice totally free of bubbles when working at 0°C. Apart from the two short periods when melting was ascertained, ensuring perfect sliding, the friction should have remained very low at a temperature close to −1°C.

• (b) driving condition: the adherence condition on the driving ring seems to be well represented if one refers to the deformation of the marker network after one and four revolutions (see Figs 14 and 15). A computation from a separate finite-element model (Reference MeyssonnierMeyssonnier, 1983) indicated that the sliding velocity over a driving cog, corresponding to a mean drag of 0.07 Mpa, is negligible (u _b > 0.1 cm a⁻¹). However, it can reach 2 m a⁻¹ in zones where the tangential stress reaches its maximum value of 1 MPa (see Fig. 8).

• (c) plane-flow hypothesis: the observations made on pieces of the ice sample collected at the end of the experiment did not show any appreciable inclination of the markers from the vertical.

Fig. 14.

Deformation of the marker network after one revolution (time t₁ in Figure 4): о observed: Δ computed (n = 3).

Fig. 15.

Deformation of the marker network after four revolutions (time t₄ in Figure 4): о observed; Δ computed (n = 3).

5.2.1. Value of exponent n in Glen’s flow law

The numerical results shown above were obtained with n= 3. This value is generally accepted for effective shear stress τ between 0.05 and 0.5 MPa (Reference Barnes, Tabor and WalkerBarnes and others, 1971; Reference DuvalDuval, 1976).

For higher values of τ, higher values of n have been reported in the literature (Reference SteinemannSteinemann, 1958; Reference Dillon and AnderslandDillon and Andersland, 1967; Reference Barnes, Tabor and WalkerBarnes and others, 1971). In general, values close to n = 5 have been quoted. Reference Goodman, Frost and AshbyGoodman and others (1981) explained these high values by the fact that at high stresses (τ > 1 MPa) the dislocation sliding velocity is no longer a linear function of stress. According to Reference Shöji and HigashiShôji and Higashi (1978) or Reference JonesJones (1982), creep at high stresses is accompanied by generalized micro-cracking, initiated by the piling up of dislocations at grain boundaries.

In the present model, τ = 1 MPa with n = 3 corresponds to a dimensionless value τ = 1.18(U/BL)^1/3 or = 1.63(U/L). The extent of the region > 1.5(U/L) covers all of the part situated up-stream of the bump (see Fig. 8). In order to take into account the results of the literature concerning high stresses, the computations were run again with values of n higher than 3.

First, n was taken equal to 4. The corresponding map of does not differ considerably from that obtained with n = 3. Therefore, the comparison is made primarily on the total deformations. Figure 17, which represents the computed deformation of the radial line θ = π/2 after four revolutions, shows that this deformation is of course more important when n = 4, but that the differences are only appreciable near the innermost stream line. For a point on this stream line, the angular difference between its positions computed with n = 3 and n = 4 after four revolutions is 5.6°. As for the deformation of the marker array, the difference in behaviour, already noticeable after one revolution, is very clear in the final position. The model with n = 3 already has a slight tendency to overestimate the deformations and this tendency becomes very pronounced with n = 4. This value must therefore be discarded.

Fig. 17.

Deformation of the radial line θ = π/2 after four revolutions, computed with n = 3 and n = 4.

In order to complete this work, the influence of a greater increase in n, but limited to the high effective shear-stress zones, was studied. To take such behaviour into account, Reference LliboutryLliboutry (1969) and Reference Colbeck and EvansColbeck and Evans (1973) have suggested adopting the polynomial relation = B ₁ τ + B ₃ τ ³ + B ₅ τ ⁵, in which the linear term represents the behaviour at low stresses. In this work we preferred to use Glen’s flow law in which n is a (step) function of т. The viscosity is still given as n = B ⁻¹ τ ¹⁻ⁿ where n = 3 and B = B ₃ for τ > τ ₅; n = 5 and B = B ₅ for τ ₅ ≤ τ; the coefficients B _j being taken to preserve the continuity of the viscosity. In practice, the exact non-dimensional value of the threshold τ ₅ which corresponds to a given stress cannot be fixed a priori (because the stress scale can only be fixed after calculation of the torque derived from the finite-element computation, and comparison to the measured one). Therefore, an approximate dimensionless value of τ ₅ is chosen, based on the results obtained with n = 3, and it is verified a posteriori that the actual value of τ ₅ corresponding to the experimental conditions is not too far from the desired threshold.

The computation was made for a threshold τ ₅ close to 1 MPa. The n = 5 region is then approximately limited by the contour = 1.5(U/L) in Figure 8. The ratio of the observed to the computed torques gives a pressure unit σ = 0.89 MPa and the corresponding actual threshold is τ ₅ = 1.05 MPa. For a marker in contact with the bump, the angular difference with its position computed with и = 3 after four revolutions lies between 2.5° and 3°, depending on its initial position. The difference between the final deformations of the marker network is less pronounced than in the case n = 4.

5.2.2. Value of coefficient B in Glen’s flow law

The torque-measuring system is very accurate but part of the measured torque is due to the friction of the mobile upper plate on its inner fixed support. This friction torque depends very much on the tightening and lubrication of the sliding pad. It can only be estimated by destroying the ice sample, allowing the ice to melt on contact with the bumps (this is done by raising the cold-room temperature and letting the machine turn without dismantling the upper plate). This procedure was not adopted here, in order to preserve the sample. A friction torque of 600 N m, known to within 400 N m, was assumed. This corresponds to the average of values obtained in other tests at the end of which the ice was melted. The resulting uncertainty in the measurement of torque is 10%. It remains acceptable for the estimation of the stresses involved, but the accuracy of S, computed for n = 3, is only of the order of 30%.

In other experiments at the end of which friction was measured, B was found to lie between 6.34 × 10_-24 and 1.27 × 10_-23 Pa⁻³S^-1 (0.2−0.4 bar_-3 a^-1).

The values of B computed with n different from 3 are:

(a) with n = 4: B = 7.6 × 10⁻³⁰Pa⁻⁴s⁻¹ (0.024 bar⁻⁴a⁻¹)

(b) with n = 3 if τ > τ ₅, n = 5 if τ ≥ τ ₅ = 1.05 Mpa;

B ₃ = 7 × 10⁻²⁴Pa⁻³s⁻¹ (= 0.22 bar⁻³a⁻¹

B ₅ = 6.3 × 10⁻³⁶Pa⁻⁵s⁻¹ (= 0.002 bar⁻⁵a⁻¹

The values corresponding to n = 3 are higher than Reference DuvalDuval’s (1976), i.e. B = 3.17 × 10⁻²⁴Pa⁻³s⁻¹ (= 0.1 bar⁻³a⁻¹) at −0.05°C, for a low water content (w = 0.03%), and in a secondary creep state. This is due to the fact that at present deformation level tertiary creep is reached: re-crystallization and fabric formation contribute to increase the value of B. Comparison with Reference DuvalDuval’s (1981) value for tertiary creep and zero water content, B = 6.3 × 10⁻²⁴Pa⁻³ s⁻¹ (= 0.2 bar⁻³a⁻¹ indicates that w remained very low during our experiment.

5.2.3. Steady-creep hypothesis

During the experiment the ice is obviously not in a state of steady creep: stresses at a given material point are very different up-stream and down-stream of the bump (see Fig. 8).

Figures 18 and 19 show the mean angular differences θ (as defined in section 4.2.4) observed for the markers which travelled close to the ice-bump interface. These averages are calculated as follows:

• (a) for each marker i, the observation allows θ at time t_n = nΔt to be calculated where Δt is the interval between. two photographs, Ω is the machine angular velocity, and θ ⁱ (t_n ) is the polar angle of marker ; measured at time t_n;

• (b) a value of θ ⁱ (θ_m) for is linearly interpolated through

• (c) assuming that θ ⁱ (θ₀) at θ₀ = π (see Fig. 18) or θ₀ = π/4 (see Fig. 19) for all the observed markers i, it is then possible to calculate an average ϴ(θ) corresponding to an Eulerian description of the flow.

Fig. 18.

Mean angular difference ϴ(θ) for the markers along the cavity roof, compared to those computed with n = 1, 3, or 4 (θ = 0 is assumed at the bump summit for each marker).

Fig. 19.

Mean angular difference ϴ(θ) for the markers along the bump contact, compared to those computed with n = 1, 3, or 4 (θ = 0 is assumed at point θ = π/4 for each marker).

Figure 18 shows ϴ(θ) observed for markers which have passed the bump summit and are moving along the cavity roof. A comparison is made with the values calculated for ice obeying Glen’s flow law with n = 1, 3, or 4. Figure 19 gives ϴ(θ) for markers which have passed through the rejoining point and are travelling along the ice-bump interface.

First of all, it must be noted that a purely Newtonian model is not suitable for simulating the ice behaviour. Along the first three-quarters of the cavity roof, the points calculated with n = 3 or n = 4 are ahead of the observed average (see Fig. 18). This tendency is reversed when approaching the next bump. For a marker in contact with the bump, the accordance with the computation for n = 3 is excellent, except on approaching the bump summit (see Fig. 19).

The difference in behaviour observed between a theoretical Glen body and actual ice, down-stream of the bump, could be explained by the influence of transient creep occurring on unloading (Reference DuvalDuval, 1976; Reference Duval, Ashby and AndermanDuval and others, 1983). As for (re-)loading, no appreciable influence of Andrade creep (which appears upon loading of a virgin or annealed sample) can be noted. This result is in accordance with Reference LliboutryLliboutry’s (1975) estimations.