Coupled simulation of landslide, tsunami, and ground displacement for the 2017 Nuugaatsiaq event in Greenland

We investigated the entire sequence of the tsunami event caused by a large landslide that occurred on 17 June 2017, in Karrat Fjord near Nuugaatsiaq village in western Greenland, and clarified the sequential processes of this cascading phenomenon. The volume of the landslide was estimated about 44 × 10⁶ m³ based on seismological analysis at seven observation stations across Greenland and an empirical relation. We conducted sequential simulations, consisting of (1) the landslide descent into the fjord based on topography, (2) tsunami generation and large-scale propagation, and (3) ground motion caused by water surface elevation due to the tsunami, considering both static and elasto-dynamic solutions. A 1 m water surface elevation causes ground displacement of up to 0.1–1.0 mm in coastal areas, and this can be detected by seismometers. This event provided a rare chance to verify the integrated model using only local seismic records in the case of no coastal tidal measurement. The timing of simulated processes agrees well with observations, but uncertainties in landslide volume remain a major factor affecting tsunami amplitude and coastal impact. The detailed seismic signals captured both near and far from the epicentre shed light on the multi-stage dynamics of such cascading phenomena. To monitor landslide occurrence and estimate potential tsunamis early, such analyses should be systematically conducted in real time and could be utilized in early warning systems.

1 Introduction

As the climate is warming, coastal and marine-terminated glaciers in Greenland are retreating (e.g., Meredith et al., 2019; Constable et al., 2022). The retreat of marine-terminating glaciers in Greenland has various consequences for coastal biochemistry and ecosystems, including a reduction of nutrients and productivity once these glaciers have retreated inland (Meredith et al., 2019). Moreover, glacier retreat in general has exacerbated the hazards associated with tsunamis: the emergence of new water areas, the release of constraints on nearby unstable slopes, the degradation of permafrost, and reduced freezing conditions all contribute to increased tsunami hazards in fjords depending on the geological context (Higman et al., 2015; Svennevig et al., 2020). In Polar Regions such as Greenland, such hazards can have important social consequences, such as the temporary or permanent relocation of communities (Matti et al., 2023).

The 17 June 2017 landslide took place in the Karrat landslide complex, located on Ummiammakku Mountain in Karrat Fjord, which includes three landslide-prone areas (Schiermeier, 2017; Bessette-Kirton et al., 2018; see also map in Svennevig et al., 2020). The fjord includes several marine-terminating glaciers, but they are located upstream and did not play a role in this specific event. The 2017 event involved 58 million m³ of material, including 45 million m³ reaching the fjord and generating the tsunami (Gauthier et al., 2018). It is important to note that other events took place earlier; at least 3 rock avalanches were identified in 2009, 2016, and 2017 on the Karrat landslide complex using Sentinel 1 & 2, and Landsat imagery (Svennevig et al., 2020). Svennevig et al. (2020) suggest that permafrost thawing due to climate change favours landslides on this type of steep and unstable slope. This type of event represents a threat to human life and has important social consequences. The landslide of 17 June 2017 led to a tsunami, which devastated the village of Nuugaatsiaq 30 km away and subsequently led to the decision to relocate people after the event (Matti et al., 2023).

Greenland's fjords, shaped by glacial erosion, feature steep walls and deep basins, making them particularly vulnerable to landslides and subsequent tsunamis. Typical fjords in Greenland vary in width from a few hundred meters to several kilometres, with depths often exceeding 500 m and lengths of tens to hundreds of kilometres (e.g., Batchelor et al., 2019). The steep topography combined with permafrost degradation increases the likelihood of large-scale mass movements, as seen in the Karrat Fjord landslide of 2017. While detailed statistics on landslides in Greenland remain limited, the risk of similar events has been studied in the Uummannaq fjord system with retreating glaciers and unstable slopes (NGI, 2021). The first recorded tsunami triggered by a landslide is the 1952 Niiortuut landslide-tsunami event (Svennevig et al., 2023), which also had an associated single fatality and is attributed to permafrost degradation in western Greenland. More recently, a tsunamigenic landslide on 16 September 2023 in an uninhabited fjord in East Greenland was also recorded on seismic networks globally (Svennevig et al., 2024; Carrillo-Ponce et al., 2024), while more local records have been found in lake sediments (Korsgaard et al., 2024). Glaciers and permafrost in Greenland have been experiencing a growing mass deficit in response to warming temperatures (Otosaka et al., 2023). There are consequently growing risks associated with climate-induced changes in Greenland's coastal and fjord environments.

This paper aims to demonstrate the potential of seismic monitoring to enhance an alert system on tsunamis along coastal area including fjords. Seismic monitoring has a long history in Greenland (e.g., Dahl-Jensen et al., 2010; Clinton et al., 2014), and the technique has been particularly used for glaciology (Veitch and Nettles, 2012; Walter et al., 2013; Röösli et al., 2014). The calving of icebergs can be recorded and analysed through seismograms (e.g., Sergeant et al., 2016) and geodetic observations, for example (Nettles et al., 2008). Glacial earthquakes have also been detected seismically associated with both the abrupt sliding of fast-moving ice streams as well as icebergs (Joughin et al., 2008). A widespread seismic monitoring system can be used for the early detection of seismic waves brought by an earthquake or a tsunami if the data can be processed and correctly interpreted in near real-time, as the velocity of these seismic waves is faster than the tsunami wave propagation. However, the coastal measurements of seawater level change remain very few. To monitor the tsunami propagation, it is requested to make the most of the seismic data. In this paper, we first analyse the seismic data to estimate the volume of the landslide. Then we perform numerical modelling of the landslide and tsunami using the VolcFlow and FUNWAVE codes. Finally, we model the ground displacement induced by the modelled tsunami in order to compare with the seismic observation recorded at Nuugaatsiaq village (station NUUG). The whole forward modelling process from landslide, tsunami to ground motion will clarify the probable and detectable ground displacement level near the coastline according to the ongoing tsunami propagation. Overall, this paper aims to demonstrate the concept, yet there are limitations discussed in Sect. 4, and further research and developments would be needed before this concept can be effectively applied.

2 Method and model setting for the 2017 Landslide event

2.1 Seismic station at Nuugaatsiaq (NUUG)

The 17 June 2007 landslide was detected as an earthquake (71.640° N, 52.344° W, 0 km depth) equivalent to a magnitude of 4.2 and the origin time of the event was 23:39:12 UTC (USGS, https://earthquake.usgs.gov/earthquakes/eventpage/us20009nlg/executive, last access: 18 June 2026) (Fig. 1a). The source mechanism has been seismologically studied using the waveform inversions (Poli, 2017; Chao et al., 2018; Xie et al., 2020), and tsunami simulation has been carried out (Chao et al., 2018; Paris et al., 2019). Although there is no direct measurement of water surface elevation during this tsunami event, several videos filmed by the inhabitants are available on YouTube (Underwood, 2017). It is estimated that the tsunami reached 1–1.5 m in height at Nuugaatsiaq and runup flooded up to 9 m in height (Strzelecki and Jaskólski, 2020). It is also reported from the field survey that tsunamis reached as high as 90 m along the coastline on the same side as the landslide and 50 m across the Karrat Fjord near the landslide point (Georgia Institute of Technology, 2017). The seismic station NUUG in Nuugaatisaq (https://ds.iris.edu/ds/nodes/dmc/specialevents/2017/06/22/nuugaatsiaq-greenland-landslide-and-tsunami/, last access: 18 June 2026) recorded the ground motions due to the landslide and probably the following tsunami (Fig. 1). We are particularly interested in the late oscillation (Fig. 1b) seen on a long period range. Chao et al. (2018) considered that this oscillation might have been generated due to the tsunami wave push applied to the coastline in the middle between the tsunami source and NUUG station. Paris et al. (2019) considered that this might have been caused by the quasi-static water surface elevation near NUUG station. In general, the seismographs are useful to detect distant events to determine the source parameters rapidly, as the seismic waves generally propagate with a velocity of 3 km s⁻¹ or higher for S-waves in the crust. On the other hand, tsunami waves may propagate with a velocity of tens to hundreds of m s⁻¹ according to the sea depth, briefly one-tenth of the elastic wave velocity. Thus, this difference in travel time is used to give early warning of tsunami propagation at other locations globally.

Figure 1 Map of the seismic stations and the UD displacement (HHZ component in centimetres) at the NUUG station. The central time zero is set to the origin time (23:39:12 UTC). (a) Map of stations used in this study. The open star shows the epicentre location of the landslide event determined by USGS, and the grey triangles show the seismic stations. (b) high-pass filtered at 0.005 Hz, (c) high-pass filtered at 0.1 Hz, and (d) band-pass filtered at 0.02–0.1 Hz.

Figure 1b, c, and d show the UD displacement at the NUUG station with different frequency bands. We removed the DC offset and linear trend, and applied a taper to the segmented waveform, which included sufficient pre-event and post-event portions. We then corrected for the instrumental response so that the amplitude response becomes flat between 0.005 and 45 Hz. The sensor at the NUUG station is STS-2, which has a flat response between $\mathrm{1}/\mathrm{120}$ and 50 Hz, so we think that the frequency response after the instrumental correction is reliable. We applied a non-causal 4th order Butterworth filter with different frequency ranges. The top panel shows the high-pass filtered record at 0.005 Hz. It shows a harmonic signal with a period of 150 s starting about 400 s after the origin time of the landslide, which is interpreted as the main tsunami-induced response. The middle panel shows the high-pass filtered record at 0.1 Hz. The high-frequency ground motion succeeding the landslide block movement is commonly observed, but the source of the high-frequency waveforms is not well resolved. Yamada et al. (2013) suggested the internal collision with the returning mass and other complicated processes on the slipping surface, and Doi and Maeda (2020) as well as Yamada et al. (2020) suggested the slide and debris flow. Thus, we interpret this high-frequency waves may be generated by the interaction between the water and landslide mass dropping into the water. The bottom panel shows the band-pass filtered record at 0.02–0.1 Hz that includes the main movement of the landslide and a chain of vibrations of smaller amplitude, which is related to the tsunami propagation.

2.2 Simulation Procedure

The phenomena are complex, including the landslide, tsunami propagation, and seismic wave propagation. However, since each process occurs in distinct domains and time scales, these processes can be modelled sequentially. Figure 2 shows our simulation strategy for the whole phenomenon.

Figure 2 Flowchart of the simulation strategy for the landslide-tsunami-ground motion phenomenon.

2.2.1 Single force inversion from seismograms

We start from a purely seismological approach of the single-force inversion with long-period seismic waveforms (e.g., Ekström and Stark, 2013; Yamada et al., 2013). We use the far-field six regional stations in Greenland (Fig. 1a) and then estimate the source time function and mass of the landslide. We removed the closest NUUG station from this inversion since it is too close to apply the point-source approximation. The inversion was performed in the frequency domain with a limited frequency window (Nakano et al., 2008). The source is assumed to be a point source at the location of the landslide (52.34° W, 71.64° N, Depth 0 km). We used the AK135 velocity structure to compute the Green's function (Kennett et al., 1995). The instrumental response was removed, and a 4th-order Butterworth filter with corner frequencies of 0.02–0.1 Hz was applied to the seismic records.

2.2.2 Landslide and tsunami generation and propagation

We perform a landslide tsunami simulation according to Mulia et al. (2020) by coupling the pyroclastic flow model, VolcFlow, with the Boussinesq wave model, FUNWAVE-TVD, using a 2-D numerical grid. We utilize the IBCAO v4.0 digital elevation model from the GEBCO database for both bathymetry and topography. The data are provided in the Polar Stereographic projection coordinates (EPSG:3996, units in meters) with a true scale at 75° N, referenced to the WGS 84 horizontal datum and Mean Sea Level vertical datums. The original 200 m spatial resolution is resampled to 100 m to suit our simulation requirements. We model the generation phase of the landslide tsunami applying VolcFlow, a 2-D depth-averaged numerical model. VolcFlow was originally developed to compute the dynamics and emplacement morphology of avalanches (Kelfoun and Druitt, 2005) and was extended to simulate the dynamics of both the landslide and water (Kelfoun et al., 2010). The model considers the landslide mass and water as two fluids and simulates the dynamics by solving the shallow water equations. This provides a more realistic landslide-induced tsunami initiation than prescribing a static initial water displacement or a rigid body slide. At each simulation time step, tsunami generation is simulated by accounting for the interaction between the landslide and water. The landslide displaces the water through momentum transfer and the elevation of the seafloor by its mass thickness, while the water affects the landslide by reducing the effective mass density underwater and exerting drag on the mass. VolcFlow has been widely applied in various landslide tsunami studies (Giachetti et al., 2011, 2012; Paris et al., 2011; Mulia et al., 2020).

VolcFlow is employed over a duration of 120 s to cover the entire process of landslide and initial tsunami generation. We assume an avalanche density of 2500 kg m⁻³ and a water density of 1027 kg m⁻³ to provide realistic environmental parameterization. The landslide source is horizontally covered by an area of 10 × 9 pixels on the 100 m grid, corresponding to a net surface of approximately 1200 m × 580 m along the major and minor axes. In total, 54 grid cells are assigned to the landslide. The numerical precision is influenced not only by the model resolution and grid spacing but also by the landslide source representation and local bathymetry. Considering the long-period nature of tsunami waves and the smoothly resampled topography, the inclusion of a finer grid is expected to have a negligible impact on the simulations.

To incorporate frequency dispersion effects during tsunami propagation, FUNWAVE-TVD is employed to model wave dynamics starting from 120 s onward. This transition time is selected to ensure the landslide-generated waves are fully developed, as VolcFlow results indicate the avalanche ceased at about 120 s. The output velocity field and water elevation from VolcFlow are implemented as boundary conditions in FUNWAVE-TVD. Unlike seismic tsunamis, nearshore landslide tsunamis have wavelengths in the intermediate water depth regime (e.g., wavelength-to-depth ratio ≈ 10). Therefore, frequency dispersion and non-hydrostatic effects should be considered. To resolve these dynamics, FUNWAVE-TVD (Shi et al., 2012) is a fully nonlinear Boussinesq wave model developed from the original FUNWAVE (Wei et al., 1995; Kirby et al., 1998). It resolves strongly nonlinear wave interactions (large wave-height-to-depth ratios) and improves linear dispersion in intermediate water depths by evaluating horizontal velocity at an optimized reference depth. Furthermore, the model employs robust numerical schemes, including wetting-drying treatment and shock-capturing (Shi et al., 2012), that enable stable and efficient simulations. The model has been extensively validated (Wei et al., 1995; Kirby et al., 1998; Shi et al., 2012) and is widely applied to wave propagation, nearshore wave transformation, tsunami inundation, and meteorological tsunamis (Abadie et al., 2012, 2020; Grilli et al., 2013; Kirby et al., 2013, 2016; Schambach et al., 2019; Paris and Ulvrova, 2019; Ho et al., 2023).

2.2.3 Ground motion due to seawater level change

Finally, we estimate the ground displacement at the NUUG station in Nuugaatsiaq. As shown in Fig. 2, we adopt two approaches. One is based on the Boussinesq problem (see Sects. S2 and S3 in the Supplement), which provides the analytical solution in a semi-finite elastic medium due to the vertical charge on the surface (Boussinesq, 1885). This solution or similar analytical solutions in the elastic medium are used for various geoscience applications to estimate the elastic displacement field due to the surface charge and discharge (for example, from ice sheets and glaciers, surface water reservoirs, mining exploitation, etc.) (e.g., Pinel et al., 2007; Bettinelli et al., 2008). On the other hand, we adopt a finite difference method (FDM) for calculating the ground motion in space and time in the elastodynamic equation, where the seismic waves are propagating (e.g., Aochi and Madariaga, 2003). The seismic waves usually propagate at a speed of a few kilometres per second, and the static displacement remains after the passage of the seismic waves. For both approaches, we assume a homogeneous elastic medium with a rigidity μ = 34.1 GPa, corresponding to the S-wave velocity of v_s = 3530 m s⁻¹, which is a typical value for crustal bedrock from CRUST1 Global Model of Earth's Crust (Laske et al., 2013). The Boussinesq solution is calculated once per second using the tsunami height at the same time, while the FDM simulation is carried out with a time step of 0.005 s continuously from the beginning of the landslide simulation. The synthetics obtained from FDM are integrated once with respect to time to obtain the displacement.

3 Simulation Results

3.1 Single-force inversion from seismograms

Figure 3 compares the synthetic and observed waveforms at NUUG and JIG3 stations from the single force inversion from seismograms. The synthetic waveforms are computed from the source time function estimated from the waveform inversion, while the NUUG station is not used to constrain the model parameters during the inversion. The detailed result is shown in Sect. S1. The waveform agreement is good at station JIG3 and the other stations. The obtained source time function (top panels of Fig. 3) shows that the horizontal particle motion is dominant in the north-south direction, and the vertical component is larger than the horizontal ones. Based on the vertical component, the duration of the event is about 80 s. The maximum amplitude of the source time function is 0.2 × 10¹² N. According to the scaling law in Ekström and Stark (2013), the mass is estimated as 0.11 × 10¹² kg. Assuming an average rock density of 2.5 × 10³ kg m⁻³, the total volume is roughly 44 × 10⁶ m³. According to Ekström and Stark (2013), the empirical relationships have a standard deviation of ±50 % in estimation between the maximum force and the mass.

Figure 3 Source time function and waveform fitting. The x-axis shows the time after the origin of the landslide. Columns show the EW, NS, and UD components from the left. Top: Source time function obtained from the waveform inversion with the frequency of 0.02–0.1 Hz. Middle: Comparison between the observed (black) and synthetic (red) waveforms at NUUG station, not included in the inversion. Note that the horizontal components at NUUG station were oriented to 79 and 169° from North, so they were rotated correctly. Bottom: Comparison for the JIG3 station used in the inversion.

The synthetic waveforms at the NUUG station also show a good fit. We note a long-period signal at 100–250 s in the NS component, similar to the obtained source time function. The landslide may have been finished in about 100 s, while the signal force inversion detected the beginning of the tsunami generation (100–250 s), which reached up to 90 m locally (Georgia Institute of Technology, 2017). The obtained source time function here is more complex than the previous analyses in the literature (Chao et al., 2018; Xie et al., 2020), inferring the complex process of landslide and tsunami generation.

3.2 Landslide and tsunami generation and propagation

We first used the bathymetry data used in Paris et al. (2019) to independently estimate the landslide volume, and we obtained an estimate of 49.7 × 10⁶ m³, hereafter called the AP model. This is consistent with our seismological estimate (44 × 10⁶ m³) and the originally reported one (∼ 50 × 10⁶ m³) by Paris et al. (2019). However, this is subject to uncertainty between 33.4 and 76 × 10⁶ m³ (Bassette-Kirton et al., 2017; Chao et al., 2018; Gauthier et al., 2018; Paris et al., 2019). Figures 4 and 5 show snapshots of the landslide-induced tsunami generation until 120 s and the tsunami propagation after 120 s for the reference landslide volume (AP model). We observe that the first tsunami wave front arrives in front of NUUG station at about 500 s, and the long-period harmonic wave oscillates until 1000 s. The wave height in front of the NUUG station reaches about 0.7–0.8 m. Later, the tsunami wave field becomes more dissipative due to the multiple reflections, and the tsunami wavelength and width become smaller and more varied.

Figure 4 Snapshot of the simulation of landslide and tsunami generation. The black area indicates an ongoing landslide on the target slope. The landslide arrives at the sea in around 30 s. The generated tsunami height is shown in the colour scale in meters. The map uses the Polar Stereographic projection coordinates (EPSG:3996, true scale set at 75° N). The map scale is in kilometres.

Figure 5 Snapshot of tsunami generation and propagation. The white area corresponds to the land above sea level. The triangle shows the position of NUUG. The map coordinates are the same as in Fig. 4.

3.3 Ground displacement due to seawater level change

According to the water surface level change obtained in Fig. 5, we compute the displacement field in an elastic medium. Figure 6 compares the two approaches at the NUUG station position. Additional comparisons for simpler cases are given in Sect. S4. Without a filter, the two estimations are very similar in displacement, allowing to analyse the causality better (Fig. 6a). The travel time of the seismic body waves is about 5–10 s at the distance of 30 km, so the seismic response appears quickly. The first movement appears about 30–40 s, reflecting the timing of the tsunami generation due to the landslide. In the observation record (Fig. 6b), the long-period ground motion generated by the tsunami is much larger than the body waves associated to landslide mass movement. On the filtered ground motions (Fig. 6c), we observe an oscillation amplitude of the static solution (Boussinesq solution) is larger than the dynamic one (FDM), although the observed amplitude is much larger (Fig. 6d). The impact of the tsunami wave approaching the NUUG station becomes apparent after 400 s and periodic (mono frequency) oscillations are observed between 500 and 1200 s. These timings correspond to the passage of the tsunami wave near the station inferred from the snapshot of the tsunami propagation (Fig. 5). The phases between Fig. 6c and d may not be perfectly matched. This is because the time series are aligned at the origin time of the landslide, and the origin of the tsunami generation may not be sufficiently constrained.

Figure 6 Comparison of the surface displacement calculated at Nuugaatsiaq, (x,y) = (−1636.7 km, −1223.4 km). (a) Raw synthetic ground displacement in cm as calculated by the analytical Boussinesq solution and the FDM simulation. (b) Observed ground displacement in cm, lowpass filtered for 0.005 Hz. (c) Filtered synthetic ground displacement in cm from panel (a) at the bandpass window between 0.005 and 0.01 Hz. (d) Observed ground displacement in cm, filtered between 0.005 and 0.01 Hz. Note that the filtered records do not start from zero since we used a non-causal Butterworth filter.

Figure 7 Spatial distribution of peak ground displacement (PGD). The synthetics, calculated at the centre of each mesh, are filtered between 0.005 and 0.01 Hz. A star represents the landslide position, and a triangle indicates the position of NUUG. The left and right panels show the horizontal and vertical components, respectively.

Figure 7 presents the spatial distribution of the maximum displacement of the ground surface simulated by FDM. Here we apply a bandpass filter between 0.005 and 0.01 Hz (100–200 s). We find that the vertical component of the ground motion is dominant and larger than the horizontal components. We observe a displacement of up to 0.01 cm along the coastline, attenuating with distance from the coast. We demonstrated here that the tsunami can deform the surrounding solid earth and therefore, the tsunami movement is detectable from a seismic station near the coastline. It should be visible in such a fjord context in which the seawater level changes by an amplitude of 1 m with a wavelength of a few kilometres. The observed ground motions during the passage of the tsunami at NUUG were larger than the simulated ones. This may be due to frequency limitation in the numerical simulations and uncertainty in the model parameters (see next section).

4 Discussion

4.1 Uncertainty in amplitude

We have conducted a chain of simulations that connect the occurrence of the landslide and the subsequent generation and propagation of the tsunami to the elastic ground displacement at the site of a seismic station (NUUG), located at a distance of 30 km from the landslide. While the simulated ground motions show good agreement in time with those recorded by the seismic station, the amplitude of the displacements is underestimated.

There are several possible reasons for this underestimation of ground motion amplitude. First, we assumed a quite high rigidity μ of the medium of 34.1 GPa. A lower rigidity is plausible, starting from 10 GPa, which corresponds to v_s = 2000 m s⁻¹. As the amplitude is proportional to $\mathrm{1}/\mathit{\mu }$ (see Sect. S2), an up to three times larger response could be expected.

Figure 8 Comparison of the simulated coastal water levels during the tsunami passage at Nuugaatsiaq for different model parameters. The parameters and main results of the simulations are given in Table 1.

Table 1 Comparison of landslide models and main effects on the resulting simulated tsunamis at Nuugaatsiaq. It is worth noting that our seismological analysis gave an approximate estimation of 44 × 10⁶ m³ for the landslide volume. Smooth 3 and smooth 5 are the averaged values over 3 × 3 or 5 × 5 grids over the original topography. Note that the Water-Reach Time indicates the time when the landslide mass reached the water, and the First-wave time indicates the time of the first negative wave amplitude appearing in the water. Since the tsunami wave height is first positive as seen in Fig. 4, the First wave time is the sum of the Water reach time and half-period of the tsunami.

Secondly, the estimation of landslide volume may remain uncertain up to a factor of 2. Let us carry out parameter studies by changing the volume of the landslide (up to twice the AP model) and the topography where the mass slides. Figure 8 compares the seawater level change at Nuugaatsiaq during the tsunami passage for various landslide models (Table 1). The amplitude of water surface elevation at Nuugaatsiaq becomes approximately double (0.8 m in the AP model to 1.6 m in the AP × 2 model). On the other hand, the phase of the water level time series does not change significantly, since the timing of the tsunami generation process (the time when the landslide reaches the water and the first wave is generated) is unchanged. Since ground displacement is proportional to the instantaneous change in water level, an uncertainty factor of two in the water level amplitude implies an equal uncertainty factor in the ground displacement estimation. We also tested the smoothing effect on the topography model (smooth 3 and smooth 5), but found no significant differences in the resulting tsunami generation and propagation.

4.2 Comparison with other tsunami seismograms

There have been several cases in which the terrestrial geophysical instruments could detect the tsunami propagation. Table 2 lists a few examples of tsunami (or seawater change) observations from the seismograms and tiltmeters. The tiltmeter is very sensitive to a small deformation, and the seawater change of tens of centimetres can be detected at a far distance (50 km). For example, Nawa et al. (2007) analysed the records from pressure gauges and broadband seismometers in Antarctica for signals from the 2004 Indian Ocean tsunami generated by an earthquake in Sumatra. A tilt effect of several µGal (10⁻⁸ m s⁻²) for about 0.2 m of seawater level change was observed at a frequency range between 0.3 and 0.6 mHz. During the 2010 M_w 8.8 earthquake in Maule, Chile, tsunami propagation was detected by tiltmeters along the Chilean coastline (Boudin et al., 2013). The tilt response was about 0.05–0.01 µm at 7 km from the coastline for a water surface elevation of about 10 cm. The same tsunami was observed even along the Japanese coastline (Kimura et al., 2013; Kubota et al., 2020), with water surface level anomalies of about 20–40 cm. This could be observed up to 50 km away from the coastline with about 5 × 10⁻³ µrad. Compared to the past observations, the estimated tsunami height of the 2017 Karrat Fjord event was locally higher, and the seismic station was closer to the coastline (< 1 km). Therefore, the ground oscillation could be clearly observed in broadband seismograms in velocity and displacement without the help of tiltmeters and gravimeters. Our parameter study (Sects. S3–S4) shows that the seawater level should change by 1 m in the surrounding of a few-kilometre area in order to measure a ground displacement of 1 mm or larger. Such a configuration is plausible in fjord environments.

Table 2 Summary of tsunami observation by seismograms and tiltmeters.

4.3 Implications for risk management

Cascading risks induced by climate change are explicitly considered in the 6th assessment report (AR6) of the Intergovernmental Panel on Climate Change (IPCC), specifically in the Polar cross-chapter paper of the report of the working group II (Constable et al., 2022). However, the specific issue of tsunamis triggered by increasingly unstable slopes in a context of retreating glaciers is only implicitly considered, as part of a broad range of cascading impacts from climate change. We argue that it is important to recognize and assess this risk due to its potential to become a substantial threat to human life and key infrastructure, especially considering that these cascading risks are also well identified as a gap of knowledge in AR6. In Greenland, mapping of existing landslide deposits in Western Greenland and inhabited parts of eastern Greenland is already underway (Svennevig, 2019). In addition, climate model projections – particularly those identifying areas vulnerable to permafrost thaw, whether currently or in the future – and the location of large calving glaciers provide valuable insights for assessing future hazard risks.

Our results from the parameter studies and the case application demonstrate that the propagation of tsunamis in fjords can be effectively monitored in near real time using seismic data. Analysing the source mechanism of seismic wave radiation allows us to identify landslides and earthquakes (i.e., moment tensor and single force), even if the earthquake location includes uncertainty. This is relevant for Arctic communities living close to fjords as it paves the way for the development of tsunami alert systems. To realize this potential, the concepts presented here should evolve toward an operational system, which requires advanced demonstration in real environments and validation. For example, future applied research and development projects could consist of demonstrating the concept on local sites, including by deploying adequate optical-fibre cables along coastlines. Concurrently, a systematic identification of coastal settlements concerned by this hazard will be essential to assess its importance and the need for alert system deployments. Recent research has made substantial progress in the area of human settlements at risk from permafrost thaw and sea-level rise in the Arctic (Tanguy et al., 2024). Similar efforts could be undertaken to assess the potential of the threat from tsunamis to the settlements identified by Tanguy et al. (2024).

5 Conclusion

We conducted a comprehensive simulation chain of the cascading tsunami event that occurred on 17 June 2017 in Karrat Fjord, western Greenland. Using seismic records from seven stations across Greenland, we derived a source time function based on a single-force model and estimated the landslide volume (approximately 44 × 10⁶ m³) using empirical relationships. This estimate lies within the range of previously reported values. We then simulated the descent process of the landslide and its interaction with the fjord water, and subsequently reproduced the generation and propagation of the resulting tsunami using available topography and bathymetry data. Ground displacement associated with sea-level changes caused by the tsunami was also modelled and compared with the seismic observations recorded at the NUUG station located along the fjord. This event provided a rare opportunity to validate an integrated modelling framework using only local seismic records, in the absence of direct coastal sea-level measurements. The results show good agreement between the timing of the simulated signals and the observed waveforms, indicating that both the landslide initiation and the subsequent tsunami propagation speed were accurately represented. These findings are consistent with the results of Paris et al. (2019), who analysed the same tsunami event. On the other hand, although the predominant tsunami period was successfully reproduced, the subsequent wave trains lost consistency due to multiple reflections within the fjord. Tsunami amplitudes remain sensitive to uncertainties in the initial landslide volume. Ground displacement calculated using both Boussinesq theory and the 3D finite-difference method showed consistent results, suggesting that the ground displacement response to the sea-level change was largely quasi-static. Simulated vertical displacements reached the millimetre scale for a water surface elevation of 1 m. It should also be noted that filtering seismic records is sensitive at low frequencies. To our best, the filtered observation shown in Fig. 6 seems correctly extracted. Therefore, there remains an important discrepancy, and this is likely because of uncertainties in subsurface rigidity and landslide parameters. Finally, this study demonstrated that tsunami propagation can be traced along the coastline via seismic observations. Coupled modelling – from landslide dynamics to tsunami propagation and seismic response – offers detailed insights into the spatiotemporal evolution of such events. The dominant period between 100 and 200 s is similarly found both in the observation and simulation. The direct observation implies propagating tsunami without inverting any source parameters, so that the target-oriented placement of the captor is reliable for this period. Such an approach can potentially enhance tsunami hazard assessment in fjord environments and contribute to early warning capabilities in the future. Seismic data can constrain not only the landslide parameters but also the tsunami progression in near real-time. Looking ahead, the deployment of the new seismic observation techniques, such as optical-fibre sensing along coastlines or the seafloor, could significantly improve our ability to monitor and understand similar cascading hazards, including glacier-related seismic activity.