Evidence of tsunami forces shaping clifftop boulder fields from integrated morphometric and numerical modeling approaches

TL;DR

This study shows that tsunamis, not storms, likely moved large boulders on Mallorca's cliffs, using modeling to improve coastal hazard assessments.

Contribution

A multidisciplinary approach combining morphometric analysis and numerical modeling to identify tsunami-induced boulder transport mechanisms.

Findings

Velocity thresholds, not wave height alone, are key to explaining boulder transport by tsunamis.

Only long-period, high-magnitude southern waves can displace observed boulders in Mallorca.

Storm wave simulations failed to displace boulders, supporting the tsunami hypothesis.

Abstract

Cliff-top boulder deposits are rare but striking geomorphic features found along high-energy coastlines. Their formation mechanisms—whether from storm waves or tsunamis–remain a subject of scientific debate. This study examines tsunami-induced boulder transport on the steep southern coast of Mallorca Island (Western Mediterranean), integrating detailed morphometric analysis with high-resolution simulations using a non-linear numerical model. Boulder characteristics, including volume, elevation, distance from cliff edges, and orientation, were mapped and used to compute mobilization thresholds under subaerial, submerged, and joint-bounded block (JBB) conditions. Our simulations tested 360 tsunami scenarios varying in wave height, period, and direction. Results show that wave height alone is insufficient to explain boulder transport; instead, velocity thresholds–particularly for saltation and JBB lifting—provide more accurate indicators. The model successfully replicates tsunami propagation, cliff impact, and inland flooding, revealing that only long-period, high-magnitude waves from southern directions exceed the thresholds needed to displace the observed boulders. Boulder deposits differ across the study area: smaller, higher-elevation blocks in the western sector are only mobilized by a narrow range of extreme tsunami conditions, while larger, lower-elevation eastern boulders respond to a broader spectrum. Storm wave simulations failed to reach or displace the boulders, strengthening the tsunami hypothesis. These findings highlight the critical role of numerical modeling in tsunami hazard assessment and call for re-evaluation of other cliff-top deposits globally. The methodology presented here demonstrates a robust, multidisciplinary approach to interpreting boulder emplacement and contributes to refining coastal hazard mapping in tsunami-prone regions.

Figures

Tsunami (red) and wind wave (blue) $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document}$ and $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document}$ mobilization threshold for each block considering (A) sub

Maximum wave height *H* (A) and total velocity *u* (B) at the impact point of each boulder for all the tsunami waves simulations. Results are grouped in colours by offshore direction $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _0$$\end{document}$ , the brighter the colour the higher the period *T*. Background lines stand for wave height (A) and velocity (B) thresholds. Total velocity magnitude *u* is co

Mean velocity at the boulders position at the impact instant respect to all the offshore wave heights and periods and directions for both sectors. Black lines stand for: sliding velocity threshold, point lines; rolling, dashed lines; saltation, dot-dashed line; JBB saltation/lifting, solid line.

Tsunami waves propagation phases: (a) drawback, (b) wave approach, (c) impact, (d) flooding and (e) ebb-flood. Upper panels show the tsunami wave bidimensional approach at Els Bancals. Red dots stand for the positions of the boulders. Lower panels show the cross-shore profile of the ES (black line), the sea-surface elevation (blue line, left axis) and the cross-shore velocity component $\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt

Sea surface elevation at Els Bancals on the Impact time-step for: (A) a 10 m tsunami wave (T=1800 s) and (B) a 10 m storm wave (T=10 s). Red dots stand for the position of the boulders.

Full text

Introduction

Tsunamis are long gravity waves that are generated in the ocean by abrupt geophysical events such as undersea earthquakes, submarine landslides, or volcanic eruptions. Among them, those caused by undersea earthquakes are most frequent^1^. Sometimes these processes work in conjunction with each other to enhance the amplitudes of the tsunamis. The energy carried by tsunami waves can cause significantly greater damage than typical storm surge waves due to their longer wavelengths, higher velocities, and ability to propagate across entire ocean basins with minimal energy loss. The primary effects of tsunamis in coastal areas include extensive flooding, high-energy wave impacts, coastal erosion, strong currents, and their ability to transport large volumes of debris.

In recent decades, significant advancements in monitoring systems, data acquisition technologies, and numerical modeling techniques have facilitated precise recording, analysis, and investigation of the complex processes governing tsunami generation, propagation, and subsequent inundation impacts. These developments have been particularly evident in recent major events such as the 2004 Indian Ocean tsunami^2,3^, the 2016 Fukushima tsunami^4,5^, and the recent 2022 Hunga Tonga-Hunga Ha’apai eruption^6,7^. The high-resolution data obtained from these events have substantially improved the accuracy of tsunami hazard quantification and risk assessment. However, such comprehensive datasets remain largely confined to these recent occurrences.

Although there are no data sources comparable to those available for recent events, tsunamis have left behind geological evidence in the form of tsunamites, sedimentary deposits that provide invaluable insights into past seismic activity and its impacts on coastal regions. These deposits are typically characterized by unique sediment layers that can include a mix of marine and terrestrial materials, anomalous boulder placements, and microfossils that indicate sudden marine incursion onto land. Studies of tsunamites around the world, from the Mediterranean to the Pacific rims, have enabled scientists to reconstruct the frequency, magnitude, and effects of historical tsunamis, extending our understanding of these phenomena beyond the limits of written records. Ultimately, the analysis of tsunamites not only helps in the historical understanding of tsunamis but also aids in refining current models used for predicting potential future tsunami risks.

The formation of cliff-top boulder deposits is a subject of ongoing debate among geologists and oceanographers, primarily concerning the differentiation between deposits left by tsunamis and those created by other high-energy events like storms. Critics argue that certain geological features attributed to tsunamis deposits, such as the size, the sorting of sediments or the imbrication of larger clasts, can similarly result from storm surges, making it challenging to definitively attribute such layers to tsunami activity without additional contextual evidence. Focusing on tsunami boulder deposits, large clasts displaced from coastal or nearshore environments by high-energy tsunami waves, it is evident that their size, orientation, and spatial distribution are influenced by hydrodynamic forces, substrate characteristics, and coastal geomorphology. This type of deposit is found worldwide (^8–11^, among others).^12^, along with subsequent studies by^13,14^, developed parametric approaches to elucidate the formation mechanisms of tsunami boulder deposits. These approaches are based on boulder morphometrics and provide wave height and velocity thresholds for boulder mobilization, considering various mobilization conditions, including submerged, subaerial, and boulder lifting scenarios. These parameterizations enable the estimation of the energy required for the transport of boulders. Several studies focused on the identification and genesis of these deposits rely on these formulations to ascertain their origins (^15–19^, among others).

However, the processes of wave propagation, both for tsunamis and storm waves, are highly nonlinear and complex, particularly as they approach shallow waters^20^. Therefore, parametric or simplified approaches to these processes are insufficient for determining the specific phenomena that have caused the displacement of boulders at the wave-rock impact point. This complexity underscores the necessity for a multidisciplinary approach that integrates geomorphological data with historical tsunami records and contemporary modeling to accurately identify and attribute the origins of cliff-top boulders associated with tsunamis.

In this study, we investigate the potential mechanisms underlying the formation of cliff-top boulder deposits by integrating detailed geomorphometric measurements with numerical modeling, specifically focusing on boulder deposits along the southern coast of Mallorca Island (Balearic Islands, Western Mediterranean).Fig. 1(A) Location of the Balearic Islands (Minorca, Mallorca, Ibiza) and USGS earthquake data since 1700 at the Western Mediterranean (through Qquake QGIS software,^21^). (B) Boulders location at the study site. Colours stand for boulder volume. (C) Boulders axes definition. (D) Western boulder field at Es Bancals. Boulders are located at 37 m above sea level. (E) Eastern boulder field at Es Bancals. Numbers refer to boulders. A person in a circle as scale.

Es Bancals cliff-top boulders

A flat lying Upper Miocene Reef complex builds up Es Bancals study site (location on Fig. 1A,B), in the south of Mallorca Island (Spain). On top of the Reef Complex, Middle Pleistocene carbonate aeolianites and colluvial deposits exist^22^. A single cliff, 34 m high, formed by the two units, exists at the western sector (Fig. 1C). To the best of the authors’ knowledge, this sector hosts the highest-positioned coastal boulders reported in the Western Mediterranean: 37 m above sea level and at a 34 m distance from the cliff edge, with volumes of up to 0.99 m \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^3$$\end{document} (Fig. 1D). A wave-cut platform at 10 m a.s.l. separates both units in the eastern sector of the study area (Fig. 1E). In the Eastern Sector the three eolianites levels are more resistant than the colluvial (mostly marls) deposits and form two minor (3–4 m high) topographic steps, separated by two flat areas located at 15–18 and 22 m above sea level. At these upper terraces, large imbricate coastal boulders appear (Fig. 1E). Further details of boulder position and size are provided at^23^.

The southern part of Mallorca Island, where Es Bancals study site is located, is the continental area closest to the main tsunamigenic sources of the Western Mediterranean which, according to^24^, are located offshore Algeria (Fig. 1A). At least six tsunamis have been catalogued over the last centuries affecting the Balearic Islands: in 1365, 1756, 1802, 1856, 1980 and 2003^25^.^26^ described in his chronicles a very large tsunami in 1756 that flooded more than 2.4 km inland in Cala Santanyì (SE Mallorca).

The most recent tsunami reaching the Balearic Islands shores took place on May 21, 2003 and was generated by the Zemmouri earthquake (offshore Algeria). This tsunami affected parts of Spain, France and Italy, but it mostly damaged harbour facilities south and east of the main Balearic Islands, including 3 m high waves on Ibiza Island. Tsunami simulations of the 1856^27^ and 2003^28^ events predicted tsunami wave impact in the exact places where cliff-top boulders on the Balearic coasts appear (Fig. 1 B). Despite storm wave height and energy being much larger from the north and west at the Balearic Islands^29^, these boulders are located in the eastern and southern areas of the islands^18^.

Morphometric analysis and mobilization thresholds

Cliff-top boulders at the study area are primarily accumulated in two distinct groups: the western sector (WS) and the eastern sector (ES). The WS includes 11 studied boulders located at an average height of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$36.46\pm 1.22\,\hbox {m}$$\end{document} above sea level and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$37.14\pm 6.88\,\hbox {m}$$\end{document} from the cliff edge. These boulders exhibit an orientation of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$220\pm 39.94^{\circ }$$\end{document} relative to the north. The boulders are relatively homogeneous, with average dimensions of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$1.08\times 0.70\times 0.34\,\hbox {m}$$\end{document} and an average volume of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$0.31\pm 0.31\,\hbox {m}^{3}$$\end{document} . The ES contains 29 studied boulders at an average height of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$20.24\pm 3.45\,\hbox {m}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$49.09\pm 16.11\,\hbox {m}$$\end{document} . In contrast to the WS, the boulder dimensions in the ES are more heterogeneous, with volumes ranging from 0.23 to 6.72 m³. Notably, they are spatially distributed from west to east within the sector, showing a gradual increase in volume along the East cliff ledge (Fig. 1 C). Their average orientation aligns with that of the WS at \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$220.53\pm 22.81^{\circ }$$\end{document} , with a dip angle towards the sea of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$26.95\pm 13.73^{\circ }$$\end{document} . Boulders characteristics for both sectors are detailed in the Table 1 and Table 2, respectively.Fig. 2. Tsunami (red) and wind wave (blue) \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} mobilization threshold for each block considering (A) subaerial, (B) submerged and (C) JBB conditions. Right panels show the boxplot of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} for both sectors. In (D) Velocity thresholds for each block considering sliding, rolling, saltation and JBB lifting conditions. Vertical dashed line separates western and eastern sectors. Right panels show the boxplot of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{slid}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{roll}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{salt}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{lift}$$\end{document} for both sectors.

Three initial conditions are considered to determine the boulders mobilization threshold: (1) subaerial boulders, which are fully exposed above mean sea level (MSL); (2) submerged boulders, which remain consistently underwater; and (3) Joint-Bounded Boulders (JBB), whose size and shape are primarily controlled by pre-existing joints, fractures, or planes of weakness in the bedrock, rather than by surface processes such as abrasion or weathering. JBBs are typically located at or near the transition zone between the subaerial and submerged environments. The modified^12^ equations, detailed in Table 3, are applied to each boulder in terms of tsunami wave height \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} and sea wave heights \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} and of velocity u. These thresholds provide the wave height at the impact point able to mobilize the boulder. Figure 2A–C show the resulting values of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} for each condition, along with boxplots for each threshold heights and sector. Results are also summarized in Tables 4 and 5. Notable differences between the western sector (WS) and the eastern sector (ES) are evident for both subaerial and submerged conditions, with the WS consistently exhibiting lower mobilization thresholds for both \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} across all conditions.

For subaerial conditions, the mean \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} threshold values are 0.40±0.16 m for the WS and 0.79±0.31 m for the ES, while the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} thresholds are 1.62±0.64 m and 3.15±1.25 m , respectively. In submerged conditions, tsunami and wave height thresholds increase for both sectors, with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} values of 0.55±0.64 m (WS) and 1.87±1.01 m (ES) and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} values of 2.20±1.33 m (WS) and 7.49±4.02 m (ES). Under JBB conditions, both tsunami and wave height thresholds increase for both sectors (WS: 2.19±0.93 m and 8.75±3.74 m; ES: 3.39±1.39 m and 13.58±5.56 m), remaining consistent the trends between sectors. These variations can be directly attributed to differences in boulder dimensions between the two sectors.

Regarding velocity thresholds for boulder mobilization type (Fig. 2d) in subaerial/submerged conditions, as expected, the sliding velocity threshold ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{sl}$$\end{document} ) is smaller than the rolling ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{rol}$$\end{document} ) and saltation ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_{slt}$$\end{document} ) thresholds. For the WS, the thresholds range from 3.17±0.64 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} for sliding, 4.31±1.25 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} for rolling, and 6.80±2.12 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} for saltation. In the ES, the respective values are 4.33±0.77 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} , 6.19±1.63 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} and 8.88±2.26 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} . Similar to wave height thresholds, velocity thresholds follow a comparable pattern across mobilization types, with lower velocities observed in the western sector (WS) compared to the eastern sector (ES). The only exception is for saltation, where velocities are relatively consistent across both sectors, although slightly higher in the ES. As expected, JBB velocity thresholds, which represent the saltation/lifting mechanisms of mobilization, are higher in both sectors. These thresholds are 7.77±1.62 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} for the western sector (WS) and 9.40±1.92 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathrm {ms^{-1}}$$\end{document} for the eastern sector (ES).

Modelling tsunami waves and boulder transport capacity

Morphometric analysis and analytical equations provide the tsunami and wind waves thresholds in terms of height and velocity at the initial position of the blocks to be mobilized. However, wave propagation processes are nonlinear and complex, involving various mechanisms such as shoaling, breaking, and energy dissipation. Therefore, determining which tsunami wave conditions may impact the cliff edge with sufficient height and velocity requires the use of advanced modeling tools.Fig. 3. Maximum wave height H (A) and total velocity u (B) at the impact point of each boulder for all the tsunami waves simulations. Results are grouped in colours by offshore direction \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _0$$\end{document} , the brighter the colour the higher the period T. Background lines stand for wave height (A) and velocity (B) thresholds. Total velocity magnitude u is computed considering the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_x$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_y$$\end{document} components of the velocity field.

To this end, and to assess the transport capacity of cliff-top boulders, the propagation of different tsunami waves, with varying heights, periods, and directions, has been simulated from offshore depths towards the cliff using a Fully Nonlinear Boussinesq Wave model^30^, totaling 360 simulations of tsunami waves. Further details of the numerical model as well as model setup and assumptions can be found in Section Methods.

For all simulations, the point of contact between the incoming wave and the position of each boulder is extracted in terms of H, the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_x$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_y$$\end{document} components, and the corresponding velocity magnitude u is calculated. These variables have been compared with the wave height and velocity thresholds provided in the previous section.

Figure 3 presents the maximum water column height (H) and the total flow velocity (u) at the boulder impact location, as obtained from the numerical simulations. Each point along the vertical axis corresponds to one of the initial tsunami wave conditions tested in the numerical simulations. Note that if the water does not reach the impact position, no result is shown. The results are color-coded by approximation angles, with color intensity representing the wave period. In 54.4% of the simulations, the water never reached the clifftop at any boulder position. Of the simulations reaching the cliff-top simulations, 100% reach ES and the 4% reach the WS. Regarding the direction of incoming waves reaching the cliff-top, 12.6% correspond to southern incoming waves, progressively increasing to 34.75% for waves from the southwest (SW) direction. Notably, all incoming waves reaching the western sector (WS) originate from the south (S) and south-southwest (SSW) directions. In contrast, wave directions reaching the eastern sector (ES) range from S to SW.

Regardless of the case, all simulations where water reached the clifftop surpassed the three wave height thresholds considered for each boulder (Fig. 3A). However, when analyzing the results in terms of total velocity at each boulder position, it is evident that, unlike the wave height thresholds, not all simulations exceeded the mobilization velocity threshold (Fig. 3B). In the eastern sector, among all simulations where an impact has occurred, 70% exceed the sliding threshold, 60% surpass the rolling threshold, 47% exceed the saltation threshold, and 42% surpass the JBB lifting threshold. In the western sector, among the simulations that result in an impact, 99% exceed the sliding threshold, 86% surpass the rolling threshold, 47% exceed the saltation threshold, and 32% surpass the JBB lifting threshold.

To ease the interpretation and analysis of model results, the mean velocity per sector, modeled at the position of the boulders for the different values of H, T, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta$$\end{document} considered, has been computed (Fig. 4). Contour lines represent the different thresholds associated with mobilization types (sliding, rotation, saltation and JBB lifting). Model results indicate that boulder mobilization at the cliff-top in the WS sector is restricted to particularly high and long waves. Under these conditions, the primary transport mechanisms are associated with subaerial and submerged scenarios, most likely through sliding and rolling transport. Notably, such conditions in the WS sector are only feasible with southward incoming tsunami waves. In contrast, the ES sector exhibits a broader range of offshore conditions capable of mobilizing boulders, regardless of the specific transport type or wave direction. Notably, incoming wave angles around \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$200^0$$\end{document} present the most likely conditions for boulder displacement across all scenarios.Fig. 4. Mean velocity at the boulders position at the impact instant respect to all the offshore wave heights and periods and directions for both sectors. Black lines stand for: sliding velocity threshold, point lines; rolling, dashed lines; saltation, dot-dashed line; JBB saltation/lifting, solid line.

Formation of Es Bancals boulder deposits

Morphometric thresholds and modelling results point out that tsunami forcing is a plausible mechanism of formation of these deposits. The morphometric analysis reveals that smaller blocks with lower mobilization thresholds are concentrated in the WS at higher elevations, whereas larger blocks, transported further inland, are found in the eastern sector, which sits at a lower elevation. The numerical results for impact cases, considering wave height and velocity fields, align with these lateral variations along Els Bancals. This increased potential for transport is attributed to the lower elevation of the ES sector, which results in higher wave velocities, further facilitating the transport of boulders.

Regarding block mobilization conditions, these are limited to waves with offshore parameters of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t>10$$\end{document} m and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$T>800$$\end{document} s, falling within the tsunami spectrum. Under these conditions, various block mobilization mechanisms can occur across different sectors, ranging from lower-energy processes (sliding) to more energetic ones (JBB lifting). The most likely mobilization mechanism is sliding, consistent with flume experiments on tsunami boulder transport^31^. Notably, a wider range of mobilization mechanisms is plausible in the eastern sector (ES), with probabilities ranging from 70 to 49%, covering both low- and high-energy processes. In contrast, in the western sector (WS), mobilization is primarily restricted to sliding and rolling, with minimal occurrence of saltation or JBB lifting. It is important to highlight the difference in results when considering wave heights as a threshold versus using the velocity field. In this regard, based on the numerical results, it is evident that wave height thresholds must be treated with caution, as they apply to the wave height at the exact location of boulder displacement, rather than offshore wave conditions. In coastal areas, dissipation and breaking processes play a crucial role, often resulting in a significant reduction in wave height. Furthermore, due to the dynamics of wave motion and the unique morphology of the study area, exceeding the wave height threshold at the impact point does not necessarily imply that the velocity thresholds have also been surpassed. The velocity threshold is then the key variable to determine whether the incoming tsunami waves can transport boulders. Numerical modeling provides valuable insights into wave propagation and the associated velocity fields, underscoring the importance of a multidisciplinary approach in understanding the formation of these deposits.

If we consider the tsunami directions most likely responsible for generating the tsunamite deposits at Els Bancals, these align with previous studies on tsunami generation in the southern Western Mediterranean, at the coast of Algeria (^27,32,33^, among others) as well as with observations of tsunamite deposits along the Western Mediterranean coast (^18,19,34,35^, among others). Based on the results, directions within the SSW-S range could potentially generate these deposits in the eastern sector, while the western sector is limited to SSW directions. Wave propagation processes, including refraction and shoaling across the platform, play a crucial role in facilitating boulder transport for this specific tsunami wave direction.

While the current data do not allow us to determine whether the tsunamite formations in both sectors occurred simultaneously or were caused by the same tsunami event, numerical results suggest that S-component events could potentially generate both deposits. This direction aligns with events associated with the North Algerian seismic sources identified by^32^, particularly the S-3 source. However, events with SSW components are also plausible, corresponding to the S-1 and S-2 seismic sources from North Algeria.

Since the initial position of the blocks before mobilization is unknown, verifying the model results at that position is not possible. However, the current position of the boulders is considered a reliable indicator, as wave heights and velocities would have been greater at their original offshore location compared to their final transported position. In addition, although the precise source of the boulders along the cliff cannot be identified due to the absence of clear erosional or depositional features, their origin is considered local. The c-axis dimensions of the boulders are of the same order of magnitude as the calcoarenite layers underlying them, strongly suggesting that these strata are their source. In addition, the presence of well-developed marine terraces and benchs in the area further supports a local origin for the boulders.

Numerical modelling of tsunami waves in shallow waters

This study focuses more on the processes of tsunami wave propagation rather than tsunami wave generation. While there are numerical studies that address both the origin and propagation of tsunami waves^32,36^, the resolution of those models may overlook critical dissipation and shoaling processes (transition from kinetic to potential energy), which are essential for assessing tsunami impacts in very shallow areas. By forcing the high-resolution model with long waves that cover the tsunami wave spectrum^37^, the aim is to simulate tsunami waves that have already propagated from their source, having experienced significant shoaling effects as they transition from deep waters to the continental shelf.

Simulations successfully capture the phases of tsunami dynamics in shallow waters. Figure 5 illustrates the different stages of wave propagation for one of the cases considered: drawback flow, approaching wave, impact, flooding and ebb-flood. During the drawback state, sea-surface elevation decreases due to negative (off-shore directed) fluxes in front of the tsunami wave (water retirement functioning similarly to set-down of waves), which results in a lifting of tsunami wave height (Fig. 5A). This process facilitates the transition from vertical to horizontal flow^20^. After this transition, the tsunami wave propagates fast towards the shore, with high onshore directed velocities (Fig. 5B). Wave shoaling results in the amplification of the wave height on approaching the cliff and moving in very shallow waters. Unlike the run-up of a tsunami on beaches or over a constant slope^38,39^, in this case, the tsunami wave climbs the cliff^40,41^ reaching the impact point (boulder position), with significant onshore-directed velocities (Fig. 5C). Once the wave surpasses the cliff-top, it continues propagating inland as a progressive bore^42^, with lower but still considerable velocities (Fig. 5D). Finally, the tsunami retreat, or ebb-flood, occurs with offshore-directed velocities (Fig. 5E). There are numerous observational, experimental, analytical, and numerical studies on tsunami transformation in coastal areas, particularly for constant sloping coasts^43–45^. Special attention has been given to research focusing on tsunami amplification factors caused by shoaling processes (^46–48^, among others). These studies emphasize the importance of tsunami amplification in shallow coastal areas and its significant impact on the total tsunami run-up, which is crucial for accurate tsunami hazard assessment and coastal planning. These works estimate amplification factors ranging from 3 to 6 in regions with gentle slopes. However, tsunami wave amplification due to cliffs has been explored in only a few studies, primarily based on wave tank experiments and numerical modeling of small-scale experiments^40,41^.

Following the definition of amplification factor^40^ as \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H=H_m/H_0$$\end{document} and setting \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_m$$\end{document} at the impact point, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_0$$\end{document} as the off-shore tsunami wave height, modelling results show a mean amplification factor of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2.28\pm 0.39$$\end{document} for all the simulations and boulders position. This factor decreases with increasing wave height, ranging from 3.08 to 1.86 from \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_0=10$$\end{document} to 20 m. This trend is consistent with existing laboratory experiments and numerical studies on tsunami wave interactions with cliffs, which show that smaller incoming waves exhibit higher amplification factors^40,41^. Noteworthy, this consideration should be taken with caution, since these works are based on small-scale experiments with forcing wave heights two orders of magnitude smaller than the current analysis. Clearly, further data is needed to validate these results.Fig. 5. Tsunami waves propagation phases: (a) drawback, (b) wave approach, (c) impact, (d) flooding and (e) ebb-flood. Upper panels show the tsunami wave bidimensional approach at Els Bancals. Red dots stand for the positions of the boulders. Lower panels show the cross-shore profile of the ES (black line), the sea-surface elevation (blue line, left axis) and the cross-shore velocity component \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$u_x$$\end{document} (red line, left axis). Black squares stand for the mean position of the boulders.

Tsunami waves or storm waves?

Cliff-top boulders, found in elevated coastal areas, present questions about their transport mechanisms, with both tsunamis and storm waves considered as potential forces^19,43^. Tsunami waves are characterized by long wavelengths and slow, sustained energy over vast distances, while storm waves are shorter, more frequent, and primarily driven by wind, resulting in distinct transport capabilities and coastal impacts.

To assess whether extreme storm wave conditions can reach the velocity thresholds required for different transport mechanisms at the study site, irregular waves corresponding to these extreme conditions are simulated. The offshore storm wave conditions are modeled using a wave spectrum with wave heights ranging from 5 to 10 m, and a peak period of 10 s (minimum frequency \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f_{min} = 0.07 \mathrm {s^{-1}}$$\end{document} and maximum frequency \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f_{max} = 0.125 \mathrm {s^{-1}}$$\end{document} ), representing the most energetic scenarios in the area (retrieved from Puertos del Estado). To prevent wave dissipation across the computational domain -since tsunami waves are forced in the offshore boundary at 6 km from cliff toe- the storm waves are initiated at a depth of 20 m, approximately 1000 m from the cliff base. Note that \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} threshold from the applied parameterizations to these boulders result in storm wave height values never recorded in 60 years of records.

None of the simulations involving short wind waves reached the position of the blocks, except for sporadic time steps where waves reached the cliff-top, capturing splash processes from their impact on the cliff. Figure 6 illustrates the impact time-steps of a 10 m tsunami wave with a duration of 1800 s and a 10 m storm wave with a period of 10 s, highlighting the differences in wave propagation and impact processes. These results reveal that the formation of Es Bancals boulder deposits is more plausible linked to tsunami waves rather than wind waves.

As aforementioned, most works related to tsunami boulders study formation mechanism of these deposits on slope coast. In addition, most approaches are based on morphometric data and analytic solutions of wave propagation processes. However, wave propagation is a high non-linear process. Complex bathymetries, wave and currents interaction, wave refraction and reflection and positive coupling should be addressed with numerical models that are able to capture these processes. In this sense, to discern whether storm waves or tsunami waves are the most plausible mechanism for the formation of such deposits, requires a transdisciplinary approach, such as the one employed in this study, to enhance our understanding of these processes.Fig. 6. Sea surface elevation at Els Bancals on the Impact time-step for: (A) a 10 m tsunami wave (T=1800 s) and (B) a 10 m storm wave (T=10 s). Red dots stand for the position of the boulders.

Conclusions

This study explores the mechanisms behind the formation of cliff-top boulder deposits by combining detailed geomorphometric measurements with numerical modeling, focusing on boulder deposits along the southern coast of Mallorca Island (Balearic Islands, Western Mediterranean).

Numerical results provide strong evidence that tsunami waves are the primary mechanism responsible for the transport and deposition of cliff-top deposits at the study area. In addition, the most probable tsunami directions responsible for these deposits align with seismic sources found in the southern Western Mediterranean (^17,27,32,33^, among others).

The combined analysis of morphometric data and numerical simulations of tsunami wave propagation at impact points, considering wave height and velocity fields, reveals that impact velocities are the key parameter in determining whether boulders can be transported by the flow and the type of mobilization they undergo. Although analytical solutions for wave propagation can approximate wave height and velocity at impact points^19,43^, they do not account for critical processes such as wave interaction with complex bathymetries, propagation, refraction, and other two-dimensional effects. These processes are essential for capturing lateral variations within a given area.

The current approach prioritizes tsunami wave propagation and transformation over generation, highlighting the importance of high-resolution simulations in capturing key wave shoaling and dissipation effects. By simulating tsunami waves that have already traveled from their source, the study effectively represents critical processes such as drawback flow, wave amplification, impact, and flooding–essential for understanding tsunami behavior in shallow coastal areas.

Furthermore, the numerical results align with previous research on tsunami impacts in cliffed environments, reinforcing our understanding of wave amplification and boulder transport in steep coastal settings. However, such studies remain limited, underscoring the need for further research. A comprehensive review of similar deposits in high cliffs is necessary, as some may have been misclassified as wave-driven rather than tsunami-driven, or vice versa.

The high-resolution numerical modeling approach used in this study is particularly relevant for tsunami risk assessment, as it captures key wave transformation processes, including shoaling, dissipation, and impact dynamics in both sloped and cliffed coastal environments. By accurately simulating tsunami wave propagation and its interaction with complex coastal topographies, this methodology provides valuable insights into potential inundation zones, wave amplification effects, and boulder transport mechanisms. Such detailed analyses enhance hazard mapping and coastal planning efforts, enabling more precise risk evaluations and the development of targeted mitigation strategies to protect vulnerable coastal communities from future tsunami events.

Methods

Morphometric Analysis and mobilization thresholds

Morphometric analysis has been conducted by identifying the largest boulders along the platform at the study site. For each boulder, the maximum axis (a), the minimum axis (b), and the block thickness (c), representing the largest horizontal dimension, the largest vertical dimension, and the width of the boulder, respectively, have been measured, as shown in Figure 1C, together with their orientation \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\rho _{bl}$$\end{document} , bed-slope angle \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{sl}$$\end{document} , height above sea level, and distance from the cliff edge. Orientation has been determined by measuring the azimuth of the long axis relative to the true north. The morphometric values of each boulder have been calibrated using the triangulation method^49^. Morphometric data for both sectors is detailed in Tables 1 and 2.Table 1. Morphometric data of western sector Boulders at Es Bancals.Boulder numbera (m)b (m)c (m)Direction \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\rho _{bl}$$\end{document} (deg.) \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{sl}$$\end{document} (rad.)Height (m)Distance to the sea (m)Boulder volume ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$m^{3}$$\end{document} )11.100.700.350.000.0036.3027.200.2721.601.370.45227.000.3636.5028.000.9931.070.400.350.000.0037.5039.000.1540.620.590.32202.001.4137.2039.300.1251.170.520.460.000.0037.8045.400.2860.800.520.21241.000.5637.5044.200.0971.600.800.25148.000.4237.7047.400.328a0.720.680.21286.000.7335.2034.800.108b1.230.530.37239.000.6335.6036.300.248c0.670.560.15224.000.3535.8037.100.0691.301.000.65199.000.2134.0029.800.85Average1.080.700.34220.750.5836.4637.140.31Std.0.350.280.1439.940.371.226.880.31Table 2Morphometric data of eastern sector boulders at Es Bancals.Boulder numbera (m)b (m)c (m)Direction \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\rho _{bl}$$\end{document} (deg.) \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{sl}$$\end{document} (rad.)Height (m)Distance to the sea (m)Boulder volume ( \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$m^{3}$$\end{document} )112.001.050.650.00.0024.5023.101.36121.301.201.100.000.0024.8023.901.72131.030.750.560.000.0024.8023.400.43141.200.800.43246.000.5425.3026.200.41152.001.200.24248.000.3824.8026.900.58162.301.800.300.000.0025.7034.801.24172.101.200.40233.000.2125.9038.801.01181.801.200.20249.000.2626.2043.600.43191.900.800.500.000.0023.7061.200.76200.501.200.700.000.0019.4051.900.42211.800.800.40243.000.6817.6038.300.58221.801.100.80230.000.3719.1052.801.58231.300.700.30234.000.6520.5062.600.27241.301.100.800.000.0018.4065.001.14251.400.800.700.000.0020.3070.400.78261.901.900.300.000.0018.4065.801.08272.401.600.80200.000.2818.7068.203.07283.002.200.800.000.0017.0061.505.28292.001.700.60189.000.5218.8069.802.04301.801.600.80230.000.1718.3072.602.30313.101.700.60235.000.6819.1073.903.16323.202.400.50243.000.5215.3054.703.84333.402.700.600.000.0014.9051.705.51342.701.900.40219.000.3516.2063.002.05352.101.300.30180.001.0118.6042.300.82361.300.900.40204.000.9318.6043.400.47373.101.400.40210.000.2618.5043.001.74384.002.100.80218.000.2117.4035.906.72391.201.000.20194.000.3518.1040.400.24401.100.700.30185.000.5618.3043.700.23Average2.001.360.53220.530.4720.2449.091.71Std.0.820.540.2322.810.243.4516.111.69

Mobilization thresholds considering the subaerial, submerged and joint bounded block (JBB) conditions in terms of wave height and velocity have been computed by considering the modified^12^ equations^13,14^ as summarized in Table 3.Table 3. Mobilization thresholds equations for subaerial, submerged and joint bounded block conditions, where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\rho _b$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\rho _w$$\end{document} stand for boulder and water density and set to 2.6 and 1.025 g cm \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$^{-3}$$\end{document} , respectively, V is the corrected boulder volume, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mu$$\end{document} is the static friction coefficient (set to 0.65), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{sl}$$\end{document} is the bed slope angle of the pre-transport setting, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C_d$$\end{document} is the drag coefficient (set to 1.95,^18^), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C_l$$\end{document} is the lift coefficient (set to 0.178,^18^), q is the empirical aspect ratio of each boulder set to 0.73^18^ and g is the acceleration due to gravity.SubaerialSubmergedJoint bounded blockTsunami wave height \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t$$\end{document} (m) \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t = 0.5 \mu \dfrac{\rho _b V C_d a c q}{\rho _w}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t = 0.25 \left( \dfrac{\rho _b - \rho _w}{\rho _w}\right) ^2 a C_d \left( \dfrac{ac}{b^2}\right) + C_l$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_t = \dfrac{(\rho _b - \rho _w) V (\cos \theta _{sl} + \mu \sin \theta _{sl})^2}{\rho _w C_l a b q}$$\end{document} Wind waves height \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s$$\end{document} (m) \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s = 2 \mu \dfrac{\rho _b V C_d a c q}{\rho _w}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s = \left( \dfrac{\rho _b - \rho _w}{\rho _w}\right) ^2 a C_d \left( \dfrac{ac}{b^2}\right) + C_l$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_s = \dfrac{(\rho _b - \rho _w) V (\cos \theta _{sl} + \mu \sin \theta _{sl})^{0.5}}{\rho _w C_l a b q}$$\end{document} Velocity thresholds u (m/s) \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {Sliding: } u^2 = 2 \dfrac{(\rho _b/\rho _w - 1) g c (\cos \theta _{sl} + \mu _s \sin \theta _{sl})}{C_d (c/b) + \mu _s C_l}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {Rolling: } u^2 = 2 \dfrac{(\rho _b/\rho _w - 1) g c (\cos \theta _{sl} + \mu _s \sin \theta _{sl})}{C_d (c^2/b^2) + \mu _s C_l}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {Saltation: } u^2 = 2 \dfrac{(\rho _b/\rho _w - 1) g c \cos \theta _{sl}}{C_l}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text {Saltation/lifting: } u^2 = 2 \dfrac{(\rho _b/\rho _w - 1) g c (\cos \theta _{sl} + \mu _s \sin \theta _{sl})}{C_l}$$\end{document}

Numerical modeling

The FUNWAVE-TVD (Fully Nonlinear Boussinesq Wave model with Total Variation Diminishing scheme), developed by^30^, is a numerical model designed for the comprehensive simulation of fully nonlinear and dispersive water wave phenomena. The model accurately represents a range of wave conditions, including interactions with complex coastal bathymetry, shoaling, and breaking. Widely utilized in coastal engineering and research, FUNWAVE-TVD is particularly effective in studying nearshore processes, providing insights into coastal morphodynamics, sediment transport, and wave-induced changes in coastal areas. Notably, its utilization extends to tsunami research, where its robust capabilities contribute to understanding and modeling the complexities of tsunami wave propagation and inundation dynamics. As developed in^50^, the computational domain embraces 6075x6000 m with an equispaced grid of 10x10 m. Topobathymetry has been generated by considering the high-resolution 2.5x2.5m LIDAR digital terrain model (DTM) topography^51^ and the bathymetry from the Balearic Islands Territorial Information System.

Since we focus on the high-resolution solution of the tsunami waves propagation towards the coast, accounting for shallow water conditions and physical processes (i.e. dissipation, shoaling, refraction, diffraction, reflection and wave breaking), the model is forced at the off-shore boundary of the computational domain (40 m depth), following preliminary results shown in^50^, considering different off-shore wave conditions. A total of 360 simulations were conducted, covering off-shore wave heights \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_0$$\end{document} from 5 to 20 m, periods T from 300 to 1800 s, and directions \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _0$$\end{document} ranging from \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$180^{\circ }$$\end{document} (S) to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$225^{\circ }$$\end{document} (SW), matching the direction of potential tsunamis forcing origin^32,36^, an assuming that the MSL is the current reference position. Each simulation represents a single offshore wave condition defined by the tuple [ \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H_0$$\end{document} , T, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _0$$\end{document} ] imposed at the wavemaker position. Periodic lateral boundary conditions are considered. Note that, to avoid numerical artefacts or instabilities in the results, the lateral extent of the computational domain is sufficiently long on both sides of the area of interest.Table 4. Mobilization velocity thresholds for western sector boulders at Es Bancals.Boulder numberSubaerialSubmergedSubaerial/submergedJBBHt (m)Hs (m)Ht (m)Hs (m)Sliding u (m/s)Rolling u (m/s)Saltation u (m/s)Ht (m)Hs (m)Salt./lift u (m/s)10.411.620.491.983.114.187.702.078.287.7020.793.171.335.314.576.978.443.1012.429.4330.230.930.170.692.412.567.702.078.287.7040.341.370.371.482.563.342.911.516.046.5850.301.200.220.892.742.918.832.7210.888.8360.301.200.441.772.894.175.491.485.926.5170.461.850.883.523.505.396.221.746.967.068a0.391.580.682.713.224.975.141.465.856.478b0.311.230.281.133.003.537.122.6110.428.648c0.321.300.642.562.874.534.891.034.125.4390.582.320.552.193.974.8210.384.2817.1111.07Average0.401.620.552.203.174.316.802.198.757.77Std.0.160.640.331.330.641.252.120.933.741.62Table 5Mobilization velocity thresholds for eastern sector boulders at Es Bancals.Boulder numberSubaerialSubmergedSubaerial/submergedJBBHt (m)Hs (m)Ht (m)Hs (m)Sliding u (m/s)Rolling u (m/s)Saltation u (m/s)Ht (m)Hs (m)Salt./lift u (m/s)110.612.430.622.483.854.7710.493.8415.3710.49120.692.780.471.894.174.3513.656.5026.0213.65130.431.740.361.453.283.759.743.3113.249.74140.461.850.532.113.644.777.903.0312.129.32150.692.781.867.424.096.616.141.666.646.90161.044.172.9811.914.537.307.131.777.107.13170.692.781.234.914.196.368.142.6310.538.68180.692.782.088.313.946.345.721.345.366.20190.461.850.481.903.364.149.202.9611.839.20200.692.780.592.364.105.2010.894.1416.5610.89210.461.850.582.333.624.877.262.8111.228.96220.642.550.552.214.284.9511.255.5222.0712.57230.411.620.582.313.364.776.372.118.447.78240.642.550.542.153.974.5911.644.7318.9211.64250.461.850.341.363.403.6210.894.1416.5610.89261.104.403.0012.024.627.427.131.777.107.13270.933.711.124.505.026.7511.415.3921.5812.43281.275.102.018.055.418.0411.644.7318.9211.64290.983.941.566.225.187.759.384.2316.9011.00300.933.711.084.344.936.6211.555.1920.7712.20310.983.941.666.655.177.748.894.2116.8310.98321.395.563.4213.665.869.478.563.5214.0810.04331.566.253.6114.445.749.2410.083.5514.1910.08341.104.402.7210.905.168.337.982.7510.998.87350.753.011.787.134.166.675.191.927.677.41360.522.080.702.793.715.196.392.6510.608.71370.813.241.696.754.517.058.092.6810.738.77381.224.861.937.725.598.2111.515.2721.0612.28390.582.321.435.713.726.015.641.375.506.27400.411.620.572.273.374.776.562.128.467.78Average0.793.151.877.494.336.198.883.3913.589.40Std.0.311.251.014.020.771.632.261.395.561.92