Notation

Introduction

Large-scale numerical models of ice-ocean dynamics (e.g. Hibler, Reference Hibler III1986) include no parameterization for incoming ocean waves, yet the influence of waves and swell appears to be significant in comparison with other terms in the stress balance equation. The effect of waves is felt through radiation stress, an excess flux in momentum due to the presence of the waves which acts both into and along the ice edge. At an ice edge, indeed at any floating or fixed body, the radiation stress depends directly on the proportion of energy reflected by the body together with the proportion allowed through into the open seas beyond the body. Hence, open water reflection and transmission coefficients are crucial.

In a recent series of papers, Fox and Squire (Reference Fox and Squire1990, Reference Fox and Squire1991a,b) and Squire and Fox (1990) have reworked the classical problem of ocean waves propagating into an ice sheet, obtaining precise solutions where only approximate solutions were hitherto available. The new analysis includes a mathematically complete set of wave potentials on both sides of the ice edge, and the problem is solved variationally as will be discussed briefly later. The precision of the solution is established by an independent asymptotic energy flux calculation. It is this work which is used to provide the reflection coefficient at the edge of the marginal ice zone, enabling the contribution of radiation stress in the ice edge stress balance to be assessed. In modelling the ice cover as a continuum we are certainly idealizing the real world. However, for numerical practicality the continuum hypothesis underpins all of the large scale models with which we are concerned, and it is commonly used in connection with wave propagation problems in ice-infested seas (e.g. Weber (Reference Weber1987), who models the sea ice as a viscous Newtonian fluid; or Liu and Mollo-Christensen (Reference Liu and Mollo-Christensen1988), who assume a thin elastic plate under compression). The continuum hypothesis will be most valid when the entire ice cover moves horizontally as what is effectively a rigid or near rigid lid. This will occur at high sea-ice concentrations such as are observed at the Greenland Sea ice edge where high waves pommel ice floes into a slurry of 100% concentration. A similar ice cover exists in the Weddell Sea in winter for hundreds of kilometres (Wadhams and others, 1989). Radiation stress ideas have been invoked by Martin and others (Reference Martin, Kauffman and Parkinson.1983) and Wadhams (Reference Wadhams1983) in connection with the formation and subsequent dynamics of ice-edge bands, and by Squire (Reference Squire1989) for the special case of obliquely incident waves at super-critical angles. In each of these papers either perfect or zero absorption was assumed so that the issue of the accurate calculation of the reflection coefficient did not arise.

The reflection and transmission coefficients

The new theoretical solution of Fox and Squire (Reference Fox and Squire1990, Reference Fox and Squire1991a,Reference Fox and Squire1991b) and Squire and Fox (Reference Fox and Squire1990) provides accurate estimates for the reflection and transmission coefficients at the edge of an ice sheet as a function of wave period, ice thickness and water depth. For completeness, we include here a brief account of the method of solution for normally incident ocean waves; the more complicated case of the reflection of obliquely incident waves has now been solved and will be reported in subsequent publications.

With normal incidence, the problem is essentially two-dimensional and we may consider an open water region x < 0, 0 < y < Ή with the usual free surface boundary conditions, and an ice-covered region, x > 0, 0 < y < H with suitably modified boundary conditions. The ice sheet is assumed to behave as a continuous thin elastic plate of uniform thickness, with zero bending moment and shear at its free end. Thus, assuming linearity and irrotational flow, a wave-like velocity potential ϕ(x,y)e^iwt is to be found such that ϕ(x, y) satisfies Laplace’s equation ∇²ϕ = 0 for −∞ < x < ∞, 0 < y < H, with all boundary conditions satisfied. It is assumed that ϕ(x, y) may be expressed as the sum of separated solutions (modes) of Laplace’s equation, where the modes are found as the products of eigenfunctions of the separate x-variable and y-variable differential operators in that equation, chosen to satisfy the appropriate boundary conditions. In the open water region, the modes have the form e^kx cos ky where k satisfies

This equation has two imaginary roots ±ik_T corresponding to travelling waves and infinitely many real roots ±k_n (n = 1,2,…), where

corresponding to evanescent modes. These modes form a complete set and are orthogonal. In the ice-covered region the modes have the same form (e^kx cos ky), and k satisfies

This dispersion equation also has two imaginary roots ±ik_T which correspond to travelling waves, four complex roots

where * indicates complex conjugate, corresponding to damped travelling modes, and infinitely many real roots ±k_n (n = 1,2,…), where (n − 1)π/Η < k_n < nπ/Η, each giving an evanescent mode. These modes also form a complete set but are not in general orthogonal.

Including only the physically significant roots of Equations (1) and (2), i.e. those satisfying boundedness at large x, we may write down the potential for the open water region:

and for the ice-covered region:

To determine the reflection and transmission coefficients a ϕ(x, y) must be found such that ϕ_k = ϕ_i and ϕ_w = ϕ_i on x = 0,0 < y < H, noting that bending moment and shear vanish at x = 0⁺, y = H. The problem is solved variationally, after terminating the series in Equations (3) and (4) at a finite number of terms, by minimising the error term

on x = 0 subject to unit incident amplitude and the edge boundary conditions. The weight μ is adjusted to improve convergence. The error term ε is a quadratic form in the real and imaginary parts of the coefficients I, R, Τ, {a_n }, b₊ , b₋ and {b_n }. The minimisation is performed, with respect to these coefficients, using a preconditioned linear conjugate gradient algorithm (Gill and others, Reference Gill, Murray and Wright1981) restricted to the subspace in which the incident amplitude is unity and the edge conditions are satisfied. For more details see Fox and Squire (Reference Fox and Squire1990, Reference Fox and Squire1991a, Reference Fox and Squire1991b).

Radiation stress

The phenomenon of radiation stress in water waves is now well understood and the reader is referred to Longuet-Higgins and Stewart (Reference Longuet-Higgins and Stewart1964) for a complete treatment of the theory with some applications. The results of a later work (Longuet-Higgins, Reference Longuet-Higgins1977) are also relevant to the present discussion. That paper begins with the following few sentences which we repeat here as a concise description of the radiation stress effect. “Water waves transport both energy and momentum, and any solid body which absorbs or reflects wave energy must absorb or reflect horizontal momentum also. Hence the body is subject to a mean horizontal force. In low waves, the force may be calculated immediately when the incident, reflected and transmitted wave amplitudes are known.” Using the method of Fox and Squire (Reference Fox and Squire1990) to obtain accurate values for the transmission and reflection coefficients at an ice edge, Longuet-Higgins (Reference Longuet-Higgins1977) enabled the radiation stress to be found for water of finite or infinite depth. Longuet-Higgins (Reference Longuet-Higgins1977) considered the effect of two-dimensional, irrotational waves on any number of floating or submerged bodies of finite horizontal range; the bodies may absorb or generate energy at the same frequency. He showed that the horizontal force acting normally on the bodies is

per unit length, where the terms on the right are respectively the momentum fluxes in the incident, reflected and transmitted waves. On water of finite depth H, this expression becomes

if we assume unit incident amplitude. Note that the R and Τ coefficients both refer to waves travelling in open water, R is the amplitude coefficient of waves reflected from the floating bodies, whereas Τ is the amplitude transmission coefficient of waves in the open water beyond the floating bodies. The factor (1 + 2k_TH(sinh 2k_TH)⁻¹), which introduces the effect of water depth, varies between 1 for infinitely deep water and 2 for shallow water.

We suppose now that the ice cover may be represented as a continuum and that an elastic constitutive relation may be used to determine the energy balance at the ice edge itself. Then three types of wave, as described above, will exist in the ice: an undamped travelling wave; a pair of damped travelling waves; and an infinite set of evanescent modes. Only the undamped travelling wave will propagate any distance into the sea ice. The others will be rapidly attenuated away close to the ice edge. We also suppose that the travelling wave is gradually damped as it propagates through the ice cover so that at some penetration it too will have negligible amplitude. Thus, the transmission coefficient Τ in Equation (7) will be zero, as it represents the input of energy into open water beyond the ice sheet, and the equation will become

We have made a slightly different assumption to those of Martin and others (Reference Martin, Kauffman and Parkinson.1983) and Wadhams (Reference Wadhams1983), which were both concerned with bands of ice. In Martin and others (Reference Martin, Kauffman and Parkinson.1983) the band was assumed to absorb all the wave energy with no transmission or reflection, and in Wadhams (Reference Wadhams1983) zero absorption was assumed.

It is Equation (8) which is used to compute the force per unit length due to radiation stress at an ice edge.

The Incoming Wave Field

Monochromatic waves

We consider first the effect of monochromatic waves of various frequencies on ice covers of different thicknesses. The force per unit length for 1 m amplitude incident waves, as given by Equation (8), is found at three water depths; 10 m in Figure 1, 100 m in Figure 2, and 1000 m in Figure 3. In each case the curves are labelled with the appropriate ice thickness, and it is clear that the effect of radiation stress increases significantly as h is increased. This is due principally to the reflection coefficient R which, as might intuitively be expected, increases with increased ice thickness (sec Fig. 6, Fox and Squire, Reference Fox and Squire1990). Furthermore, Figures 1, 2 and 3 have some common features because F contains: the factor (1 + 2k_TH/ sinh 2k_TH) which always lies between 1 and 2, tending to 1 for short waves and 2 for long waves; and a factor (1+R ²) which tends to 2 for short waves and tends to 1 for long waves because the latter are relatively unaffected by the presence of the ice sheet. Thus as period → 0 or period

. Between these two extremes, F approaches ρg although it may not reach this value. Consequently, Figures 1, 2 and 3 are superficially similar, with each curve decreasing to a minimum from its asymptotic short period limit (as period is increased) and then increasing to its long period asymptote. It is only for thick ice on shallow water (i.e. for 5 m and 10 m in Fig. 1) that more structure is evident where, at intermediate periods, the curves overshoot the ρg asymptote

Fig. 1. The radiation stress in Ν m⁻¹ plotted against wave period (on a logarithmic, scale) for ice of thickness 0.5m, 1m, 2m, 5m and 10 m on water of 10 m depth. The incoming wave amplitude is 1 m.

Fig. 2. The radiation stress in N m⁻¹ plotted against wave period or various thicknesses of ice on water of 100m depth. The incoming wave amplitude, is 1 m.

Fig. 3. The radiation stress in N m⁻¹ plotted against wave period for various thicknesses of ice on water of 1000 m depth. The incoming wave amplitude is 1 m.

Fig. 4. The Pierson-Moskowitz spectrum assumed to act on the ice cover. The peak period is at 10 s (0.1 Hz) and significant energy is limited to the range 4–20s, i.e. 0.05 to 0.25 Hz.

before tending to this value from above. In Figure 3 all curves drop rapidly to their minimum and remain at this value until the incoming waves have sufficient length to feel the sea floor.

These calculations can, of course, be done at any depth with the caveat that a larger number of evanescent modes are required to reach a predefined accuracy for very deep water, but that there is little point in introducing additional depth for a wave which effectively already sees the water as deep.

Spectral input

Wc now consider a more realistic ocean wave forcing, namely that of an incident sea made up of waves at various frequencies. In line with earlier work, the spectral model we shall use is the Pierson-Moskowitz spectrum described by Phillips (Reference Phillips1977). This spectrum, which models a fully developed sea, is represented by the equation

It is assumed that a Pierson-Moskowitz spectrum with a peak frequency f_m = 1/10 Hz (Fig. 4) acts upon an ice cover floating on water of depth 100 m. Five thicknesses of sea ice are considered: 0.5 m, 1 m, 2 m, 5 m and 10 m. The radiation stress spectral density functions corresponding to each thickness at water depths of 10 m and 100 m are plotted in 5 and 6, respectively. The equivalent plot at 1000 m is not included, being very similar to that at 100 m. Because of the convex-up nature of Figure 4, as opposed to the concave-up nature of Figures 1–3, the radiation stress spectrum tends to be slightly less strongly dependent on frequency than the forcing spectrum, as the two effects counteract one another. The bell-shaped form is preserved in all cases, but there is a slight change in peak frequency. The significantly greater radiation stress spectral densities of Figure 5 as compared to Figure 6 are due to the (1 + 2k_TH/sinh2k_TH) factor in Equation (8). This approaches 2 for most waves in this shallow water.

Fig. 5. Force per unit length spectra for 10 m of ice resulting from the Pierson-Moskowitz spectrum of Figure 4.

Fig. 6. Force per unit length spectra for 100 m of ice resulting from the Pierson-Moskowitz spectrum of Figure 4.

An indication of the magnitude of the radiation stress effect is provided by finding the significant force per unit length, F_s , defined by analogy with significant wave height as

where m_o is the variance or area beneath the spectrum. This is done in Table 1.

Conclusions

The significant force per unit length estimates given in Table 1 totally dominate those which could arise due to wind stresses, current stresses, sea surface tilt or Coriolis force. Suppose, for example, we consider the effect of the normally incident fully developed Pierson-Moskowitz sea as described above on a 10 km² area of 1 m ice, and we then compare this with the force due to the action of a 10 m s⁻¹ wind, assuming an air-ice drag coefficient of 3 × 10⁻³. The force at the ice edge due to the incoming waves is approximately 1.3 × 10⁸ Ν for deep water, and nearly double this for shallow water. The force imparted on the ice due to the wind is 3.9 × 10⁷ Ν. Likewise, if we consider the effect of a current beneath the ice cover with, say, a water-ice drag coefficient of 7 × 10⁻³ and a current velocity of 0.2 ms⁻¹, we obtain a force of about 2.9 × 10⁷ N. In part, the significance of the wave radiation stress is due to the input forcing spectrum chosen, which has a large significant wave height (~ 4 m) as it represents a fully developed sea. However, even in more moderate seas, the radiation stress is a significant contribution which has been neglected hitherto.

Acknowledgements

The authors are indebted to the New Zealand Foundation for Research, Science and Technology, the University of Otago, and the New Zealand University Grants Committee for their continued support. This work was completed while Colin Fox was a University Grants Committee Post-Doctoral Research Fellow.