Introduction and Background
Theoretically, a large ice sheet resting on a flat base will flow in a pseudo-radial fashion from its ice divide. In reality, however, the regular flow of ice sheets is disturbed by the effect of basal topography on ice dynamics. Variations in ice flow can result in changes to the surface morphology of an ice sheet. In central regions of East Antarctica, some of the most obvious surface morphological features result from ice flow over subglacial lakes. These basal water bodies, often hundreds to thousands of square kilometres in area, manifest themselves on the ice-sheet surface as distinct smooth, flat regions (Reference Cudlip and McIntyreCudlip and McIntyre, 1987; Reference Ridley, Cudlip and LaxonRidley and others, 1993; Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996; Reference Siegert and RidleySiegert and Ridley, 1998). However, the influence of subglacial lakes on ice dynamics has yet to be quantified fully.
The largest subglacial lake, by an order of magnitude, exists beneath, and up to 230 km north of, Vostok station, central East Antarctica. Vostok lake has an area of 14 000 km2 and a volume of approximately 1800 km3 (Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996). European Remote-sensing Satellite (ERS-1) radar altimetry shows that the margins of the flat region above the Vostok lake are located within <1 km of the lake shore, as determined from radio-echo sounding (RES) (Fig. 1). Subglacial lakes of a considerable size are easily distinguishable from the surrounding ice-sheet base on RES data due to their bright, flat radar reflection when compared to that of bedrock (Reference Oswald and RobinOswald and Robin, 1973; Reference McIntyreMcIntyre, 1983; Reference Siegert, Dowdeswell, Gorman and McIntyreSiegert and others, 1996).
Fig. 1. ERS-1 altimeter-derived surface for the ice sheet over Vostok lake; contour spacing 10 m (after Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996). The lake area is indicated by shading Our flowline starts at position X (A) and ends at position ϒ.
The surface and subsurface slopes of the ice sheet along the 230 km length of Vostok lake indicate that the ice sheet is effectively floating on the lake water (Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996). Moreover, due to the lack of friction at the ice–water boundary, the basal shear stress above subglacial lakes is negligible. Thus, the subglacial lake shore represents an ice-sheet transition zone in which grounded ice, with considerable basal friction, changes to floating ice with zero basal drag. We conclude that this grounded/floating transition at the edge of large subglacial lakes is comparable to the grounded/ floating transition at the ice-sheet/ice-shelf boundary.
In this paper, we use a numerical model of the ice-sheet/ice-shelf transition (Reference MayerMayer, 1996) to calculate the ice-sheet velocity, temperature, stress and strain fields across Vostok lake. Measurements of the ice-sheet surface and basal topography are used to construct a conceptual flowline over which the model runs. The results are then compared to the pattern of internal RES layering in order to derive quantitative information about (a) ice flow over Vostok lake, and (b) rates of basal melting and freezing caused by the lake.
Glaciological Setting
The regional flow of the ice sheet around Vostok lake is from west to east. The ice sheet has a surface slope of ∼0.08° in this area. Over the subglacial lake, however, the ice-sheet slope is 0.01° from north to south (Fig. 1). Ice thickness is about 4200 m over the northern end of the lake and steadily decreases to about 3700 m at the southern end, close to Vostok station. RES data indicate subglacial hills surrounding the lake (Reference Robin, Drewry and MeldrumRobin and others, 1977). These data show that the thickness of ice surrounding the lake varies by as much as 500 m in 5 km (Fig. 2b). The surface elevation of the Vostok lake region (Fig. 1) shows interesting morphological features along the western and eastern margins of the lake that may be related to the ice flow across the lake. The ice sheet over the western lake margin is characterized by a small (approximately 2 m) trough in the ice-sheet surface, about 10 km wide and 200 km long. In contrast, a few kilometres east of the eastern margin of the lake, a small rise occurs in the ice surface which is about 10 m high, about 10 km wide and >200 km long (Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996; Reference Siegert and RidleySiegert and Ridley, 1998).
Fig. 2. RES profiles across the lake (from Reference Robin, Drewry and MeldrumRobin and others, 1977) as indicated in Figure 1: (a) 60 MHz RES profile along line AB; (b) 60 MHz RES profile along line CD. In (a), internal structures are clearly detected across the left (western) lake margin. The bending of the right part of all the reflections originates from an increase in flight altitude.
No systematic mapping of Vostok lake has been done. The only RES data available are from the 1970s survey of Antarctica, held in archive at the Scott Polar Research Institute, University of Cambridge, U.K. Two flight-lines from this database cross the western margin of the lake. Across this lake boundary, AB (Fig. 1), RES data show that the ice-sheet internal layers bend downwards towards the base of the ice sheet by 500 m (Fig. 2a). The location of this apparent subsidence corresponds with the position of the shallow surface trough. In contrast, the internal ice-sheet layering across the eastern margin of the lake, CD (Fig. 2b), does not display any recognizable bending. Thus, the topographic rise in the ice-sheet surface over the eastern shore of the lake does not correlate with any recognized internal structure of the ice sheet.
Model Flowline and Boundary Conditions
There are only eight RES flight-lines flown across the lake that yield information on ice thickness (Reference Robin, Drewry and MeldrumRobin and others, 1977). The areal coverage of these data is not dense enough to allow a gridded interpolation. However, assuming hydrostatic equilibrium for the floating part of the ice sheet, it is possible to calculate an areal ice-thickness distribution for this region from the observed surface elevation (Fig. 1). The RES ice-thickness profiles can be used to calibrate the distribution of floating ice over the lake.
In addition, RES profiles aligned along the ice-flow direction can be used to formulate model boundary conditions. In general, ice flow is dominated by the driving stress generated by the surface slope (Reference PatersonPaterson, 1994), which determines the direction of the ice flow. One of the RES profiles (flight 123, 1974/75, A–B) crosses the contour lines of the surface elevation in the west of the lake almost perpendicularly in a general direction of 72°. This RES line is used to provide ice-thickness data for our flowline calculations. The ice-sheet surface gradient to the west of the lake is almost constant along the entire lake margin, so that parallel flow can be expected, justifying the application of a two-dimensional flowline model.
A comparison of the surface slopes on and around the lake area reveals a ten times higher gradient on the adjacent grounded region than on the floating ice. The direction of the slope over the lake is from north to south, with a gradient of −1.6 × 10−4. If the ice flow were in the direction of the surface slope over the lake, the entire inflow of ice across the western lake margin would leave the lake over its southeastern margin. Assuming mean flow velocities of about 3 m a−1, an ice flux of approximately 3 × 109 m3 a−1 occurs across the 280 km long western margin of the lake. In addition, the accumulation of 0.03 m a−1 ice equivalent contributes a mass flux of 4.2 × 108 m3 a−1 over the 14 000 km2 area of the lake (Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996). If all the ice left the lake across the southern 50 km wide margin, a mean ice velocity of 18 m a−1 would result at Vostok station, which is in disagreement with the observed value of 3.7 ± 0.7 m a−1 (Reference Liebert and LeonhardtLiebert and Leonhardt, 1973).
This indicates that the flow direction over the lake is controlled by the grounded-ice flow regime, from west to east. Therefore, a flowline was constructed across the lake from west to east, with only a slight southwards deviation over the centre of the lake (Fig. 3). Unfortunately, none of the RES flights shows ice-bottom reflections across the entire lake. Therefore, part of the flowline had to be derived from the surface elevation data using the assumption of hydrostatic equilibrium over the lake. Since the surface topography allows the clear determination of the grounding line, the ice thickness prior to regrounding was calculated at 3950 m. East of the lake, there are very few measurements of ice thickness. However, we need to account for grounded ice east of the lake for our flowline to be complete. We therefore use an estimate of the grounded-ice thickness based on a simplified relation between ice thickness, surface slope and a yield stress of 55 kPa (Reference PatersonPaterson, 1994, p. 240):
Fig. 3. Surface topography in the central part of Vostok lake and its surroundings. The thick dashed line is the constructed flowline for the model calculations. The thin dashed line indicates the part of the RES profile where the bed reflection cannot be detected.
We note that the value of ice thickness east of the lake does not adversely affect model results over the lake itself. The geometry of our model flowline is shown in Figure 4 (where the ice flow is directed from X to Y).
Fig. 4. Hypothetical flowline geometry across Vostok lake constructed from available RES data (A-B) and ERS-1 data (entire flowline). The ice flow is directed from X to ϒ. The ice thickness between B and T is estimated as explained in the text. The lake is located at 27–84 km. The upper panel shows the detailed surface topography in high vertical resolution; the lower panel shows the entire ice-sheet geometry.
The thermal boundary conditions are governed by an annual mean temperature of −55°C at the surface, a geothermal heat flux of 54.6 mW m−2 and temperature at the pressure-melting point at the floating-ice subsurface. The accumulation rate in the centre of Antarctica is very low. In our model we used the accumulation rate at Vostok station (about 80 km from our flowline) of 0.027 kg m−2 a−1 w.e. (Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996). The thermal conditions were calculated from an initially linear temperature gradient between surface temperature and pressure-melting point at the bottom.
To summarize, the model requires boundary conditions in the form of surface and subglacial elevation along an ice flowline, surface temperature, basal heat flux and surface accumulation. Model results are compared with internal RES layering along the flowline. As discussed, these boundary conditions are established from available RES and field measurements. The largest source of potential error is the location of our ice flowline. However, recent interferometric synthetic-aperture radar data of ice flow over Vostok lake indicate that our flowline is suitably positioned over the lake (personal communication from R. Kwok, 1999). We are therefore confident of the validity of our model boundary conditions.
Description of the Model
Ice-sheet flow is governed by shear stresses, whereas in floating ice longitudinal deviatoric stresses dominate. In an ice shelf, longitudinal deviatoric stress gradients are counteracted by lateral drag along the sides of the shelf. The influence of lateral drag on a central flowline reduces with increase in the ratio of width to length. Our flowline is situated almost in the centre of the lake, with a relatively high width-to-length ratio. Thus, numerical investigations of the flow pattern over a subglacial lake must include longitudinal deviatoric stresses. It is expected that a clear transition of the stress regime along the flowline exists across the grounding/floating transition, similar to that measured in an ice-shelf transition zone.
The numerical model used to determine the velocity field on a vertical plane along the flowline was developed from a three-dimensional, isothermal version described in Reference MayerMayer (1996). This model is able to calculate diagnostic and prognostic solutions of ice dynamics across the transition zone between ice sheets and ice shelves. It was reduced to a two-dimensional version, which was coupled with the thermodynamical temperature solution by the flow parameter (A) of Glen’s flow law for ice.
The numerical formulation is based on the two-dimensional Stokes equations:
All physical quantities used for the description of the model are explained in Table 1.
Table 1. Physical quantities and constants used in the formulation of the model
The full stresses were separated according toReference Van der Veen and WhillansVan der Veen and Whillans (1989) into resistive stresses and litho- (cryo) static stresses:
where H 0 is the overburden ice thickness. The resulting set of equations can be expressed in terms of deviatoric stresses τij , the cryostatic load gradient −ρg∂ (H + h − z)/∂x and the resistive stress in the vertical Rzz (using Rxx = 2rxx + Rzz (Reference Van der Veen and WhillansVan der Veen and Whillans, 1989)):
where H is the ice thickness, h is the ice-base elevation and z is the vertical coordinate, positive upward.
Using the Leibniz rule for changing the order of differentiation and integration and the upper boundary condition according to Reference Van der Veen and WhillansVan der Veen and Whillans (1989), the vertical resistive stress can then be expressed as:
For closing the problem, the incompressibility equation
is used to derive the vertical velocity in combination with the appropriate boundary condition.
The relation between the stresses and the resulting strain rate is described by Glen’s flow law (Reference GlenGlen, 1955),
, which is reformulated to:
where
is the effective deformation rate and A is the Arrhenius-type flow parameter:
It is convenient to scale the mathematical formulation in order to allow a simpler numerical treatment. The ice bottom is defined as s = 0, whereas the ice surface equals s = 1, which gives a scaled vertical coordinate:
The transformed system of equations is then solved in a finite-difference scheme, using central differences and a point relaxation method for iteratively obtaining the velocities u, w in horizontal and vertical direction.
As a boundary condition there exists no coupling between ice surface and atmosphere, so that both stress components (normal and shear stress) equal zero:
where
is the resulting stress vector and
is the normal unity vector of the surface (Reference MacAyealMacAyeal, 1997). At the bottom of the floating ice, the hydrostatic pressure necessary to float the ice is equal to the isotropic pressure exerted by the lake water almost everywhere. The deviation from hydrostatic conditions is expressed by the vertical resistive stress Rzz
(Equation (6)). Because water cannot exert shear stress on the ice, the shear stress must vanish at the ice–water interface:
In grounded ice, the velocity at the ice-sheet base depends on the thermal conditions. For a cold-based ice column the velocities at the lower boundary are zero. If the basal temperature is at the pressure-melting point, a sliding relation must be introduced. A basic empirical relation between the sliding velocity u s, the basal shear stress τ b and the effective pressure N is:
where p and q are positive integers (here p = 3 and q = 1). k depends on thermal and mechanical properties (Reference PatersonPaterson, 1994, p. 151) and is very small, of the order of 0–1.0 × 10−15 m s−1 Pa−2. For assumed hydrostatic conditions over the lake (which is realistic for almost the entire lake), the effective pressure can be calculated in terms of a local height above buoyancy h * for a theoretical lake surface h 1 (h * = H + ρ w/ρ(h − h 1)):
To avoid a mathematical conflict at the grounding line, where N becomes zero and therefore the basal sliding velocity would be infinite, a boundary value for N close to the grounding line was introduced. This basic relation between basal sliding velocity, basal shear stress and effective pressure is used in several ice-sheet models (Reference Budd, Jenssen and SmithBudd and others, 1984; Reference HuybrechtsHuybrechts, 1990), but is only a rough parameterization for largely unknown conditions at the lower boundary (Reference PatersonPaterson, 1994).
The temperature field in the ice depends on boundary conditions of (1) the mean annual temperature at the surface, and (2) the geothermal heat flux at the bottom of the grounded ice sheet. The lower boundary of the floating ice is assumed to be at the pressure-melting point (Reference Siegert and DowdeswellSiegert and Dowdeswell, 1996). The time-dependent temperature distribution can then be calculated (Reference Huybrechts and OerlemansHuybrechts and Oerlemans, 1988):
which is governed by vertical diffusion, horizontal and vertical advection, and internal deformation, caused by horizontal shear strain rates. Equation (16) is in the vertically scaled form and is used for coupling with the temperature-dependent flow parameter A. Units used in Equations (1–16) are summarized in Table 1.
Model Results
Modelling the dynamics of a viscose material, such as ice, is very sensitive to boundary conditions. For the Vostok lake area, we have identified an ice flowline with appropriate boundary conditions of surface elevation (from the ERS-1 altimeter), ice thickness and subglacial topography (from RES). RES data also provide good information regarding particle flow within the ice sheet. The RES profile (Fig. 2a) shows clear internal reflections in the ice which are interpreted as isochrons (Reference SmithSmith, 1996; Reference Siegert and RidleySiegert and Ridley, 1998). In the lower regions of ice sheets, ice flow is parallel to the base of the mobile part of the ice sheet. Our RES data indicate that isochrons are also parallel to the ice-sheet base (Fig. 2). Thus, we can infer particle flow paths from the pattern of isochronous surfaces in the lower regions of the ice column (e.g. Reference SmithSmith, 1996). We contend that a model detailing ice flow across Vostok lake must be capable of matching particle paths with these internal structures.
An inverse-type modelling procedure was employed where particle-path results were forced to be compatible with internal RES layering through adjustments to model parameters. The model parameters varied in this procedure included (a) the basal sliding algorithm, (b) ice rheology and ice temperature, and (c) the addition and subtraction of mass at the ice-sheet base.
One problem encountered in our investigation of model boundary conditions was the lack of ice-velocity information along our flowline. However, ice velocity at the inflowing end of the flowline is required as a model boundary condition. Assuming that the velocity is of the order of that measured at Vostok station (3.7 m a−1), we performed a series of preliminary model runs in which a variety of inflow ice velocities less than 3.7 m a−1 were used. Model results were then examined with respect to mass conservation across the flowline. A reasonable steady-state solution was obtained for an inflow velocity of 2.2 m a−1 and this was used in all model runs presented.
Mass conservation and velocity fields
Figure 5a shows the horizontal velocity distribution along our constructed flowline. The ice flow is relatively undisturbed over the lake and does not show the typical characteristic of a floating ice shelf with free expansion towards the open sea. The velocity even decreases towards the lake from 2.1 m a−1 to 1.8 m a−1 and remains almost constant across it. The subglacial bedrock perturbations upstream of the lake lead to higher velocities above them to satisfy continuity of ice flux. On the lee side of the bumps (in relation to the flow direction) the basal velocity reduces to almost zero. The reduced basal velocity close to the western grounding line indicates a downwelling process over the bedrock bump nearest to the lake. The ice thickness across the lake reduces by 26 m over the 57 km length of the lake (20 m of this reduction happens within the first 2 km; Fig. 4). Thus, the horizontal deformation rate
is close to zero over the lake.
Fig. 5. Calculated velocity distributions (m a−1) across the flowline according to the boundary conditions described in the text: (a) horizontal, (b) vertical.
East of the lake, the horizontal velocity increases in response to the steeper ice-sheet surface gradient. It should be reiterated that the flowline east of the lake was constructed theoretically and has a relatively smooth ice–bed-rock interface. Consequently there is no variation in the pattern of basal friction to the east of the lake. It seems likely that the real flow of ice on this side of the lake is more similar to the pattern on the western lake margin, with a strongly undulated bedrock profile. The increase of the mean ice velocity from that over the lake is governed by the surface slope derived from ERS-1 data, not by the constructed bedrock topography. Thus, the bedrock geometry across the eastern side of the lake does not affect the main results of this paper.
The vertical velocity distribution (Fig. 5b) displays a pattern typical of the transition between grounded and floating ice. The largest gradients occur close to the western grounding line. Across the transition zone, the vertical ice velocity changes from −0.27 m a−1 at the ice-sheet base to very small values within 5 km. Here the ice slides down the steep bedrock bump (Fig. 4) and thickens considerably at the western edge of the lake.
If the glaciological situation was that of a parallel-sided ice shelf, open to the ocean at one side, negligible basal friction would lead to an immediate thinning of the ice column. In the case of Vostok lake, however, the ice regrounds after 57 km and thus is not freely floating. In ice shelves the horizontal gradient of τxx controls the thinning rate. In the case of Vostok lake the horizontal gradient of τxx vanishes across the floating section of the ice sheet. Therefore, ice flow across the lake is controlled by the slow increase in ice velocity past the eastern ice margin. The ice sheet over the lake acts rather passively in terms of ice dynamics, but controls the influx of ice across the western margin due to the thickness of ice in this region.
Stress distribution
The comparison of longitudinal stress and shear stress within the ice (Fig. 6) reveals the different stress regimes between grounded and floating ice. The stress regime deviates considerably from that of freely floating ice shelves. There is no distinct horizontal gradient in the deviatoric stress until the ice leaves the area of the lake. The floating part of the ice sheet acts as a region of stress conversion from compression in the west to extension in the east, with very low stress values over the lake. A more typical behaviour of grounded/floating ice is revealed by the shear stress. The western part of the flowline shows an alternating pattern of high and low shear-stress values. The floating part itself must of course be shear- stress free, and the continuation on the grounded part (constructed from a basic model, described in the flowline section) behaves more like a “shelfy stream”, with almost no basal friction (Reference HindmarshHindmarsh, 1993).
Fig. 6. Stress distribution (kPa) along the flowline: (a) longitudinal deviatoric stress, (b) horizontal shear stress.
Based on these model velocity distributions, the western grounding/floating transition zone is approximately 5–10 km long, corresponding to 2–3 ice thicknesses. This value fits well into the range resulting from theoretical work (Reference Herterich, Van der Veen and OerlemansHerterich, 1987; Reference MayerMayer, 1996).
The undulating bedrock surface, together with the associated differences in basal sliding, initiates strong shear-stress variations, which are highest in the lower 25% of the ice column but are also present in the upper regions. These values might be large enough to create fabric changes in the ice, which then could lead to difficulties in the analysis of ice cores.
Basal melting conditions from particle paths
Our inverse modelling procedure allowed us to experiment with a variety of basal-sliding, ice-rheology and basal mass-balance configurations. We found that the pattern of particle paths is most sensitive to the subtraction and addition of material at the ice base. Further, under one configuration, particle paths were matched with RES layers. The resulting mass-flux and basal mass-balance pattern along the flowline is shown in Figure 7. The theoretical mass flux is calculated from the inflow over the eastern model boundary and the mass balance at the surface and ice base assuming steady-state conditions and mass conservation (e.g. the continuity equation). The modelled mass flux shows that the model matches these requirements everywhere along the flowline. Although there are several possible reasons why ice might be subtracted or added at the ice-sheet base, the measurements available are more compatible with “melting and freezing” events than other causes (see discussion below).
Fig. 7. (a) Mass flux along the flowline, calculated from the model results (solid line) and the continuity equation (dashed line with stars); (b) melt rates at the ice base.
Along the flowline, melting and freezing occurs only over the first 18 km of the floating section (Fig. 7). The melting acts to reduce the ice thickness over 2 km by 145 m, which is partly counteracted by enhanced inflow resulting in a net reduction of 20 m. Refreezing then compensates the dynamic thinning, and the ice thickness remains constant over the next 16 km. After this point the ice thins slightly from 3796 m to 3790 m before becoming grounded again.
Several internal radio-echo layers can be resolved across the lake (Fig. 2a). We use one of the strongest internal reflections, which can be traced across the western grounding line, as the standard of comparison for our calculated particle paths. Figure 8 displays the flowline geometry, including the internal reflection derived from the RES profile, and particle paths for a number of ice depths. The particle path which enters our flowline model at 900 m depth matches very well with the observed internal reflection. Across the bedrock bump west of the lake the particle path and the radar reflection follow the bedrock geometry as expected. Assuming no basal melting, the slope of the internal reflection decreases until it becomes flat over the lake. However, the RES data show a “dipping” of the reflector close to the grounding line, before it flattens over the lake. The only model run which matched particle paths with the observed dipping was when a reduction of basal material as the ice becomes afloat was included.
Fig. 8. Flowline geometry and calculated particle paths. The thick dotted line displays the internal reflection found in the RES profile.
Melting of 145 m of ice across the first 2 km of the lake requires a melt rate of 150 mm a−1. The basal conditions then switch to a freezing rate of 80 mm a−1 which decreases to zero after 10 km (Fig. 7b). We suggest that the mechanism responsible for this exchange of basal mass is similar to that which occurs beneath large ice shelves (Reference Lewis and PerkinLewis and Perkin, 1986; Reference ThyssenThyssen, 1988; Reference Jenkins and DoakeJenkins and Doake, 1991). Since our model calculates only ice dynamics, it cannot specify the melting/freezing mechanism. Nevertheless, the subtraction and addition of basal ice along our transect may be associated with circulation of the lake water. Further, our melting/freezing assumption has implications for the residence time, and thereby the age, of the lake water. If basal melting occurs along the entire western coast in the same order as in the flowline, roughly 8 × 107 m3 of ice will be melted every year. If 20% of the meltwater mixes with the total volume of the lake water (estimated to be 1800 km3 (Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others, 1996)), an equal amount of the lake water must refreeze according to our model. In this case, the entire water reservoir will be exchanged within 112 500 years. This result matches very well with the theoretical residence times given by Reference Kapitsa, Ridley, Robin, Siegert and ZotikovKapitsa and others (1996). If there is complete mixing the residence time will reduce to 22 500 years.
The internal RES layer shows a sharper depression than the modelled particle path, which implies a somewhat narrower region of basal melting/freezing than we calculate in our model. Still, the order of magnitude in area and the rates at which the processes operate agree very well between the model and the observed flow features. The melting of 150 mm of ice per year on a length of 2 km requires 9.14 × 107 J of energy, which must be provided by cooling the water at the ice–lake interface. The rather cold ice comes in contact with the lake water which has a temperature of at least 270.32 K, corresponding to the pressure-melting point at this depth. If melting occurs, the water mass adjacent to the ice base cools and becomes buoyant because fresh water has its maximum density at higher temperatures than the pressure-melting-point temperature (Reference GillGill, 1982, p. 599f). The cold water rises along the inclined ice subsurface towards the centre of the lake. During this ascent, the cold water reaches an area with a higher pressure-melting point, whereupon refreezing occurs.
Other explanations for internal layer patterns at the floating ice-sheet margin
We attribute the addition and subtraction of basal ice across our two-dimensional flowline to melting and freezing processes. However, a number of additional explanations require consideration. Similar RES layer patterns to those observed on the western margin of the Vostok lake (Fig. 2a) have been measured at the grounded-floating transition of Carlson Inlet, Fletcher Promontory, West Antarctica (Reference Vaughan, Corr, Doake and WaddingtonVaughan and others, 1999). In both cases, the bending of internal layers increases with ice depth, and the location of the layer pattern is associated with a change in the ice-surface slope. The cause of this layering pattern at the margin of Carlson Inlet is thought to be anomalous rates of accumulation at the grounded-floating transition. This may also be the case for the ice surface at the western margin of Vostok lake, and is a plausible explanation of how ice lost at the base is replaced on the surface. However, accumulation patterns alone cannot explain the pattern of layering at Vostok lake. This is because even though the deepest internal layer shows the bending characteristic, the ice-sheet base does not. Therefore, as ice flows across the margin, mass must be lost along our flowline between the deepest layer and the ice-sheet base.
Internal RES layers have also been shown to bend due to the three-dimensional flow of ice (Reference Robin and MillarRobin and Millar, 1982). This may represent a cause for the internal layer patterns observed across the western margin of Vostok lake. The idea behind this explanation is as follows. If the subglacial topography is not flat, ice may converge within subglacial valleys and diverge around subglacial hills. This may cause the thickness between internal layers to be large over valleys and small over hills. RES layers do therefore show bending in response to ice flow across bedrock obstacles. For the case of layering over the western margin of Vostok lake, it is conceivable that RES layer bending is caused by the following two conditions: (a) a subglacial valley exists at the margin of the lake where internal layering folds downwards, and (b) our RES transect is located slightly off the line of true ice flow. In this case, part of the transect may be derived from ice flow over a subglacial hill (causing internal layering relatively high in the ice column), and another part may originate from within a valley (resulting in internal layers relatively deep in the ice column). Thus, the internal layers would appear to fold from one part of the transect to another. Although we cannot discount this possibility, the limited evidence available tends to support our melting/ freezing proposal more than the three-dimensional ice-flow explanation. This is because if three-dimensional ice flow did cause folding of internal layers at the western margin of the lake, there would be no reason to expect the layers to recover to their original ice depth after folding had occurred. However, we have two examples of radar layering across the western margin of the lake which demonstrate that the layers actually do return to their original depths after folding. For our melting and freezing theory, this feature is necessary if the volume of melted water equals the volume of water refrozen. Thus, we think that internal layer bending across the western margin of Vostok lake is more likely to be caused by subglacial melting than by the three-dimensional ice flow upstream.
Discussion
Our model results are linked with RES data from the central part of the western shoreline of Vostok lake. The RES data show a distinct bending in the internal layering that can be modelled when zones of melting and freezing (in the order of several centimetres per year) are introduced. Internal layering features similar to those in Figure 2 occur across other parts of the western lake margin to the north of our flowline (Reference Siegert and RidleySiegert and Ridley, 1998). If these RES features are also caused by zones of melting and freezing, then a continuous zone of basal melting may exist from our investigated flowline, northwards, along the western margin of the lake. After the ice sheet has moved across the melting zone, the model predicts refreezing of lake water onto the ice-sheet base. Such refreezing would result in the accretion of several tens of metres of basal ice by the time ice regrounds at the eastern margin. Our two-dimensional model suggests that Vostok lake causes up to 145 m of ice to be converted from normal deposited glacial ice to refrozen basal ice across the centre of the lake. Thus, across at least the southern end of the lake and our flowline, the lower 4% of the ice column is converted from glacier ice to basal ice as a result of ice flow over Vostok lake. Because basal ice has physical properties different to those of normal glacier ice, its rheology may also be different. Thus, the actual dynamics of ice flow downstream of Vostok lake may be subtly different to those upstream.
Our two-dimensional model does not allow a prediction of the current “state of health” of Vostok lake. However, it does indicate that the lake is a relatively dynamic system which causes basal melting and freezing over 10 000 km2 of the Antarctic ice sheet. In our transect the net volume of meltwater matches with the net volume of accreted ice. Therefore, the lake exists in relative equilibrium state.
Our model indicates that Vostok lake is an active feature beneath the central East Antarctic ice sheet, rather than an inert object. Our results suggest the lake influences ice-sheet mass balance and the physical properties of up to 5% of the overlying ice column. However, under modern environmental and glaciological conditions, the relative flow of ice west of Vostok lake is similar to that to the east. Therefore, the overall net effect of this large subglacial lake on large-scale ice-sheet motion appears to be small despite the dynamism that our model suggests it induces.
Conclusions
A model of the ice-shelf/ice-sheet transition zone, which incorporates (a) a complete two-dimensional solution of the stress–strain field along a flowline, including longitudinal stresses within the ice sheet, (b) full coupling between the basal boundary conditions (expressed as relative drag and sliding velocity) experienced under the ice sheet and the internal ice rheology, and (c) full coupling between the solutions to the two-dimensional mechanical and thermal regimes of the ice sheet, is applied on an ice-sheet flowline across the subglacial lake near Vostok station. The model requires inputs of surface and subglacial topography, thermal boundary conditions and accumulation rates. Model results indicate a number of interesting features associated with the dynamics of the ice sheet as ice flows across the lake.
As expected, the ice flow across Vostok lake is restricted spatially by regrounding, and therefore does not resemble the flow characteristics of an ice shelf spreading into the open sea. Stretching of the ice flow does not occur when floating starts. However, as ice is transferred across the grounding/floating transition, the horizontal stress components are reduced from around 100 kPa to very low values, whereas the shear stress diminishes as expected within floating ice. We calculate that there is very little dynamism within the ice sheet as it flows across the lake, as is evident from the regular, steady velocity and stress distributions.
Our results indicate that the western transition between grounded and floating ice is between 5 and 10 km (2–3 ice thicknesses), whilst across the regrounding region the transition occurs very gradually over several tens of kilometres due to the idealized geometry and the basal sliding conditions. Within the transition zone there exist large stress gradients, which might be observed as fault structures or even crevassed zones at the ice-sheet surface.
Model results replicate particle flowpaths established from radio-echo layering when basal melting is included along a 2 km section of the flowline as soon as ice crosses the lake margin. After this zone of melting, we expect basal refreezing to occur. This should excite circulation and mixing within at least the upper layers of the lake. If 20% of the meltwater mixes with lake water before refreezing, we calculate the average residence time of liquid water to be around 100 000 years.
By analyzing RES data from other regions of the lake, we conclude that a zone of basal melting may occur from our flowline northwards along the western lake margin. In other regions, net freezing is predicted to occur. By the time ice regrounds across the eastern margin, up to 200 m of ice may have been converted from normal glacial ice to basal ice.
Acknowledgements
Funding from the European Science Foundation (EISMINT Exchange Grant 97–3) for collaboration during a visit to Aberystwyth is gratefully acknowledged. M.J. Siegert acknowledges funding from a College Research Grant, University of Wales, Aberystwyth. This is Alfred Wegener Institute contribution No. 1758.
