Introduction
Recent field studies have identified major changes occurring in the flow dynamics in the vicinity of Ice Streams B and C in West Antarctica (Reference Van der veen and OerlemansVan der Veen and Oerlemans, 1987; Reference Stephenson and BindschadlerStephenson andBindschadler, 1988). Mass-balance calculations on these ice streams indicate Ice Stream B is in negative balance by about 40% while Ice Stream C is in positive balance (Reference Shabtaie and BentleyShabtaie and Bentley, 1988; Reference Shabtaie and BentleyShabtaie and others, 1988; Reference Whillans and BindschadlerWhillans and Bindschadler, 1988). The discharge of these ice streams into the ice shelf affects the thickness distribution of the ice shelf and its pattern of flow (Reference MacAyeal and LangeMacAyeal and Lange, 1988). These changes in the ice shelf, in turn, affect the ice-stream flow by altering the longitudinal resistance exerted on the ice stream.Reference Thomas Thomas (1979) referred to this effect as “back-pressure”. Ice rises are one of the major sources of this back-pressure and, as such, are integral parts of the dynamics of the ice-shelf–ice-stream system.
Crary Ice Rise (CIR) lies directly in the path of Ice Stream B (Fig. 1).Reference MacAyeal and Lange MacAyeal and others (1987) have calculated the force exerted on Ice Stream B by CIR. They also have calculated (see MacAyeal and others, 1989) that the region around CIR is thickening at a rate of 0.13 ± 0.09 m/a. Other authors have identified other features of CIR which further suggest that CIR is not in equilibrium: thickness anomalies in the ice shelf down-stream of CIR (Reference JezekJezek, 1984); a non-equilibrium transverse surface profile (MacAyeal and Thomas, 1980); and separation of a large piece of the ice rise (a “raftD) (Reference Bindschadler, Stephenson, Roberts, MacAyeal and LindstromBindschadler and others, 1988b). In this paper we calculate the regional pattern of net mass balance in the vicinity of CIR to identify where rapid changes are taking place and discuss possible causes of the changes and their ramifications on ice flow in the area.
Method of Calculation
Net mass balance refers to the summation of all gains and losses in ice volume for a region of space. If this region is enclosed by the curve Γ, then the net mass balance, M, can be expressed as
(Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others, 1987, Equation (15)) where ρ̄ is the depth-averaged ice densjty, H is the ice thickness, H is the horizontal ice velocity, u is the unit vector normal to Γ, A is unit length along Γ, A is the accumulation (surface and base), and S is the horizontal area contained within Γ. The first term in Equation (1) is the mass change due to advection while the second term is the mass change due to accumulation or ablation of ice. As inReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others (1987), we calculate the average density as
where H is the ice thickness, ρ ice is the density of pure ice (917 kg/m3), and α and β are parameters used to fit Equation (2) to firn densities near the surface. For the measured density profile at J9, α = 608 kg/m3 and ß = 0.043 m−1; we use these same values.
Application of Equation (1) to actual data involved defining Γ and approximating the first integral of Equation (1) by a series of summations along line segments. The line segments usually connected occupied field stations where accurate geographical coordinates and ice velocity were known. Ice thickness was measured at intermediate points along most of the line segments and velocity data were interpolated to these points to permit calculation of mass flux across sub-segments.
While similar in concept, the method used in this paper to approximate Equation (1) differs in a number of ways from the method applied in Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others (1987). There are four main differences: scale, projection, interpolation of velocity along line segments, and mass-flux calculation for sub-segments. The effect of the differences is often negligible, but not always. Many times the difference is sensitive to the particular geometry of line azimuth and/or velocity gradient. It is a general requirement of net mass-balance calculations that calculation errors of mass flux be minimized because the net mass balance is usually a small difference between large, but very nearly equal, mass fluxes. Therefore, it is worth detailing the differences for the sake of comparing our results with Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others (1987) and for future calculations of mass balance.
Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others (1987) used an Earth radius of 6370997.0 m and a polar-stereographic projection of geographic coordinates. In this paper, the local surface is assumed spherical with a radius of 6357020.8 m. This radius is derived from the radius of the WGS-72 ellipsoid at lat. 84 ° and corrected for the geoid height at lat. 84 °S., long. 171 °W. using the Rapp Set A geoid. The difference between the values of Earth radii is 0.2%; all distances are different fromReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others (1987) by this amount.
InReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others (1987), the distortion caused by the polar-stereographic projection was minimized by choosing the pole of projection near the center of CIR, but, as with any polar projection, scale distortion increased with distance from the projection pole. This paper uses no grid projection. All distances are great-circle distances rather than straight-line distances on a plane of projection. Quantitatively, this effect is small near the pole (only 0.02% at a distance of 100 km from the pole).
Each line segment is divided into a varying number of sub-segments to allow an adequate representation of the ice-thickness variation. Velocity data were then interpolated to these same points. In this paper, the velocity magnitudes and azimuths are linearly interpolated along a given line segment. At each of these interior points, the exact bearing of the great-circle line segment was used to calculate the corresponding normal component of velocity. By contrast, in Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others (1987) an average bearing was used for each line segment, the velocity component normal to this bearing at each end of the line segment was calculated, and it was these two end values which were linearly interpolated to other points on the line segment. The latter method is not a true linear interpolation of the velocity field; its non-linear character increases as the angle between the line segment and velocity azimuth decreases. The use of a constant average bearing also introduces a non-linearity to the interpolation which depends on both the rate at which the true bearing changes along the line segment and the velocity gradient.
The final difference between methods concerns the calculation of mass flux for a sub-section of a line segment. With velocity, line bearing, and ice thickness all specified at each end of the sub-segment, we assume in this paper that the velocity, thickness, and density each vary linearly along the sub-segment,
where U i is the velocity magnitude, V i is the velocity component normal to the sub-segment, θ1 is the angle between the velocity azimuth and the bearing of the sub-segment, and γ is the fractional distance along the sub-segment λi to λi+1 The mass flux due to advection is then taken from the first integral in Equation (1),
Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others (1987) used a trapezoidal approximation which computes the mass flux per unit width at each end of the sub-segment and multiplied the average by the sub-segment length. In our notation, that approximation is written
Equation (10) ignores the cross terms which Equation (9) contains. The error in using Equation (10) will increase with the gradients of velocity and ice thickness.
Due to the complexities, and in some cases the peculiarities, of applying each method to specific data, an analytical treatment of the differences is not given here. Rather, because the same line segments used inReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others (1987) appear in this paper, a comparison of the mass fluxes derived by the two methods appears later in Table II.
The error in the calculated values of advective mass flux follows from Equation (9),
where δU, δθ, and δH are the errors in velocity magnitude, velocity azimuth, and ice thickness, respectively. In practice, there is also an error due to position errors: the position of some of the flight lines over which ice thicknesses were measured did not match precisely with the station coordinates. We incorporate this error in position as an error in ice thickness.
The net mass balance is then the sum of mass fluxes across line segments comprising any closed curve plus the accumulation collected within that curve.
Data
Figure 2 indicates the area of interest, the station names and areas over which the net mass balances are calculated. Detailed features in Figure 2 represent all crevasses and undulations seen on aerial photography of the area (Reference Bindschadler, Vornberger, Stephenson, Roberts, Shabtaie and MacAyealBindschadler and others, 1988b, figs 1 and 2). The location of the ice rise in the center of this area is evident from the lack of crevassing or undulations. Line segments comprising the perimeter of the area are identical with those used by Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others (1987) and connect sites occupied by surface field parties. Sub-divisions of this area are defined by connections to surface stations on the ice rise.
Table I lists the locations, velocities, and ice thicknesses at each station. This list includes stations on the ice rise, around the perimeter of the area, and some additional stations generated as described below. Positions of most stations indicated on Figure 2 are from occupation by satellite-positioning receiver. Velocities are from repeated geoceiver occupation or, where parenthesized, from interpolated methods described below. Ice thicknesses are from airborne radar soundings.
Positions and velocities of most of the perimeter stations were measured sometime between 1983 and 1985 (Reference Bindschadler, Stephenson, Roberts, MacAyeal and LindstromBindschadler and others, 1988a). For stations J101 J11, and I11, data were acquired during the Ross Ice Shelf Geophysical and Glaciological Survey (RIGGS) (Reference Thomas, MacAyeal, Eilers and GaylordThomas and others, 1984). The velocities of J11 and I11 were interpolated from adjacent RIGGS station velocities.
The perimeter-line segments avoid crevassed areas where high velocity gradients can be expected. However, this is not true of the interior segments which cross on to the ice rise. On these line segments, additional stations were introduced to permit a more variable gradient of velocity. The character of the crevassing is distinctly different on opposite sides of CIR and a different scheme was used on each side. Figure 3 shows the velocity profile along each perimeter-line segment.
For the line segments on the south-west side of CIR (E3-J1, E4-J2, and LP1-J3). a velocity profile was constructed which reproduced the velocity profile measured at the station L1 stake network (Reference Bindschadler, Vornberger, Stephenson, Roberts, Shabtaie and MacAyealBindschadler and others, 1988b, fig. 7). This network, also on the south-west side of CIR, extends from a point on CIR 2.3 km from the ice-rise–ice-shelf boundary to stakes on the ice shelf 10km from the boundary (see Fig. 2). The measured profile of velocity can be closely approximated by linear interpolation between four points: zero velocity at the beginning of the line on the ice rise; 10 m/a at the ice-rise–ice-shelf boundary in a direction parallel to the boundary; 145 m/a at a point 4.4 km out on the ice shelf; and 266 m/a at a point 9.5 km out on the ice shelf. The velocity azimuths for these last two points were linearly interpolated between the azimuth of the boundary and the azimuth of the station velocity at the ice-shelf end of the line segment. This profile is assumed at each of the three line segments on this side of CIR, and Table I and Figure 3 show the intermediate stations which were generated to provide this approximation of the velocity profile. The errors assigned to the velocity magnitudes are the velocity errors of the ice-rise stations and the azimuth errors are set conservatively to ±5°, larger than any measured error. Motion on the ice rise is so small that satellite positioning produced velocities which were within one standard deviation of 0 m/a and azimuths which were unreliable.
On the north-east side of CIR, the velocity profiles along the lines E3–K1 and LP1–K3 were defined by adding two points each: one at the edge of the ice rise, where the velocity was set to zero, and the other point at the edge of the intense crevassing (see Fig. 2). The velocities at each of the latter points was calculated by extrapolating the component of velocity normal to the line segment using the appropriate shear strain-rate component from the measured strain-rate at the corresponding ice-shelf station (Reference Bindschadler, Stephenson, Roberts, MacAyeal and LindstromBindschadler and others, 1988a). The magnitude error was set to the error at the ice-shelf station. The azimuths of these velocities were set normal to the line segment. Although this method only estimates one component of the full velocity, it is only this component which is needed for the mass-balance calculations.
The final line segment extends from the point of the ice rise farthest down-stream along a “suture line” separating crevassed ice from rifted ice. This line appears to represent a natural flow boundary where ice flowing around the ice rise converges at the down-stream end. The ice-rise point, IR, was assigned a velocity of zero, while the velocity at the point J10′ (at the intersection of this line segment and J10–J11) was calculated from the linear interpolation of velocity and azimuth along J10-J11 and using the larger errors of J11. The azimuths of the IR-J10′ line segment and the J10′ velocity differed by only 4°, indicating that very little mass flowed across this line segment and that the interpolation scheme was accurate within the stated errors.
Ice thicknesses at each station and along most line segments were measured by radar sounding. These data were kindly provided to us by S. Shabtaie (personal communication). In 1984, airborne radar profiles were flown along all segments up-stream of and including K3-LP1-J3. Errors for these segments were estimated by Shabtaie to be ±10 m. Ice thicknesses at I11, J10, and J11 were measured during RIGGS but no radar flights between these stations were made. The ice-thickness profiles between these stations (including J3-I11 and K3-J10) were supplied by Shabtaie from a contour map he prepared which included all available ice-thickness data in this area. His data sources included some ground-based soundings as well as airborne measurements from RIGGS and earlier field studies. He estimated the uncertainty of ice thickness along these profiles as ±20 m.
Accumulation data were collected during RIGGS both on and around CIR. From dating of stratigraphic layers of increased radioactivity deposited during 1954 and 1964–65, Reference Clausen, Dansgaard, Nielsen and CloughClausen and others (1979) deduced a pattern of accumulation rate over the Ross Ice Shelf. The values around CIR ranged from 0.09 to 0.12 m/a (ice equivalent) with the single station on CIR somewhat lower (0.07 m/a). In Reference Thomas, MacAyeal, Eilers and GaylordThomas and others (1984), measurements on exposed stakes were analyzed and gave a similar pattern over the ice shelf but somewhat larger values: 0.09–0.16 m/a (ice equivalent) with 0.10 m/a on the ice rise. The general pattern is of increasing accumulation to the south-west but there is a large local variability. The lower value on CIR is likely due to relatively higher winds on the elevated ice rise. The relative transport of mass by the wind would cause the stations down-wind to exhibit enhanced accumulation. Rather than attempt to account for local variability which probably has not been adequately sampled with the RIGGS stations, we use a constant accumulation rate of 0.12 m/a and assign a probable error of ±0.03 m/a to cover the range of values measured by both groups.
Mass Fluxes
Table II lists the results of the mass-flux calculations over all line segments and compares our calculations with those ofReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others (1987). The only difference which exceeds one standard deviation occurs for line segment C2–0; most line segments agree to well within one standard deviation. The difference for the C2–0 line segment is caused by different position coordinates for station C2 used in this paper and MacAyeal and others (1987). We feel the coordinates in this paper are more accurate and are consistent with the satellite-receiver positions used to derive the velocity used in this paper. The only other difference worth noting is for line segment G2–KI. This difference is also caused by a slight difference in coordinates for station G2.
As expected, most mass enters through the line between G2 and C2 (408 842 kg/s) and most departs between J10 and III (292 901 kg/s). On the sides of the perimeter, the general flow is outward up-stream of CIR and inward on the down-stream side. The net results of the flow around the sides are a loss of 51 511 kg/s to the south-west and a loss of 41 467 kg/s to the north-east.
To ease the comparison of mass flux along line segments of different length, Table II includes the mass flux/km. Across the up-stream line segments, the two largest values occur. Also, there is a strong gradient in mass flux/km increasing from G2 to C2. The mass flux leaving the perimeter between the down-stream stations is less than enters up-stream. This is due to a combination of thinner ice and lower velocities (see Fig. 3). Along the sides of the perimeter, the mass flux is considerably less than either the up-stream or down-stream ends. For the interior lines, the sharp gradient of flow in the boundary layer around the ice rise is evident.
Net Mass Balance
Table III and Figure 4 present the net mass balance (Equation (1)) for a number of areas defined by various combinations of line segments. Our calculations agree quite closely withReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others (1989) that the advective mass fluxes nearly balance for the complete area (22 962 ± 24 639 kg/s from Table III; 8800 ± 26 981 kg/s inReference MacAyeal, Bindschadler, Shabtaie, Stephenson and Bentley MacAyeal and others, 1989) but the table makes it clear that within the complete area there are significant redistributions of mass taking place. Mass is being added to areas 1, 3, and 4 at high rates while mass is being removed rapidly from areas 5 and 6.
The advective mass fluxes are large enough that adding the accumulation flux to obtain net mass flux modifies the pattern only slightly. Data on the basal conditions are not available so the results of Table III for net mass balance assume negligible basal melting or freezing. It is possible that some ice accumulation may be added by freezing of water exposed to the air within the large rifts which occur in area 5 down-stream of CIR, but the areal extent of these rifts is limited and unlikely to alter the strongly negative net mass balance. Similarly, any gain or loss of mass at the base of the ice shelf would have to be large to alter significantly the strong regional differences shown by Table III and Figure 4.
Discussion
Although the net mass balance of the combined areas is nearly zero, the regional differences do not represent a redistribution of mass already within the perimeter. This would require the regions of mass gain to be down-stream of the regions of mass loss. Instead, the situation is one of mass being added in the up-stream and south-west areas from the ice streams and mass being removed at the downstream end by the flow of the ice shelf. This gradient will encourage ClR to migrate up-stream as thickening up-stream promotes grounding and thinning down-stream promotes ungrounding. From the radar data, the bedrock underneath the up-stream part of CIR dips up-stream at an average grade of 2 × 10−3. Using the average thickening rate of 0.76 m/a for area 1, the predicted rate of grounding-line movement up-stream is 380 m/a. Reference JezekJezek (1984) had an earlier hypothesis of grounding-line retreat based on deviations of the velocity directions from the trace of debris tracks but made no estimate of retreat rate.
The retreat of the grounding line down-stream appears to be more episodic, related to the creation of rifts and the removal of rafts of ice. This zone of rifts differs markedly from the adjacent ice across the “suture zone” (line segment IR-J10′; Fig. 4). Table III shows that there is little, if any, mass flow across this boundary. The natural boundary appears to be wave-like and suggests a periodic variation in rift-formation activity. Each full cycle measures about 11 km. Using an ice velocity of 444 m/a at J10′, a period of 25 years per cycle results. It is not known what mechanism would be responsible for driving rift formation at this period.
There may be a second, much longer cycle of rift- (or raft-) formation activity.Reference Jezek Jezek (1984) noted “dome-hollow” pairs in the ice-thickness map of the Ross Ice Shelf down-stream of CIR. From their distance down-stream and the current velocity field, he estimated their dates since leaving CIR to be 400 and 550 years. These features may be evidence for historic activity similar to what is occurring at present.
There is also a significant transverse gradient of net mass balance: net gain of mass on the south-west side of CIR and net loss on the north-east side. This pattern is consistent with the calculation that the net force exerted by the ice shelf on CIR is directed at an azimuth of 0.7° true (Reference MacAyeal, Bindschadler, Shabtaie, Stephenson and BentleyMacAyeal and others, 1987, table IV). This orientation would promote mass deposition on the southern side of CIR and mass removal on the northern side. It also would help to explain the difference in the character of the crevassing on the two sides of CIR as well as the raft separation and rift formation in the lee of CIR.
The orientation of the net force is a result of the ice-thickness distribution on the ice shelf. Reference Bentley, Clough, Jezek and ShabtaieBentley and others (1979) showed that the Ross Ice Shelf was 100–150 m thicker on the south-west side of CIR than on the northeast side. This is probably a result of the thicker ice discharged into the ice shelf by Ice Stream A which flows in a much deeper channel than Ice Stream B (Reference Shabtaie and BentleyShabtaie and Bentley, 1988).
With this new knowledge of the pattern of net mass balance around CIR, certain features of CIR can now be explained. The general surface topography of CIR consists of a long ridge oriented along the line IR–LP1–E4, a dome north-east of the ridge roughly centered on the line segment LP1–K3, and a broad depression between the ridge and the dome (Reference Shabtaie and BentleyShabtaie and Bentley, 1988; Reference Bindschadler, Koci and IkenBindschadler and others, in press). The major peculiarity with this shape is the absence of any high spots underlying the long ridge, in contrast to the dome which is underlain by a submarine rise. The ridge rests on bedrock which dips to the south-west at a slope of 6 × 10−3. The existence of the ridge can now be explained by the excess mass accumulating on the south-west side of CIR. Also, the depressed topography between the ridge and the dome is the result of the removal of mass in the embayment between the ridge and the dome.
The fact that the present regional pattern of net mass balance provides an explanation of the current shape of CIR implies that this pattern has persisted for a substantial length of time. At its highest point, the ice thickness exceeds the flotation thickness by 80 m. At the current maximum rate of thickening given by Table III (0.88 m/a for area 4), this would take nearly 100 years. This age estimate should be viewed as a lower bound because local thickening rates would probably not be sustained at this maximum rate. On the other hand, measurements of the temperature profile through CIR indicate that the initiation of the ice rise occurred within the last few hunded years (Reference Bindschadler, Koci and IkenBindschadler and others, in press). These two techniques, mass-balance calculation and temperature measurement, may provide a powerful method to bracket the age of CIR within narrow limits.
Acknowledgements
Our thanks are given to the many field workers who assisted in the collection of these data. The description of this work was aided by comments from C.R. Bentley, C.S. Lingle, S.N. Stephenson, P.L. Vornberger, I.M. Whillans, H.J. Zwally, and two anonymous reviewers. The work was supported by U.S. National Science Foundation grants DPP-8514543 and DPP-8614407.