Evidence for a more extensive Greenland Ice Sheet in southwestern Greenland during the Last Glacial Maximum

Reconstructions of past ice-sheet dimensions serve as important prescribed features in transient climatic model simulations (He et al., 2013; Zweck and Huybrechts, 2005) and are also used to independently assess paleoice–sheet model performance (Downs et al., 2020; Lecavalier et al., 2014; Tarasov and Peltier, 2004). Yet, such reconstructions reveal conflicting estimates of the maximum size of the Greenland Ice Sheet during the Last Glacial Maximum (LGM, 26–19.5 ka) around some parts of Greenland. For example, while most reconstructions of Greenland Ice Sheet extent during the LGM depict a terminus located somewhere on the continental shelf, its exact position is uncertain (Fig. 1; Bradley et al., 2018; Fleming and Lambeck, 2004; Funder et al., 2011; Lecavalier et al., 2014). Placing direct constraints on the LGM ice extent requires offshore geophysical surveys, seafloor mapping, and sediment coring; however, only a few areas have yielded detailed knowledge about past ice extent (Dowdeswell et al., 2014; Hogan et al., 2016; Jennings et al., 2017; Newton et al., 2017; Ó Cofaigh et al., 2013; Ó Cofaigh et al., 2004).

Offshore LGM ice extent can also be inferred using terrestrial evidence. Glacial sediments and landforms can be used to constrain former ice extent, but these archives rarely extend back into the LGM on Greenland (Håkansson et al., 2007; Larsen et al., 2018; Young et al., 2020). Relative sea-level curves have been used in combination with glacially weathered trimlines to infer the size and lateral extent of ice sheets, but these estimates are typically used as a minimum ice-sheet estimate (Funder and Hansen, 1996). A growing number of Greenland's coastal regions have been studied using cosmogenic nuclide exposure dating, which provides absolute age constraints on the timing of ice retreat onto land (e.g., Winsor et al., 2015; Roberts et al., 2009; Søndergaard et al., 2020). However, deglaciation chronologies are less straightforward when used in areas where frozen bedded conditions dominated during the LGM, such as in eastern and northern Greenland, or elsewhere at higher elevations (Briner et al., 2013; Corbett et al., 2013; Graham et al., 2019; Kelly, 1980).

Southwestern Greenland has an extensive terrestrial record of the Greenland Ice Sheet during the Holocene (Young et al., 2020), but like other regions in Greenland, its glacial extent during the LGM is unclear. A commonly cited reconstruction places the southwestern Greenland Ice Sheet LGM terminus a mere 40 km offshore of the modernday coastline (Funder et al., 2011), yet studies on the wide continental shelf off Disko Bugt support an LGM ice limit at the shelf break (Ó Cofaigh et al., 2013). Terrain in southwestern Greenland above ~650 m above sea level (asl) is often heavily weathered and has apparent ¹⁰Be exposure ages that far exceed the LGM period (Kelly and Bennike, 1985; Rinterknecht et al., 2009; Roberts et al., 2009; Winsor et al., 2015). These results have been interpreted to indicate that the surface of the LGM may not have been more than 650 m asl, which would mean that some mountain summits survived as nunataks during the LGM.

The exposure history of summits interpreted to be potential nunataks could be complicated by non-erosive LGM ice cover. For example, weathered uplands, even those with pre-LGM cosmogenicnuclide inventories in bedrock, and sometimes even in erratics, have been shown to be at least intermittently covered by ice during the LGM (Briner et al., 2006; Corbett et al., 2013; Young et al., 2018). In south-western Greenland, if coastal uplands were covered by ice during the LGM, it would point to a terminus position far onto the continental shelf and possibly at the shelf break. These issues collectively have led to an uncertain reconstruction of the lateral extent and maximum elevation of the Greenland Ice Sheet during the LGM in some areas around Greenland.

Here, we provide new evidence to support an extensive Greenland Ice Sheet configuration in southwestern Greenland during the LGM. We investigate the relationship between the elevation of marine limits and their locations relative to LGM ice extent. Next, we develop a 1-D crustal rebound model for a model-data comparison of relative sea-level history. Finally, we constrain ice thickness with cosmogenic nuclide measurements from erratic and bedrock surfaces atop coastal summits in south-western Greenland. Combined, our data suggest that the southwestern Greenland Ice Sheet terminus extended to the shelf break during the LGM.

We explore the relationship between marine limit elevation and its distance to ice-sheet margins by compiling marine limit data from around the Laurentide, Fennoscandia, Barents–Kara, and Greenland Ice Sheets. We use LGM ice-sheet outlines from published compilations (Funder et al., 2011; Dalton et al., 2020; Hughes et al., 2016) and include 62 marine limit observations from northeastern Canada (Dyke et al., 2005), 55 from Greenland, 22 from Svalbard (Forman, 1990), and 42 from Ireland, the UK, and Scandinavia (Hogan et al., 2016) (Fig. 2). We measure the distance from each marine limit location to the mapped LGM ice-sheet margin. In places where the mapped LGM extent did not reach the shelf break, we made a second measurement from the marine limit location to the shelf break.

We modeled two versions of the ice sheet, one with an LGM terminus at the shelf break (“Max Ice”) and another with the LGM terminus at the location at the Fiskebanke moraines on the inner shelf (“Min Ice” in Fig. 3). We enforced these configurations from 30 ka until 15 ka, which was long enough to relax the underlying crust into equilibrium with the overlying ice load. We applied the same retreat history to both the Max Ice and Min Ice configurations. We initiated deglaciation at 15 ka based on the timing of retreat from the shelf break farther north (e.g., Ó Cofaigh et al., 2013) and then used five mapped Greenland Ice Sheet margin positions for Holocene time steps at 11.6 ka, 10.4 ka, 9.1 ka, 8.1 ka, and 7.3 ka based on data from Young et al. (2020) and Lesnek et al. (2020), with linear retreat rates between each ice margin position. There is evidence that the Greenland Ice Sheet retreated inland of its current position during the middle Holocene. Although the exact distance remains unknown, it likely was no more than a few tens of kilometers from its current position (Young et al., 2021; Lesnek et al., 2020; Ruskeeniemi et al., 2018; Young and Briner, 2015). We represented this by positioning the margin at the modern terminus at 5.5 ka, then 20 km inland of the modern terminus from 4.0 ka to 1.2 ka, then at its modern position at 0 ka, with linear rates of margin migration between these defined times.

To validate this simple model for visco-elastic deformation of the lithosphere, we ran the Benchmark A comparison test described in Martinec et al. (2018) and Spada et al. (2011). Our modeled vertical displacement under the outer 70% of the glaciated area agreed with the Benchmark A results (Martinec et al., 2018) to within 10 m (~10%).

Three previous studies with cosmogenic nuclide measurements are relevant to our work because they include ¹⁰Be ages in the Sisimiut area. Here, we recalculated published ages using the methods described below. Samples collected from bedrock above 650 m asl yielded some ages in excess of 100 ka (Roberts et al., 2009), and samples from lower elevations have ¹⁰Be ages that are >15 ka and thus could be interpreted to indicate ice-free LGM conditions or early deglaciation (Rinterknecht et al., 2009; Roberts et al., 2009). Roberts et al. (2009) interpreted a boundary at ~650 m asl as a possible maximum LGM ice-sheet elevation for the region and supported a limited maximum lateral extent with the LGM terminus placed just offshore (Kelly and Bennike, 1985; Weidick and Bennike, 2007). Winsor et al. (2015) reported ages of ca. 35 ka from boulders at ~600 m asl, which could be interpreted to indicate an even lower LGM ice surface. Winsor et al. (2015) also reported ¹⁰Be ages of boulders that average ca. 14.5 ± 1.0 ka from ~150 m asl and suggested that deglaciation of the Sisimiut area took place at this time. The ca. 14.5 ka age of deglaciation in Winsor et al. (2015) contrasts with a single ¹⁰Be age of ca. 20 ka from a boulder at ~400 m asl only ~10 km away reported by Rinterknecht et al. (2009). Finally, Young et al. (2020) reported three ¹⁰Be ages from ~230 m asl that average 11.5 ± 0.3 ka, which is notably younger than the nearby ages from Winsor et al. (2015). All of these previously dated sites lie to the west of the Taserqat Moraine, which Young et al. (2020) dated to 11.6 ± 0.4 ka (n = 5).

To further explore the unresolved nature of upland glaciation and former Greenland Ice Sheet thickness in the Sisimiut area, we use multiple nuclides, ¹⁰Be, ²⁶Al, and in situ ¹⁴C, due to their combined utility in elucidating glacier history in areas with complex exposure histories (e.g., Beel et al., 2016; Bierman et al., 1999; Graham et al., 2019; Young et al., 2018). We collected three bedrock and six boulder samples from a summit adjacent to the ocean with sample elevations ranging from 355 m asl to 397 m asl (Figs. 4–5). The summit at this “western site” displays evidence of past ice cover due to the presence of glacial erratics. Signs of advanced subaerial weathering, such as protruding quartz veins with ~1–2 cm of relief above the surrounding bedrock, are also visible at this site, which indicates a potentially complex exposure history. We also collected three bedrock samples from a summit at 798 m asl located 32 km inland from the western site (Figs. 4–5). The summit bed-rock surfaces at this “eastern site” exhibit signs of advanced weathering with deep weathering pits; no boulders suitable for sampling were present.

We collected samples using a rock hammer, chisel, and angle grinder from the top ~2 cm of each surface. The boulders were perched in positions that imply stability since deposition. Bedrock samples were collected from local topographic high points to minimize the effect of shielding from seasonal snow cover and surrounding topography.

We processed samples using a series of physical and chemical treatments. At the University at Buffalo (Buffalo, New York, USA), the rock samples were crushed and then sieved to 250–500 µm. Quartz was isolated using froth flotation in an acetic acid solution followed by hydrochloric and hydrofluoric/nitric acid baths. We dissolved clean quartz with beryllium carrier with a concentration of 1074 ppm of ⁹Be and added ²⁶Al carrier (1001 ppm) to the bedrock samples (Corbett et al., 2016). Once beryllium and aluminum were isolated and purified, they were oxidized and loaded into cathodes for accelerator mass spectrometry (AMS) measurement at PRIME Lab at Purdue University (West Lafayette, Indiana, USA). ¹⁰Be and ²⁶Al ages were calculated with the CRONUS v3 online age calculator (https://hess.ess.washington.edu), and ¹⁰Be and ²⁶Al ages were calculated using the Baffin Bay production rate of Young et al. (2013) using the Lm scaling scheme (Balco et al., 2008; Lal, 1991; Stone, 2000). We also processed six quartz separates for in situ ¹⁴C analysis at Lamont-Doherty Earth Observatory following standard procedures detailed in Lamp et al. (2019). To avoid potential carbon contamination from the froth flotation procedure at the University at Buffalo, we leached samples in concentrated nitric acid for several days before processing them further. The measured ¹⁴C concentrations were blank-corrected using a long-term laboratory blank (Lamp et al., 2019). In situ ¹⁴C ages were calculated using the western Greenland (¹⁴C) production rate calibration data set (Young et al., 2014) and Lm scaling (Balco et al., 2008; Lal, 1991; Stone, 2000).

Our marine limit analysis shows that marine limit elevation is significantly related to distance from the hypothesized LGM margin (Fig. 6). We fit data for all four ice sheets to Equation 1 and found r² = 0.73 with C = 4.7 and n = 0.50 in agreement with the theoretical value of n = 0.5 and C of order 1 (Andrews, 1968a). Some of the scatter, we assume, is due to the added influence that ice margin history during deglaciation (variable retreat rate, stillstands, and re-advances) has on the relationship. Of particular interest are the data from southwestern Greenland using the Greenland Ice Sheet ice limit from Funder et al. (2011), which fall outside the trend of data from all other ice sheets (Fig. 6). On the other hand, when using an LGM margin at the shelf break, the best-fit equation resembles that of other ice sheets much more closely.

The LGM reconstruction of Greenland Ice Sheet extent in the Sisimiut region proposed by Dyke (2004) and Funder et al. (2011) results in a marine limit versus distance relationship that differs from that of other Northern Hemisphere ice sheets (Fig. 6) and requires a relationship different from the theoretical parabolic profile (Equation 1). On the other hand, our “Max Ice” scenario, in which the ice-sheet terminus reaches the shelf break, the marine limit versus distance plot becomes similar to the results from other ice sheets (Fig. 6) and fits well with a parabolic profile. These findings support an LGM extent at the shelf break off Sisimiut.

The results from our crustal rebound simulations (see animations in Supplemental Material¹) show relative sea-level histories similar to those from independent observations. Our simulated relative sea-level curve using “Max Ice” and a deglaciation age of 12 ka yields a marine limit for Sisimiut of 136 m asl. This is close to the observed marine limit of 140 m asl and generally compatible with field constraints of relative sea level (Bennike et al., 2011; Funder and Hansen, 1996; Long et al., 2011). The simulated marine limit using “Min Ice” is 47 m asl and has a poorer fit with observations (Fig. 7). Thus, the overall better fit between our simulated “Max Ice” relative sea-level history and observations favors a “Max Ice” scenario.

Despite the agreement between geological observations and the results from our “Max Ice” crustal rebound simulations, the ice-sheet history we use in our model is not well constrained prior to 11.6 ka. In our “Max Ice” scenario, we force the ice margin to retreat from the continental shelf break at 15 ka (see Section 2.2) and then prescribe a linear retreat rate over the next 4.4 k.y. until the ice margin reaches the next recorded position ~115 km inland at 11.6 ka. In our “Min Ice” scenario, we also begin retreat at 15 ka at a constant rate until 11.6 ka, but the ice margin begins much farther inland (Fig. 3). If the Greenland Ice Sheet had a complex ice-margin history with cycles of repeated advance and retreat (e.g., retreat from the LGM limit during the Bølling followed by re-advance and final retreat after the Younger Dryas), such an ice-margin history could lead to a different loading history. Similarly, our assumption that the crust achieved equilibrium with the LGM ice load may not be valid if the Greenland Ice Sheet was at its maximum extent for a relatively short period, such as arriving to its LGM limit not long before it receded from that limit, a result that is supported by our cosmogenic nuclide results (see Section 3.3). If the maximum LGM ice extent were reached relatively late and occupied for only a brief period of time, then the Greenland Ice Sheet would need to have been thicker than in our “Max Ice” simulation to yield a marine limit that matches independent observations. In any case, despite these uncertainties, our simplified ice-sheet history for “Max Ice” nevertheless produced a reasonably close fit with observations (Fig. 7).

We also explored how thick the Greenland Ice Sheet would have to be in our “Min Ice” scenario to reproduce the observed 140 m asl marine limit near Sisimiut. This requires an ice thickness of ~1400 m at the modern coastline, which would require an increase in the yield stress from 110 kPa to 310 kPa to maintain an inner-shelf terminus position. In addition to being a factor of three higher than the typically accepted value (van der Veen, 2013), a yield stress of 310 kPa results in an ice-sheet thickness over 5 km at the ice divide. We find it unlikely that a minimum ice configuration could result in a reasonable data-model fit.

From the western site, three bedrock samples have ¹⁰Be ages of 24.0 ± 0.6 ka, 18.6 ± 0.5 ka, and 15.8 ± 0.4 ka (Tables 1–2). The five samples from erratics at the western site have ¹⁰Be ages of 14.7 ± 0.5 ka, 11.8 ± 0.4 ka, 11.7 ± 0.4 ka, 11.0 ± 0.4 ka, and 10.3 ± 0.3 ka. The ²⁶Al ages at the western site are in good agreement with the ¹⁰Be ages, and the ²⁶Al/¹⁰Be ratios of the samples range from 6.77 ± 0.26–7.78 ± 0.31 (Table 1). We measured in situ ¹⁴C in four samples from the western site (one bedrock and three erratics) to further constrain the exposure history of the site (Table 2). In situ ¹⁴C ages from two erratics are 11.1 ± 1.1 ka and 11.5 ± 1.4 ka, and the bedrock sample has an age of 11.7 ± 1.2 ka. One of the in situ ¹⁴C measurements from an erratic resulted in a concentration of ~210 k atoms/g, which is higher than the site-specific saturation value of ~190 k atoms/g. We suspect this sample was contaminated during the extraction process and disregard it from further consideration.

Two of the erratic samples (18GRO-6 and 18GRO-11) at the western site have ¹⁰Be ages that are beyond the 1σ age range of the others (one older, one younger; Fig. 4), and the average age of the remaining three erratics is 11.5 ± 0.4 ka. We interpret this age of 11.5 ± 0.4 ka as the timing of deglaciation at the western site.

The in situ ¹⁴C ages for the samples from the western site strengthen the case for the 11.5 ± 0.4 ka age best representing the timing of deglaciation. Sample 18GRO-9, an erratic with a ¹⁰Be age of 11.0 ± 0.4 ka, has an in situ ¹⁴C age of 11.1 ± 1.1 ka. Both 18GRO-4 (bedrock) and 18GRO-11 (erratic) have ¹⁰Be and ²⁶Al ages (ca. 15–17 ka) that are older than our interpreted age of deglaciation (ca. 11.5 ka), which we suspect is due to inheritance. The in situ ¹⁴C ages for these same samples (11.7 ka and 11.5 ka, respectively) are within 1σ of our erratic ¹⁰Be-based age of deglaciation of 11.5 ± 0.4 ka. The average of all three in situ ¹⁴C ages from the western site is 11.5 ± 0.3 ka, the same as our average age of deglaciation based on ¹⁰Be ages of erratics (11.5 ± 0.4 ka). Overall, these ages support a deglaciation of the western site around the Younger Dryas/Holocene boundary and imply ice cover prior to that time.

At the eastern site, the three bedrock samples (18GRO-1A, 18GRO-, and 18GRO-2) have ¹⁰Be ages of 79.4 ± 1.3 ka, 63.3 ± 1.0 ka, and 69.0 ± 1.7 ka, with corresponding ²⁶Al ages of 66.3 ± 2.9 ka, 56.1 ± 1.1 ka, and 55.7 ± 1.4 ka (Table 1; Fig. 4). The three samples have ²⁶Al/¹⁰Be ratios ranging from 5.85 ± 0.20–6.40 ± 0.16 ka (Table 1). We also measured in situ ¹⁴C in two of the bedrock samples from the eastern site, which are 17.1 ± 2.5 ka and 19.0 ± 3.2 ka (Table 2; Fig. 4).

The cosmogenic nuclide data from the eastern site exhibit evidence for a complex burial and exposure history. The ²⁶Al and ¹⁰Be ages from the three bedrock surfaces are far older than the regional timing of deglaciation (e.g., Young et al., 2020). The samples have ²⁶Al/¹⁰Be ratios that are 6.10 ± 0.13, significantly lower than 7.3, which has been proposed as an appropriate constant exposure value for high latitudes (Corbett et al., 2017; Halsted et al., 2021). This indicates that these samples experienced episodes of burial, upwards of 250 k.y., likely by overriding ice (Bierman et al., 1999; Beel et al., 2016; Strunk et al., 2017). These durations increase when using a ratio of 7.3. These ²⁶Al/¹⁰Be ratios, which are similar to values that others have found in western Greenland uplands (Roberts et al., 2009; Beel et al., 2016; Corbett et al., 2017), do not allow us to definitively say when this burial occurred, but it is possible that a portion of the burial occurred during the LGM (Beel et al., 2016).

Without further information, whether the eastern site was covered during the LGM remains unknown. However, the apparent in situ ¹⁴C ages provide additional insight. Samples 18GRO-1a and 18GRO-2 are much younger than their ¹⁰Be and ²⁶Al ages and average 18.1 ± 1.3 ka. There are two ways that these results could be interpreted. First, the in situ ¹⁴C ages could represent the timing of deglaciation of the eastern site. Second, the ages could contain some amount of inherited ¹⁴C due to only a brief, minimally erosive overriding event during the LGM, in which case the site was deglaciated sometime after ca. 18.1 ka. In the next section, we explore these interpretations in light of the larger data set from both our data from the western site and prior work.

The cosmogenic nuclide ages at the western site agree with ¹⁴C-dated sediments collected from lakes in the region. Our ages are older than the basal lake ¹⁴C sediment age of 10.5–11.0 cal ka B.P. (re-calibrated) that provides a minimum limiting deglaciation age for the region (Bennike et al., 2011). The deglaciation of the western site also agrees with the average of three ¹⁰Be ages from erratic boulders (11.5 ± 0.3 ka) collected ~10 km south of Sisimiut reported by Young et al. (2020) (Fig. 4). The ²⁶Al/¹⁰Be ratios from the western site samples reveal insignificant burial even considering the two pairs from bedrock with low levels of inheritance (Fig. 8). Based on this, we interpret the age of deglaciation for the Sisimiut coast to be ca. 11.5 ka.

We suggest that both the western and eastern sites were deglaciated at around the same time, in which case the in situ ¹⁴C concentrations are likely due to (1) low-magnitude erosion (if any) of the bed due to non-erosive ice conditions during the LGM and/or (2) an insufficient duration of ice cover to result in complete ¹⁴C decay. The 11.6 ka Taserqat Moraine is draped across the massif of the eastern site (Fig. 4), a mere ~100 m lower in elevation than the summit and only ~1 km to the SE. Given this geometry, we find it unlikely that the eastern site would have remained a nunatak during the LGM. Moreover, considering our exposure ages from the western site and across the general field area, we suggest our eastern site likely deglaciated just prior to Taserqat moraine deposition. Both the western site and the island site of Young et al. (2020) deglaciated at ca. 11.5 ka, the same age as the Taserqat Moraine (11.6 ± 0.4 ka), which implies rapid deglaciation of the region at this time. Thus, the Greenland Ice Sheet likely retreated fast enough to negate any difference in ice thickness at our two study sites, causing them to become ice-free at roughly the same time. Finally, although higher in elevation, the eastern site is farther inland, and generalized ice-surface profiles indicate that both sites would deglaciate at similar times. LGM ice probably would have been sufficiently thick to prevent continued cosmogenic nuclide accumulation through ice during the LGM, a phenomenon that has been observed at thin ice caps elsewhere in the Arctic (Pendleton et al., 2019).

Our exposure ages agree with those of Young et al. (2020) and, combined, suggest that the Sisimiut region deglaciated at ca. 11.5 ka. We note that one of the boulders that we dated (18GRO-9) yields a ¹⁰Be exposure age of 11.0 ± 0.4 ka, but a previous study reported a ¹⁰Be exposure age of 19.6 ± 5.2 ka from the same boulder (GRE-1; Rinterknecht et al., 2009); however, the two ages overlap at 2 sigma. The cause of the age discrepancy may be due to improving laboratory methods like minimizing background ¹⁰B isobaric interference and using improved ¹⁰Be ion exchange chromatography for increased sample purity (Corbett et al., 2016). Winsor et al. (2015) collected nine samples from boulders for ¹⁰Be exposure dating at two sites close to Sisimiut. Four of five samples collected at ~585 m asl yield an average age of 36.2 ± 2.1 ka with one younger age of 17.9 ± 0.8 ka. Winsor et al. (2015) raised the possibility that the site remained above the LGM ice surface; in light of our results, it seems that their alternative interpretation is more likely, which is that the samples contain inheritance and the site was covered by LGM ice. Winsor et al.'s (2015) four additional ¹⁰Be ages, each on a boulder from ~150 m asl, average 14.4 ± 0.9 ka after excluding one young age (10.2 ± 0.6 ka), which they interpreted as an outlier. These ages cannot be easily explained because they are at odds with our results and those of Young et al. (2020), but our results and the relatively younger ages from Young et al. (2020) suggest that small amounts of inherited nuclides may explain these slightly older ages.

Finally, we use likely timing of deglaciation (11.5 ka) to invert the in situ ¹⁴C concentration for a plausible duration of LGM ice coverage of the eastern site. Generally, in situ ¹⁴C becomes either fully saturated when constantly exposed or has decayed to a non-detectable level when fully buried at ~30 k.y. (Graham et al., 2019). Offshore evidence supports a glaciation in southwestern Greenland during Marine Isotope Stage 4 (MIS 4, 75–58 ka) that was more extensive than current LGM reconstructions suggest (Seidenkrantz et al., 2019). Thereafter, the Greenland Ice Sheet is generally thought to have been about as extensive as its current size throughout MIS 3 (58–27 ka) (England, 1999; Larsen et al., 2018). Thus, a plausible glacial history of the eastern site is that it became deglaciated at ca. 58 ka from an MIS 4 overriding event and subsequently remained ice free until the LGM. Following Graham et al. (2019), we model the in situ ¹⁴C concentration in the samples from the eastern site assuming they were saturated with respect to ¹⁴C upon arrival of the ice sheet during the LGM. Next, by assuming no erosion during ice cover and re-exposure of the surfaces at 11.5 ka, we calculate the onset of ice-margin advance and surface burial to be between 19.6 ka and 14.6 ka (Fig. 9). Regardless of the exact timing and duration of LGM ice cover, our in situ ¹⁴C measurements that yield finite exposure ages below the saturation limit mandate that these bedrock surfaces were not exposed throughout the LGM.

The new results presented here point to a more extensive LGM configuration in southwestern Greenland than is currently depicted in most literature. Our comparison of marine limit values from southwestern Greenland to those from other sectors of the ice sheet suggests that an inner-shelf LGM position is incompatible with the relatively high marine limits observed at coastal sites. Our simplified crustal rebound modeling also reveals that a shelf-break LGM configuration is more compatible with existing relative sea-level data than an ice sheet that terminated on the inner shelf. Finally, our new series of cosmogenic nuclide exposure ages supports ice cover of coastal summits and a final deglaciation age for the coast that occurred no earlier than the late Pleistocene/Holocene boundary. Thus, we favor the view that the Greenland Ice Sheet likely terminated at the continental shelf break during the LGM, with coastal deglaciation occurring just prior to ca. 11.6 ka.

Our conclusion is consistent with recent studies that place the LGM terminus position at the shelf break in the Disko Bugt and Uummannaq Fjord regions (e.g., Ó Cofaigh et al., 2013). Furthermore, the results highlight the need for direct sediment dating of the continental shelves around Greenland. Specifically, if the depositional age of the Helle-fiske or Fiskebanke Moraines were known, then it would be possible to constrain the age and location of the terminus of the Greenland Ice Sheet throughout southwestern Greenland more accurately. Our evidence from southwestern Greenland joins other work suggesting that shelf break termini may have been common throughout Greenland (Fig. 1A; England, 1999; Evans et al., 2009; Larsen et al., 2018; Ó Cofaigh et al., 2013; Roberts et al., 2009; Vasskog et al., 2015). Additional offshore studies and further applications of the modeling approaches taken here could test this hypothesis. If the Greenland Ice Sheet terminated at the shelf break throughout Greenland, this maximum extent would be consistent with the LGM extents of other major coastal ice sheets (Hughes et al., 2016; Mangerud et al., 1998; Margold et al., 2018).

Whereas our findings indicate that the coastal mountain summits in southwestern Greenland near Sisimiut were covered during the last glaciation, simple modeling of in situ ¹⁴C inventories suggests the duration of ice overriding was only 5.5 ± 2.5 k.y. starting at 17.1 ± 2.5 ka. Such a late maximum in ice volume is consistent with some ice-sheet modeling studies (e.g., Lecavalier et al., 2014; Buizert et al., 2018) that suggest a late LGM ice phase, perhaps correlative with Heinrich Stadial 1 (ca. 18.5–15 ka). Nevertheless, our cosmogenic nuclide measurements, modeled relative sea-level curves, and modeled ice-sheet profiles suggest that the Greenland Ice Sheet terminated at the shelf break at some point during the LGM or shortly thereafter in south-western Greenland.

We acknowledge Joseph Tulenko for assistance with generating maps of the study location. We thank T. Tuna and E. Bard for careful gas-source measurements at Centre de Recherche et d'Enseignement de Géosciences de l'Environnement (CEREGE). We acknowledge support by National Science Foundation Arctic System Sciences grants ARC-1504267 to J.P. Briner and E.K. Thomas and ARC-1503959 to N.E. Young. We extend thanks to Polar Field Services for field logistical support, and to two anonymous reviewers who took time to critically review this work.

All MATLAB code and input data used to model the viscoelastic relaxation of the lithosphere for comparison to relative sealevel data is available at http://hdl.handle.net/10477/82607.