1. Introduction
Marine infiltration has played a significant role in the dynamics of the Northern Hemisphere ice sheets since the Last Glacial Maximum (LGM). Glacio-isostatic modelling suggests that one of the earliest features of the retreat was the deglaciation of the Barents Sea ice sheet at 18 kyr cal BP (Reference Knies, Kleiber, Matthiessen, Muller and NowaczykLambeck, 1995), in which sea-level induced marine calving and corresponding retreat of the ice margin is thought to have played a substantial role (Reference Lambeck and ChappellLandvik and others, 1998). For the Laurentide ice sheet the Hudson Bay region is thought to have become ice-free between 9.4 and 9 kyr cal BP (Reference Dyke, Vincent, Andrews, Dredge, Cowan and FultonDyke and others, 1989). This event is one of the last major ice-margin retreats associated with deglaciation since the LGM. Both of these events are significant for the history and stability of these ice sheets. Moreover, it is our view that they represent good observational constraints on the marine-retreat history of Northern Hemisphere ice sheets, a process which is poorly understood and appears to exhibit variability on both temporal and spatial scales (Reference Brown, Sikonia, Post, Rasmussen and MeierBrown and others, 1983; Hughes, 2002). Additionally, refined measurements of eustatic sea-level change since the LGM have identified periods of fast and significant sea-level rise (Reference FairbanksFairbanks, 1989; Reference Vieli, Funk and BlatterYokoyama and others, 2000), which, although not correlated in time with known large-scale marine ice retreat, are certainly possible results of large-scale marine calving. Calving events related to the retreat of the marine margin of Northern Hemisphere ice sheets have been linked to the production of meltwater pulses affecting oceanographic circulation and producing rapid climatic change (Reference Lindstrom and MacayealMacAyeal, 1993; Reference Clarke, Marshall, C., Bilodeau, Veiga-Pires, Clark, Webb and KeigwinClarke and others, 1999). An accurate representation of the change in the marine extent in ice-sheet models is therefore necessary in order both to generate realistic ice-sheet dynamics and to infer the role of ice sheets in the larger climate system.
In general, however, ice-sheet models applied to the Northern Hemisphere ice sheets have focused on land-based ice and have seldom treated the interaction with the ocean in a physically coherent way. One avenue to overcome the problem is to consider rigorous coupling with an ice shelf over the Arctic Ocean in analogy with the situation in Antarctica (Reference Le Meur and HuybrechtsLindstrom and MacAyeal, 1989). In that case, the marine extent coincides with the grounding line, but it is questionable whether a thick ice shelf actually existed over the Arctic basin, whether it was continuous in space and time, and whether it was connected to a circumarctic grounded ice sheet that fed it. Anything else would make such a coupling scheme very complex to implement in a numerical model. Therefore most treatments of the marine margin have not considered floating ice dynamics explicitly. Instead the assumption has been made that marine calving is the dominant process limiting the extent of marine-based ice sheets. Because the spatial resolution of ice-sheet models is too low to allow explicit physical treatment of marine calving, and the process by itself is too poorly understood, parameterizations are used. The most common element of the varied parameterizations of marine calving is that water depth (and/or ice thickness) at the margin is considered to most strongly control the process (Reference Clarke, Marshall, C., Bilodeau, Veiga-Pires, Clark, Webb and KeigwinClarke and others, 1999; Reference Polyak, Forman, Herlihy, Ivanov and KrinitskySiegert and others, 2001a). In the more sophisticated modelling schemes, the parameterized calving flux is considered as an ablation term at the margin (Reference PeltierPfeffer and others, 1997; Reference Marshall and ClarkeMarshall and others, 2000), while most other studies in effect derive a depth contour beyond which calving removes all ice (Reference HuybrechtsHuybrechts and T’siobbel, 1997; Reference SvendsenTarasov and Peltier, 1997; Reference Charbit, Ritz and RamsteinCharbit and others, 2002).
The formulation developed in this paper is essentially of the second type, as our main interest is in producing the geomorphologically observed marine extent of Northern Hemisphere ice sheets in the simplest possible way. To this end, we follow the work of Reference Brown, Sikonia, Post, Rasmussen and MeierBrown and others (1983) and Reference Marshall, Tarasov, Clarke and PeltierMeier and Post (1987) on the observed marine-calving induced change in grounded ice extent for tidewater glaciers, and assume that water depth is the most important physical property which controls the extent of the marine margin. By assuming that, below a certain critical water depth, marine calving is complete, we empirically parameterize the allowable grounded marine extent of the ice sheets as a function of water depth, using a water-depth dependence on eustatic sea level. By relating water depth to eustatic sea-level changes, we reproduce the first-order characteristics of marine ice behaviour (i.e. advance over the continental shelf during periods of low sea level, and retreat during periods of rising sea level) without introducing so many model parameters as to under-determine the magnitude of important processes.
2. Marine-Extent Parameterization
The relation used here to determine the marine extent follows that of Reference HughesHuybrechts (2002) in an application for the Greenland ice sheet. It is similar in several respects to that of Tarasov and Peltier (1997), Reference HuybrechtsHuybrechts and T’siobbel (1997) and Reference Charbit, Ritz and RamsteinCharbit and others (2002). In the Tarasov and Peltier ice-sheet model, the marine extent is set at the 400 m marine bathymetric limit of the present-day continental shelves. This effectively assumes an Arctic continental shelf devoid of water where the surface mass balance solely controls the ice-sheet extent until it extends to the fixed bathymetric limit. Reference HuybrechtsHuybrechts and T’siobbel (1997) also use a constant 500 m value for the marine limit, but found that this formulation was insufficient to model marine extent, as it allowed marine advance but provided no effective mechanism for marine retreat. In the Reference Charbit, Ritz and RamsteinCharbit and others (2002) model, the marine extent follows a flotation criterion. Here we also use a water-depth criterion, but, instead of invoking the marine extent for a constant water depth, we consider the water depth to be a function of eustatic sea level. Any grounded ice that extends to a location with marine bathymetry deeper then this critical water depth is considered to be beyond the ice-sheet margin, and calving here is complete.
We find that by qualitatively tuning the eustatic-sea-level–marine-extent relationship, the best fit to the (sparse) geomorphological observations since the LG M is given by:
where Hc is the maximum water depth to which ice can ground, and ΔHsl is the time-dependent eustatic sea-level change. The definition of water depth used here is the contemporaneous marine bathymetry, which changes with time due to hydro-isostatic loading from eustatic sea-level change and glacio-isostatic loading from the advance and retreat of the ice sheets. This generates a feedback between the marine extent and marine grounding from the time-evolving marine bathymetry. We do not include changes in eustatic sea level in the calculation of water depth here, as this complication was found to have a second-order effect compared to Equation (1), which is already a function of eustatic sea level. Given the irregular topography of the continental shelves of oceanic regions, we define a marine margin for the ice sheets where Equation (1) applies. This area is a band around the edge of the ice sheet where the ice front is adjacent to the open ocean and is 5 1200 m thick. Ignoring this constraint essentially allows marine calving to occur in the interior of the ice sheets wherever the marine bathymetry is less than Hc. The marine-extent relationship is formulated so that during periods of low eustatic sea level the ice sheets expand over the continental shelf, and during periods of high sea level the ice sheets retreat towards the present-day coastline.
The marine-extent determination given in Equation (1) is the simplest relationship we have found which satisfies the geomorphological observations. Essential characteristics are the hybridization (at –100 m here) and choice of gradients. The first-order features of the relationship are that at the present day (ΔHsl = 0 m) the marine extent is close to the present-day coastal marine margin. At the LGM (ΔHsl ≈ –130 m), the marine extent occurs at the L GM marine bathymetric contour, which accounting for glacial isostasy and lowered eustatic sea level is near the 600 m depth contour of the continental shelves of present-day ocean bathymetry.
3. The Ice-Sheet Model
The marine margin formulation given in Equation (1) was implemented in an upgraded version of the three-dimensional thermomechanical ice-sheet model of Reference HuybrechtsHuybrechts and T’siobbel (1997). This model has a horizontal mesh size of 50 km, 17 layers in the vertical, and covers all of the Northern Hemisphere where widespread glaciation is believed to have taken place. The main inputs to the model are bed topography (derived from the ETOP05 dataset), mean annual surface temperature and the mass balance. The latter two represent the climatic forcing. They are based on the present-day mean monthly surface temperature at sea level and the mean monthly precipitation rate, to which perturbations are applied in the anomaly mode to account for different climates. We retained the precipitation climatology from the Jaeger maps (Reference Imbrie, A., Imbrie, J. Hays, Kukla and SaltzmanJaeger, 1976), but exchanged the ECHAM-1 temperature climatology for the U.S. National Centers for Environmental Prediction (NCEP) re-analysis of long-term monthly mean temperature interpolated to the 1000 mbar level provided by the U.S. National Oceanic and Atmospheric Administration–Cooperative Institute for Research in Environmental Sciences (NOAA–CIRES) Climate Diagnostics Center. The mass-balance model distinguishes between snowfall, rainfall, ice- and snowmelting, meltwater retention and subsequent runoff, the latter of which is based on a degree-day model with positive degree-day (PDD) factors of 3 mma–1PDD–1 for snow, and 8 mm a–1PDD–1 of water equivalent for ice.
The treatment of glacial isostasy is important for the time-dependent change of marine bathymetry that is used as a control on the extent of the marine margin. A lithosphere with flexural rigidity 1025 Nm is used to calculate the steady-state deflection caused by changing ice and water loads over the model domain. Isostatic equilibrium is attained by a schematic mantle adjustment towards this steady-state deflection using a point-wise relaxation with a decay time of 3 kyr. Although simple, this formulation is superior to the local glacio-hydrostatic lithosphere formulation adopted previously, and sufficient to reproduce the major features of more sophisticated viscoelastic formulations (Reference LandvikLe Meur and Huybrechts, 1996). Another important parameter in the model is the multiplier in the rate factor of Glen’s flow law (enhancement factor), which was set to 26. Although this value is too high for the Greenland ice sheet (Reference Tarasov and PeltierTarasov and Peltier, 2000), it is chosen so that the model results are comparable to the Northern Hemisphere model results of Reference HuybrechtsHuybrechts and T’siobbel (1997).
4. Experimental Design
We ran the ice-sheet model over the last glacial cycle, with initial conditions of a glaciological steady-state ice configuration forced by present-day climate. The time-dependent climate forcing is built around two U.K. Meteorological Office (UKMO) V.3.2 Paleoclimate Modeling Intercomparison Project (PMIP) (Reference Hewitt and MitchellHewitt and Mitchell, 1996) precipitation and temperature time slices, from an LG M simulation and from a present-day simulation. The monthly temperature data were reduced to the same topographical level using a lapse rate of –0.008°C m–1, and required some corrections to eliminate artifacts arising from land/sea mask changes between the LGM and today. To use these fields over an entire glacial cycle, we consider them as “extremes” in the glacial/interglacial climate contrast, and use a normalized “glacial index” derived from a synthesized Greenland Icecore Project (GRIP) δ18O/Vostok deuterium record (Reference HughesHuybrechts, 2002) to rescale the fields over an entire glacial cycle (Fig. 1a). The anomalies in climate between the LGM and present-day UKMO PMIP fields are then applied as temperature differences and precipitation ratios to the observed present-day climate to generate time-dependent changes in temperature and precipitation over an entire glacial cycle.
Fig. 1. Time evolution of ice-sheet model forcing and output. (a) Glacial index used to force the precipitation and temperature fields derived from a synthesized Vostok–GRIP climate record. (b) Eustatic sea level based on SPECMAP, but with corrections from Reference LambeckLambeck and Chappell (2001) for 21 kyr cal BP to present day. (c) Time-dependent change in total ice volume. (d) Surface area evolution of the grounded ice sheet. (e) Ice-volume growth and decay converted into metres of mean sea-level equivalent. The conversion shown takes into account the volume of ice displacing ocean water (HI, hydroisostatic) and disregards ice volume outside the ice-sheet masks implied by the work of Reference Dyke, Vincent, Andrews, Dredge, Cowan and FultonDyke and others (1989) and Reference StuiverSvendsen and others (1999) (HI and GE, geomorphological extent).
The eustatic sea-level forcing is based on the spectral-mapping project (SPECMAP; Reference Huybrechts and T’siobbelImbrie and others, 1984) but with modifications from Reference LambeckLambeck and Chappell (2001) so that the minimum in eustatic sea level occurs at 21 kyr cal BP and not at 18 kyr cal BP. This modified sea-level record is shown as a function of time in Figure 1b.
5. Results
5.1. Glacial cycle
Figure 1e shows our estimate of the far-field eustatic sea-level contribution of the Northern Hemisphere ice sheets over a glacial cycle. It is based on total ice-sheet volume (Fig. 1c), but has been modified to show a more realistic equivalent sea-level change (see below for explicit details of the modification). At the LGM, the corresponding eustatic sea-level change is slightly above 100 m. This value is about 75% of the estimated global eustatic sea-level change at maximum glaciation (Reference ClarkClark and Mix , 2000). Taken together with the Reference HughesHuybrechts (2002) estimate of a 14–18 m sea-level contribution from the Antarctic ice sheet, this gives a total estimate of eustatic sea-level change of 110–120 m. Given the climate-related uncertainties related to ice-sheet model input, we believe that the modelled Northern Hemisphere ice-volume changes over a glacial cycle reported here are therefore quite reasonable.
To evaluate the marine-extent formulation, we focus on changes in the marine margin since the LGM , because these allow us to compare the model results with geomorphological observations in the area. Attempts have been made to reconstruct the extent of the Northern Hemisphere ice sheets prior to the LGM (Reference Clark and MixClark and others, 1993; Reference MacAyealMangerud and others, 1998), but much of the evidence has been overridden by subsequent ice-sheet advance. In our comparison, we convert cited 14C years to calendar years using the INTCAL98 time-scale (Reference M., Kromer, Becker and FergusonStuiver and others, 1998).
5.2. Last Glacial Maximum
Time slices of modelled ice-sheet elevation at various moments are displayed in Figure 2. At 19 kyr cal BP the model reproduces the LGM ice extent in accord with the geomorphological observations (Reference StuiverSvendsen and others, 1999; Reference DykeDyke and others, 2002) except that it includes excessive glaciation over the Taimyr Peninsula and the Kara Sea in northern Russia (see Fig. 3 for locations). Although there is some evidence for regional glaciation in Severnaya Zemlya (Reference JaegerKnies and others, 2001), the ice extent generated here is probably too great. Sensitivity studies indicate that this excessive glaciation is likely a result of the UKMO PMIP climate forcing. To generate LGM climate fields, boundary conditions for topography and surface albedo were used from the Ice4G reconstruction of Reference Miller and KaufmanPeltier (1994). This reconstruction contains excessive glaciation over the Taimyr Peninsula compared to the work of Reference StuiverSvendsen and others (1999), and as a result this pattern is imprinted in the climate forcing of UKMO which the ice-sheet model uses to generate the ice cover. Glaciation also occurs in this region for other models that use PMIP climate input (Reference PfefferPollard and Thompson, 1997; Reference Charbit, Ritz and RamsteinCharbit and others, 2002), providing additional evidence that this is a result of the climate formulation.
Fig. 2. Evolution of ice-sheet elevation (in metres) since the LGM using the marine extent relationship given in Equation (1). Times are given in each panel and are expressed in kyr cal BP.
Fig. 3. Index map showing geographic locations mentioned in the text. Colour contours give the present-day bathymetry of the continental shelf.
The existence of ice over the Kara Sea in our model, a region thought to be ice-free at the LGM (Reference Knies, Kleiber, Matthiessen, Muller and NowaczykLambeck, 1995; Reference MangerudMangerud and others, 2002), is also thought to result from inaccurate climatic input, as mass-balance reconstructions which generate no ice sheet over the Kara Sea require special climate modification, in particular to mimic very dry conditions (Reference Siegert and MarsiatSiegert and others, 1999; Reference Siegert, Dowdeswell and MellesSiegert and Marsiat, 2001). In most other areas, however, the ice sheet modelled for 19 kyr cal BP matches the inferred LGM ice extent quite well. Siberia and Alaska are predominantly ice-free, and the southern margins of the Laurentide and Fennoscandian ice sheets are reasonably well matched with the geomorphological reconstructions of Reference DykeDyke and others (2002) and Reference Andersen, Denton and HughesAndersen (1981). We show the 19 kyr cal BP elevation as LGM because that is close to when the modelled icevolume reaches a maximum (at 18.4 kyr cal BP; Fig. 1c). Because of the problems associated with glaciation over the Taimyr Peninsula, we have modified the ice-volume conversion in Figure 1e to only include ice that occurs within the geomorphologically inferred LGM ice-sheet extent of Reference Dyke, Vincent, Andrews, Dredge, Cowan and FultonDyke and others (1989) and Reference StuiverSvendsen and others (1999). The conversion from ice-sheet volume to eustatic sea-level contribution also accounts for changes in the volumetric capacity of the Northern Hemisphere ocean basins as a result of hydro-isostatic adjustment. However, these modifications do not affect the timings of maximum ice-sheet volume and the corresponding minimum contribution to sea level, which occur simultaneously.
5.3. Early marine retreat
The maximum ice-sheet volume for our model occurs about 3 kyr later than is usually inferred from the observed global
eustatic sea-level record (19–22kyr cal BP; Reference Vieli, Funk and BlatterYokoyama and others, 2000). That is because of the phase lag before the model ice dynamics reacts to the climate index forcing. Reference Vieli, Funk and BlatterYokoyama and others (2000) suggest that sea level began to rise rapidly at 19 kyr cal BP, and this time is also thought to mark the beginning of the retreat of the Barents Sea ice sheet (Reference Lambeck and ChappellLandvik and others, 1998). However, despite the climate/ice-volume lag, the marine-extent formulation produces a retreat in the Barents Sea at this time because the sea-level forcing already started to rise at 21 kyrcal BP (Fig. 1b). At 18 kyr cal BP, the central ice domes of the Laurentide and Fennoscandian ice sheets are higher than at 19kyr cal BP, but the St Anna Trough is already ice-free as a result of a retreat in the marine extent (Fig. 2). The modelled retreat history at this location is complicated by the excessive glaciation over the Kara Sea, but the St Anna Trough was thought to have been glaciated during the LGM (Reference PollardPolyak and others, 1997). The marine bathymetry here is also the deepest of any part of the continental shelf to which the Fennoscandian ice sheet extended, with a present-day water depth of 600–700m. Therefore with our water-depth-dependent marine extent relationship it is the first region of the Barents Sea ice sheet to become ice-free. By 17 kyr cal BP (Fig. 2), the Barents Sea is open, with glaciation confined to Svalbard, Franz Josef Land and Novaya Zemlya and regions with shallow marine bathymetry close to the coast. In the ice-sheet model, the Barents Sea opens to the sea so that Svalbard and Franz Josef Land are isolated at 17.4 kyr cal BP. This isolation occurs at least 1.6 kyr later in the model than suggested by Reference Andersen, Denton and HughesAndersen (1981), although Reference Lambeck and ChappellLandvik and others (1998) suggest that it occurred between 19 and 18 kyr cal BP. However, considering that our climate forcing generates a somewhat late LGM ice sheet, the uncertainty in the timing of the event and the errors associated with dating this event (radiocarbon time-scale calibration (Reference Siegert, Dowdeswell, Hald and SvendsenStuiver and others, 1986) and variations in 14C production (Reference Bard, Arnold, Hamelin, Tisnerat-Laborde and CabiochBard and others, 1998)), we conclude that the marine retreat here generates the observed infiltration of the Barents Sea quite well.
5.4. Late marine retreat
The modelled marine retreat for the Laurentide and Innuitian ice sheets begins late in the last deglaciation at about 13 kyr cal BP. By this time, the Innuitian ice sheet begins to recede in the northwest, starting with the outer islands (Fig. 2). The other significant marine retreat during the later stage of deglaciation is the marine infiltration of Hudson Bay by 9 kyrcal BP. Figure 2 shows the modelled ice-sheet elevation at 10 kyr cal BP. A major retreat of the Fennoscandian ice sheet has already occurred, although the Taimyr Peninsula remains glaciated. Themodel does not predict the separation, at this point, of the Laurentide and Cordilleran ice sheets, which is thought to have occurred by 15 kyr cal BP (Reference Dyke, Vincent, Andrews, Dredge, Cowan and FultonDyke and others, 1989). Hudson Bay is still glaciated, although Hudson Strait is ice-free. At 9 kyrcal BP (Fig. 2), Hudson Bay is ice-free, as is Foxe Basin. After 9 kyr cal BP, changes in the marine limit occur only close to Svalbard and Franz Josef Land and in northern Canada, and are of negligible magnitude. Finally, it is worth noting that an advantage of parameterizing marine extent in terms of eustatic sea level is that the model is able to reproduce a full return to interglacial ice-free conditions in the Holocene, without marine ice sheets remaining on the Arctic continental shelf as in some previous modelling studies (e.g. Reference Greve, Wyrwoll and EisenhauerGreve and others, 1999).
6. Discussion
6.1. Marine retreat in Eurasia
Using the simple marine-extent parameterization expressed in Equation (1), we are able to reproduce the large-scale marine infiltration events associated with the Fenno-scandian ice sheet since the LGM. There is a delay in the retreat of the Barents Sea ice sheet, but this difference in timing is small compared to the observational errors and can be corrected for using an alternate climate forcing. The essential component of the marine-extent model which generates the retreat is that we define a maximum water depth for the marine margin, or equivalently a“minimum” depth for complete calving. Our required water depth which determines marine extent is thus an overestimate compared to explicit calving-law determinations of the marine limit in which the rate of calving is directly proportional to water depth (e.g. Reference Mangerud, Astakov and SvendsenMarshall and others, 2000; Reference Polyak, Forman, Herlihy, Ivanov and KrinitskySiegert and others, 2001). Regardless of the details of the formulation, the marine margin changes most where the marine bathymetry is deepest. For the Barents Sea ice sheet, this is at the St Anna Trough and Bjørnøyrenna (Fig. 3), and in all model runs used to tune the marine-limit relationship these are the locations where ice loss occurs first. The relative timing between marine retreat at Bjørnøyrenna and retreat at the St Anna Trough is of interest in examining the importance of ice thickness as a control on marine calving. The present-day bathymetry of the deepest (and thus most sensitive area in our formulation) part of the St Anna Trough is approximately 100 m deeper than that of Bjørnøyrenna. The isostatic depression of these areas in our model is approximately the same for both areas, so that, using our marine-extent relation, the bathymetry is deeper for the St Anna Trough than for Bjørnøyrenna. It is a characteristic of our marine-margin determination that retreat of the marine margin near the St Anna Trough always precedes, or occurs coincidentally with, retreat at Bjørnøyrenna, as the water depth at Bjørnøyrenna is always shallower. However, if ice thickness is a crucial parameter in marine calving (e.g. Reference PeltierPfeffer and others, 1997; Marshall and others, 2000) and the ice sheet was thicker over Bjørnøyrenna than over the St Anna Trough, marine calving near Bjørnøyrenna could have been greater, resulting in the marine margin retreating here first. Large-scale marine retreat at Bjørnøyrenna preceding retreat at the St Anna Trough is assumed in some Fenno-scandian deglaciation chronologies (e.g. Ice4G (Reference Miller and KaufmanPeltier, 1994); maximum model reconstruction of Lambeck (1995)). However, given the uncertainty of the LGM extent of the Barents Sea ice sheet near the St Anna Trough (Reference JaegerKnies and others, 2001; Reference MangerudMangerud and others, 2002) and the uncertainties regarding LGM climate in this region (Reference Pollard and S. L. ThompsonPollard and others, 2000; Reference Polyak, Forman, Herlihy, Ivanov and KrinitskySiegert and Marsiat, 2001), we make no attempt here to discriminate between models governing marine-margin changes on this basis.
6.2. Marine retreat in North America
For the North American region, the marine-extent relationship also reproduces the geomorphologically inferred retreat reasonably well. The marine retreat of the Innuitian ice sheet starts in the northwest at about 13 kyr cal BP, which is broadly consistent with geomorphological inference (Reference DykeDyke and others, 2002). For the Laurentide ice sheet, the marine retreat of Hudson Bay begins at the easternmost part of Hudson Strait at 12 kyr cal BP. By this time, approximately half of Hudson Strait is ice-free as is Ungava Bay. Our model contains no ice-stream physics, so the extent of ice in Hudson Strait is probably underestimated (e.g. Reference Mangerud, Astakov and SvendsenMarshall and Reference Clark and MixClarke, 1997). The easternmost part of Hudson Strait has the deepest marine bathymetry of any area of the Hudson Bay/Hudson Strait basin system (present-day depth of 600 m), so the modelled marine retreat occurs relatively early here. The ice-sheet model does not produce the northward advance of ice across Hudson Strait associated with the Younger Dryas cooling event reported by Reference Meier and PostMiller and Kaufman (1990) and modelled by Reference PeltierPfeffer and others (1997). The Younger Dryas period is the only time since the LGM that the ice-sheet model generates an advance in the marine limit, and the advance occurs only in coastal regions of Norway, British Columbia and Alaska. Although evidence for marine advance in Norway has been found during this period (Reference BaumannBaumann and others, 1995), in the model the advance is generated by climate forcing and not from the marine–grounding feedback which is included in the model by relating the marine extent with the isostatically adjusting marine bathymetry. Over the continental shelf, the magnitude and rate of isostatic adjustment is small, as any part of the ice margin which becomes ice-free in the model changes from glacio-isostatic loading to hydro-isostatic loading. Thus we conclude that this feedback is of second-order importance compared to the eustatic sea-level dependent forcing (Equation (1)) of the marine extent.
The marine infiltration of Hudson Bay by 9 kyr cal BP is reproduced reasonably well by the model, although Foxe Basin is also ice-free at this time. The reconstruction of Reference Dyke, Vincent, Andrews, Dredge, Cowan and FultonDyke and others (1989) suggests that Foxe Basin was still glaciated at 9 kyr cal BP, becoming ice-free only by 8 kyr cal BP. The minimum marine bathymetry of this basin is shallower than that of Hudson Bay, so that, according to our marine extent, relationship retreat would begin in Foxe Basin only after retreat in Hudson Bay. However, as Hudson Bay is thought to be close to the maximum elevation of the Laurentide ice sheet at the LGM, it is probably the area in our model that is most sensitive to changes in the marine extent caused by glacio-isostatic adjustment. The simple treatment of isostatic adjustment used here contains the assumption that the Earth is in isostatic equilibrium at the present day. Accounting for the observed present-day residual isostatic adjustment in the model would generate a shallower marine bathymetry in Hudson Bay and effect changes in the marine extent accordingly. The fact that both Ungava Bay and Foxe Basin record an early deglaciation in the model could also be related to feedbacks in ice dynamics which result from changes in the marine extent. On our model resolution of 50 km these areas are represented by only a few gridpoints, and the change in ice velocity accompanied by a large marine retreat is therefore quite sufficient to drain these basins. The quantification of the relative contributions of these processes to the modelled change in ice extent is, however, not insightful given the absence of physically relevant processes in the model such as fast ice streams, ice-shelf dynamics and enhanced basal flow. As a physically realistic model of marine-calving induced changes in marine extent is beyond the resolution of 50 km, we focus here on larger-scale marine infiltration events only.
6.3. Marine-extent formulation
The marine-extent relationship is a hybrid equation (separated at ΔHsl = –100 m), with different gradients for different segments of eustatic sea-level change. Experience from tuning the marine-extent formulation has shown that the hybridization is required to produce the correct timing for deglaciation of the major Northern Hemisphere ice-sheet marine margins. For example, using a simple sea-level equation of the form Hc = 6ΔHsl produces a correct Laurentide and Fennoscandian ice-sheet marine extent, but the retreat of the Barents Sea marine limit occurs only at 13 kyr cal BP and that of Hudson Bay not at all. For Hc = 4ΔHsl , the LGM ice-sheet configuration has regions of Hudson Bay and the Barents Sea ice-free. For Hc = 5ΔHsl , Hudson Bay is glaciated, but, contrary to the geomorphological observations, the St Anna Trough in the Barents Sea is ice-free (Reference PollardPolyak and others, 1997). Therefore a simple linear relationship between marine extent and eustatic sea level is insufficient to reproduce the observed marine infiltration of the Northern Hemisphere ice sheets since the LGM.
It is possible that different marine-extent relationships are appropriate to different ice sheets (as has been suggested for ice-sheet rheology (Reference Tarasov and PeltierTarasov and Peltier, 2000)). Reference Tarasov and PeltierVan der Veen (2002) suggests that observational evidence showing that the calving rate for tidewater glaciers is an order of magnitude faster than for lacustrine glaciers is attributable to differences in basal sliding velocity, resulting from different water densities. He also argues that the particular form of glacier calving is related to whether the glacier is polar or temperate, and whether the glacier extent is grounded or floating. Although this suggests that it would be reasonable to assume regional heterogeneity in the treatment of the marine extent in ice-sheet models, it is our intent here to generate the simplest marine-extent formulation that explains most of the geomorphologically inferred marine-retreat history of the Northern Hemisphere ice sheets.
We suggest that our results regarding a non-linear relationship between marine extent and eustatic sea level are analogous to those of Reference Clarke, Marshall, C., Bilodeau, Veiga-Pires, Clark, Webb and KeigwinClarke and others (1999) and Reference Marshall and ClarkeMarshall and others (2000) relating to tidewater calving. They discuss the desirability of a“calving switch”, between fast and slow periods of marine calving with corresponding changes to the marine extent. The correspondence between their results and those presented here suggests that the requirement for Equation (1) to be of hybrid form is not a result of the neglect of ice thickness in the parameterization used here. The fit of Equation (1) to the geomorphologically inferred observations would suggest that the role of marine-extent changes in the retreat of the Northern Hemisphere ice sheets has two components. The first (for Hsl <–100m) predominantly relates to large-scale calving of predominantly marine-based ice sheets (e.g. the Barents Sea ice sheet) and appears to be relevant to the instability of glacial dynamics when subjected to large marine-infiltration events. The second (for Hsl = –100m) relates more to marine calving produced by eustatic sea-level changes and appears to be relevant to the process of marine infiltration of narrow waterways and coastal regions.
The results presented here show an interesting correspondence to those focused on the rate of marine calving of tidewater glaciers and the role of water depth. For example, using a flotation criterion in a numerical model, Reference Van der VeenVieli and others (2001) conclude that a linear relationship between calving rate and water depth is only valid during periods of slow changes to the grounded marine extent. They agree with the conclusions of Reference Marshall, Tarasov, Clarke and PeltierMeier and Post (1987) that during periods of fast changes, changes in marine extent are dominated by the amount by which the thickness at the ice edge exceeds its flotation height. If these results are applicable here, this would suggest that the non-linearity of Equation (1) may result from the neglect of imposing a flotation criterion.
Another possible explanation for the non-linearity of Equation (1) is that it reflects a climate signal that is not encapsulated in the ice-sheet model. Figure 1 shows that, since the LGM, eustatic sea level has changed inversely with temperature change as reflected by the GRIP ice core. Therefore low sea level is correlated with colder ice sheets which would have the capacity to expand in deeper water. As the Barents Sea ice sheet was much further to the north than the Hudson Bay region, this leaves open the possibility of an entirely different form of calving for the two ice sheets (Reference Tarasov and PeltierVan der Veen, 2002). As our fit to the geomorphological retreat history is empirical in nature, it is these observations which constrain our understanding of the relevant processes in marine calving and changes in the marine margin.
7. Conclusions
Using a simple marine-extent relationship, we are able to reproduce the first-order processes associated with changes in the marine margin of the Northern Hemisphere ice sheets since the LGM. The timing errors associated with our marine-extent relationship appear to be related to the climate forcing used in the ice-sheet model, and this forcing must be improved before we are able to state that this formulation is sufficient to reproduce all of the observed marine retreat. However, at this stage we are able to conclude that, before explicit model physics of calving, more appropriate margin dynamics, and ice-stream physics are included in ice-sheet models, a simple parameterization of the form used here serves as a useful and physically reasonable tool with which to model the marine retreat of the Northern Hemisphere ice sheets since the LGM.
Acknowledgements
This work was performed within the HGF-Strategiefonds Projekt 2000/13 Seal (Sea Level Change). We especially thank C. Hewitt, J. Gregory and G. Munhoven for their help in obtaining and processing the UKMO PMIP climate fields over our model domain, and S. J. Marshall and L. Tarasov for thoughtful reviews.
