1. Introduction

The ideas of genericity and of structural stability, which originated as concepts within mathematics (Reference ThomThom, 1975; Reference ArnoldArnold, 1983, Reference Arnold1986; Reference Guckcnheimcr and HolmesGuckenheimer and Holmes, 1983), have much to teach us in our attempts to model the natural world. In mathematics they have, of course, the rigour that is appropriate to that subject; in science it is natural that they should be applied in a somewhat looser way. In particular, they are of wide applicability in geophysics, where problems characteristically involve scalar, vector or tensor fields that have no special symmetry. An example in glaciology is the time-dependent behaviour of the edge of an irregular ice sheet (Reference NyeNye, 1990). This paper studies the related problem of the behaviour of a flow centre in an irregular ice sheet.

The point of view to be adopted is unconventional in glaciology, and needs to be stated explicitly. The conventional way of treating ice sheets theoretically is first to set up a mathematical model, which is an idealization of the real ice sheet and which is chosen so as to be analytically or numerically tractable. The model will commonly contain assumptions about the flow law for ice, the isotropy of the ice, its thermal properties, the form and nature of the bed, the various boundary conditions, and so on. Here, however, we proceed differently.

Consider the surface of a real ice sheet at a particular time. A topographic map made at a certain scale represents a smoothed version of the real surface, and it is this smoothed surface that we deal with. Consider also the real two-dimensional horizontal velocity field of the ice at the ice-sheet surface, and let this now be smoothed at the same map scale. These smoothing operations completely, define our model, and they represent the correspondence between model and reality. Thus, we are not dealing with an ice-sheet model in the usual sense, but with an ice sheet which is simply a smoothed version of a real one, such as the Antarctic or the Greenland ice sheet, in all its geographical complexity.

With such a general approach, it might be thought that little of value could be said, but this is not the case. Certain conclusions follow purely from topological considerations, which are absolute for generic smooth surfaces and functions. The word generic, as we use it in this paper, simply means typical. It refers to something that will occur when there are no special conditions. Likewise, non-generic means "of vanishingly small probability". This is admittedly a loose definition; in mathematics, as one would expect, generic has a strictly defined meaning, but such rigour would not be appropriate here.

The importance of the topological conclusions that follow from genericity derives precisely from the fact that they are independent of specific mechanisms of flow, and of specific models, in the usual sense of that term. Specific models, however good, are necessarily approximate; topological conclusions are absolute.

2. Structural Stability: Maxima, Minima and Saddles

As an introductory example, let us consider how the surface of an ice sheet can change with time. The smoothed surface that we have constructed will generically show three kinds of points where the surface is level: maxima(summits), saddles and minima (hollows), although the latter will occur rarely. At a slightly later time, as the ice sheet evolves, all these points will be found to have moved a little, but they will all retain their identities: they are said to be structurally stable. Generically, there are only two ways in which they can destroy themselves: (1) by moving to the edge of the ice sheet and disappearing, or (2) by a maximum or minimum colliding with a saddle to leave no level place (Fig. 1). If the lefthand diagram of Figure 1 represents a maximum and a saddle, the transition to the righthand diagram could be accomplished by progressively tilting the whole surface until eventually it all slopes downwards to the right.

Fig. 1. Height contours of an ice sheet. A saddle S (index -1) and a maximum M (index +1) annihilate one another, the total singularity index remaining zero throughout.

3. Directions of Ice Flow

Consider now, at a given instant, the smoothed two-dimensional distribution of horizontal velocity at the surface of the ice sheet, represented by the two-dimensional vector field u(r), where r = (x, y) denotes position. As in any such field there are, in general, points where u(r) = 0. An ice centre is a point where u(r) = 0 and where the two principal horizontal strain rates e ₁ and e ₂ are both positive. (It is non-generic for either of them to be zero precisely at a point where u(r) is also zero. That would be a structurally unstable situation; for example, a change in the scale of smoothing would destroy the coincidence.)

Although an ice centre will normally occur near a summit, it is non-generic that it should coincide with a summit exactly; this would entail conditions so special that they would not be found in Nature. Similarly, it is non-generic that the rotation rate at an ice centre should be exactly zero, although it may very well be quite small. Most detailed ice-sheet models are based on a flow rule of the type

where ϕr is the height of the surface, and s(r) is some scalar function of position (which would include, for example, a dependence on bed roughness, as well as on powers of

and of ice thickness). This is the most general rule that ensures that u(r) is perpendicular to the surface contours (Reference NyeNye, 1952; Reference HutterHutter, 1983, p. 459), as is commonly observed to be approximately the case. In such a model, an ice centre would necessarily occur exactly at a summit (grad ϕ=0). It follows from rule (1) that, since curlgrad= 0,

Thus, with this rule, not only would there be no velocity at the summit, but no rotation either. However, there would be rotation wherever grad s(r) is not parallel or anti-parallel to grad ϕr. It is useful to note these consequences of rule (1) but, since the rule is only approximate and our aim is to avoid all approximation, we shall not, in fact, apply any such constraint in what follows. Our point of view is that there is nothing to prevent a rotation at an ice centre, and therefore it will occur in Nature, even if only to a small degree. A possible physical mechanism is provided by the nature of the ice-sheet bed: this will necessarily favour cither clockwise or anticlockwise motion; it will never favour both equally, for that would be non-generic. Of course, this is only one of a large number of possible physical mechanisms.

Taking an ice centre as origin O and x, y axes parallel to the principal horizontal strain rates e ₁ and e ₂ at O, the velocity components u(x,y) and v(x,y), at a given instant, are given to lowest order by

where ω is the instantaneous rate of rotation at O. To find the possible patterns of instantaneous flowlines, we write

, thus

If the velocity distribution were steady with time, these equations would give the actual particle motions. However, generically, the velocity distribution changes with time, and so the equations merely describe the instantaneous motion of the particles. We now follow a standard procedure from dynamical systems analysis (Reference PippardPip-pard, 1985, p. 22). First eliminate y to give

where T is the trace, e ₁ + e ₂, and ∆ the determinant, e₁e₂ + ω², of the matrix of coefficients in Equations (2). Since Equation (3) is the equation for damped harmonic motion, and T is positive, the solution grows exponentially:

The resulting flowline patterns are drawn in Figure 2a, b, c and d for e ₁/e ₂ = 4. When

is below the critical value given by , there are two straight flowlines, generically not at right-angles, one of them being a common tangent at 0 to all the other flowlines. When is above the critical value, the flowlines cease to possess a common tangent. We expect the rotation rate to be small at an ice centre, as already discussed, and therefore a pattern between those of Figure 2a and b is to be expected. As time progresses and the ice sheet evolves, an ice centre will move, but because it is structurally stable the general features of its local flowline pattern will be preserved.

Fig. 2. Flowline’patterns at an ice centre for varying amounts of rotation. The ratio 0F the principal strains e₂/e₁ is taken as 4. The straight flowlines come together as w increases.(a)

, the critical value, (d) w = 2.22e₁.

If rule (1) were obeyed, the ice centre would be exactly at a summit and there would be no instantaneous rotation at this point. It is readily shown that in this case the contour lines for the surface would be concentric ellipses of the same shape as the representation quadric for strain rate, namely

These are the ellipses shown in Figure 2a. The instantaneous flowlines would be the one-parameter (A) set of curves

which are drawn in Figure 2a for e ₂/ e ₁ = 4. They would fan out from the direction of lesser principal strain rate.

A similar treatment can be given for the flowlines around a zero-velocity point with e ₁ and e ₂ having opposite signs, as would be expected near a surface saddle. If rule (1) applied, u would be exactly zero at a surface saddle, and so would ω. Therefore, we expect a u = 0 point near a surface saddle and that ω will be small at this point. Since

is likely to be small compared with , we expect the determinant A to be negative, and therefore, in dynamical systems terminology, a saddle-point instability. However, if, exceptionally, were large enough , four other possibilities would appear as follows (see Reference PippardPippard, 1985, p. 22 for terminology):

where "fi" means focal instability, "ni" means nodal instability, "ns" means nodal stability and "fs" means focal stability.

4. Patterns of Strain Rate

The pattern of horizontal strain-rate directions near an ice centre behaves quite differently from the pattern of flowlines. To explain this, it is convenient to begin with an ice sheet that is circularly symmetric in plan view in all its features (topography, accumulation rate, temperature, etc.). Of course, this is quite non-generic, but we shall in a moment perturb the ice sheet in a general way. Before perturbation, the ice spreads out radially from the centre in all directions and the horizontal strain rate eij (i,j= x,y) at the centre is isotropic:

say, the constant K being positive. (At an isotropic point for strain rate, in general, e has this form, but K need not be positive.)

Away from the centre the strain will not be isotropic, and in the simplest case the trajectories of principal strain rate would be radial lines and circles (Fig. 3a). However, circular symmetry would still be retained if the trajectories had a spiral character, as in Figure 3b. It can be shown that this implies that the fiowlines are spirals also. In addition, there is a rotation rate that varies with radial distance. (These effects involve terms two orders higher than were considered in section 3, that is, cubic terms in Equation (2).)

Fig. 3. Trajectories of principal strain rate for a circularly symmetric ice sheet, (a) without, and (b) with rotation.

The singularity index of the patterns of orthogonal strain-rate directions in Figure 3a and b is +1. However, the singularity index for the trajectory pattern of an isotropic point, in general, is ±½ as can be seen by inspection of the patterns in Figure 4 (Reference NyeNye, 1983, Reference Nye1986a). (The reason for the

, as compared with 1 for the contour patterns, is that rotation of the strain-rate tensor at a point by π brings it into self-coincidence, whereas a contour, which has a high side and a low side, has to be rotated by 2π to be brought into self-coincidence.) Therefore, we can expect that, on perturbation, our +1 isotropic point will split into two, each of index . In other words, the +1 isotropic point is structurally unstable. (A similar splitting occurs with the isotropic points for curvature situated at each end of the axis of revolution of a spheroid; when the rotational symmetry is broken, each of these splits into two. In crystal optics the transition from a uniaxial crystal to a biaxial crystal is somewhat similar. However, these analogies are not perfect, as we shall see.) After perturbation, the local flowline pattern will be like those of Figure 2; unlike the central isotropic point, the flow centre is not split.

Fig. 4. The pattern of principal strain-rate trajectories around an isotropic point can, in general, he of three different kinds. The lemon and monstar patterns have index

while the star pattern has index

. On the other hand, the star and monstar have three straight trajectories, while the lemon has only one. The nomenclature is due to Reference Berry and HannayBerry and Hannay (1977).

In principle, the two isotropic points for strain rate could be lemon or monstar (Fig. 4), but not star, because this has index -½. Isotropic points are also characterized by the contour property. Around each isotropic point one can draw contours of equal principal strain rate (lines along which a principal strain rate e ₁ ord e ₂ takes a fixed value). There are two possibilities: the contours for both e ₁ and e ₂ are either two families of ellipses or two families of hyperbolas. Thus, in principle, not only could our two isotropic points be lemon or monstar, but in addition they could be elliptic or hyperbolic. We now show that they are, in fact, always lemon and nearly always hyperbolic.

To represent the unperturbed e_ij within a small neighbourhood of the central degenerate isotropic point, taken as the origin, write the quadratic expansion

(the strain-rate distribution is considered smooth for the same reason as the surface topography). Linear terms have been omitted, because the pattern is centro-symmetric and therefore the tensor must be invariant under the transformation

. Circular symmetry of the kind shown in Figure 3a requires, further, that (a) the directions of principal strain rate are radial and circumferential, and (b) the two strain-rate invariants are functions of r and not 0 (r, 0 being polar coordinates). A short analysis shows that these conditions reduce e to the form

If, however, we admit a rotation, as in the generic case of Figure 3b, there are additional terms involving c:

where

Working with this more general form, we now perturb it by a small amount in lowest order to make the origin non-isotropic, thus:

where

where ϵ1 and ϵ2 are small constants. "Small" means that the separation of the new isotropic points is small compared with a characteristic radial dimension of the original strain-rate distribution. (A more general perturbation in q1, by ϵ₃ say, rather than by –ϵ1, would have no effect other than to change the mean normal strain rate at the origin.) It is necessary to impose one further condition: because the strain-rate tensor is the gradient of a single-valued velocity function, its components must satisfy what in elasticity would be called the equations of compatibility. In two dimensions, as here, they reduce to the single condition:

Applied to Equation (7) this gives a = 3b, so that the tensor now has the form

where

The conditions for an isotropic point, e_xx = e_yy, e_xy = 0, reveal that there are now two such points symmetrically placed and, without loss of generality, we can make them occur on the x, y axes (different from those of section 3 because ϵ₂ ≠ 0) by setting

Taking a suitable sign for ϵ₁, they now occur at x =

say, and y = 0. Shifting the origin to x-+ϵ, y=0, and retaining only the linear terms in the expansion about this point, we obtain

where K' is a constant,

. To characterize the isotropic point, we use the discriminants given in the Appendix. Thus, to check the singularity index, we evaluate the discriminant D ₁ to find

Since this is always positive, the index is

as was anticipated.

To decide whether there are one or three straight trajectories passing through the point (straight meaning zero curvature at the point itself), that is, whether the pattern is lemon or monstar, we evaluate the line discriminant D _L, (see Appendix). Thus

Since D ₁ is always negative, the line number is 1 and the pattern is lemon rather than monstar.

For the contour discriminant D _c (see Appendix) that decides whether the contours of constant principal strain rate are elliptic or hyperbolic, we find

Thus, if c ² < 12b ², the contours are hyperbolic and, if c ² > 12b ², they are elliptic. Recalling that c is the parameter that controls whether the pattern is spiral or not, we sec that, unless the spiral quality is pronounced, the contours are hyperbolic; the elliptic case would be rare. It may be surprising that the contours at the isotropic points can be hyperbolas when the unperturbed contours were simply concentric circles; the details of the transition are illustrated in Nye (1991).

Although the tensor in Equation (9) with condition (10) contains four parameters, namely K,ϵ₁, b, c, the pattern of principal strain-rate directions it describes is, in fact, entirely controlled by a single parameter, apart from an overall scaling. To see this, first subtract off the isotropic strain rate K and divide the remaining strain rate by ϵ₁, neither of which operations affects the pattern. Then apply the isotropic scaling

to give, in place of Equation (9), the strain-rate tensor

where c’ = c/2b. In this form the only free parameter is c’, and the isotropic points occur at X = ±1, Y = 0. Figure 5a, b, c and d shows the patterns of principal strain-rate directions for c’ = 0, 0.5, 1, 10. Notice that it is only for c' = 0 (i.e. c = 0) that the two lemons share a common trajectory. Such sharing would not be typical in Nature. Generically, the lines from the two points miss one another, and the whole pattern has a spiral character, however slight this may be in practice. The change-over from hyperbolic to elliptic takes place at

Fig. 5. Patterns of principal strain-rate directions (Fig. 3b) after perturbation. The values of c’ are (a) 0, (b) 0.5, (e) 1, (d) 10. The trajectories mere integrated by a low-order Runge-Kutta method with an interval proportional. to the difference of the two principal values.

We have not been specific about the nature of the small perturbation that acted on the original degenerate isotropic point at the centre of a circular ice sheet; this is precisely because it does not matter, so long as the perturbation is such as to break the original circular symmetry, leaving either none at all, or, at most, two-fold symmetry (an example of what happens when four-fold symmetry is preserved by the perturbation is shown in Figure 13 of Reference NyeNye (1986a)). The parameter ϵ, which is called a control parameter, could represent time, so that we should then be watching the evolution of the ice sheet. More interestingly, it could represent the scale of resolution on which we choose to observe, that is, the map scale. Viewed from a distance, the patterns in Figure 5a, b, c and d are those of the unperturbed strain rate. Thus, generically, a pattern which appears, on a certain map scale, to be as in Figure 3b will inevitably be found, on higher resolution, to be broken up in the way shown by figure 5b, c and d. (Likewise, a pattern which, on a certain scale, shows no isotropic points may be found, as the resolution is increased, to contain one or more close pairs of isotropic points of opposite index.) As another example of the interpretation of ϵ, suppose we are considering a circularly symmetric ice sheet modelled numerically, and suppose r is some constant which breaks the symmetry of the model when it is changed (it might be a constant governing the temperature distribution, for example). Then the central degenerate isotropic point will break up in the way we have calculated.

Conclusions

We have shown, on purely topological grounds without reference to the mechanism of flow, that when the perfect symmetry of a circular ice sheet is perturbed in a general way the isotropic point for strain rate at its centre, being structurally unstable, breaks up into two structurally stable components. Around them, the strain-rate trajectories always have the lemon pattern, and the contours of equal principal strain rate are usually hyperbolic. This splitting is in contrast to the behaviour of the ice centre itself, which remains single and intact.

As the ice sheet becomes progressively less circular, we expect the two isotropic points to separate further, and, since their indices cannot change discontinuously, they will remain

, Meanwhile, new isotropic points will be born together as twins of opposite index, or will enter the ice sheet at its boundary. Thus, the original isotropic points will continue to exist until either they interact with one of the new ones or escape across the boundary. Although any one isotropic point retains its index (between interactions), nevertheless its line and contour character can change.

Of course, the distribution of the horizontal component of velocity in an ice sheet is only one example of a two-dimensional flow field which exhibits phenomena similar to those described here; many other examples could be envisaged. We also remark, without proof (to be published), that there is a perfect analogy between the breakup of the central isotropic point and the behaviour of the focus of a slightly imperfect lens used with circularly polarized light. The strain-rate ellipse is analogous to the polarization ellipse. However, the analogy with the breakup of a circularly symmetric degenerate umbilic point of a surface is not exact because of the presence of the parameter c.