Evidence of Middle Holocene landslide-generated tsunamis recorded in lake sediments from Saqqaq, West Greenland

The Vaigat strait (Sullorsuaq) in West Greenland is well known for its susceptibility to landslides and historical landslide-generated tsunamis. Recent mapping of the seabed in the Vaigat strait has revealed several prehistoric giga-scale (volumes of 10⁹ m³) tsunamigenic landslides. However, the timing of these giga-scale tsunamis is largely unconstrained, but they are assumed to have occurred after the last deglaciation. Here, we report on lake sediment core records from four coastal lakes located between 19 and 91 m above sea level (a.s.l.) on the Saqqaq foreland at the eastern end of the Vaigat strait. We use a multiproxy approach including X-ray fluorescence (XRF) and magnetic susceptibility core scanning along with a screening for marine diatoms to identify at least two tsunami deposits in two of the four sediment cores. Radiocarbon dating of aquatic macrofossils and bulk samples suggest that the tsunami events occurred at ca. 7.6 and 7.3 ka cal BP. Using a previously published relative sea level curve from Vaskebugt, Arveprinsen Ejland (Alluttoq), located 40 km southeast of Saqqaq, we infer wave runup heights of ca. 41–66 and 45–70 m, respectively, for the two tsunami events. These runup heights from prehistoric tsunamis are 1 to 2 orders of magnitude higher than the historic landslide-tsunami runup heights at Saqqaq which only reached an elevation of ca. 3 m in November 2000. While we found deposits from two tsunami events in the lake sediments, landforms from at least nine giga-scale landslides are present on the seafloor of Vaigat. We infer that these deposits probably represent the two most recent tsunamis identified in the Vaigat strait and that the older tsunamis must have happened between the last deglaciation and the oldest sediment in the lakes, i.e., between ca. 10.0 and 8.5 ka cal BP.

1 Introduction

Landslide tsunamis are among the most devastating natural disasters in fjord settings with several recent examples from Norway, Alaska, and Chile (e.g., Blikra et al., 2006; Sepulveda et al., 2010; Higman et al., 2018; Svennevig et al., 2020, 2023a). A landslide tsunami occurs when large rock masses suddenly fail and slide into a waterbody, usually in the form of a rock avalanche (sensu Hungr et al., 2014, and Hermanns et al., 2021), displacing large amounts of water that travel through the fjord rising to a tsunami wave that can inundate nearby coastal areas (Hermanns et al., 2013). As such, landslides with tsunamigenic potential pose a direct threat to coastal communities, and under a warming climate they are expected to occur more frequently (Patton et al., 2019).

The Vaigat strait (Sullorsuaq; Fig. 1) in West Greenland experienced two landslide-triggered tsunamis in 1952 and 2000 (Pedersen et al., 2002; Dahl-Jensen et al., 2004; Svennevig et al., 2023a). In 1952, the tsunami caused damage to infrastructure in the settlement of Qullissat and the loss of one human life (Svennevig et al., 2023a). In 2000, the landslide at Paatuut generated large waves which hit the settlement of Saqqaq, while the now-abandoned settlement of Qullissat experienced severe damage to buildings at lower elevations (Pedersen et al., 2002). In 2021, a third landslide occurred at Assapaat in between the Niiortuut and the Paatuut landslides (Fig. 1). Although ca. 4×10⁶ m³ (0.004 km³) of material entered the Vaigat strait, no tsunami was observed, which is probably due to local topographical control and the specific dynamics of this landslide (Svennevig et al., 2022). It has recently been shown that the historical landslides in the Vaigat strait were preconditioned by permafrost degradation (Svennevig et al., 2022, 2023a). Farther north at Karrat Isfjord, a large rock avalanche triggered a tsunami in June 2017 that inundated the nearby settlement of Nuugaatsiaq and claimed four human lives (Paris et al., 2019; Svennevig et al., 2020), showing that landslides are also occurring outside of the Vaigat strait. Svennevig et al. (2020) mapped three areas of continued slope deformation near the 2017 landslide, and today the two settlements of Illorsuit and Nuugaatsiaq with a combined 170 inhabitants are closed due to continued risk of landslide-induced tsunamis in the area.

Figure 1 Landslide mapping by Svennevig (2019) and Svennevig et al. (2023b). Bathymetry from BedMachine Greenland v4 (Morlighem et al., 2021) and other map data from SDFI (2018) and ArcticDEM (Porter et al., 2018). A seismic profile line from Vaigat available from Marcussen et al. (2001) is shown, where the sediments along the profile were interpreted as either postglacial slumps or glacial moraines, surrounded by horizontal Holocene sediments.

Svennevig (2019) compiled a preliminary inventory of landslides in Greenland based on available digital elevation models (DEMs) and bathymetric data. In the Vaigat strait, nine large underwater landslide deposits were mapped. They are recognizable by their hummocky and blocky topography in the otherwise smooth U-shaped glacial trough with a Holocene sediment cover comprising the bottom of the Vaigat strait (Fig. 1) (Svennevig et al., 2023b). Seismic profiles from the Vaigat strait show the presence of localized chaotic accumulations of sediment interpreted as either old submarine slides or submarine aggradations from subaerial landslides, with a general thickness of 50–100 m that can locally exceed 200 m (Marcussen et al., 2001; Pedersen et al., 2002). Based on the long runout, large volumes, and giant displaced blocks Svennevig et al. (2023b) suggest that these rock avalanches had a significant tsunamigenic potential, which should leave onshore evidence of tsunamis. However, nothing is known of when the landslides occurred after the last deglaciation or what magnitude the tsunamis could have generated.

Coastal lake sediments are excellent archives for capturing and preserving evidence of past tsunami events as these events have a distinct lithological signature that is different than the normal lacustrine sedimentation (e.g., Bondevik et al., 1997a; Long et al., 2015). If several lakes are cored in the same area at different elevations, i.e., the staircase approach, the sedimentological records can furthermore be used to estimate the runup height of a tsunami and thus its magnitude (e.g., Bondevik et al., 1997b). Evidence of tsunamis is in general rare in Greenland and has only been encountered by coincidence three times in relation to paleoclimatic studies or relative sea level (RSL) reconstructions using lake sediment records. In Loon Lake in East Greenland and in the Ammassalik area in southeast Greenland, deposits attributed to the Storegga submarine slide have been described (Wagner et al., 2006; Long et al., 2008). In Disko Bugt (Qeqertarsuup Tunua) in West Greenland, tsunami deposits attributed to rolling icebergs have been identified (Long et al., 2015).

Until now no systematic study of lake sediment cores has been undertaken to reconstruct the Holocene tsunami history in the Vaigat strait. In this study we use six lake sediment cores from the coastal foreland at Saqqaq in the Vaigat strait to constrain the timing and magnitude of prehistoric tsunamis.

2 Study site: Saqqaq foreland and surroundings

The settlement of Saqqaq (population 132 in 2020) in West Greenland is situated on the southern coast of the Saqqaq foreland protruding into Vaigat (Fig. 1). The Vaigat strait is generally between 500 and 600 m deep with depths of up to 620 m south of Saqqaq, except where landslide deposits, dikes, and sills (Paleocene intrusions) are present. The elevation of the foreland gradually increases from sea level up to ca. 200 m about 5 km from the coast at the foot of the ca. 900 m high mountains (Fig. 2). West of the foreland lies Saqqaqdalen, which is a valley with a flat floor consisting of a fluvial plain that continues into the Umîvik delta.

Figure 2 Lake sediment coring sites at the Saqqaq foreland and surroundings. The foreland elevation ranges from 0–200 m a.s.l. and has multiple lakes in this height interval. The elevation color scale has been partitioned into intervals below and above 200 m to emphasize the morphology. Inferred directional ranges of the tsunamis are indicated with wavy blue arrows. Catchments of the cored lakes are local, and water entering the lakes is from rainfall or snowmelt and is sheltered from the rivers draining the mountains. The 70 m contour represents the coastline at ca. 10.0 ka cal BP as inferred by a radiocarbon dated marine shell (K-994 from Weidick, 1972). DEM and map data are from SDFI (2018) except for the streams which were derived from the DEM using GRASS (2022). Geological map data are from Pedersen et al. (2013).

This exposed coastal configuration and the multiple lakes situated at different elevations, also known as a staircase, makes the Saqqaq foreland ideal for retaining evidence of past tsunamis to determine the timing and elevation of tsunami runup. Weidick and Bennike (2007) provide a minimum age of the last deglaciation of around 10 ka cal BP from dating marine shells at 70 m a.s.l. (Fig. 2). The susceptibility of the lakes to inundation from tsunami waves would have been higher in the Early Holocene considering Saqqaqdalen would have been submerged more than 70 m due to glacial isostatic subsidence, and in the Middle Holocene the glacial isostatic subsidence would depress the landscape ca. 6–25 m (5.8–7.6 ka cal BP) relative to today.

The landscape surrounding the lakes consists of undulating bedrock surfaces with a thin and discontinuous sediment cover. The bedrock of the Saqqaq foreland generally consists of Precambrian granodioritic augen gneiss (Garde and Steenfelt, 1994). An outcrop of the Kingittoq Member, a Late Cretaceous (Cenomanian) fluvial and deltaic deposit, is located on the mountain slopes (Fig. 2) and consists of successions of mudstones, heteroliths, and well-sorted sandstones (Pulvertaft, 1989). At the type section of this member 25 km northwest of Saqqaq, sandstone sheets are often overlain by carbonaceous heterolithic mudstones with thin coal beds (Dam et al., 2009). At Saqqaq, the Kingittoq Member is more coarse-grained and consists of cross-bedded pebble conglomerates, with no coal beds reported (Pulvertaft, 1989).

2.1 The November 2000 tsunami in Saqqaq

Detailed observations of runup heights for the tsunami in Saqqaq in November 2000 have not previously been published. To compare these to the paleotsunamis reported below we have mapped the runup height of this event using recognizable buildings in photographs taken after the event (Fig. 3a). As the land was covered in snow at the time of the event and the tsunami occurred at mean tide, the maximum runup is readily recognizable. West of the peninsula, the maximum runup height was ca. 3 m, whereas east of the peninsula, no visible tsunami runup was observed (Fig. 3a). Saqqaq is not in direct line of sight to the Paatuut landslide 55 km away (Fig. 1); thus the waves that were observed have been reflected and refracted.

Figure 3 (a) Photograph from the Illulissat police the day after the tsunami that reached Saqqaq after the 21 November 2000 landslide at Paatuut. The maximum runup is seen as discoloration, and the removal of the snow is indicated with X, Y, and Z. (b) Topographic map of Saqqaq from ASIAQ Greenland Survey with the viewing direction and field of view of the photo in (a). The maximum tsunami runup is indicated with a red line, and X, Y, and Z refer to the same marked localities in (a).

3 Methods

Four lakes were selected using high-resolution satellite images and a 2 m resolution DEM (SDFI, 2018). Lakes with small, local catchments were chosen, and lakes with large tributaries or visible influx of sediment were avoided. Their sill elevations range from 19 to 91 m a.s.l. (Fig. 2). The lake depth was measured using a GARMIN echo sounder, and the sediment cores were obtained at the deepest part of the lake by hammering a 60 mm piston corer into the lake floor from a Zodiac inflatable boat.

After recovering the sediment cores, they were kept upright and drained before being shipped to Denmark where they were stored under refrigerated conditions at 2–4 °C. In the laboratory, the cores were split, and the stratigraphy was described. The cores were scanned using an ITRAX core scanner to obtain X-ray fluorescence (XRF), X-ray imagery, magnetic susceptibility (MS), and an optical line scan image. The cores were scanned with a 1 mm resolution using a standard setting of 30 kV and 50 mA for the X-ray fluorescence scanning and a dwell time of 30 s. The X-ray settings were 3 kV and 50 mA with a dwell time of 1 s and a step size of 0.2 mm. Both XRF and X-ray used a Rh tube. The MS has a step size of 4 mm, and the optical scan has a resolution of 0.047 mm. The ITRAX core scanner provides counts of ca. 40 elements, and it has been used to identify tsunami deposits onshore and in lake cores (e.g., Chagué-Goff et al., 2017; Shinozaki, 2021). In this study, we use the element titanium (Ti) as a proxy for minerogenic input to the lake (i.e., tsunami deposits) and the element ratio Ca $/$ Fe as a proxy for the marine influence in the sediment. Titanium is suitable to differentiate between marine deposits with a high minerogenic content and lacustrine gyttja deposits with low minerogenic content. As Ti is a common constituent of rocks such as gneisses or schists, it primarily indicates a terrigenous continental source (Rothwell and Croudace, 2015). In our setting, high Ti values are interpreted as minerogenic deposits flushed into the lake from either shoreline or onshore surroundings as the tsunami wave loses energy and the water runs back to the sea. The Ca $/$ Fe ratio is used as a proxy for the input of biogenic carbonates of marine origin (Ca) relative to detrital clay (Fe) of terrigenous origin into the lake and can be used to indicate tsunami deposits (Chagué-Goff et al., 2017).

We describe the tsunami deposits in the cores using typical stratigraphic signatures identified by Bondevik et al. (1997a) and summarized by Long et al. (2015). These are stratigraphical signatures in lake sediments that have been found in lakes inundated by the tsunami waves from the Storegga submarine landslide off the coast of Norway. Here, the base of the tsunami deposit consists of massive or graded sand, suggesting a rapid deposition from suspension overlying an erosional unconformity. The erosional unconformity was caused by the turbulence from the tsunami wave as it crossed the lake threshold. Above the sand, a facies termed “organic conglomerate” may be present. This facies consists of coarse organic detritus with large rip-up clasts of gyttja, peat, and silt. The large rip-up clasts are of irregular form and sizes and occur with twigs and plant fragments in a matrix of gyttja and silt. Lenses of sand may also be present (Bondevik et al., 1997a). One criterion from Bondevik et al. (1997a) we are not able to test is the lateral extent of the tsunami deposits as we did not survey the lake stratigraphy with multiple cores during fieldwork.

The chronology is based on accelerator mass spectrometry (AMS) radiocarbon (¹⁴C) dating of aquatic plant macrofossils (mosses) or bulk sediment samples, as we were unable to identify any terrestrial macrofossils. The local bedrock consists of non-calcareous rocks; hence this specific problem can be excluded when using aquatic macrofossils or bulk sediment samples for dating the cores. The macrofossil samples are unspecified aquatic mosses (bryophytes), and this is not without implications for how the obtained ages can be used. Aquatic mosses will fix reservoir-aged dissolved CO₂ from the water column or sediment pore water in quantities large enough to introduce old-carbon reservoir effects to the radiocarbon dating of the moss (Marty and Myrbo, 2014). Radiocarbon ages above the tsunami unit also have complications. In principle, the sediment just above the tsunami unit should show the minimum age of the tsunami. However, Bondevik et al. (1997b) found that for lakes in western Norway radiocarbon ages above the tsunami deposits are commonly older than the actual age of the tsunami inundating the lakes.

In total, we dated 12 samples to constrain the age of the cores and provide age control of tsunami events (Table 1). We used three sampling strategies to determine the age of a tsunami deposit. To provide a maximum age of the events we dated samples below the tsunami deposit under an erosional unconformity. Bulk sampling of rip-up clasts or macrofossils embedded in the tsunami deposit was also dated to provide a maximum age, while a minimum age is obtained by dating lacustrine samples above the tsunami deposit where lacustrine sediments re-occurred. The ¹⁴C ages were calibrated into calendar years using OxCal v4.4 (Bronk Ramsey, 2009) and the IntCal20 calibration curve (Reimer et al., 2020).

Table 1 All macrofossil samples in this study are aquatic mosses (bryophytes). All radiocarbon ages were calibrated using OxCal v4.4 (Bronk Ramsey, 2009) and the IntCal20 curve (Reimer et al., 2020), except for sample K-994, which has not been normalized for isotopic fractionation to a δ¹³C value of −25 ‰; i.e., we have added 400 years before calibrating the date into calendar years using Marine20 (Heaton et al., 2020) using a local d$R=-\mathrm{49}\pm \mathrm{59}$ from West Greenland (Pearce et al., 2023). Errors in sill heights in our study are calculated from Coeurdevey and Soubirane (2013).

There is no local RSL curve for Saqqaq, so we recalibrate existing isolation basin radiocarbon dates from Long et al. (1999) from Vaskebugt (Kangerluarsuk), Arveprinsen Ejland (Alluttoq) (see Fig. 6). This location is situated 40 km to the southeast of Saqqaq (Fig. 8) and would have deglaciated a few hundred years later than Saqqaq, which would have seen slightly higher uplift rates (Long et al., 1999). Consequently, when using the Arveprinsen Ejland RSL curve on Saqqaq, it will underestimate how much uplift has occurred at any given time.

In addition to the radiocarbon dates used for the RSL curve from Long et al. (1999) we also use radiocarbon dates from Weidick (1968, 1972), where shells provide a minimum age of the last deglaciation at Saqqaq (Weidick, 1968) and a minimum height of the RSL at that time (Weidick, 1972). We have recalibrated all samples so that they are standardized with the radiocarbon ages used in this study (Table 1). The ¹⁴C ages for samples from this study and Long et al. (1999) were calibrated into calendar years using OxCal v4.4 (Bronk Ramsey, 2009) and the IntCal20 calibration curve (Reimer et al., 2020). K-994 has not been normalized for isotopic fractionation to a delta ¹³C value of −25 ‰; i.e., we have added 400 years before calibrating the date into calendar years using Marine20 (Heaton et al., 2020) using a local d$R=-\mathrm{49}\pm \mathrm{59}$ from West Greenland (Pearce et al., 2023).

We screened sediment core SAQ21-06 from relevant depths for diatom assemblages to find marine species that could be used as an indicator of marine water intrusion to lakes. In total five diatom microscopy slides were prepared according to Battarbee et al. (2001). We treated ca. 0.2–0.5 g of fresh sediment with 10 % hydrochloric acid (HCl) and 30 % hydrogen peroxide (H₂O₂) in a water bath for 4 h to remove carbonates and organic material from the sediment. Samples were rinsed with distilled water, and slides were mounted using Naphrax^® with a refractive index of 1.73. Diatoms were identified with an Olympus BX51 light microscope using phase contrast optics at 1000× magnification.

4 Results

We retrieved four lake cores of which two contain tsunami deposits and another two help constrain the age of the tsunami deposits. Lake coring site characteristics are summarized in Table 2, and cores are plotted with sediment proxies in Figs. 4 and 5. All core lengths reflect the thickness of the sediment package until the corer encountered bedrock or an impenetrable substrate. The lakes are situated below the marine limit, and marine sediments should be expected in the bottom of the cores. However, in several of our attempts to penetrate the deepest part of the succession, the core tubes fractured probably because of the stiff nature of the sediments. Accordingly, we were unable to penetrate the entire succession and retrieve marine sediments.

Figure 4 Sediment proxies (optical core and X-ray images, stratigraphy, and XRF) and calibrated radiocarbon ages from the four lakes. T1 and T2 show tsunami sediment units.

Figure 5 Correlation of tsunami facies T1 and T2 in lake cores SAQ21-06 and SAQ21-09. A–D show details of the tsunami units T1 and T2.

Table 2 Lake and coring site characteristics.

4.1 Lacustrine and glaciolacustrine deposits

4.1.1 Gyttja facies

The organic-rich gyttja is the main lithology in all cores. It is brown to brownish-grey and varies from being finely laminated to massive in composition (Figs. 4, 5). Both the laminated and massive gyttjas have a very low minerogenic content characterized by low Ti and MS values.

4.1.2 Interpretation of gyttja facies

The gyttja facies is interpreted as lacustrine sediments deposited under normal conditions in the lakes with input from glaciers. The variations in gyttja composition between the lake records probably reflect different lake sizes and depths, as well as catchment properties.

4.2 Tsunami deposits

We have identified deposits from two tsunami events, T1 and T2, which are both present in cores SAQ21-06 and SAQ21-09 (Figs. 4 and 5).

4.2.1 Tsunami unit 1 (T1)

The younger (T1) tsunami deposit is 28 and 41 cm thick in SAQ21-06 and SAQ21-09, respectively. In both cores, it has a distinct erosional base followed by a layer of massive coarse to very coarse sand (grain size 1–2 mm). In SAQ21-09, T1 has eroded into the sediments of T2, while in SAQ21-06 there is a laminated gyttja facies in between the two tsunami deposits, making it possible to separate T1 from T2. In SAQ21-09 there is a 3.5 cm massive sand layer above the erosional unconformity. The massive sand is followed by 38 cm of organic conglomerate consisting of irregular clasts composed of silt or sand and with plant fragments. Some of the rip-up clasts consist of laminated gyttja, and the size of the clasts decreases upwards in T1 with several lenses of coarse sand visible within the organic conglomerate. The contact to the massive gyttja above is transitional for 1–2 cm.

In SAQ21-06 the corresponding massive sand layer is 0.5 cm thick with an erosional unconformity at the base. It is overlain by 27 cm of organic conglomerate consisting of large rip-up clasts of irregular form and size and with plant fragments, silts, and sands. There are at least two erosional unconformities inside the organic conglomerate with layers of sand. A thin layer of sand constitutes the sharp contact between the organic conglomerate and 20 cm of lacustrine laminated gyttja above.

Our sediment proxies (Fig. 4) – MS, Ca $/$ Fe, and Ti – are all very high at the massive sand at the base of the tsunami sediments and where sand lenses are located in the organic conglomerates, making the sediment proxies much more variable relative to the gyttja.

4.2.2 Tsunami unit (T2)

In core SAQ21-09 the oldest tsunami unit (T2) makes up the bottom of the core. Here, the top of the 26 cm thick T2 sediments has been eroded by the overlying T1. In the other core, SAQ21-06, T2 is 7 cm thick and enclosed by laminated gyttja above and below.

The base of T2 in core SAQ21-09 consists of a 0.2 cm layer of massive coarse sand and very coarse sand. Most of the coarse sand flushed to the bottom of the tube when it was sealed and had to be cut off from the main core. It is not shown here as primary sedimentary structures were not preserved. Above the massive sand is a 3 cm thick layer of organic detritus with sand grading upwards, which is followed by 6 cm of coarse organic detritus with small gyttja clasts and plant fragments in a matrix of fine–medium sand. Above is a 8 cm block of laminated gyttja, where the laminations of the gyttja are cross-cut by the unit below, indicating that the gyttja was deposited prior to the unit below and has been displaced subsequently. Then follows a 0.1 cm layer of coarse sand grading to medium–fine for 1 cm. The contact between the coarse sand and the laminated gyttja below is sharp. Overlying this is 8 cm of coarse organic detritus with gyttja clasts and plant fragments and fine to medium sand.

The 7 cm thick T2 unit in SAQ21-06 consists of coarse organic detritus with clasts with a thin layer of sand at its base and a sharp contact without visible erosion of the laminated gyttja below. There are small clasts of gyttja and more plant fragments than in SAQ21-09. T2 has a sharp contact to the laminated gyttja above.

Our sediment proxies (Fig. 4) – MS, Ca $/$ Fe, and Ti – are all very high in the organic detritus facies and low in the tilted block of laminated gyttja in SAQ21-09. In the upper organic detritus facies in SAQ21-09 the variability is greater due to larger clasts.

4.2.3 Interpretation of tsunami facies (T1 and T2)

We interpret T1 and T2 as tsunami deposits as they contain many of the diagnostic characteristics described by Bondevik et al. (1997a). The lowest part of T1 and T2 consists of coarse to very coarse sands with clear erosional unconformities at the base, and the sand is angular and similar to the sand found in the terrain surrounding the lake. Thus, it is likely that the massive, coarse sand has been flushed into the lakes from the surroundings when the tsunami wave withdrew. In both deposits, this is followed by an organic conglomerate with rip-up clasts, sands, and plant fragments. The organic-rich conglomerates have erosional unconformities and lenses of sand which we interpret as being caused by subsequent waves from the same landslide-tsunami wave train. We interpret the deposits in the two cores as tsunami deposits. While identification of the tsunami deposits is based on visual description of sediments and structures, as well as sedimentological proxies, the correlation of the units between the two cores with tsunami deposits is primarily based on the laminated gyttja separating the two tsunami units T1 and T2 in core SAQ21-06, constrained by the age control. This correlation is supported by the visual appearance and sedimentological proxies of T1 and T2 in the two cores. Since we did not survey the entire lake stratigraphy, we cannot exclude the small possibility that the ca. 42 cm unit of laminated gyttja could be a large rip-up clast and that T1 and T2 may be the result of one tsunami event.

The T2 tsunami facies in the two lakes have the same signatures, and we correlate them based on visual appearance. Although T2 has the same general appearance in the two lakes, there are some notable differences. These reflect that lake SAQ21-09 is much more exposed to tsunamis and has experienced much more turbulence and erosion when invaded by the tsunami compared to SAQ21-06. SAQ21-09 has a coarser sand fraction at the base and no gyttja below the base of the tsunami facies which is laying directly on top of an impenetrable substrate. The sediment core also has sand lenses from distinct waves from the same landslide-tsunami wave train, indicating that the aspect to the strait and exposure to the incoming tsunamis are more important than lake height for capturing tsunamis in the lake sediments. This may also explain why there is no erosional unconformity in SAQ21-06.

4.3 Diatoms

All prepared diatom microscope slides were screened to find marine diatom species that do not normally live in freshwater systems and thus could be used as an indicator for marine water invasion. Diatom assemblages in the samples included species that typically inhabit freshwater and weakly brackish systems. The most common species that were found in almost all the samples were Navicula vulpina, Navicula similis, Neidium iridis, Stauroneis phoenicenteron, and Tabellaria flocculosa. Only one sample (from 61–62 cm depth) in core SAQ21-06 included diatoms that have a marine origin. This sample was taken from the massive sand layer in the base of the T1 deposit and included a diatom species, Cocconeis scutellum, which is a common species in marine and brackish coastal waters (Cremer, 1998; Witkowski et al., 2000; Pearce et al., 2014; Oksman et al., 2022).

4.4 Radiocarbon chronology

The tsunami units (T1 and T2) are present in two cores (SAQ21-06 and SAQ21-09), and combined with the absence of tsunami deposits in cores from two other lakes (SAQ21-07 and SAQ21-11), we constrain the age of the two tsunami using the radiocarbon dating of macrofossils and bulk sediment (Table 1).

Samples for radiocarbon dating of T2 were taken below and above the tsunami sediments in SAQ21-06. The bulk dates below are 8.1 and 7.9 ka cal BP, and the bulk date above is 7.9 ka cal BP, suggesting an age of T2 of around 7.9 ka cal BP. However, a macrofossil (aquatic moss) sampled at the top of T2 yields an age of 7.6 ka cal BP, indicating T2 is younger than the 7.9 ka cal BP that the bulk sediment samples bracketing T2 yield. We suggest that the small discrepancy between the macrofossil (aquatic moss) date at the top of T2 and the bulk date above T2 is due to the re-sedimentation of older carbon in the bulk dates that makes them a few hundred years too old. This is attributed to erosion of older organic material by the tsunami and the subsequent redeposition into the lake as the older organic material is washed back into the lake after the tsunami event (Bondevik et al., 1997b). In core SAQ21-06 the 42 cm long sequence of laminated gyttja between T1 and T2 dates to 7.9 ka cal BP at the bottom and 7.2 ka cal BP at the top, respectively above T2 and below T1. This combined with the close grouping of ages in and around the T2 unit in cores SAQ21-06 and SAQ21-09, indicates that some of the moss and bulk sediment samples are a few hundred years too old probably because of a local reservoir effect (Marty and Myrbo, 2014; Strunk et al., 2020) and re-mobilization of old organic material by the tsunami (cf. Bondevik et al., 1997b). Accordingly, we conclude that T2 occurred around 7.6 ka cal BP as this is the youngest age related to T2.

A bulk sediment sample below the erosional unconformity of T1 in SAQ21-06 yields a maximum age of 7.2 ka cal BP. Above T1 are two dates yielding ages of 5.8 and 6.7 ka cal BP; i.e., T1 seems to have occurred between 7.2 and 6.7 ka cal BP. The oldest of the two dates, 6.7 ka cal BP, is taken immediately above the T1 unit, while the 5.8 ka cal BP date is sampled above T1 and may be less affected by re-sedimentation of older sediments washing into the lakes from the surroundings after the tsunami occurred (Bondevik et al., 1997b). However, in a nearby lake at lower elevation (SAQ21-07, 19 m a.s.l.) the sedimentary record shows continuous deposition of lacustrine sediments extending back to 7.3 ka cal BP without any evidence of T1. To explain this discrepancy we see two possibilities. Either T1 is not a tsunami deposit but a local phenomenon or something is wrong with the dating of T1. As T1 is recorded in two lakes (SAQ21-09 and SAQ21-06) and is stratigraphically consistent, it is not likely that it represents a slumping event or another local sedimentological phenomenon, and we interpret it as a tsunami deposit. Instead, we suggest that T1 must have occurred at ca. 7.3 ka cal BP (within the uncertainty of the ¹⁴C ages) to satisfy the chronology in both SAQ21-06 and SAQ21-07.

5 Discussion

We consider that the topographical setting, the relative sea level history, and the age of the events, in addition to the signatures for tsunami sediment facies in lakes used by Bondevik et al. (1997a), are helpful in evaluating possible depositional processes. For example, the topographical setting of the lake basins makes it impossible for flash flooding from rainstorms and glacial meltwater to reach the lakes (Fig. 2.), and these processes cannot explain the lake stratigraphy. Earthquakes can cause slumping in the lake sediments; however, West Greenland is considered tectonically stable with minimal active faulting (Voss et al., 2007), and the geometry and succession of the sediment facies do not suggest that slumping from the lake sides has occurred. Lake inundation by tsunami waves is the remaining depositional process that can explain the presence of sediments with the stratigraphical signatures of the deposits identified here. In the following, we use the minimum runup heights as a tool to assess the source and magnitude of the tsunami triggering mechanism.

5.1 Tsunami runup heights and magnitude

Our results show that at least two tsunamis, T1 and T2, invaded the lakes at Saqqaq at ca. 7.3 and 7.6 ka cal BP, respectively (Fig. 5). We use the recalibrated RSL curve from Long et al. (1999) to plot the two events on curves showing the elevation of the lakes relative to the sea surface over time since deglaciation (Fig. 6). The older tsunami (T2) entered lake SAQ21-09 in 7.6 ka cal BP when the lake had a minimum sill height of 41 m. Lake SAQ21-11, which is located 700 m further inland, has no evidence of tsunami deposits since 8.5 ka cal BP. It then had a minimum sill height of 66 m. The bracket of minimum and maximum runup height of tsunami T2 is thus 41–66 m (Fig 7). T2 is also recorded in lake SAQ21-06, which at the time had a minimum sill height of 21 m.

Figure 6 Recalibrated relative sea level curve from Vaskebugt, Arveprinsen Ejland (Long et al., 1999). Sample V8 was excluded from the polynomial fit.

When the younger tsunami (T1) entered lake SAQ21-09 at 7.3 ka cal BP, the lake had a minimum sill height of 45 m marking the minimum runup height. Like for T2, the maximum runup height for T1 is also constrained by lake SAQ21-11, which at the time had a minimum sill height of 70 m. This brackets the runup height of T1 to 45–70 m. T1 is also recorded in lake SAQ21-06, which at the time had a minimum sill height of 25 m.

The absence of traces of tsunamis in the recorded lake stratigraphy in lake SAQ21-11 does not preclude that the lake was invaded by T1 and T2, as a tsunami wave loses energy as it moves inland and may leave no trace in the sedimentological record. It does, however, provide an upper bound on elevation where the tsunami was able to leave its signature in the stratigraphy.

5.2 Triggering mechanism of the tsunamis

Figure 7 illustrates that 41 m is the minimum runup height for any tsunami found in this study, and we use this elevation to evaluate potential tsunami triggering mechanisms. Several mechanisms exist that can produce tsunamis in a setting such as Saqqaq, such as calving from glaciers, rolling icebergs, wave runup during storms, earthquakes, and marine-terminating landslides.

Figure 7 Lake elevations of four lakes extended back in time using present-day elevation and the recalibrated relative sea level curve (Fig. 6). Tsunamis T1 and T2 have been plotted, and the intersection with the lake elevation curves shows the lake elevations when lakes SAQ21-06 and SAQ21-09 were invaded by tsunami waves. This provides the minimum runup height of T1 and T2. There is no sediment older than 8.5 ka cal BP in our core record, and after ca. 6 ka cal BP a tsunami would have had to have a minimum runup of up to 19 m to reach a lake. These elevations and ages define fields where tsunamis are possible but not recorded by our data (light blue “tsunami possibility fields”). The maximum runup heights of two historical tsunami waves from 1952 (Svennevig et al., 2023a) and 2000 from the present study are shown for comparison.

Calving from glacier fronts or capsizing icebergs can generate tsunamis. MacAyeal et al. (2011) suggested that as a rule, the open-water upper limit to tsunami wave height from such a calving or capsizing event is 0.01H, where H is the initial vertical dimension of the iceberg. Between 7 and 8 ka cal BP ago the calving fronts of the marine-terminating outlet glaciers draining the ice sheet into Disko Bugt (Qeqertarsuaq) would have been located roughly where they are today (Weidick and Bennike, 2007), and the icebergs would have had to travel 45 km through the fjords of Kangia (Jakobshavn Isfjord) or Torsukattak before entering Disko Bugt itself (Fig. 8). The maximum water depth at the mouths of the two fjords would have been 280 m (Kangia) and 380 m (Torsukattak) at 8 ka cal BP, restricting the vertical height of the icebergs to less than 400 m and consequently a theoretical upper limit of the open-water tsunami wave height to ca. 4 m. A typical depth at the bottom of the Torsukattak fjord is 650 m; thus a tsunami initiated from calving would have an upper limit to wave height of ca. 7 m, which would dissipate during the 65 km travel distance to Saqqaq (Fig. 8). Even though we do not consider the amplification of the wave heights as they near the shoreline, ice calving and capsizing icebergs do not have the ability to create waves with a runup height of 41 m.

Figure 8 Bathymetry of Disko Bugt, Vaigat, and neighboring fjords (Pedersen et al., 2013; Morlighem et al., 2021). Areas with a high risk of landslides (Dahl-Jensen et al., 2004) are shown as purple shadings along the north and south coasts of Disko Island and Nuussuaq. Mapped postglacial landslides (Svennevig, 2019) are shown as red polygons, and historical landslides in Vaigat in 1952 (Niiortuut, yellow), 2000 (Paatuut, blue), and 2021 (Assapaat, green) are also indicated.

The wave climate in Disko Bugt is driven by a combination of locally generated wind waves and swells from the Davis Strait. The fetch for wave generation is limited, except towards the Davis Strait in a southwesterly direction. Storms mostly occur during winter when sea ice at the coast limits almost all wave activity. During the sea-ice-free period, katabatic winds blow off the ice sheet from the east where fetch is limited. Consequently, wave runup height is much less than the minimum found in any lake.

Earthquakes also have the potential to produce tsunamis. However, West Greenland is situated on a tectonically passive continental margin and is considered tectonically stable (Voss et al., 2007).

Subaerial marine-terminating landslides have the potential to produce very large runup tsunamis in confined fjord settings (e.g., Dahl-Jensen et al., 2004; Sepulveda et al., 2010; Higman et al., 2018; Svennevig et al., 2020). The 1958 Lituya Bay landslide tsunami is the largest such tsunami recorded in historical times with a near-field (ca. 1.5 km from the source) maximum runup height of 524 m (Franco et al., 2020), demonstrating the extreme tsunamigenic potential of landslides. Svennevig (2019) mapped several landslides in the Vaigat strait, and Pedersen et al. (2002) and Dahl-Jensen et al. (2004) marked the Vaigat strait as an area where the bedrock configuration preconditions landslides (Fig. 8). Svennevig et al. (2023b) found that at least nine very large postglacial rock avalanche deposits are present on the seafloor of Vaigat as little as 15 km from Saqqaq and argue, based on morphology along with the large volume and long runout, that these had a significant tsunamigenic potential.

Based on the high runup heights and the presence of at least nine giant rock avalanches in Vaigat, we conclude that tsunamis T1 and T2 recorded in our study were initiated by subaerial landslides entering the sea in the Vaigat strait.

5.3 Magnitude estimates of the tsunamigenic landslides

We cannot, based on the available data, pinpoint which specific landslides triggered the tsunamis that we have reconstructed from the lake sediments. However, based on our estimated minimum and maximum runup heights and the minimum and maximum distance to the largest landslides in the Vaigat strait (from Svennevig et al., 2023b), we can assess the landslide volumes involved using the semi-empirical SPLASH equation (Oppikofer et al., 2018). The SPLASH equation is a parametrization of a power-law relationship between landslide volume, distance, and vertical runup of the tsunami wave. Parameters were found by least-square-fitting data from a limited number of large rock slope failures and the observed runup at different distances from each landslide. It provides a semi-empirical prediction of runup height from preliminary regional studies. We use it here since it provides an estimate of the magnitude without using more advanced tsunami models, which is beyond the scope of the present paper.

The distance to the giant landslides in the Vaigat strait (where they entered the sea) is 90 to 25 km (Fig. 8), and the runup for T1 is 45–70 m and for T2 41–66 m, as presented above (Fig. 7). Applying these numbers to the SPLASH equation we get that a landslide initiated 90 km away would have had a volume of 1.8 to 3.9 km³ to produce tsunami T1 and 1.5 to 3.3 km³ to produce tsunami T2. A landslide at the closest giant landslide deposit, just across the Vaigat strait, 25 km to the southwest of Saqqaq (Fig. 8), would have had a volume of 0.34 to 0.74 km³ to produce tsunami T1 and 0.27 to 0.63 km³ to produce tsunami T2.

Svennevig et al. (2023b) calculate volumes of the three largest giant postglacial landslides in the Vaigat strait (that may represent several events) to be ca. 1.7 to 8.4 km³. The largest of these, the Ujarrassusuk rock avalanche complex, is situated 25 km from Saqqaq just across the strait. Had this rock avalanche occurred in a single event, a runup height of 280 m would be expected in the Vaigat strait (also using SPLASH; Oppikofer et al., 2018).

This is significantly different from the up to 0.74 km³ we infer above for the highest possible tsunami (T1). From the nine postglacial rock avalanche complexes described by Svennevig et al. (2023b), we would expect more evidence of tsunami activity than the two events (T1 and T2) reported above and significantly higher waves. As our study does not cover the age interval from deglaciation at 10 to 8.5 ka cal BP, it is possible that deposits from the inferred multiple larger tsunamis are present in the tsunami probability field bracketed by these ages (Fig. 7).

Two recent landslide-triggered tsunamis were observed in the Vaigat strait in 1952 and in 2000. Both landslides originated from the slopes at the Nuussuaq Peninsula, the same coast Saqqaq is located on. In 1952 after the Niiortuut landslide, a volume between 1.8 and 4.5×10⁶ m³ (0.0018–0.0045 km³) entered the sea and caused an estimated runup in Saqqaq between 1 and 1.8 m (Svennevig et al., 2023a). During the 2000 Paatuut landslide 30×10⁶ m³ (0.03 km³) entered the sea (Dahl-Jensen et al., 2004) and triggered a tsunami that had a runup of 3 m in Saqqaq (Fig. 3). Residents reported “strange waves” for 2.5 h. The location of the landslide and the duration of the wave runup indicate multiple reflections off the steep slopes of the Vaigat strait. We conclude the T1 and T2 tsunamis also had an area of origin in the Vaigat strait. However, the historical landslides in Vaigat in 1952 and 2000 were at least an order of magnitude smaller than those that produced the T1 and T2 tsunamis we find in the lake sediment cores. We also note that there is a tsunami probability field of relatively large young tsunamis with a runup of up to 19 m younger than 6 ka cal BP (Fig. 7). Constraining events in the two tsunami probability fields should be the priority of future studies.

With an increasing elevation and distance to the shore, landward lakes are expected to capture fewer and thinner sand beds, with a decreasing grain size (Bondevik et al., 1997a). The sampled lakes are at different elevations with varying exposure to reflected tsunami waves. Lakes SAQ21-09 and SAQ21-11 are far more exposed to direct waves coming from the Vaigat strait than SAQ21-06 and SAQ21-07, and a landslide tsunami reaching the latter two lakes would either wash over the Saqqaq foreland or reflect off the northern coast of Disko Island before reaching Saqqaq (Figs. 2 and 8). However, eyewitness accounts and videos from the 2017 Karrat fjord landslide tsunami show the complexity of waves reflected, refracted, and diffracted in a fjord, especially when sea ice is present (Paris et al., 2019). The first wave does not necessarily have the highest amplitude, and waves can come from multiple directions. Therefore, at the Saqqaq foreland, the exposure of lakes to tsunami waves coming from the Vaigat strait can be controlled either by the landward distance and elevation or by the exposure of lakes on the western side of the foreland. Potentially, for every tsunami, waves are captured by each lake at the Saqqaq foreland.

5.4 Implications for other records of tsunamigenic lake sediments in Disko Bugt

Other tsunami deposits have been found in lakes surrounding Disko Bugt. At Vaskebugt, Arveprinsen Ejland (Fig. 8), Long et al. (2015) re-appraise two lake cores used for basin isolation dating in a previous study (Long et al., 1999) and propose that tsunamis generated by rolling icebergs could have deposited sediments in some of the lakes. The lakes have an ideal setting for capturing tsunamis, and cores V1 and V8 may contain tsunami deposits that previously have been interpreted as internal slumping in the lakes (Long et al., 2015). Of particular interest is the tsunami unit in core V1, which is a 0.5–2 cm thick unit of coarse sand and angular gravel with an erosive base. It is enclosed in freshwater gyttja; however, in some instances it directly overlies the marine sediments, resembling the tsunami sediments described by Bondevik et al. (1997a). Presently, the sill height of lake V1 is 61 m a.s.l. and would have been located at 0 m a.s.l. at 9.5 ka cal BP (Fig. 6). The amount of erosion is unknown, but the tsunami must have occurred after 9.5 ka cal BP. The lack of absolute age control on the tsunami facies limits the use of sill heights for discrimination between tsunamigenic mechanisms, but if for example, the tsunami had occurred 500 years after basin isolation – at 9 ka cal BP – the lake would have risen to a sill height of 11 m a.s.l., which would exclude rolling icebergs as the tsunamigenic mechanism based on the arguments discussed above in Sect. 5.2. Lake V8 is also located at Vaskebugt and presently has a sill height of ca. 13 m and a modeled age of isolation of 6.7 ka cal BP. It is suggested by Long et al. (2015) that the potential tsunami unit in this lake could represent a marine flooding of the lake basin and adjoining surfaces. The V8 radiocarbon sample (5.2 ka cal BP) appears not to have been sampled at the marine contact but above the tsunami facies, and it is excluded from the RSL reconstruction in both Long et al. (1999) and the present study (Fig. 4). Thus, the tsunami responsible for the tsunami facies in lake V8 should have occurred between 6.7 and 5.2 ka cal BP when the sill heights were between 0 and 12 m a.s.l. This is inside the younger tsunami possibility field (Fig. 7), and both landslides and rolling icebergs could have been the mechanism initiating the tsunami.

Long et al. (2015) also describe a likely tsunami deposit in lake IV4 at Tasiusarsuit in the southern part of Disko Bugt (Fig. 8) and attribute it to an iceberg-generated tsunami occurring at 6.1 ka cal BP (6242–5910 ka cal BP). This single event is recorded in lake IV4 in the 200-year window where the estimated > 3 m runup height from an iceberg-generated wave was able to reach the lake. The interpretation that icebergs generated the waves that entered the lakes at Tasiusarsuit and Vaskebugt relies on a topographical analysis and the then-current knowledge of landslide- and iceberg-generated tsunamis in Disko Bugt (Long et al., 2015). An important assumption made by Long et al. (2015) is that the 2000 Paatuut landslide is the largest tsunamigenic landslide that is likely to have occurred in the Vaigat strait, and due to large distances from areas with a high frequency of landslides in the Vaigat strait and the south coast of Disko Island, the wave energy would have dissipated before reaching the lakes at Vaskebugt or Tasiusarsuit. Landslide tsunamis may also have reached lake IV at Tasiusarsuit from the south coast of Disko Island where Svennevig (2019) mapped giant landslides (Fig. 8). However, the seafloor off the south coast of Disko Island has only a few depth soundings (Morlighem et al., 2021), and no presence of submarine landslide landforms has been reported.

The present study shows that there is onshore evidence of giant landslide-generated tsunamis in the Vaigat strait, and these could have reached Vaskebugt and explain the presence of the tsunami units in lakes V1 and V8. The age control of these units is not of sufficient quality to unambiguously explain if these deposits were generated by waves from landslides or rolling icebergs.

Identifying additional tsunamis and constraining their timing are outstanding research questions, and the scientific literature on lake deposits (i.e., Long et al., 2015) and gaps in the bathymetric data indicate that this will be revealed by conducting additional offshore and onshore research in the Vaigat strait and Disko Bugt areas. Additional fieldwork in the lakes around Saqqaq could reveal more tsunami deposits and provide better identification and improved age control over tsunamis that are reported here.

6 Conclusion

We analyzed lake sediment cores from the Saqqaq foreland in West Greenland and found sediments that we interpret to be deposited by tsunamis. We can distinguish between at least two tsunami events (T1: 7.3 ka cal BP; T2: 7.6 ka cal BP), each with multiple waves, and determine minimum runup heights of 45 and 41 m, respectively. We find no tsunami deposits in a higher elevation, which helps us constrain the maximum possible runup height to 70 and 66 m. Due to the large runup heights and the presence of deposits from giant rock avalanches in the Vaigat strait 25 to 90 km from Saqqaq, we conclude that the tsunamis were likely triggered by giant landslides.

Two historic landslide-triggered tsunamis in 1952 and 2000 had runup heights of 1–1.8 and 3 m in Saqqaq, respectively, which shows that historic events are at least an order of magnitude smaller than T1 and T2. However, tsunamis with a runup height of up to 19 m could have occurred in the past 6000 years without being captured in the cored lakes as this interval is not covered by our survey. Likewise, the time between our oldest dated sediment core bottom (8.5 ka cal BP) and deglaciation (10 ka cal BP) defines a possibility field where older giant landslide tsunamis could have occurred. Nine giga-scale landslide deposits are described from the Vaigat strait, and as we only identify deposits from two resulting tsunamis in our study, we infer that the deposits from these other landslides could be encountered in this possibility field.