Selective erosion beneath the Antarctic Peninsula Ice Sheet during LGM retreat

Abstract

In mountainous terrain, the relationship between ice sheet dynamics and basal topography is complex, with each component influencing the other. This paper investigates how the last glacial maximum Antarctic Peninsula Ice Sheet might have modified its bed both at maximum extent and during progressive grounding line retreat. Focussing on the Marguerite Trough Ice Stream we then examine the degree to which basal topographical conditions affected the rate of ocean-forced recession. Zones of peak subglacial erosion are preferentially located in areas of convergent flow and where horizontal strain rates are highest. During ice sheet retreat, potential erosion rates increase in these areas, but the foci remain fixed. This leads to selective and progressive deepening of subglacial basins. As grounding lines migrate landward, faster retreat tends to occur over subglacial basins, especially if flow is divergent, whereas slower retreat takes place on sloping beds and where the geometry of the outlet allows convergent flow and a non-negative flux balance. In conclusion the Antarctic Peninsula Ice Sheet selectively erodes its bed beneath linear outlets and, over successive glacial cycles, progressive deepening of subglacial basins may bring about non-linear retreat of the ice sheet margin.

Introduction

The notion that ice sheets can differentially modify their beds through the course of repeated glaciations has been around since the earliest observations of terrestrial glaciations (Agassiz 1841), but systematic mapping, interpretation and synthesis of the geomorphological evidence indicative of such behaviour has advanced most notably through a comprehensive exploration of the subject spanning almost half a century (Sugden 1968, 1978, Jamieson et al. 2005, 2010, Sugden et al. 2005, Phillips et al. 2006). Collectively this body of work put forward, developed and tested the concept of ‘selective linear erosion’ by glaciers, in which successive glacial cycles repeatedly focused erosion beneath wet-based ice occupying deep troughs. Over time, this concentrates erosion along valley floors, whilst neighbouring upland areas remain relatively unscathed due to their cover of thinner, colder and less erosive ice. Following directly from the earliest of these insights and building on later expansions of the theory, many recent studies now report observations of differential subglacial erosion (Fabel et al. 2002, Briner et al. 2006, Stroeven et al. 2006, Licht & Palmer 2013) or describe numerical simulations that attempt to either quantify or explain the process (Alley et al. 2003, Bougamont & Tulaczyk 2003, Bougamont et al. 2003, Pollard & DeConto 2003). Furthermore, it is becoming increasingly apparent that subtleties of basal topography may in turn condition the behaviour of an overlying ice sheet, for example through stabilizing feedbacks, such as sedimentation at glacier grounding lines (Alley et al. 2003, 2007). In these scenarios, continued subglacial transport of sediment to the grounding line buffers the glacier from minor changes in sea level, building a sedimentary wedge that ‘fills in’ the gap beneath the ice-shelf base and the sea floor. This essentially makes the ocean shallower at such locations, and thus a greater sea level rise is required to destabilize the glacier margin. The ice sheet–topography system is therefore complex, and influences operate in both directions.

Numerical modelling presents an ideal approach with which to investigate ice sheet–topography interactions. Yet, until recently, computational demands have limited the spatial resolution of continental-scale ice sheet simulations to grids that allow tractable experimentation which, in the case of continents such as Antarctica, has typically been 20–40 km (Huybrechts 2002, Pollard & DeConto 2009). At this scale, major outlet glaciers are simulated, but finer details are lost. Thus numerical replication of selective linear erosion at a scale observable in the field has remained to some extent unachievable. Computational advances now make it possible to simulate continental ice sheets at a scale that can be compared to field observations (e.g. Golledge et al. 2013), meaning that geomorphological questions can now be addressed more easily.

The principal aim of this paper is to use high-resolution (5 km), whole-continent Antarctic ice sheet simulations to explore the changing geometry and dynamics of the Antarctic Peninsula Ice Sheet (APIS) through the last glacial termination. The specific objectives are to: i) identify the pattern of ice flow at the Last Glacial Maximum (LGM, c. 20 ka) and its consequence in terms of probable subglacial erosion, and ii) determine how the speed of ice sheet retreat following the LGM may have been modulated by topographical feedbacks.

Study area

The Antarctic Peninsula is characterized by a linear mountain chain and numerous outlying islands (Fig. 1a). Summits in Palmer Land rise to over 3000 m above sea level (a.s.l.), and accumulate ice preferentially on their western flanks as a consequence of a relatively high precipitation rate (Arthern et al. 2006, Lenaerts et al. 2012). Eastern slopes accommodate more modest glaciers due to the drier regional atmospheric conditions. The western continental shelf is more deeply incised than in the east, with well-developed linear troughs (e.g. Marguerite Trough) extending from present-day coastal areas to the continental shelf-break. Evidence from trough-mouth fans suggests that these troughs have been repeatedly occupied by ice streams throughout the Cenozoic (Livingstone et al. 2012).

Fig. 1 a. Bed topography of the Antarctic Peninsula and surrounding continental shelf, from Fretwell et al. (2013). Thin black lines show modern coast and grounding lines (Haran et al. 2005), dashed black line indicates modelled Last Glacial Maximum (LGM) margin position (Golledge et al. 2012a). AI=Alexander Island, BC=Bryan Coast, EIS=Evans Ice Stream, GVIS=George VI Sound, LIS=Larsen Ice Shelf, MB=Marguerite Bay, MT=Marguerite Trough. b. Modelled surface elevation and flow-lines of grounded ice (ice shelves shown in black). c. Modelled horizontal ice flux at the LGM (Golledge et al. 2012a). d. Transient grounding line positions during modelled retreat, driven by increasing oceanic temperature and rising sea level

Ice sheet model and experimental design

We use the Parallel Ice Sheet Model (PISM; Bueler & Brown 2009), a 3D Eulerian finite difference model that is purposely designed to enable high-resolution and long-term simulations of ice sheet flow. This model has been successfully applied to a wide variety of glaciological research, including high-resolution simulations of the Greenland and Antarctic ice sheets (Golledge et al. 2012a, Bueler et al. unpublished) and of mountain ice fields (Golledge et al. 2012b). PISM combines the shallow approximations of ice sheet (SIA) and ice-shelf (SSA) flow (Bueler & Brown 2009), so that grounded ice is influenced both by gravitational and longitudinal stresses, whereas floating ice is solely governed by the latter (i.e. spreading forces). This combination approach to ice sheet stress balance allows both sheet and shelf flow to be simulated, as well as that of transitional zones such as ice streams where ice is grounded (like an ice sheet) but has very low basal drag (like an ice shelf). An enhancement factor of three was used for the SIA, based on parameter tuning experiments (Golledge et al. 2012a), with no explicit calving of ice shelves imposed except to limit their extent to within the continental shelf area. Recently, the ability of 3D ice sheet models (such as PISM) to accurately track grounding line migrations has been investigated using simplified geometry experiments (Pattyn et al. 2013). Those simulations revealed that few finite-difference ice sheet models were capable of accurately tracking grounding line movements over short time frames and at coarse resolutions. In our implementation we hope to have mitigated some of these shortcomings by running long (multi-millennial) integrations and by using a relatively fine (5 km) spatial grid, but this remains an area where continued model refinement is required.

The simulations presented and described here are from a suite of experiments that were initialized from a thermally-equilibrated and dynamically-stable 20km resolution present-day configuration Antarctic ice sheet. This spun-up model was then allowed to re-equilibrate under environmental conditions representative of the LGM. The magnitude of imposed LGM atmospheric cooling was based on the Vostok d¹⁸O time series (Petit et al. 1999), oceanic cooling from stacked global benthic d¹⁸O records (Lisiecki & Raymo 2005), and long-term sea level changes from Imbrie & McIntyre (2006). From this steady-state LGM Antarctic ice sheet simulation we then simulated its dynamic response to a range of ocean-forcing perturbations in the form of sea level and oceanic temperature changes. Since the experiments were designed simply to examine the response of the modelled ice sheet to isochronous perturbations, arbitrary sea level increases of 25 and 50 m, and increases in sub-shelf melting (representing oceanic heat flux) of 33, 67 and 100% of the LGM to present difference were used.

Calculation of potential subglacial erosion

One way to infer the probable erosional effects of the LGM APIS is to consider subglacial erosion of bedrock as a function of basal ice velocity and the shearing force exerted at the bed (Pollard & DeConto 2003):

(1) $$\eqalignno { \theta \,=\,{\rm }\left( {{\rm 1}\,{\rm -}\,\left( {S/T_{{lim}} } \right)} \right)\,{\rm x}\,{\rm 2}e^{{{\hbox - }{\rm 1}0}} \tau _{b} U_{s} (1)$$

Here, bedrock erosion rate (θ) is dependent on ice geometry and dynamics, and the thickness of eroded material (S) which evolves until a maximum limiting threshold is reached (T _lim ). Basal shear stress (τ _b ) and the sliding velocity (U _s ) that arises from it are calculated from the combined stress balance approach described above and in Bueler & Brown (2009). In the above scheme, subglacial erosion potential under an ice sheet in equilibrium will depend primarily on its sliding speed, thickness and surface slope, but the actual amount of bedrock erosion that might occur will be further controlled by the thickness and rheology of overlying sediments as well as the coherency and erodibility of the rock type being overridden. Since these latter two components are largely unknown, except where geophysical studies have allowed local interpretations to be made, we focus here on calculating the glaciological component of subglacial erosion potential. These results are presented with dimensionless units to acknowledge that alternative parameterizations would yield different results. However, it is the relative pattern of erosion potential that is most important in this study, manifest in the different rates that characterize different areas.

Landscape preservation

A primary aim of this study is to investigate the likely pattern of bedrock erosion, as defined above, and the dominant controls that influence landscape evolution. Previous attempts at the interpretation of landscape preservation beneath former ice sheets have commonly invoked favourable basal thermal conditions (i.e. a frozen bed) as the primary control (Dyke 1993, Kleman 1994, Kleman & Hättestrand 1999, Ó Cofaigh et al. 2001, Davis et al. 2006, Reinardy et al. 2009), yet numerous process studies have reported observations of subglacial motion (and thus the potential for erosion) or evidence of bed modification at subfreezing temperatures (Shreve 1984, Echelmeyer & Wang 1987, Cuffey et al. 1999, Atkins et al. 2002, Fitzsimons 2003, Atkins 2013). In the theoretical analysis of Shreve (1984), who calculated sliding speeds for ice at temperatures as low as 20°C below the pressure-melting point, it was demonstrated that sliding could occur at such temperatures because of basal water films as thin as 1 nm (1x10^-9 m). Although basal velocities of subfreezing glaciers are typically low, their erosion potential may be significant over millennial timescales, particularly in areas where the overlying ice is thick (e.g. >1 km).

Role of topography in ice-stream retreat

The Marguerite Trough Ice Stream (MTIS) has received considerable attention both in terms of field-based geological investigations (Ó Cofaigh et al. 2005, Anderson & Fretwell 2008, Kilfeather et al. 2011) as well as modelling studies (Payne et al. 1989, Jamieson et al. 2012). The consensus of this body of research is that the ice stream was sourced in the central mountains of Palmer Land, drained into George VI Sound, flowed northward into Marguerite Bay, and from there towards the continental shelf break (Fig. 3a). Weak deformation till and the presence of mega-scale glacial lineations attest to the fast flow of this outlet (Ó Cofaigh et al. 2005), whilst the occurrence of meltwater channels incised into both soft sediment and bedrock led Anderson & Fretwell (2008) to infer a well developed subglacial hydraulic network beneath the MTIS. Submarine geomorphology is also characterized by well developed grounding zone wedges, which may have stabilized the ice stream during periods of episodic retreat (Jamieson et al. 2012). Here we use our continental-scale model to investigate this particular ice stream system in the context of topographic control on ice sheet behaviour.

Results

The Antarctic Peninsula Ice Sheet at the last glacial maximum

The relatively low aspect-ratio of the simulated LGM APIS means that ice discharging from central accumulation areas is topographically controlled to a greater extent in western areas than in the east. Modelled flow-lines illustrate the effect of this topographical influence, preferentially funnelling ice through troughs (Fig. 1b). This pattern of basal control on ice flow is much less apparent on the Weddell Sea side of the APIS, where instead ice appears to flow more or less directly towards the continental shelf-break. In the east, ice flow from the southern APIS converges with ice from West Antarctica to flow northward through an extended Evans Ice Stream; in the west, glaciers draining the Bryan Coast of West Antarctica discharge into the Bellingshausen Sea.

The modelled pattern of LGM ice flow identifies many individual locations of convergent outlet glaciers and probable ice streams, yet their rates of ice discharge are highly variable. For example, the greater ice thicknesses evident in the southern APIS coupled with higher velocities give rise to a greater horizontal ice flux per unit area than outlets in the northern APIS (Fig. 1c). In terms of discharge volume, the dominant outlets are those draining into the Bellingshausen Sea from Ellsworth Land, the extended Evans Ice Stream in the western Weddell Sea, the Marguerite Trough Ice Stream (MTIS) and to a lesser extent a smaller ice stream sourced from central Alexander Island. Each of these conduits features a well-defined ‘trunk’, and is fed by a number of tributaries. In the MTIS and the two large Bellingshausen Sea outlets, peak fluxes occur in the mid-sections of the glaciers, downstream of their onset zones but upstream of their termini on the continental shelf. By contrast, fluxes in the Evans Ice Stream increase progressively towards the shelf-break.

Subglacial erosion potential

Given the formulation of Eq. (1), it is not surprising that the greatest subglacial erosion is predicted to occur beneath fast-flowing outlets of the LGM APIS (Fig. 2a). In such areas, modelled erosion rates are at least one or two orders of magnitude higher than flanking areas, suggesting that preferential deepening of cross-shelf troughs is probable. Conversely, upland areas, as well as inter-stream sectors, are most likely to be protected from any significant subglacial erosion. This pattern is maintained during retreat of the LGM ice sheet, with zones of erosion potential largely remaining unchanged (Fig. 2b). Considering the central sector of the LGM APIS in more detail, Fig. 3a illustrates that the greatest subglacial erosion potential occurs beneath the central parts of the trunks of the largest outlets, in some cases extending into inland catchments. Beneath the MTIS, erosivity most probably reached a maximum in two locations: one zone at the mouth of George VI Sound and a second in the inner shelf area of Marguerite Bay (Fig. 3a). These locations correspond precisely with peaks in longitudinal strain rates, indicating accelerating flow (Fig. 3b). These zones of peak glacier acceleration do not migrate during deglaciation; rather, their magnitudes are enhanced with respect to ‘steady-state’ conditions at the LGM but their positions remain fixed. During MTIS retreat, along-flow discharge rates increase dramatically in George VI Sound (Fig. 3c) due to the accelerations described above, whereas flux across the continental shelf from Marguerite Bay declines as both ice thickness and horizontal velocity reduces.

Fig. 2 Relative pattern of subglacial erosion potential (Eq. 1) for a. LGM conditions and b. after ocean-forced retreat. Grey shading denotes areas of the ice sheet where calculated erosion potential is more than an order of magnitude below the lowest coloured value (i.e. <0.01).

Fig. 3 a. Pattern of potential Last Glacial Maximum (LGM) subglacial erosion in the central Antarctic Peninsula. AI=Alexander Island, GVIS=George VI Sound, MB=Marguerite Bay. Present-day coastline and LGM flow lines shown in black. Red line illustrates Marguerite Trough Ice Stream (MTIS) centreline. b. Horizontal strain rates along the MTIS centreline at LGM (bold, green) and during stages of ocean-forced retreat (dashed lines). c. Ice discharge rate (flux) along the MTIS at LGM (bold, red) and during retreat (dashed lines). d. MTIS surface and bed profiles at LGM (bold, blue) and during retreat (dashed lines).

Landscape preservation

These simulations indicate that interior areas with little or no subglacial erosion at LGM will be the same as those under a more restricted ice sheet configuration, and thus the underlying landscape in these areas will probably survive unmodified.

To investigate the primary control on subglacial erosion, calculated erosion potential was first plotted at each grid cell against the horizontal ice flux (thicknessxvelocity) in that cell, and each point shaded according to the temperature of basal ice in that location, relative to the pressure melting point (Fig. 4a). There is a clear correlation between erosion potential and ice flux, implying that the greatest modification of the glacier bed is likely to arise wherever ice discharge is high, i.e. where either velocities or thicknesses are greatest. Conversely, landscapes would be protected in areas where fluxes are low, such as in the central parts of ice sheets where velocities are low, or in parts of the ice sheet that are thin, for example over mountain summits. However, the relationship between subglacial erosion potential and basal ice temperature for the simulated LGM APIS is more complex. Whilst the very highest rates of erosion (>1) only occur beneath ice at the pressure melting point, there is no clear relationship at lower temperatures and basal erosion appears to be possible (albeit in relatively few locations) at -10°C and below. This occurs as a consequence of longitudinal stresses. Where fast-flowing wet based ice transitions into slower and colder ice in upstream catchment areas it exerts a pulling force. This allows sliding even at sub-freezing temperatures, in agreement with field observations. Since flow rates are extremely low in these areas, the potential for erosion will be correspondingly small, unless compensated for by a greater thickness of ice. The critical factor appears to be the presence of water at the bed. In our model, ice with a positive sliding velocity also has a basal water layer (Fig. 4b) and the results illustrate that even very thin films of water are sufficient to allow basal motion. However, where these films approach 1 nm, basal velocities are extremely low (Fig. 4c).

Fig. 4 a. Relationship between calculated subglacial erosion potential and the rate of ice discharge, shaded according to temperature of basal ice relative to the pressure-melting point. b. Thickness of basal water layer in relation to basal ice velocity. c. The distribution of basal ice velocities for cells where the subglacial water layer thickness is less than the minimum shown in b.

In summary, it is apparent that temperature alone does not exert complete control over modelled basal erosion potential or the ability of ice sheets to modify or preserve their underlying landscapes. Subfreezing temperatures do not necessarily preclude sliding (and by inference, basal erosion), and conversely, ice at the pressure-melting point does not necessarily lead to significant modification of the glacier bed. What is important in terms of predicting likely erosion patterns is the accurate calculation of horizontal ice flux.

Role of topography in ice-stream retreat

During deglaciation, these model experiments illustrate that, under the range of applied ocean forcing scenarios, grounding line positions tend to retreat gradually in some areas, and more rapidly in others (Fig. 1d). For example, both the outer and inner Bellingshausen Sea appear to favour closely-spaced grounding line positions, whereas the intervening section of continental shelf appears to be a less favourable location for grounded ice margins. Likewise in the western Weddell Sea, the outer part of the shelf appears to allow incremental retreat, whereas the inner shelf is conducive to more abrupt back-stepping of the grounding line. The same pattern is evident seaward of the present-day Larsen Ice Shelf, and at a smaller scale elsewhere around the APIS.

To investigate this in more detail the pattern of retreat of the MTIS along a single modelled flow-line during retreat from its LGM position is considered (Fig. 3a).

These experiments illustrate a series of progressively-retreated grounding line positions that reflect non-linear recession of the grounded ice sheet in response to the imposed oceanic forcings. Many of the 48 simulated grounding line positions overlap, suggesting that some topographic situations favour slower retreat than others (Fig. 3d). The greatest density of grounding line positions occurs on subglacial topographical highs in the inner and outer continental shelf, with far fewer overlapping positions in the deeper mid-shelf. Slower retreat appears to coincide with down-flow dipping beds and the upper parts of retrograde slopes, but is less common on the flat bed and lower parts of flanking slopes of the deep basin of the mid-shelf.

Discussion

Climatic factors appear to have produced asymmetric behaviour of the APIS at the LGM. The wetter western side hosted a larger number of fast-flowing outlet glaciers and ice streams whose positions were topographically determined. The climate and physiography of the eastern APIS gave rise to fewer, but more extensive, areas of fast flow. These differences have no doubt persisted through repeated glacial cycles of the Pleistocene and thus have favoured preferential erosion of the western continental shelf and the creation of large cross-shelf troughs. At the domain scale, the modelled APIS flow-lines (Fig. 1b) are generally consistent with flow directions reconstructed from geological data (Heroy & Anderson 2005, Golledge et al. 2013). However, locally the model fails to simulate correctly the precise flow pattern of some of the smaller APIS ice streams (e.g. Biscoe Island Ice Stream, Palmer Deep Ice Stream, Gerlache and Boyd Strait Ice Stream; Canals et al. 2000, Evans et al. 2004, Amblas et al. 2006, Domack et al. 2006). This is most probably due to the use of a smoothed topographical grid, which may also impact our ability to correctly predict the finer-scale details of basal erosion. In the central sector of the APIS, geological data suggests that a local ice dome was present over Alexander Island at the LGM (Graham & Smith 2012, Johnson et al. 2012). Our simulation for the region suggests that the northern part of Alexander Island was indeed independent of (although confluent with) the APIS, but the detail contained in our smoothed topographical input data is insufficient to reproduce the complexity of the discrete glacier systems indicated by offshore data. On the eastern side of the Antarctic Peninsula the modelled ice sheet limits agree more closely with the probable former extent of the LGM ice sheet interpreted from geological evidence (Livingstone et al. 2012). However, the model struggles to resolve fine-scale flow features (and therefore areas of potential erosion) such as in the Robertson Trough (Evans et al. 2005, Reinardy et al. 2011). Overall, it is clear that whilst this model may be able to predict the ‘selectivity’ of glacial erosion at a broad scale, fine details of ‘selective linear erosion’ in the APIS are still under-represented.

Although there may be limitations to the accuracy with which our model can mimic the geological record, it does at least allow inferences to be made with respect to the dominant controls on landscape evolution. Subglacial erosion appears to be governed more by the magnitude of horizontal ice flux than by basal thermal regime, thus locations where ice is either thin or slowly-flowing may be protected to some extent, irrespective of whether they are wet-based or not. However, once formed, topographical deeps tend to promote further subglacial erosion in areas where ice is accelerating. This tends to occur where ice flow converges in narrow troughs, since mass conservation requires an increase in velocity to compensate for the reduction in width. Since these topographical conditions remain unchanged during retreat phases of the ice sheet, zones of acceleration are more or less fixed, and it is only the velocity of the ice flowing through these zones that varies. Therefore, this suggests that subglacial erosion will most probably be focused in specific zones throughout a glacial cycle, irrespective of changes in glacier extent, until these areas are no longer covered by grounded ice. In the context of APIS expansions during glacial maxima, our model-based prediction is that subglacial erosion of favourable locations on the outer continental shelf will be episodic, whereas the inner shelf might be exposed to longer periods of ice cover and thus greater erosion. These experiments lend support to the notion of preferential overdeepening within glacial troughs, rather than progressive erosion along the entire length of trough floors during advance/retreat cycles. In this sense, the landscape of the Antarctic Peninsula can truly be considered one of ‘selective’ erosion, and in some areas is both selective and linear.

During retreat of the APIS through the last termination, topographical influences on ice sheet flow were further modified by the effects of rising sea levels and warming oceans. These environmental drivers led to thinning, flotation, and subsequent retreat of grounded ice, with faster flowing outlets responding more quickly than slower flowing ice (Golledge et al. 2012a). On a perfectly flat bed, an ice-stream of uniform width will retreat linearly in response to applied oceanic forcings (Jamieson et al. 2012, exp. C). However, topographical conditions, such as variability in trough width and bathymetry, lead to a non-linear retreat rate characterized by stillstands at favourable locations (Jamieson et al. 2012, exp. A). The marine ice sheet instability hypothesis (Weertman 1974, Thomas & Bentley 1978) predicts that a retreating marine-based ice sheet will be unconditionally unstable on a bed that slopes landward, assuming no lateral variation in the flow field. Our simulations, in common with those of Jamieson et al. (2012), suggest that topographical influences in the Antarctic Peninsula may have allowed outlets, such as the MTIS, to retreat slowly even on reverse slopes, most probably as a consequence of changes in ice-stream geometry that sufficiently modified the rate of ice flux across the grounding line. These results concur with other theoretical studies (Gudmundsson et al. 2012, Gudmundsson 2013), and highlight the critical influence on retreat rates exerted by topography in deeply incised landscapes like the Antarctic Peninsula.

Conclusions

These simulations have investigated both the role of ice sheet dynamics in modification of subglacial topography, and the role of basal topography in influencing ice sheet behaviour. Using a high-resolution model of the LGM APIS, both at its maximum extent and during stages of retreat, has demonstrated that variability in subglacial erosion potential most probably reflects spatial differences in net horizontal flux, rather than basal temperature. Furthermore, the greatest basal erosion most probably coincides with zones of convergent flow and highest strain rates. During retreat, grounding line migration may be slow even on retrograde slopes if horizontal fluxes are sufficient, but deep basins tend to produce faster retreat, particularly where ice flow is divergent. Together, these results confirm the complex inter-relationship between ice sheets and topography, and provide a mechanistic basis for the palaeo-glaciological interpretation of geological data. Within the limits of model resolution, the simulated pattern of basal erosion is consistent with the geomorphologically-based concept of selective linear erosion, but more detailed topographical input data is required before a comprehensive prediction of subglacial erosion variability may be possible.