Introduction
Floating ice masses affect shipping and offshore activity, afford a base for scientific expeditions, and are a potential resource. The flow of glaciers and ice shelves, which may result in the eventual formation of icebergs, has been studied by many investigators. Reference NyeNye (1951) studied the plastic flow within an ice mass, under its own weight, assuming a constant yield stress and Lévy-Mises equations of flow and reported the two-dimensional flow within the glacier with reference to bed slope and formation of transverse crevasses and shear faults. Later Reference NyeNye (1952) assumed laminar flow in valley glaciers in which the movement is parallel to the bed, and explained the formation of transverse crevasses and shear faults by a flow in which an excess longitudinal stress was assumed. Crevasse patterns based on the nature of the longitudinal stress were also indicated. Reference GlenGlen (1952) performed compression tests on blocks of ice and found that ice does not have a constant coefficient of viscosity; a linear relationship was seen between the logarithms of strain-rate and stress, representing a motion in-between laminar and plastic flow. Reference NyeNye (1957) analysed glacier flow based on Glen’s law. Lliboutry (1957) used plasticity theory to determine the profile of a glacier, and explain the formation of crevasses, bulges passing down a glacier, glacier erosion, and the upturning at the foot of an ice fall.
Reference WeertmanWeertman (1957) analysed the creep deformation of floating ice shelves, using Glen’s flow law and Nye’s application of plasticity theory for regions away from the edge of the shelf. Reference ReehReeh (1968) treated the floating ice shelf as a beam loaded by the difference between the hydrostatic pressure in ice, the water pressure at the ice front, and the buoyancy of water; the maximum shear stress occurred at a distance of about the thickness from the ice front. The flexure of a floating ice tongue due to tidal variations was studied by Reference HoldsworthHoldsworth (1969), who idealized the ice tongue as a beam on an elastic foundation, and studied the break-up of the ice, considering the behaviour of the material as elastic, elastic-plastic, and fully plastic. The floating ice was assumed to be sufficiently long so that the ice floated freely, away from the hinge line, ensuring that the imbalance of forces at the ice front did not affect the analysis. Various mechanisms of calving were discussed by Reference Holdsworth and HusseinyHoldsworth (1978). The expansive creep in an ice shelf, together with bending due to bed geometry at the hinge line and ocean level, was considered to study the break-up of the ice. A linear elastic theory with constant properties was used. Calving at the hinge line due to tidal variation and geometry of the bed, and at the ice front due to the imbalance between the hydrostatic pressure in the ice and the water pressure, was reported. The formation of large icebergs was explained by a vibration calving mechanism in which the ice tongue was represented as a rectangular plate with varying thickness. A finite-difference scheme was used to solve the eigenvalue problem, and the locations of maximum stress and deformations required for fracture were obtained from the mode shapes. Reference SmithSmith (1977) presented a fracture mechanics approach for calving in floating ice shelves due to tides. By this means, the floating ice shelf was modelled as a cantilever beam loaded by the buoyancy of the water to determine the stresses in the ice, and the depth of cracking was determined by comparing the stress-intensity factors due to tensile stresses caused by bending and the compressive overburden pressure. Reference IkenIken (1977) presented a plane-strain finite-element solution for an ice mass breaking off a cliff. Glen’s flow law was used to calculate the stresses and movement in the ice, and to study calving in a glacier undercut by a lake. The calving mechanism was modelled as a process in which a crack opens perpendicular to the principal tensile stress, extending up to a point of zero tensile stress, and which is then followed by a build-up of tensile stress at the crack tip until a critical value is reached when the process repeats itself.
Interest in the study of the behaviour of ice under loads has increased with the use of floating ice as a structural member. Reference WilliamsWilliams (1976) presented a theory for the analysis of a viscoelastic plate on elastic foundations, considering thermal stresses and the dependence of material properties on temperature. The plate equations were obtained from the three-dimensional theory of linear viscoelasticity. Reference YakuninYakunin (1974) presented a linear viscoelastic analysis of a floating ice sheet, idealizing the ice as a combination of Maxwell and Kelvin models. A survey of the bearing capacity of floating plates was given by Reference KerrKerr (1976). Reference Vaudrey, Katona and FrankensteinVaudrey and Katona (1975) presented a viscoelastic finite-element analysis of floating ice sheets in which the bulk and shear moduli were expressed as relaxation functions given by a constant and summation of exponential terms. Different types of viscoelastic models were obtained by varying the number of exponential terms chosen. The buoyancy of the water was represented by a Winkler foundation, and the deformation in a circular floating platform subjected to central loading was illustrated. The above formulation was extended by Reference Swamidas, El-Tahan and ReddySwamidas and others (1978) to include temperature-dependent relaxation functions, varying thickness and reinforcement, and results were presented for an artificially thickened offshore platform. Reference Tinawi and MuratTinawi and Murat (1979) developed a finite-element program for the study of floating ice sheets. In this, a polynomial variation of the elastic modulus across the ice thickness was assumed, and a temperature-dependent flow law representing creep and an incremental analysis using initial strain technique were used. Reference HutterHutter (1975) developed a plate theory to analyse floating ice sheets, assuming viscoelastic material behaviour and a non-uniform temperature distribution across the depth. The theory also accounted for shear deformation and large displacement, Reference Hutter, Williams and TrydeHutter and Williams (1980) presented a brief review of analysis of floating ice sheets, discussing flow and temperature distribution. Reference Dempster and TrydeDempster (1980) studied the size, towing, and drift characteristics of icebergs based on experimental data and indicated the feasibility of viscoelastic analysis of calving glaciers to determine the probable upper bounds on iceberg size.
This paper describes a finite-element viscoelastic analysis of an ice tongue. A computer code, VISICE, developed by Reference KatonaKatona (1974), is modified for studying calving mechanisms. The forces at the ice front were studied based on temperature-dependent relaxation functions, tidal actions, and bedrock geometry.
Statement of Problem
The behaviour of a floating ice tongue due to tides, bedrock geometry, and forces at the ice front are studied. Figure 1 shows the break-up at the mouth of a glacier. The splaying and uniform break-up prompted the use of an axisymmetric code. An assumed plan view of the ice shelf is shown in Figure 2. Numerical results are presented for an ice tongue, 300 m thick, floating for a length of 1.8 km. The finite-element mesh is shown in Figure 3, with the buoyancy of the water represented by springs at the interface nodes
Three types of possible loading of the ice tongue are considered. Figure 4 shows the salient features of marginal calving. The imbalance between the hydrostatic stress in the ice and the sea-water pressure at the ice front produces a net outward force and a clockwise moment at the edge. Figure 5 shows a step in the bedrock. The sea-water level is such that the ice level decreases seawards away from the bedrock support.
The difference in the level of water from that required for full buoyancy is taken as 3 m. The bending in the ice tongue due to tidal variations in the water level is shown in Figure 6. A tidal variation of 2 m from the mean is assumed with a period of 12 h.
A linear viscoelastic theory is used for the analysis. Relaxation functions are obtained by modelling the ice as a standard linear solid (Fig. 7) based on its short- and long-term properties. The effect of temperature on the dynamic modulus (Reference MichelMichel 1978), and the variation of Poisson’s ratio across the thickness (Reference HutterHutter 1974), have been reported to be not very significant. However, under static conditions when the load application is relatively slow, ice exhibits considerable delayed elasticity giving rise to a lower apparent elastic modulus. Under these conditions, which describe the present case, the elastic modulus and Poisson’s ratio are dependent on temperature (Reference GoldGold 1958), Steady-state temperature distribution has been obtained assuming the thermal conductivity to vary linearly as the temperature. This leads to a parabolic variation in temperature. The water is assumed to be at 0°C. Mass densities of 875 kg/m3 and 1025 kg/m3 are assumed for the ice and the water respectively.
Analysis
Temperature distribution
Reference Drouin and MichelDrouin and Michel (1974) have reported the findings of various investigators on the thermal conductivity of ice as a function of temperature. It was found that the thermal conductivity increases almost linearly with decrease in temperature. The coefficient of thermal conductivity k can be expressed as
where k 0 and b are constants, and θ is the temperature.
Applying the law of thermal conduction.
where q is the rate of heat transfer, A is the sectional area, and is the temperature gradient; the temperature distribution is given by
where θ1, θ2 are the ice surface and ice-water interface temperatures respectively, d is the depth of the floating ice, and z is the distance from the ice-water interface.
Material model
Ice is modelled as a linear viscoelastic material to account for its time-dependent behaviour (Reference KatonaKatona 1974). The viscoelastic constitutive relation for an isotropic material is given by
where σ ij are the stresses, ε ij are the strains, K and G are the bulk and shear moduli, and d̄ represents a derivative operator. The bulk and shear moduli are monotonically decreasing relaxation functions, represented by the Prony series given below:
and
where t is the time and αi, βi are the relaxation times. K 0 G 0 are the long time moduli, and the instantaneous values are given by K 0 + Σ K i and G 0 + Σ G i. Different types of viscoelastic behaviour are obtained by varying the number of terms in the summation. In the analysis only one term is considered for representing a standard linear solid (Fig. 7).
Finite-element model
An axisymmetric finite-element formulation (Reference KatonaKatona 1974) is used for the solutiun of the floating ice-tongue problem. Quadrilateral isoparametric elements with two degrees of freedom per node are used. An incremental approach is used, and the governing matrix equation for the solution of displacement increments for a given time, tn , is given as
where is the incremental stiffness matrix, {u} is the displacement vector, {R} is the vector of applied loads, {H} is the history vector up to the time tn , and [Kω ] is the stiffness matrix due to the fluid foundation. A Winkler-type foundation is assumed, and for a compatible axisymmetric element with two nodes on the fluid interface, the spring stiffness is given as
where the suffix m refers to the element m, γ is the unit weight of supporting fluid, and r 1 and r 2 are the radial distances of the two ice-water interface nodes. The equations are reduced using a solver for banded, partitioned matrices.
Results and Discussion
The short- and long-term displacements at the top and bottom of the ice tongue, due to the imbalance of forces at the ice front causing creep, are shown in Figure 8. The ice front is subject to larger bending and due to relaxation there is a reduction in bending away from the ice front. The stresses near the ice front are shown in Figure 9. A maximum radial stress of 0.28 X 103 kN/m2 occurs at the top at a distance from the ice front of about three-quarters of the depth. A similar value was obtained by Reference Holdsworth and HusseinyHoldsworth (1978) about half the depth away from the ice front. Cracking would take place at the point of maximum tensile stress and travel downwards resulting in the breaking away of an ice mass whose length is about equal to its thickness. A similar situation in the third dimension would result in evenly-dimensioned ice masses when other environmental forces are absent. Table 1 gives a summary of some iceberg features on the northern Grand Banks of Newfoundland reported by Reference AllenAllen (1973). The least ratio of sides is a little more than half and compares well with the viscoelastic calving value.
The stresses at the top and bottom layers of elements due to a tidal variation of 4 m are shown in Figure 10. The ice responds with the tide and the maximum stresses occur at the hinge line. As the stresses alternate, cracks would form and close at the hinge. The cracks will deepen as the modulus of the section is reduced after cracking and may not regain the full tensile strength when it closes. A similar situation is seen for the bending due to a rock step at the hinge (Fig. 5). The stresses at the top and bottom layer of elements are plotted in Figure 11. The maximum stresses occur at the hinge line, and, unlike tidal bending, the stresses do not alternate. Cracks will form at the hinge line and progress downward resulting in the whole ice tongue breaking away. Unlike the case of tidal bending, another maximum occurs at a point 1 km away from the hinge. Similar maxima have been reported by Reference HoldsworthHoldsworth (1969, Reference Holdsworth and Husseiny1978), and equally-spaced rifts have been noted on ice shelves in Antarctica.
Simple examples have been used to demonstrate the use of a viscoelastic model. The model can be improved by extending the code to viscoplastic analysis, and the inclusion of fracture mechanics to study crack propagation.
Acknowledgements
The authors would like to thank Dr A.A. Bruneau, Vice-President of Professional Schools and Community Services, Memorial University of Newfoundland, St John’s, Newfoundland, for continued interest and encouragement. Appreciation is expressed to Professor J.S. Tennant, Chairman, Ocean Engineering Department, Florida Atlantic University, for the opportunity afforded to the first author (DVR) to continue the work as a Visiting Professor. The support of this investigation by the National Research Council of Canada, Grant No. 8119, is gratefully acknowledged.