Submarine Landslides in the West Continental Slope of the South China Sea and Their Tsunamigenic Potential

The 109 meridian fault is located in the west of the South China Sea (SCS) connecting to the offshore Red River Shear Zone. Seismic data from the central Vietnamese shelf indicates that many submarine landslides were developed along the steep continental slope in this offshore region. Here, we analyze the potential for such landslides to trigger damaging tsunamis based on the local geological background and sedimentary environment. We assess their tsunamigenic potential along the coast of Southern Central Vietnam (SCV). We point out that the evolutionary processes of the 109° meridian fault: striking-subsidence of the adjacent basin, combined with the high sediment input from numerous montane rivers of the hinterland generate conditions that likely favor the development of submarine landslides along the well-defined and steep continental slope near SCV. To estimate the impact of tsunami waves on the SCV coastline, we conducted a pilot study using two numerical models: NHWAVE and FUNWAVE-TVD to model 4 representative landslides with volumes ranging between 1.3 and 14 km3 and water depth of 300–1000 m. The submarine landslides were treated as rigid slump and deformable slide corresponding to two different sedimentary environments. Our results show that the tsunami waves generated by rigid slump can reach up to 20 m height in the landslide source area and ∼5 m when arriving at the closest coastline. Tsunami waves could arrive at the central Vietnam coast within 30 min in eight scenarios. Our initial results also suggest that seafloor topography, i.e., waveguide effects of ocean ridges, shelf resonance and the potential bay resonance cause significant variability in potential run-up. We note that ocean ridges located in the deep basin of the SCS focus the tsunami energy propagating towards the northwest coast of Luzon Island, Philippines where tsunami wave heights of ∼2.3 m wave height are modeled. These findings underscore the importance of tsunami hazard assessments that account for both earthquake generated and earthquake triggered tsunamis. Our work also highlights a continued need to examine tsunami sources in the region as mitigation and preparedness for the socio-economically important and heavily populated coastlines of the SCS is reliant on a detailed understanding of the hazard.


INTRODUCTION
Submarine landslides, in addition to submarine earthquakes, are one of the principal causes of large tsunamis . Since the 1998 Papua New Guinea event which caused more than 2000 fatalities, people have become increasingly aware of the hazard posed by tsunamis with landslide causes . The need for increased awareness was further raised by two recent events in Indonesia in 2018: the 2018 Palu tsunami (Liu et al., 2020;Nakata et al., 2020;Gatter et al., 2021;Schambach et al., 2021) and the 2018 Anak Krakatau event (Hunt et al., 2018;Grilli et al., 2019;Ye et al., 2020;Grilli et al., 2021) both involving landslide mechanisms.
Submarine landslides have been identified at various geological environments worldwide, i.e., passive, convergent and strike-slip margins as well as volcanic regions (Hampton et al., 1996;Gee et al., 1999;McAdoo et al., 2000;Urgeles and Camerlenghi, 2013;Wang W. et al., 2018;Gatter et al., 2021;Tappin and Grilli, 2021). As one of the largest marginal seas in the world, the South China Sea (SCS) contains a variety of continental margins, some of which that provide favorable environments for the development of submarine landslides with various scale (Moscardelli and Wood, 2016). The SCS has wide and passive continental margins in the south and north, a subduction belt in the east, strike-slip margin in the west and numerous isolated carbonate platforms inside the sea basin.
Over the last 15 years several submarine explorations have identified a large number of submarine landslides with various volumes in the SCS [e.g., Gee et al. (2007), Wu et al. (2011), Wang et al. (2014), Li W. et al. (2015), Sun et al. (2017), Sun et al. (2018a), Sun et al. (2018b)]. Some of these landslides are among the largest slides in the world with volumes exceeding 1,000 km 3 , such as the Brunei Slide has estimated volume of 1,200 km 3 (Gee et al., 2007) and the Baiyun complex has~1,035 km 3 (Sun et al., 2008;Sun et al., 2017;Sun et al., 2018a;Wang J. et al., 2018). Besides the giant slides observed in the open continental slopes, the majority of submarine landslides in the SCS are relatively small, typically on the order of a few km 3 , especially those developed in the submarine canyons (He et al., 2014;Chen et al., 2016;Zhou et al., 2019), steep continental slope in the west (Fyhn et al., 2009a;Fyhn et al., 2009b;Hong Nguyen et al., 2012) and carbonate islands (Wang W. et al., 2018;Huang et al., 2020). Previous studies generally focus on interpreting the characteristics of submarine landslides themselves (e.g., morphology, geometric features, evolution processes etc.) through high-resolution bathymetric and 2D/3D seismic data (Gee et al., 2007;Sun et al., 2008;Fyhn et al., 2009a;Huang et al., 2020). Whether those landslides are tsunamigenic and to what extent they can cause hazard are still poorly understood, particularly for submarine landslides with relatively small scales. A few studies investigated the tsunamigenic potential of the two giant landslides in the SCS: the Brunei slide (Tan et al., 2017) and the Baiyun complex Ren et al., 2019;Sun and Leslie, 2020). These studies highlight the potential catastrophic tsunami effect on the densely populated coastal regions (e.g., the southern China) in the SCS . Compared with landslides with giant sizes, landslides with volumes on the order of a few km 3 generally generate shorter wavelengths and its wave amplitudes decay faster, and therefore affecting relatively limited coastal areas (Janin et al., 2019). Nonetheless, when such tsunami events occur in densely populated region, small-scale landslides can also trigger non negligible tsunami disasters. Examples include the 1998 PNG tsunami (Tappin et al., 2008), the 1979 Nice airport tsunami (Assier-Rzadkieaicz et al., 2000). Therefore, quantitative assessments of tsunamigenic potential of small-scale landslides are needed.
In this study, we investigate the potential tsunami hazard posed by submarine landslides in the western continental slope of the SCS. The available seismic data reveals that numerous landslides exist along the 109°meridian fault area (Pham Nang Vu, 2009). Evidences can also be observed in other seismic profiles in Gwang andWatkins (1998), Fyhn et al. (2009b), Tan et al. (2014), Vu et al. (2017). The 109°m eridian fault located in the western SCS belongs to a giant strike-slip fault system: Red River-Ailao fault system which runs through the west part of the South China Sea. The geological evolution of this fault system shaped the rapidly changing topography in east Sunda shelf (Fyhn et al., 2009a;Fyhn et al., 2009b) and formed an abrupt shelf break with the slope more than 10°in the offshore region near Central Vietnam. The well-developed submarine landslides are located over a steep slope (10°) starting at~200 m water depth,~40 km away from city Nha Trang in central Vietnam. Similar landslides today may become potential tsunami sources threatening the nearby region, especially the central Vietnamese coast (Figure 1).
Historical documents and field surveys suggest that Vietnamese coast has been frequently attacked by coastal flooding events (Hong Nguyen et al., 2012). Among all the historical records, five of them located near the meridian fault region ( Figure 1) are thought to be tsunami events (Nguyen et al., 2007). As part of the Red River-Ailao fault system, the 109°m eridian fault is still seismically active. During 1919-2005, 18 earthquakes were recorded along the 109°meridian fault, with the maximum magnitude of Mw 6.1 that located in 109. 0°E, 10.1°N, happened in 1923(Hong Phuong, 2001Hong Nguyen et al., 2012). Seismic activity is one of the potential triggers for mass failures on the slope (e.g., Tappin et al., 2008). The frequency of a tsunami either triggered by underwater earthquake or submarine landslide is commonly at a scale of decades, hundreds of years or even longer (Talling, 2014), therefore, the tsunami hazard is often ignored or overlooked (Giuliani et al., 2019). If a submarine landslide that is either earthquake independent or a combination failure of the meridian fault and a triggered landslide was to occur it would likely cause a tsunami that would affect the Vietnamese coast. To date this has not been quantitatively addressed.
Here, we conduct a pilot series of numerical experiments to investigate the tsunamigenic potential of submarine landslides near central Vietnam. We first briefly summarize the tectonic evolution of western SCS and analyze the potential factors affecting the slope instability. Then based on the constraints from geological information and seismic data, we specify the geometric parameters of four representative landslides. We then presents the modeling procedure which includes a nonhydrostatic 3D model NHWAVE (Tehranirad et al., 2012) for tsunami generation and the fully non-linear dispersive FUNWAVE-TVD model  for tsunami wave propagation in a larger domain. We conclude by describing the potential characteristics of tsunamis generated by the modeled landslides with reference to the implications for tsunami hazard locally in Vietnam and in the broader SCS region.

Geological History of the Western SCS
The geologic structures of the western SCS are mainly related to the large-scale strike-slip systems: The Red River-Ailao fault in the north (Briais et al., 1993) connects the 109°meridian fault offshore of central Vietnam, and extends to the Lupar line in the south (Figure 1). The Red-River Ailao fault originates from the Ailao Mountain, and stops southwest of Hainan Island, likely connecting to the northern tip of the 109°m eridian fault (Fyhn et al., 2009a) which developed along the continental shelf -slope transition zone in eastern Vietnam (Morley, 2002;Clift et al., 2008). The 109°meridian is located on the shelf break, and has a total length~800 km , covering a latitude range between 17°N and 10°N. Three successive deformation phases for the Red-River Ailao fault and the 109°meridian fault have been inferred and include leftlateral movement from~30 to 16 Ma, a transitional phase between 16 and 5.5 Ma when slip velocity decreased and eventually ceased, and right-lateral movement with low rates after 5.5 Ma (Tapponnier et al., 2001;Clift and Sun, 2006;Zhu et al., 2009). The 109°meridian fault has the Zhongjiannan Basin/Phu Khanh Basin located offshore to the east. It was likely formed by intensive spreading during the Middle-Late Eocene (Hongfang and Ling, 2006). From the Late Oligocene to Early Miocene (~28-13.8 Ma B.P), the basin experienced a right-lateral transition and extension rifting stage (Lei et al., 2015). The late Neogene (~10.4 Ma B.P) uplift occurred in the south and central Vietnam, which suggested the seaward tilting of the region onshore and offshore of central Vietnam (Fyhn et al., 2009a;Wang et al., 2013). During the same period, onshore fission track studies reveal the onset of major uplift evidenced offshore by an acceleration in the rate of late Neogene cooling (9-15°C/Ma) (Carter et al., 2000). The uplift and related increase in offshore sediment supply influenced the Zhongjiannan Basin/Phu Khanh Basin, resulting in the high slope gradient continental slope along the western SCS (Hongfang and Ling, 2006;An et al., 2013;Wang et al., 2013).
The elevation difference of the steep continental slope approaches 1,500 m at many places. Although the slope has been subsequently adjusted and modified by sediment deposition, most of the slope gradient still maintains above 10°. Moreover, the east coast of Vietnam is located in a tropical monsoon zone and has numerous montane rivers that drain the 10,350 km 2 range that extends most of the way down the Vietnam and reaches highest of 2,819 m (Souvignet et al., 2014;Tan et al., 2014). Holocene sedimentation rates on the offshore shelf of the Thu Bon River mouth are very high, ranging between 10.1 and 124.5 cm/ka (Schimanski and Stattegger, 2005). The high sedimentation rates on the continental shelf off the coast of Vietnam form fundamental conditions for developing submarine landslides. For example, in the central Vietnam offshore, during the late Miocene (~5.5 Ma), a large-scale submarine landslide with an area of approximately 18,000 km 2 and a maximum thickness of 930 m formed in the deep-water region of the Qiongdongnan basin . The identification of several large-scale submarine landslides in the geological record implies that strike-slip activity on the Red-River Ailao fault and the 109°m eridian fault along with consistent and high sediment supply have an influence on the stability of the continental slope.

Oversteepening
Oversteepening is regarded as a crucial precondition prior to many slope failures. According to the mechanics and slope readjustment model of Ross et al. (1994), if the applied shear stress exceeds the shear strength of marine slope sediments, a submarine landslide will occur. The steep, sediment laden slopes of the study area is the final product of tectonic activity which offset part of depositional profile, leading to erosion and oversteepening that eventually promoted the slope instability (Fyhn et al., 2009b). In addition to tectonic factors, high sediment-deposition rates contribute to the oversteepening process. In the process of slope evolution of the central Vietnam continental slope, the excessive steepness of the basement caused by tectonic activity and the large amount of sediment brought by the central rivers would likely result in significant excess in sediment supply to the continental slope. With these contributions, a large number of submarine landslides were developed on the continental slope of eastern Vietnam ( Figure 2).

Seismic Activity
Earthquakes are identified as the most common triggers of submarine landslides in many regions (Locat and Lee, 2002;Urgeles and Camerlenghi, 2013;Huhn et al., 2019). Previous studies suggest that earthquake recurrences are remarkably correlated to the failure possibilities of submarine landslide (Urgeles and Camerlenghi, 2013). Along the 109°meridian fault, historical documents report 18 earthquakes during 1919-2005 with the maximum magnitude up to M = 6.1 (Hong Nguyen et al., 2012). Statistically, there are 71 earthquakes ranging from M = 3.0 to M = 6.1 since 1715, including volcanogenic and tectonic events (Nguyễn Hồng Phương, 2017). The former can be regarded as a direct trigger of submarine landslides, and the latter can promote the slope instability by deducing the shear strength of sediment (Posamentier and Martinsen, 2011;Huhn et al., 2019). The earthquake swarms are largely distributed near the south of the Tuy Hòa shear zone (~12°N) (Vu et al., 2008). Hong Phuong (2001) use the peak ground acceleration (a max ) to evaluate quantitatively the geodynamic regime of the central Vietnam coast, and find that there are two anomaly zones of a max . One of them is completely superimposed on the seismically active fault segment of the 109°meridian fault system, suggesting that there is a high possibility of slope mass failure triggered by earthquakes along this fault.

Relative Sea Level Change and High Sedimentation Rate
Since the late Pleistocene (~20 ka), the regional sea level offshore Vietnam has changed dramatically, especially during 9-13 ka B.P, rising from −80 m to −20 m (Schimanski and Stattegger, 2005). According to Pham Pham Nang Vu (2009), the studied submarine landslides probably occurred during 10-12 ka B.P between the transition time of the last glacial and the last interglacial when the sea level steadily rose. Some previous studies suggested that landslides occur more frequently during periods of sea level rise and lowstand (Day and Maslin, 2005;Lebreiro et al., 2009;Locat et al., 2009;Georgiopoulou et al., 2010;Smith et al., 2013), although Urgeles and Camerlenghi (2013) pointed out that no strong global correlation of landslide frequency with sea level changes. The possible mechanisms about how the climate-driven sea level change affect the landslide frequency thus still remains controversial (Tappin, 2010;Smith et al., 2013;Urgeles and Camerlenghi, 2013). Factors associated with hydrate dissociation, increased overpressure in sediments and seismic activity are proposed as probable causes of landslide during periods of sea level rise (Smith et al., 2013). In a much shorter time-scale, the sea level fluctuations caused by storm waves affect the slope stability via altering the stress regime at the seafloor, but the direct impact of changing sea level on slope is minor (Weaver and Kuijpers, 1983;Lee et al., 1996). To summarize, the change of sea level in the offshore east Vietnam coast serves more like a precondition than a direct trigger.
Rapid sedimentation is another common cause of submarine landslide through increasing the pore pressure and decreasing the effectives stress in the sediment (Harbitz et al., 2006).
Overpressure is more likely to be found in locations with rapid sedimentation. Schimanski and Stattegger (2005) analyze sediment core from Vietnam shelf and date the sediment samples with the AMS-14 C technique. Some data of the cores, locate in the narrow central shelf, showing high-sediment rates, with a MAR (mass accumulation rates) of 60-120 g/cm 2 ky. These values exceed MAR of other shelves in the SCS at least by a factor of 2-4 (Schimanski and Stattegger, 2005), occurring in early Holocene times, which greatly increase the instability of the slope.

Volcanic Activities
The volcanic activities in the 109°meridian fault were noted to have started in the late Miocene (Fyhn et al., 2009a). Volcanic activity can influence the instability of slopes in two ways. Firstly, volcanic activity beneath East Vietnam leads to an uplift of the coastal topography which was accompanied by activation of rivers (Carter et al., 2000;Clift et al., 2004;Tan et al., 2014;Maselli et al., 2020). The increased elevation strengthens the downcutting of river. With the abundant rainfall brought by the southwest monsoon, the river transport capacity is significantly enhanced, resulting boosted sedimentation rates and increased sediment thickness on East Vietnam shelf. All these factors form the precondition of slope instability. Secondly, the volcanic activities lead to the magmatic intrusions in central Vietnam coast (Clift et al., 2004;Fyhn et al., 2009a). During the process of magma cooling which produced by, the release of hydrothermal fluid and overpressure gas can rise to surface through fractures or passages in the slope, becoming one of the triggers of submarine mass failure (Fyhn et al., 2009a). A series of normal, high angle faults forming graben-like structure in the basement (Tan et al., 2014), which means that the high-angle faults of the geologic structure may become the main passages of gas transport to the shallower layer. The rising gas pass by the finer layer and then accumulates, forming local overpressure, which reduces the effective stress and the shear strength, leading to slope mass failure.

NUMERICAL MODELING OF SUBMARINE LANDSLIDES
To estimate the tsunamigenic potential of the submarine landslides on the steep continental slope near Vietnam, we conducted a pilot study using four representative landslides based on the constraint of the limited information from 2D seismic data ( Figure 2). From north to south, four submarine landslides are parameterized with different geometry, thickness, and kinematics. The submarine landslides are estimated as 5-20 km long, and 5-10 km in width with the maximum thicknesses in the range of 150-400 m. These landslides are generally distributed in water depth of 250-2,500 m between latitudes of 11.9176°and 13.4167°(see the key parameters in Table 1). Due to the limited spatial coverage of seismic data, we simplify the shape of these landslides as ellipse. Each submarine landslide has an elliptical footprint on the slope, centered as its estimated initial location. The thickness distributions of these landslides are assumed to be quasi-Gaussian. Additionally, five virtual gauges are specified near the coastline of major cities to better visualize the tsunami propagation process. The exact locations of these gauges are listed in Table 2.
We used the combination of NHWAVE-FUNWAVE models for tsunami generation and propagation. The landslide materials are treated as rigid slump and deformable slide, corresponding to rigid block and mud flow, respectively. NHWAVE can predict the instantaneous surface elevation and 3D flow field, and provides a numerical solution to the 3D Euler equations for incompressible flow in a terrain and surface-following σcoordinates . It has been validated for highly dispersive landslide tsunami generation by comparing simulated surface elevations with laboratory experimental data (Enet and Grilli, 2007;Shi et al., 2012). The dynamic landslide failure process and the resultant initial tsunami wave were modeled using fully dispersive nonhydrostatic wave model NHWAVE . Thus, by treating landslide as a rigid slump, we expect to provide a conservative (i.e., the worst-case scenario) estimate of tsunami generation and coastal impact.
For rigid slump, as in Enet and Grilli (2007), the landslide geometry is idealized as having a "Quasi-Gaussian" shape of elevation ζ(x, y), and the steepness is controlled by the ε (here ε 0.717). The elliptical footprint of downslope length b, width w, and maximum thickness T are defined as where (ξ, χ) are the local downslope and spanwise horizontal coordinates, k b 2C/b, k w 2C/w, with C acosh(1/ε). Based on these parameters, the landslide volume is given by V s bwT I 2 C 2 ( And f μ sechμ a tan sinhg μ , g μ acosh sechμ ε For the specified ε, we find C 0.8616, I 1 0.4804, I 2 0.5672 and V s 0.3508bwT. The kinematic motion of rigid slump is prescribed based on a set of semiempirical formulas given in Enet and Grilli (2007), Grilli et al. (2002). E.g. the terminal speed was calculated according to Grilli et al. (2002), U t gb π(γ−1) 2C d sin θ , where θ is slope angle, b is the length of landslide, γ is the specific gravity of the slide material, C d = 1 is the added mass coefficients. For the deformable slide, we specify the kinematic viscosity to be 0.1 and the Manning coefficient for the viscous slide is 0.1. The Cartesian computational grids were used in tsunami simulations with 500 m resolution-based grid from GBECO (General Bathymetric Chart of Oceans), which covers the whole SCS.
In this study, the tsunami waves in the first 30 min after failure were simulated by NHWAVE and serve as initial condition of FUNWAVE-TVD for propagation modeling in larger domain. The surfaces wave height and horizontal velocity of the tsunami wave are then interpolated into the whole SCS gridded data for FUNWAVE-TVD (Tehranirad et al., 2011;Shi et al., 2012;Kirby et al., 2016). FUNWAVE-TVD is developed based on the fully nonlinear Boussinesq equations and is used a "total variation diminishing" scheme to better model wave breaking dissipation and dispersion (Chiocci et al., 2008;Shi et al., 2012). The combination of NHWAVE-FUNWAVE modeling approach has been validated for the Tohoku 2011 tsunami (Tappin et al., 2014) and applied in many events Grilli et al., 2017;Grilli et al., 2019;Li et al., 2019;Schambachand Lauren, 2020;Schambach et al., 2021). According to the seismic profiles (Pham Nang Vu, 2009), we identified four representative submarine landslides (Figure 2). Due to a lack of accurately mapped bathymetric data, the initial slide is modeled as a sediment mound of quasi-Gaussian crosssections, with maximum thickness T, and an elliptical footprint of down slope length b and width w. Table 1 gives the geometric parameters and initial location estimated for each modeled landslide base on the seismic profiles. Figure 3 shows the tsunami surface elevation generated by both deformable and rigid material of LS01-LS04 at t = 10 min after the release of each slide. The spatial distributions of tsunami surface elevation are similar in shape but differ significantly in magnitudes (Figure 3). The initial waves generated by rigid slumps of all scenarios are much higher than that generated by their deformable counterparts (Figure 3). With all the other key parameters (initial water depth and slope gradient) being similar, slides with larger volumes (e.g., LS02 and LS04) generate larger initial tsunami wave heights than those with smaller volumes (LS01 and LS03). Figure 4 shows the extracted surface elevations at t = 10 min along four transects which generally follow the central axis of each slide. The heights of initial surface elevation generated by rigid slump are 1.5-3 times larger than that generated by deformable slide (Figure 4). LS02 that has the largest volume among all the slides triggers the highest tsunami waves~15 m while the initial wave generated by the smallest slide (LS03) is~5 m (Table 1; Figure 4). The initial waves generated by rigid slump and deformable slide consist similar waveform pattern but with different level of irregularity ( Figure 4). Taking LS02 as an example ( Figure 4B), after being released, in the sliding direction, the wavefront is a large crest followed by a large trough and a series of wave train with some of them even larger than the leading wave, while in the opposite direction of the slide motion, large leading depression wave followed by a crest with similar amplitude propagates backward toward the coast (Figure 4). Figure 5 show the snapshots of surface elevations generated by LS02 and LS04 at 15 and 25 min. The spatial distributions of surface elevation generated by all slides demonstrate as radical circles/ripples, with longer wave lengths in the deep sea and shorter landward wave lengths in the continental shelf ( Figure 5). When the tsunami waves propagate toward the coastlines over continental shelf, the velocity of the landward wave trains slows down due to the shoaling effect and forms a shell-shape wave pattern ( Figure 5). For all cases, the propagation directions of leading waves are similar which largely direct towards the east. Figure 6 shows the snapshots of tsunami surface elevations generated by the rigid slumps of LS02 and LS04. Surface elevations generated by the deformable slides of LS02 and LS04 are shown in Supplementary Figure S2. In all the scenarios, tsunami waves are largely featured with a leading crest followed by a wave trough (Figure 6). Only a narrow portion of tsunami waves propagating opposite to the slide direction are led by a depression wave followed by a positive crest wave (Figure 6, also refer to Figure 5). The orientations of wave fronts are determined by the seafloor bathymetry and significantly affected by the presence of seafloor relief, including Xisha, Zhongsha, Nansha islands and the seamount chains. The tsunami waves slow down when they encounter the deep-sea islands/seamounts and bend around the islands before meeting behind them due to the strong wave reflection and diffraction effect ( Figures 6B,E). The wave amplitudes increase when tsunami waves approach the continental slope and then propagate into the continental shelf due to the shoaling effect.

Characteristics of Tsunami Propagation
We investigate the wave evolution characteristics by extracting wave profiles at t = 5, 15, 25, 35, 40, 45, 50, 55, 60 min time snapshots along transects EE′ and FF′ for rigid slump LS02 and LS04 ( Figures 7A,B). In both scenarios, the tsunami waves evolve with similar pattern. The initial tsunami wave crests and troughs could reach~30 m/−40 m for LS02 and~15 m/−30 m for LS04 in the source region right above the slumps, then the waves split into two groups and depart in two opposite directions. The one propagating toward the coast experiences shortened wavelength, first increased and then decreased wave amplitude due to the combined effect of shoaling and bottom dissipation. The waves propagating towards the deep sea experiences the opposite effect with stretched wavelength and decreased wave amplitude. We observe that the wave amplitude declines significantly in both propagation directions. For LS02 ( Figure 7A), the amplitude of landward wave crest drops from 30 m in the source region to~3 m near the coastline ( Figure 7A). From longitude 109.8°E toward east where water depths of both profiles are over 1,000 m, the amplitude of seaward wave trough decreases from −40 m to −1.5 m in the deep sea at 60 min. Similarly, in the case of LS04, the amplitude of wave crest drops from 15 to 1.0 m in the landward direction with the wave trough decaying from −40 m to −0.7 m in the seaward direction ( Figure 7B).

Nearshore Impacts in Central Vietnam
The central Vietnamese coast has extremely complex bathymetry and morphology featured with many chains of islands and many half-enclosed bays. When tsunami waves propagate towards such an irregular coastline, interference including wave reflection, diffraction, refraction and resonance are expected. Figure 8 shows the snapshots of tsunami surface elevation generated by rigid slumps of LS01-LS04. Due to the difference of geographical locations, LS01 and LS02 generate tsunami waves propagating dominantly towards southwest during the shoaling stage, while the main tsunami direction generated by LS03 and LS04 is toward northwest. Consequently, tsunami waves generated by LS01 and LS02 reach cities in the north, e.g., Qui Nhon, Tuy Hoa earlier than cities in the south, e.g., Nha Trang and Cam Ranh, while tsunami waves generated by LS03 and LS04 arrive southern cities first. Significant wave refraction occurs at the coastal island chains along the complex coastline, focusing wave energy and amplify wave height behind the lee side of the islands, i.e., Hon Lon, Hon Tre and Cu lao Xuanh (Figure 8). For coastlines with halfenclosed bays, particular attention should be given to tsunami amplifications due to resonance effects (Rabinovich, 2010;Bellotti et al., 2012). When the tsunami period matches the natural   Frontiers in Earth Science | www.frontiersin.org April 2022 | Volume 10 | Article 843173 9 oscillation modes of these bays, tsunami wave heights can be enhanced, resulting exacerbated coastal damage (Lepelletier, 1981). The tsunami hazard related to resonance effect in this region deserves more detailed future study. To predict the onshore tsunami wave heights, numerical simulations based on high-resolution bathymetric data and topographic data is required.

Maximum Surface Elevation
The maximum wave height distributions of all the eight scenarios are shown in Figure 10. We also overlay the tsunami travel time contours in every 30 min on the maps. The studied landslides, especially the deformable slides ( Figures 10A-D), generate very localized tsunami impact with most of the wave energy concentrating in the central Vietnam. With other factors (initial water depth, slide material etc.) being similar, the total volume of slide is the key factor determining the tsunami generation capacity. We can see relative smaller slides (LS01 and LS03) produce much smaller spatial distributions and wave amplitudes. Similarly, besides the geological location, the only difference between LS02 and LS04 is that the LS04 is thinner. As a result, LS04 produces relatively compact and smaller wave heights. (Figure 10).
A combination of directivity from the slide orientation and the influence of bathymetry determine the final distribution of maximum wave height. In the near-source region, the shapes of the maximum surface distribution interestingly resemble "volcanic-like" (Figures 10E-H) or "flank collapsed volcaniclike" shape from a 3D view ( Figures 10A-D). The crater of these "wave volcanos" are the source areas right above the continental slope. The collapsed flanks are formed due to the continental slope reflection which will be explained in detail in the discussion section. In the relatively far-field, seafloor topography determines the propagation direction of tsunami energy. One notable phenomenon, the ocean ridge in the middle of the SCS basin guides tsunami wave energy towards the west coast of Philippines ( Figures 10F,H), resulting relatively high tsunami waves in the coastline of western Luzon.
Using LS02 as an example, Figure 10 illustrates a more quantitative comparison between the maximum surface elevation generated by deformable slide and rigid slump. While the maximum wave heights reach 4-4.6 m in both scenarios along the coastlines near the source region, the majority of remaining coastlines (i.e., southern coast of China) only experience minor or moderate tsunami wave heights (tens of   Figure 9B). We also observe one coastal location in southeast Taiwan has exceptionally high wave amplitude (~1.6 m) compared with its nearby gauges (tens of cm). When closely examining the animation of tsunami propagation in this specific region (Supplementary Video S1), we show that such typically higher amplitude is at least partially related to the unique double-arc structure in the Luzon strait (Yang et al., 1996). The role of mid-ocean ridges in guiding tsunami propagation has been clearly demonstrated in previous tsunami events, e.g., the 1833 event in South Sumatra (Okal and Synolakis, 2008), the 2004 Indian Ocean Tsunami (Titov et al., 2005), the 2006 Kuril tsunami (Kowalik et al., 2008) and the 2011 Tohoku tsunami (Song et al., 2012). Special attention should be given to coastal regions which are located in the terminus of such waveguide structures.

Tsunami Wave Arrival Time
In all simulated scenarios the first tsunami waves arrive at central Vietnam coast within 30 min. Major cities including Quy Nhơn, Tuy Hòa, Nha Trang, Cam Ranh and Phan Rang would experience the first tsunami wave between 30-60 min ( Figure 10). Depending on the landslide locations, the tsunami arrival time is~1-1.5 h for nearby archipelagos inside the SCS, like Xisha and Zhongsha. The west Luzon coast is affected by tsunami waves (~2.5 h) in the relatively far-field where the tsunami directivity effect and relatively deep bathymetry between the studied landslide cause a relatively higher projected tsunami. Elsewhere in the SCS the arrival time is commonly between 2.5-3 h for southeast Hainan, 3.5 h for southwest Taiwan Island and north Borneo Island, 5-6 h for the southern coast of mainland China and more than 10 h for east coast of Malay Peninsula and Singapore ( Figure 10).

Is the Tsunami Hazard Enhanced by the Steep Continental Slope?
The modeling suggests the potential for one particularly interesting phenomenon; that the tsunami waves pile up in front of the sharply bathymetry change at the continental slope near central Vietnam. This leads to a stripe-like wave hump parallel to the edge of continental shelf (see the maximum surface elevation in Supplementary Figure S1). To illustrate the locations of the wave hump, we extract the profiles of maximum wave height and their corresponding bathymetry along seven transects perpendicular to the orientation of the edge of continental shelf ( Figure 11). From north to south, the peak locations of the wave humps are 110 ± 0.1°which correspond to water depth 1,500-2,000 m ( Figure 11). We hypothesize that the linear wave hump feature is likely largely caused by the combined effects of wave refraction and shoaling effect on the continental slope. Due to the steepness of this specific continental slope, a large difference between the propagation speed at the top (~200 m depth) and bottom edge (2000 m depth) leads to rapid shortening of the wavelength and dramatic localized increase in wave height. At the same time, the orientations of tsunami wavefront is adjusted by strong wave refraction so that the wave trace roughly follows the orientations of bathymetry with the largest gradient difference in the deeper water. Figure 12 show time series of tsunami wave amplitude generated by LS02 ( Figure 12) at five selected gauges and their wavelet analysis ( Figure 12, see Table 2). We first compare the tsunami waveforms generated by deformable and rigid slides, which show dramatic difference in wave height and wave period at all stations. In general, we find that the surface elevation generated by rigid slump (~2-3 m) is~2-3 times larger than that generated by deformable slide (~1 m). In most scenarios, gauges closer to the source region experience larger tsunami waves ( Figure 12). One interesting wave feature at TG04 is that a noticeably larger wave train occurs between 1-1.5 h after the first group of energetic tsunami waves. A similar feature is also shown in the wavelet plots where larger wave energy appears again 1.0-1.5 h after the first energetic wave patch. Such delayed energetic waves are commonly associated with edge waves which travel alongshore with much slower speed than normal tsunami wave (Ursell, 1952).

Spectral Analysis
The subplots of wavelet analyses reveal the frequency-time content of the tsunami waveforms recorded by five synthetic gauges ( Figure 12). The wavelet plots clearly show strong tsunami energy at the period band of 5-10 min in all the waveforms generated by both deformable and rigid slides. We attribute these energy patches to the tsunami periods associated with the landslide sources themselves. The energetic patches last much longer in tsunami signals generated by rigid slump (~5 h) than deformable slide (~2 h), indicating the tsunami energy decays much faster in scenarios of deformable slide than those of rigid slump.
The wavelet plot generated by LS02 shows distinct patches of tsunami energy at the period band of 8-15 min which we attribute to the landslide source of this tsunami ( Figure 12). The expected tsunami period of LS02 can be roughly estimated using formulae T 2L gd √ , where T is tsunami period, L is source dimension, g is gravitational acceleration (9.81 m s −1 ) and d is water depth at the source location (Rabinovich, 2010;Heidarzadeh and Satake, 2014). We calculate the tsunami period with source dimensions of 20 km for LS02 and a water depth in the range of 300-1000 m ( Table 2), and the resulting period are 6.7-12.3 min for LS02. There theoretical values are close to those shown in the spectral analysis. We also observe larger bandwidths centered in 20 min, 30 min and 60-80 min in most of the wavelet plots, especially at locations of TG1, TG2, and TG5. The waves with longer periods are very persistent in TG1 and TG5 lasting~240-330 min. Such long wave periods and long-lived waves are often related to shelf resonance which is caused by the wave reflection between the coastline and the edge of continental shelf [e.g., Yamazaki et al. (2011), Melgar et al. (2018, Wang et al. (2021)]. Stronger resonance is expected in this particular shelf due to the sharp continental slope.
Harbor resonance is another factor which could significantly amplify tsunami wave amplitude and induce strong current (Rabinovich, 2010). To better understand whether the tsunami wave could generate resonance in the harbor, we choose four bays (namely Dam Thi Nai, Xuan Dai Bay, NHA Phu Bay and Vinh Cam Ranh from north to south) in the latitude range of 11.5°-14°N (see the locations of these harbors in Figure 8) to roughly calculate their eigen periods following the formula proposed by Rabinovich (2010) as T n 4L (2n+1) gh √ , for mode n = 0,1,2, . . .where L is length of the harbor; h is the average water depth in the harbor, g is gravitational acceleration, 9.81 m/s 2 ; n stands for the number of nodal lines which equals the mode number. When n = 0, we obtain the fundamental period (Rabinovich, 2010) which is the lowest mode, known as the Helmholtz mode which is of particular importance for any given harbor. The formula given above provide an approximately estimated calculation for the period of the Helmholtz and other harbor modes. The normalized periods of various modes have a specific ratio with a simple geometry of the harbor, which into a relationship of 2/(2n+1) on the basis of the period for closed-basin with constant water depth, T n 2L gh √ (Rabinovich, 2010). For rough calculations, we use the most basic pattern here, and the others can refer to (Wilson, 1972) depending on the harbor geometry and water depth.
The lengths of the bays are 9.7, 6.5, 16.7 and 4.22 km (Figure 8). By extracting the bathymetry profiles, we obtain that the average water depths of the entrance are 2.5 m, 4.4 m, 37 m, and 1 m. We take n = 0, 1, 2 for a rough calculation and the results show that the longest natural period of harbor is 130.6 min and the shortest is 11.7 min with the common range of 18-30 min (Supplementary Table S1). Except waveforms recorded at TG1 and TG4 in the scenarios of deformable slide, waveforms recorded by other tide gauges contain wave periods of 18-30 min, indicating that there is the possibility that harbor resonance could be induced. Although the tsunami affected regions are limited in many of the presented scenarios, exacerbated tsunami hazard may occur in certain bays or harbors when the natural resonance periods match the tsunami periods or fundamental period of shelf resonance (Kulikov et al., 1994;Pattiaratchi and Wijeratne, 2009). It would be very interesting to further investigate the role of bay and shelf resonance in tsunami amplification along this specific coast based on more sophisticated approach (Yamazaki et al., 2011;Wang et al., 2021), e.g., using finite element numerical models (Bellotti et al., 2012).
We should note that the scenarios we presented in this study only represent limited number of typical cases in this region. Due to the limited spatial coverage of seismic data, we reconstructed these landslides with simplified shapes. In reality, the geometric parameters of landslides could be much more diverse. To reconstruct the initial spatial distribution of landslides more accurately, more detailed geophysical and geotechnical data, including 2D/3D seismic data, multi-beam bathymetric data and sediment cores should be collected in future. Another limitation of this study is that we treat the landslide material as rigid block and viscous mudflow. The tsunamigenic capacity highly depends on the landslide dynamics (sliding, slump, debris flow and turbidity flow etc.) (Zengaffinen et al., 2020) (Harbitz et al., 2014). Future studies should consider different landslide material as well as rheological behaviors, etc. (Lovholt et al., 2015;Løvholt et al., 2017;Schambach et al., 2018;Kim et al., 2019), e.g., granular flow or multi-phase flow Si et al., 2018;Zhang et al., 2021a;Zhang et al., 2021b).

CONCLUSION
In this study, we summarize the potential for tsunamigenic landslide failures in the west continental slope of SCS along the 109°fault line. By reviewing previous studies about the geological evolution of the western SCS, we conclude that slope physiography is the final product of series of tectonic activities, seismic deformation, relative sea level changes, relatively high sediment rates and localized volcanism that directly or indirectly promote the slope instability through changing the shear stress of sediment on the slope.
Using a coupled approach: NHWAVE and FUNWAVE-TVD, we simulate the tsunami generation and propagation process of four representative landslides. By examining the hydrodynamic characteristics of the tsunami waves, we summarize our key findings below: 1) The tsunami hazard level depends heavily on the landslide materials and their volumes. Among the studied landslides with volumes ranging between 1.3-14 km 3 , the rigid slumps with larger volumes produce larger tsunami impact, e.g., initial surface elevation generated by LS02 reaches 30 m in the offshore region right above the source area and the nearshore tsunami wave heights could reach~5 m in certain coastline. The LS02 with rigid material is also the only scenario which could cause considerable tsunami impact in the relatively far-field. For all the other scenarios, the tsunami energy decays very fast, especially for those with deformable material, resulting very localized tsunami impacts. Large tsunami wave energy only concentrates in the coastline of central Vietnam and the immediately offshore region near the source. The maximum wave heights generated by these landslides differ significantly, ranging between 5 and 30 m. We roughly calculate the areas of ocean surface with wave height larger than 1 m and find the coverages are 192.3-10614.8 km 2 for deformable slides and 7,367.3 km 2 -107327.0 km 2 for rigid slumps. The strong disturbance caused by tsunami waves pose grave threat for ocean engineering in the deep sea.
2) Besides being directly influenced by landslide geometries and their kinematic features, the spatial distribution of tsunami impact is also strongly affected by the seafloor bathymetry along the path of tsunami wave propagation. We observe very interesting hydrodynamic phenomenon associated with the seamount chain, steep continental slope and the double-arc structure in Luzon strait. These seafloor structures either serve as waveguide which facilitate the focused tsunami energy in relatively far-field, or cause strong refraction along the continental slope, forming the stripe-like wave hump parallel to the 109°meridian fault.
3) The wavelet analysis of five virtual gauges in the coastal region disclose that the periods of the landslide tsunami are~5-10 min. We notice that there are tsunami waves with longer periods centered at 20 min, 30 min and 60-80 min in some locations. Such persistent energy patches are commonly associated with shelf resonance. These periods are found to be within the range of the natural resonance periods of some coastal bays (8-96 min). The potential coincidence suggest amplified tsunami waves with much longer oscillation duration may occur in this region. 4) We should note that the tsunami waves in all simulated scenarios reach the coast of central Vietnam within 30 min after the landslide failure. Loss and damage could be particularly serious for such local tsunamis due to the concentrated tsunami energy and limited evacuation time, as illustrated by the 1998 Papua New Guinea event.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding authors.

AUTHOR CONTRIBUTIONS
XP, LL, and HN designed the study. XP conducted the study and wrote the manuscript with contributions from all other coauthors. All authors contributed to the discussion and interpretation of the results.