Journal cover
Journal topic
**Natural Hazards and Earth System Sciences**
An interactive open-access journal of the European Geosciences Union

Journal topic

- About
- Editorial board
- Articles
- Special issues
- Highlight articles
- Manuscript tracking
- Book reviews
- Subscribe to alerts
- Peer review
- For authors
- For reviewers
- EGU publications
- Imprint
- Data protection

- About
- Editorial board
- Articles
- Special issues
- Highlight articles
- Manuscript tracking
- Book reviews
- Subscribe to alerts
- Peer review
- For authors
- For reviewers
- EGU publications
- Imprint
- Data protection

**Research article**
01 Aug 2018

**Research article** | 01 Aug 2018

Assessment of the peak tsunami amplitude in the region of Taiwan

^{1}Department of Earth Sciences, National Central University, Taoyuan City 32001, Taiwan, R.O.C.^{2}Earthquake-Disaster & Risk Evaluation and Management Center, National Central University, Taoyuan City 32001, Taiwan, R.O.C.^{3}Graduate Institute of Hydrological and Oceanic Sciences, National Central University, Taoyuan City 32001, Taiwan, R.O.C.^{4}Center of Excellence for Ocean Engineering, National Taiwan Ocean University, Keelung, Taiwan

^{1}Department of Earth Sciences, National Central University, Taoyuan City 32001, Taiwan, R.O.C.^{2}Earthquake-Disaster & Risk Evaluation and Management Center, National Central University, Taoyuan City 32001, Taiwan, R.O.C.^{3}Graduate Institute of Hydrological and Oceanic Sciences, National Central University, Taoyuan City 32001, Taiwan, R.O.C.^{4}Center of Excellence for Ocean Engineering, National Taiwan Ocean University, Keelung, Taiwan

Abstract

Back to toptop
The southernmost portion of the Ryukyu Trench near the island of
Taiwan potentially generates tsunamigenic earthquakes with magnitudes from
7.5 to 8.7 through shallow rupture. The fault model for this potential region
dips 10^{∘} northward with a rupture length of 120 km and a width of
70 km. An earthquake magnitude of *M*_{w} 8.15 is estimated by the fault
geometry with an average slip of 8.25 m as a constraint on the earthquake
scenario. Heterogeneous slip distributions over the rupture surface are
generated by a stochastic slip model, which represents the decaying slip spectrum
according to *k*^{−2} in the wave number domain. These synthetic slip
distributions are consistent with the abovementioned identical seismic
conditions. The results from tsunami simulations illustrate that the
propagation of tsunami waves and the peak wave heights largely vary in
response to the slip distribution. Changes in the wave phase are possible as
the waves propagate, even under the same seismic conditions. The tsunami
energy path not only follows the bathymetry but also depends on the slip
distribution. The probabilistic distributions of the peak tsunami amplitude
calculated by 100 different slip patterns from 30 recording stations reveal
that the uncertainty decreases with increasing distance from the tsunami
source. The highest wave amplitude for 30 recording points is 7.32 m at
Hualien for 100 different slips. Compared with the stochastic-slip
distributions, the uniform slip distribution will be highly
underestimated, especially in the near field. In general, the uniform slip
assumption only represents the average phenomenon and will consequently
ignore the possibility of tsunami waves. These results indicate that
considering the effects of heterogeneous slip distributions is necessary for
assessing tsunami hazards to provide additional information about tsunami
uncertainties and facilitate a more comprehensive estimation.

How to cite

Back to top
top
How to cite.

Sun, Y.-S., Chen, P.-F., Chen, C.-C., Lee, Y.-T., Ma, K.-F., and Wu, T.-R.: Assessment of the peak tsunami amplitude associated with a large earthquake occurring along the southernmost Ryukyu subduction zone in the region of Taiwan, Nat. Hazards Earth Syst. Sci., 18, 2081-2092, https://doi.org/10.5194/nhess-18-2081-2018, 2018.

1 Introduction

Back to toptop
Almost all destructive tsunamis are generated by shallow earthquakes that
occur within subduction zones. Numerous destructive tsunami events, including
the *M*_{w} 9.1 Sumatra earthquake in 2004 (Lay et al., 2005), the
*M*_{w} 8.8 Chile earthquake in 2010 (Lay et al., 2010; Fritz et al., 2011)
and the *M*_{w} 9.0 Tohoku earthquake in 2011 (Goda et al., 2015; Goda and
Song, 2016), all of which occurred in subduction zones, have occurred
recently. The island of Taiwan, which is located at the convergent boundary
between the Philippine Sea Plate and the Eurasian Plate, is constantly under
possible threat of a tsunami. The convergence rate in this area is
approximately 80–85 mm yr^{−1} (Seno et al., 1993; Yu et al., 1997; Sella et al., 2002; Hsu et
al., 2009, 2012). Thus, earthquakes occur frequently in and around Taiwan.
The shallow earthquakes that occur in the Manila Trench to the south and the
Ryukyu Trench to the northeast are particularly tsunamigenic, and earthquakes
occur more actively in the southernmost Ryukyu Trench than in the northern
Manila Trench (Wu et al., 2013). The most well-known historic tsunami events
that occurred in northeastern Taiwan are the 1867 Keelung earthquake
(*M*_{w} 7.0) (Tsai, 1985; Ma and Lee, 1997; Cheng et al., 2016; Yu et al.,
2016) and the 1771 Yaeyama (Japan) earthquake (*M*_{w} ∼ 8) (Nakamura,
2009). Accordingly, these historic recordings demonstrate that Taiwan is
under a potential tsunami threat. Furthermore, the 2011 Tohoku earthquake
induced a powerful tsunami that destroyed coastal areas and caused nuclear
accidents (Mimura et al., 2011). As there are four nuclear power plants along
the coast of Taiwan, it is necessary to carefully estimate the tsunami hazard
in addition to the hazards of compound disasters.

Probabilistic tsunami hazard analysis (PTHA) is a modification of probabilistic seismic hazard analysis (PSHA) (Cornell, 1968; SSHAC, 1997), and it is intended to forecast the probability of tsunami hazards for a given region as comprehensively as possible. The recurrence rates of earthquakes have typically been estimated using the Gutenberg–Richter relationship (Gutenberg and Richter, 1944) for a defined source region in consideration of tsunamis triggered by earthquakes. The assessment of the wave height is one of the primary differences between PTHA and PSHA. PSHA assesses the ground motion based on empirical attenuation relationships (Wang et al., 2016), while PTHA assesses tsunami wave heights using empirical approaches or tsunami simulations (Geist, 2002; Geist and Parsons, 2006, 2009). Geist and Parsons (2006) mentioned that the tsunami wave height follows a definable frequency–size distribution over a sufficiently long period of time within a given coastal region (Soloviev, 1969; Houston et al., 1977; Horikawa and Shuto, 1983; Burroughs and Tebbens, 2005). This method is of great use in establishing the tsunami probability for a region if there is an extensive catalog of observed tsunami wave heights. However, given the wide distribution of global tsunamigenic earthquakes within seafloor regions throughout subduction zones, the tsunami records obtained from coastal gauges or/and ocean buoys are too sparse to comprehensively assess the associated hazards and the recording time, since their deployment is too short to enable a study on the recurrence intervals of tsunamis and earthquakes. Consequently, because the existing tsunami catalog is limited, simulations represent an effective approach. Conventional tsunami simulation adopts a simple source approximation and applies elastic dislocation theory to calculate the deformation of the seafloor surface assuming a uniform slip over the entire fault surface (Okada, 1985; Okal, 1982). However, the complexities of earthquake rupture processes play a substantial role in the generation of tsunamis. Conventional approaches are therefore unable to capture various features of short-wavelength tsunamis in the near field (Geist, 2002; Geist and Parsons, 2009). The results of previous studies that simulated tsunamis originating from historical earthquakes around Taiwan (Ma and Lee, 1997; Wu et al., 2008) using uniform slip models agreed only with long-wavelength observations. For the purposes of hazard mitigation, it is critical to predict the amplitudes of tsunamis along various coastlines for a given earthquake as accurately as possible. To make such predictions, the effects of the rupture complexity must be taken into consideration. Recent developments in PTHA have included the adoption of stochastic slip distributions of earthquakes to determine the overall probability of particular tsunami heights (Geist and Parsons, 2006, 2009). The adoption of stochastic slip distributions is able to quantify the variations in reasonable evaluations of the probabilities of specified tsunami heights at individual locations resulting from a specific fault.

In this study, we assess the heights of tsunamis along the coastline of
Taiwan generated by the potential tsunamigenic zone at the southernmost end
of the Ryukyu subduction zone. This potential zone is located close to
Taiwan, and at least 10 earthquakes (*M*_{w} > 7) have occurred
over the past 100 years (Hsu et al., 2012), the largest of which was the
*M*_{w} 7.7 in 1920 (Theunissen et al., 2010). For this area, the plausible
magnitude of the strongest earthquake is determined within the range between 7.5
and 8.7 (*M*_{w}) (Hsu et al., 2012). The fault zone is bounded by the
Longitudinal Valley Fault to the west and the Gagua Ridge to the east (Hsu et
al., 2012). This fault geometry with a defined rupture length and width is
employed herein, and an earthquake with a magnitude of 8.15 is used in the
tsunami simulation. The stochastic slip model is invoked to describe the
uncertainty in the rupture pattern over the fault plane to enable a more
realistic assessment of the tsunami probability.

2 Earthquake scenario and tsunami simulation

Back to toptop
The estimated maximum magnitude of a possible earthquake scenario is
essential for establishing the fundamental seismic conditions of the tsunami
simulation. The scenario of a potential rupture fault extending to a depth of
13 km proposed by Hsu et al. (2012) occurs along the southernmost Ryukyu
trench with a rupture length of 120 km, a width of 70 km and a dip of
10^{∘}. Kanamori and Anderson (1975) investigated the relation between
the rupture area and moment, and revealed that most of the average stress
drops (Δ*σ*) vary between 10 and 100 bars. The average stress
drop for most interplate earthquakes is approximately 30 bars, and thus, we
set an average stress drop of 30 bars. The stress drop and seismic moment
(*M*_{0}) relation along a dip-slip fault is described as follows (Kanamori
and Anderson, 1975):

$$\begin{array}{}\text{(1)}& {M}_{\mathrm{0}}={\displaystyle \frac{\mathit{\pi}\left(\mathit{\lambda}+\mathrm{2}\mathit{\mu}\right)}{\mathrm{4}\left(\mathit{\lambda}+\mathit{\mu}\right)}}\mathrm{\Delta}\mathit{\sigma}{W}^{\mathrm{2}}L,\end{array}$$

where *W* and *L* are the width and length of the rupture plane. We can obtain the moment for this scenario under an average
stress drop of 30 bars with the assumed rupture geometry. In Eq. (1), *μ*
denotes the rigidity and *λ* is the Lamè parameter. We assume that
the crust is elastic and homogeneous. Hence, $\mathit{\mu}=\mathit{\lambda}=\mathrm{30}$ GPa (Fowler,
2004; Piombo et al., 2007). Additionally, the seismic moment can be
represented by the rupture area and average slip as follows (Lay and Wallace,
1995):

$$\begin{array}{}\text{(2)}& {M}_{\mathrm{0}}=\mathit{\mu}A\stackrel{\mathrm{\u203e}}{D}.\end{array}$$

Moreover, the seismic moment is dependent on the rupture area (*A*) and
average slip ($\stackrel{\mathrm{\u203e}}{D})$; thus, the average slip can be estimated by
Eq. (2), and it is calculated to be 8.25 m. Then, the seismic moment can be
transformed into the magnitude *M*_{w} by the following (Hanks and Kanamori,
1979):

$$\begin{array}{}\text{(3)}& {M}_{\mathrm{w}}=\left({\displaystyle \frac{\mathrm{log}{M}_{\mathrm{0}}}{\mathrm{1.5}}}\right)-\mathrm{10.73}.\end{array}$$

Therefore, the maximum possible earthquake magnitude is *M*_{w} 8.15
(${M}_{\mathrm{0}}=\mathrm{2.07}\times {\mathrm{10}}^{\mathrm{28}}$ dyne-cm).

The rupture process of an earthquake is extremely complex. Seismic inversion
results reveal that the slip distribution of a rupture has a heterogeneous
spatio-temporal development. Consequently, using a simplified uniform slip
distribution to simulate a tsunami captures only the long-wavelength portion
of the tsunami field (Geist and Dmowska, 1999). In addition, the temporal
description of the seismic rupture process can be ignored because the
propagation velocity of the tsunami wave is substantially slower than the
seismic rupture velocity (Dean and Dalrymple, 1991; Ma et al., 1991; Wang and
Liu, 2006). Andrews (1980) showed that the static slip distribution is
directly related to stress changes and that the spectrum of the slip
distribution is proportional to *k*^{−2} decay in the wave number domain:

$$\begin{array}{}\text{(4)}& \left|{F}_{s,t}\left[{D}_{x,y}\right]\right|\propto {k}^{-\mathrm{2}}\end{array}$$

where *D*_{x,y} is the slip distribution over a 2-D lattice, *F*_{s,t} is the
2-D Fourier transform, and $k=\sqrt{{k}_{x}^{\mathrm{2}}+{k}_{y}^{\mathrm{2}}}$ is the radial
wave number. The *k*^{−2} power law indicates that the slip distribution has
self-similar characteristics; moreover, this characteristic can also be
demonstrated from a fractal perspective (Tsai, 1997). Based on
self-similarity, Herrero and Bernard (1994) introduced the *k*-square model,
which leads to the *ω*-square model (Aki, 1967). The slip spectrum
follows *k*^{−2} decay beyond the corner radial wave number (*k*_{c}),
which is proportional to 1/*L*_{c}. The *L*_{c} depends on the
characteristic rupture dimension (Geist, 2002).

The heterogeneous slip distribution is proportional to *k*^{−2} and is
similar to fractional Brownian motion as a stochastic process (Tsai, 1997).
The stochastic slip distribution can be described by multiplication in the
Fourier domain:

$$\begin{array}{}\text{(5)}& {D}_{x,y}\propto {F}_{x,y}^{-\mathrm{1}}\left[{F}_{s,t}\left[{X}_{x,y}\right]\times {k}^{-\mathrm{2}}\right],\end{array}$$

where *X*_{x,y} is a random variable for the spatial distribution that
randomizes the phase, and ${F}_{x,y}^{-\mathrm{1}}$ is the inverse 2-D Fourier
transform. The random distribution of *X*, which is best described by a
non-Gaussian distribution, especially by a Lèvy distribution, can be
calculated by reversing Eq. (5) (Lavallée and Archuleta, 2003;
Lavallée et al., 2006). The Lèvy distribution can be described by
four parameters, namely *α*, *β*, *γ* and *μ*_{L}, as
follows:

$$\begin{array}{ll}{\displaystyle}\mathit{\phi}& {\displaystyle}\left(t\right)=\\ \text{(6)}& {\displaystyle}& {\displaystyle}\left\{\begin{array}{ll}\mathrm{exp}\left(-{\mathit{\gamma}}^{\mathit{\alpha}}{\left|t\right|}^{\mathit{\alpha}}\left[\mathrm{1}+i\mathit{\beta}\mathrm{sign}\left(t\right)\mathrm{tan}{\displaystyle \frac{\mathit{\pi}\mathit{\alpha}}{\mathrm{2}}}\right.\right.& \\ \left.\left.\left({\left|\mathit{\gamma}t\right|}^{\mathrm{1}-\mathit{\alpha}}-\mathrm{1}\right)\right]+i{\mathit{\mu}}_{L}t\right),& \mathit{\alpha}\ne \mathrm{1}\\ \mathrm{exp}\left(-\mathit{\gamma}\left|t\right|\left[\mathrm{1}+i\mathit{\beta}{\displaystyle \frac{\mathrm{2}}{\mathit{\pi}}}\mathrm{sign}\left(t\right)\right.\right.& \\ \left.\left.\left(\mathrm{ln}\left|t\right|+\mathrm{ln}\mathit{\gamma}\right)\right]+i{\mathit{\mu}}_{L}t\right),& \mathit{\alpha}=\mathrm{1}\end{array}.\right.\end{array}$$

The parameter *α*, 0 < *α*≤2, affects the falloff
rate of the probability density function (PDF) for the tail. The parameter
*β*, $-\mathrm{1}\le \mathit{\beta}\le \mathrm{1}$ controls the skewness of the PDF, and the
parameter *γ*, *γ* > 0 controls the width of the PDF.
The parameter *μ*_{L}, −∞ < *μ*_{L} < ∞, is related to the location of
the PDF. The Lèvy distribution is effective at describing the
distribution of a random variable, i.e., *X*, from real earthquake events,
implying that the slip distribution without self-similarity has a heavy tail
behavior (Lavallée et al., 2006). Based on experiments of generating
stochastic slip distributions, this heavy tail behavior affects the intensity
of an extreme value (Lavallée and Archuleta, 2003).

The stochastic slip distribution is generated by a 2-D spatially random
distribution by imposing a self-similar characteristic beyond the corner
radial wave number, which is constrained by the rupture dimension, in the
wave number domain. In this study, the potential rupture fault is divided into
5 × 5 km^{2} subfaults. The grid is composed of 24×14
meshes along the strike and dip directions. The produced
variable with a spatially random distribution adopts the Lèvy
distribution (*α*=1.51, *β*=0.2, *γ*=28.3, ${\mathit{\mu}}_{L}=-\mathrm{0.9}$),
which is the dip-slip result from Lavallée et al. (2006) shown in
Fig. 1a. In Lavallée et al. (2006), the slip distribution of the
Northridge earthquake was divided into the dip-slip and strike-slip
directions, and they were calculated by an inverse 2-D stochastic model to
obtain the values of the Lèvy PDF. The values of the Lèvy PDF, which
are mentioned above, are indicative of the result of the dip-slip direction. The
Northridge earthquake is a thrust earthquake (Davis, 1994), and thus, it has
a faulting mechanism that is approximately similar to our scenario fault
model. There are no inverted slip models of past earthquakes in the study
area that can be used to conduct an analysis of the Levy PDF parameters; therefore, the
Lèvy distribution in Lavallée et al. (2006) is adopted in this study.
From the perspective of mathematical operations, the slip distribution in
Eq. (5) represents a filtered random distribution. However, for consistency
with the physical behavior over the rupture surface suggested by the results
of the inverse modeling, truncation of the Lèvy distribution must be
performed to constrain the extreme slip value. The synthetic slip
distribution (Fig. 1b) produced by the spatially random distribution in
Fig. 1a is heterogeneous, and its power spectrum obeys a *k*-square model at
high wave numbers (Fig. 1c). The average slip of this synthetic slip
distribution is 8.25 m, indicating that the earthquake energy is constant as
estimated above, and the maximum slip is 31.02 m. One hundred different slip
distributions are produced for the tsunami simulation representing the
uncertainty in the results associated with complex rupture processes. In the
100 sets of results, the maximum slip range is between 20.17 and 37.97 m.
Smooth processes are not included, nor are additional regional constraints
for the slip distribution. There are two reasons for this application. The
first is that we do not have information regarding where the plate interface
is locked or the locations of asperities often repeat in historical events.
The second is that some studies reported that the asperities extend to the
boundary of the fault model (Ide et al., 2011; Lay et al., 2011; Shao et al.,
2011; Yue and Lay, 2011). According to these reasons, we do not apply any
additional constraints for stochastic slip distributions. Similarly, the
uniform slip case constitutes a complete uniform slip distribution over the
whole fault plane. Figure 1b and d demonstrate the stochastic distribution of
the scenario source models causing the maximum and minimum wave heights
at recording station 26 (Hualien) (Fig. 2). Both patterns
affecting the propagation will be discussed in Sect. 3.1.

Figure 2 shows the computational domain, recording stations and fault model.
The potential rupture fault is divided into 5 × 5 km^{2}
subfaults, and the stochastic slip distribution model is applied to determine
the amount of discrete slip on each subfault. Vertical seafloor displacements
caused by slip along the rupture plane are calculated using elastic
dislocation theory (Okada, 1985). The Cornell Multi-grid Coupled Tsunami
model (COMCOT) is used to run the tsunami simulations. COMCOT is capable
of efficiently studying the entire life span of a tsunami, including its
generation, propagation, runup and inundation (Wang, 2009), and it has been
widely used in studying many historical tsunami events, such as the 1960
Chilean tsunami (Liu et al., 1995), 1992 Flores Island tsunami (Liu et al.,
1995), 2003 Algeria tsunami (Wang and Liu, 2005), 2004 Indian Ocean tsunami
(Wang and Liu, 2006, 2007) and 2006 Pingtung tsunami, Taiwan (Wu et al.,
2008; Chen et al., 2008). COMCOT solves linear or nonlinear shallow water
equations for spherical or Cartesian coordinates using the finite difference
method. With a flexible nested grid system, it can properly guarantee both
the efficiency and the accuracy from the near-coastal region to the far-field
region. Two grid layers are used to simulate the propagation of tsunamis. The
Manning coefficient is 0.013 in this study, which assumes a sandy sea bottom (Wu
et al., 2008). The bathymetry-adopted open data from the National Oceanic and
Atmospheric Administration (NOAA) that can be download from
https://maps.ngdc.noaa.gov/viewers/wcs-client/, last access: 26 March 2018 (Amante and Eakins, 2009). The resolution of the outer
layer is 4 min for the solution of the linear shallow water equation, and
the resolution of the inner layer is 1 min for the solution of the nonlinear
form of the shallow water equation. There are 30 recording stations referring
to the positions of tidal gauges maintained by the Central Weather Bureau
(CWB) along the coastlines of Taiwan and the outlying islands. The CWB
website presents the locations of the tide stations
(http://e-service.cwb.gov.tw/HistoryDataQuery/index.jsp, last access: 26 March 2018 and
https://www.cwb.gov.tw/V7e/climate/marine_stat/tide.htm, last access: 26 March 2018). These locations are shifted slightly to the grid
nodes to accurately record the data. Table 1 presents the locations and water
depths of the recording stations in the computational mesh.

3 The effect of heterogeneous slip on tsunamis

Back to toptop
The stochastic slip model produces different slip distributions with the same
fault geometry in addition to a constant average slip and a constant seismic
moment. The model is used to describe the heterogeneous slip pattern of an
earthquake and to further examine its effect on the tsunamis originating from
the southernmost end of the Ryukyu subduction zone adjacent to Taiwan.
According to the previous sections, the maximum possible earthquake magnitude
is determined to be *M*_{w} 8.15 with an average slip of 8.25 m.
Furthermore, the uniform slip distribution on the rupture plane is also used
to simulate a tsunami to facilitate a discussion of the difference between the
effects of uniform and heterogeneous slip on tsunamis.

The static vertical displacement of the ocean floor is modeled using elastic dislocation theory (Okada, 1985) with a static slip distribution. The vertical seafloor displacement is modeled as the initial water level, and the horizontal component of the seabed displacement is not included in the simulation. Figure 3a shows the initial water elevations produced by a uniform slip distribution, and Fig. 3b exhibits the maximum free-surface elevation during the propagation. Figure 3c and e demonstrate the initial water elevations produced by the stochastic slip distributions (Fig. 1b and d). The initial water elevation with a uniform slip distribution is simple and smooth, but those with stochastic slip models are more complex and relatively heterogeneous. Nonuniform slip causes an apparent change in the wavelength distribution of the initial free-surface elevation (i.e., the potential energy distribution), which affects the path of energy propagation. In the uniform slip scenario, the maximum free-surface elevation pattern is straightforward and clearly controlled by the topography. However, many strong and seemingly chaotic paths of wave energy appear in the nonuniform slip scenarios, and the free-surface field exhibits additional uncertainties in terms of the flow. In Fig. 3b, the maximum free-surface elevation mainly propagates toward two places where the seafloor bathymetry becomes shallow relative to the deep areas northeast of Taiwan, as shown in Fig. 2. Although the propagation paths due to the nonuniform slip distributions (Fig. 3d and f) also have the same characteristics, it is notable that the paths followed by the wave energy differ depending on the rupture pattern. To the northeast of Taiwan in Fig. 3f, there is a strong wave path connecting the two higher-elevation areas of bathymetry. However, this behavior is not observed in Fig. 3b and d. In addition, the maximum elevation on the footwall in Fig. 3d is higher than that in Fig. 3f. In Fig. 3b, the high elevation appears only along the coast on the footwall side. These results indicate that the wave energy variation depends on the rupture pattern, thereby causing differences in the wave paths and leading to completely different tsunami amplitudes.

Thirty stations located along the coastlines are available for recording the amplitude of the tsunami wave height. Relative to the other stations, stations 25 (Shihti), 26 (Hualien) and 27 (Suao) are situated near the potential rupture fault, and they have high wave amplitudes and enormous variations in the tsunami simulations of 100 different slip distributions; consequently, the time series of the wave heights at these stations are shown as an example (Fig. 4). The time series of the wave heights at the other stations are shown in the Supplement. The variability in the distribution of the initial free-surface elevation results in substantial phase changes and different wave heights. It is worth noting that the average of the disordered and chaotic time series produced by the 100 different slip distributions is almost identical to the results of the time series produced by the uniform case. This implies that the uniform slip distribution simply represents an average result and that it cannot represent all of the possible situations.

According to the statistical results from 100 different slip patterns (Table 1) for 30 stations, Hualien station has a maximum wave amplitude of 7.32 m, and its maximum wave amplitude interval ranges from 1.87 to 7.32 m, which constitutes the widest interval for any recording site, and the standard deviation of this distribution is 1.024 m. These findings indicate that Hualien station has a high uncertainty in this scenario. However, the maximum wave amplitudes from the uniform slip distribution are relatively lower than those from the stochastic results. Following the above findings, we need to consider whether the estimations from the uniform slip case are appropriate for hazard analysis by focusing on the maximum wave amplitude issue.

According to the results of our simulations, we calculate the probability of the peak tsunami amplitude (PTA) at each recording station as shown in the histogram of Fig. 5. To verify the representativeness of the PTA probability distributions, another 100 sets of different slip distributions are produced and simulated under the same seismic conditions. In Fig. 5, the shapes of the PTA distributions from another 100 sets (black lines) are similar to the shapes of the histograms from the first 100 sets. These results verify the representativeness of the PTA probability distributions produced from 100 sets of slip distributions. This test also reinforces the reproducibility of our simulations and demonstrates that the number of simulations is roughly satisfactory for statistical analysis. Of course, the more slip distributions we use, the more comprehensive and stable the range we obtain.

In Fig. 5, the PTA distributions for the stations in eastern Taiwan (red
markers) have much higher values than those in western Taiwan (blue
markers) due to the specified location of the tsunami source. The shapes of
the PTA distributions in eastern Taiwan resemble lognormal distributions,
while those in western Taiwan resemble normal distributions. We suppose that
the attenuation of the wave propagation causes the lognormal distributions to
degenerate into normal distributions. The PTAs produced by a uniform slip
distribution are generally located in the middle of the PTA distributions.
Both PTA values (i.e., the value of the PTA from the uniform slip
distribution and those from the stochastic slip distribution models) decrease
with distance from the potential fault due to the attenuation of the wave
propagation (Fig. 5 shows the results for all stations, and Fig. 6 shows the
results for stations 20 through 30 in eastern Taiwan). However, some
stations, e.g., stations 17, 19 and 21, do not precisely follow this trend;
this could be the result of the coastal topography and the presence of an
energy channel. From Fig. 3d, in comparison with the adjacent coastline,
station 21 is located exactly where the wave energy gathers. In addition,
broad distributions are frequently observed at promontories along the
coastline and are caused by complex propagation path effects between the
source region and the recording locations (Geist, 2002). There are many
compound factors that affect the tsunami propagation and maximum wave height.
Figure 6 presents the relation between the distance and wave height and shows
the PTA distributions following Fig. 5. The *x* axis presents the shortest
distance between the stations and fault plane. On the footwall side, stations
20 and 22, which do not directly face the energy propagation path (Fig. 3f),
are located on islands off the coast of Taiwan; consequently, their PTA
distributions are lower than those of stations 21 and 23, even though the
distances from the potential fault are similar. On the hanging wall, station
29 is farther from the coastline of Taiwan than other stations; however,
because of the real location of the station and its numerical grid setting,
its PTA distribution is lower than that of station 30 (Fig. 3b). The ranges
of the PTA distributions converge with increasing distance on both sides
of the fault. Moreover, the PTA distributions and their average values
roughly exhibit a linear decrease with increasing distance except for
stations 25 and 26. In contrast, these two stations in the near field are
directly affected by initial water elevation, and thus, the PTAs caused by
uniform slip are quite low.

Although the seismic parameters have already been defined as constants in our experiment, there exists an uncertainty in the PTA, which is not a constant value. Hence, the uniform case cannot provide this uncertainty, and thus, the PTA could be underestimated. The results give specific PTA ranges, which represent the wave height uncertainties for the scenario of earthquakes originating from the Ryukyu Trench. It is therefore necessary to consider the effects of a heterogeneous slip distribution to comprehensively assess the tsunami hazard.

4 Discussion

Back to toptop
Most coastlines threatened by near-field tsunamis, such as the coasts of
Chile, Japan and Indonesia, are parallel to the trench axis of the associated
subduction zones. Many tsunami events, including the *M*_{w} 8.8 Chile
earthquake in 2010 (Lay et al., 2010; Fritz et al., 2011), the *M*_{w} 9.0
Tohoku earthquake in 2011 (Goda et al., 2015; Goda and Song, 2016), the
*M*_{w} 9.1 Sumatra earthquake in 2004 (Lay et al., 2005), and the
*M*_{w} 8.1 Mentawai earthquake in 2010 (Satake et al., 2013), have occurred
along these regions. However, the potential rupture fault in this study along
the southernmost Ryukyu subduction zone is perpendicular to the coast of
Taiwan, which directly affects the first wave motion. The first motion on the
footwall is up; conversely the first motion on the hanging wall is down. As a
result, the coastline retreats from the land to the sea as the first tsunami
wave approaches, allowing people additional time to leave the seafront.

The effect of a heterogeneous slip distribution is important and necessary
for near-field estimations (Geist, 2002 and Ruiz et al., 2015).
Figure 5 shows that the PTA distributions in the near field are broad, and
they narrow with increasing distance from the potential fault. The
uncertainty in the near field is higher than in the far field. At most
of the eastern stations, the values of the average PTA approach uniform
results, but the uniform slip results at stations 25 and 26 are close to the
minimum PTA (Table 1). Geist (2002) presented the average and extreme
nearshore PTA calculated for 100 different slip distributions and compared
them with the uniform slip results (Fig. 6a in Geist, 2002). The range of the
PTA also becomes narrower with increasing distance. The values from the
uniform slip distribution and the average PTA are similar, but some of the
average values are close to the minimum PTA between approximately 19 and
19.5^{∘} N. Similar characteristics of the average PTA and the results
from the uniform case are observed in different regions. The average PTA is
equal to the uniform slip result in the nearshore region, but this could be
caused by other factors (e.g., distance to the tsunami source, propagation
path, initial water elevation) that shift the average PTA toward the
minimum PTA.

Four nuclear power plants (NPPs) are located on the island of Taiwan. According to the numerical results, we infer that the mean PTA in the coastal area of NPP4 ranges from approximately 2 to 3 m. The distribution at this plant may be wider than those at other nuclear power plants due to its position relative to the tsunami source. Moreover, NPP4 is located on the shore of a bay with a curved shape; the magnification effect from the geometrical shape of the bay may serve to enhance the PTA therein. NPP3 also exhibits this condition insomuch that the energy is concentrated at the location of the plant (Fig. 3b, d and f). For the coastal areas around NPP1 and NPP2, the PTA distributions are between 1 and 2 m. The coastlines of these two nuclear power plants slightly face the direction of tsunami propagation, and thus, their PTAs should be higher than those along adjacent coastlines (Fig. 3b, d and f). In general, under this scenario, the coastline at NPP4 has the largest threat. Although NPP3 is far from the tsunami source, it faces a wave height of approximately 1.5 m on average with a ±0.5 m range of uncertainty. However, NPP3 is closer to the Manila subduction zone, and thus, it could be threatened by a tsunami originating from the Manila Trench. In contrast, the coastlines of NPP1 and NPP2 are relatively safe and have fewer uncertainties with regard to the PTA.

The use of heterogeneous slip patterns clearly delineates the range of possible waveforms and provides more information on latent uncertainties in the wave height. The 95 % confidence intervals for the wave height from 100 sets in each time series provide us with a specific range for the amplitude of the tsunami wave (Fig. 4). According to these time series, we are aware of the periods of tsunami runup and runoff and can prepare supporting policies to reduce associated disasters. For example, a nuclear power plant uses a trench from the ocean to intake water to cool the reactor; thus, if the sea level is too low to take in water, the temperature of the reactor will rise excessively, causing a nuclear disaster. Based on the results of simulations, we can estimate how much water should be stored for tsunami runoff. This issue requires more attention in Taiwan because four nuclear power plants are located near the coast.

The results of the tsunami simulations illustrate that the effect of the slip
distribution on the rupture plane has significant effects on the wave
propagation and wave height. The correctness of this slip distribution
determines whether the wave height calculations represent a useful reference.
However, some parameters of the stochastic models could influence the
synthetic slip distributions. For instance, the exponent of the slip spectrum
is associated with the roughness of the slip distribution. Higher exponential
values inhibit the powers of high wave numbers, leading to smoother slip
distributions; conversely, lower values lead to rougher slip distributions.
In general, the *k*-square model needs to be followed. Furthermore, the
interpolation of the slip distribution for a given geometry will affect the
exponent of *k* (Tsai, 1997). Interpolation will smooth the original pattern.
The powers of short wave numbers will be depressed and the powers of long
wave numbers will be enhanced. Moreover, the random spatial variability of the
slip distribution is relatively critical. According to Lavallée and
Archuleta (2003) and Lavallée et al. (2006), we adopt the truncated
non-Gaussian distribution for the spatial variability. This truncation limits
the non-Gaussian distribution to a particular range. However, extreme
truncation will cause the heavy-tailed characteristic of this distribution to
become less pronounced or even disappear, similarly to a Gaussian distribution.
In mathematics, the synthetic slip distribution is a filtering process
insomuch that the characteristics of a heavy-tailed distribution affect the
extrema of the slip distribution. The maximum slip will be greater as the
truncated range increases, and the maximum slip may exceed reasonable values
if the truncated range is excessively wide. Therefore, the parameters must be
chosen carefully to match the observations acquired by inversion.

5 Conclusions

Back to toptop
The maximum possible earthquake magnitude is *M*_{w} 8.15 with an average
slip of 8.25 m in the southernmost portion of the Ryukyu Trench. One hundred
slip distributions of the seismic rupture surface were generated by a
stochastic slip model. The maximum slip range is between 20.17 and 37.97 m,
and the average slip of each model is consistent with 8.25 m. A
heterogeneous slip distribution induces variability in the tsunami wave
heights and the associated paths of propagation. The simulated results
demonstrate that the complexity of the rupture plane has a significant
influence on the near field for local tsunamis. The PTA distribution provides
a specific range for the wave height and the probability of occurrence in
this scenario. These distributions and their average values exhibit an
approximately linear decrease with increasing distance. The coastline, which
is situated very close to or even atop the tsunami source, is directly
affected by the rupture slip distribution. Then, the range of the PTA
distribution will converge with increasing distance from the tsunami source.
In this study, Hualien station, which is located directly above the source,
has the widest PTA interval (1.87–7.32 m) and the highest wave amplitude.
The statistical summary reveals that this station has a standard deviation
of 1.63 m, which is larger than those of the other stations, and the largest
uncertainty. However, the PTA caused by the uniform slip distribution is only
1.63 m, which is much lower and even below the average (3.36 m) at this
station. This finding implies that a simplified earthquake source cannot
completely represent the tsunami amplitudes in reality. If we adopt a uniform
slip distribution to assess tsunami hazards, those hazards will be critically
underestimated. Furthermore, the variances of tsunami amplitudes are
imperative for assessing tsunami hazards, and the quantitative technique
employed is also important.

Data availability

Back to toptop
Data availability.

The data are available upon request from the authors.

Supplement

Back to toptop
Supplement.

The supplement related to this article is available online at: https://doi.org/10.5194/nhess-18-2081-2018-supplement.

Competing interests

Back to toptop
Competing interests.

The authors declare that they have no conflict of interest.

Acknowledgements

Back to toptop
Acknowledgements.

The authors
are grateful for research support from the Ministry of Science and Technology (ROC) and
the Department of Earth Science, National Central University (ROC). This work is supported
by “Earthquake-Disaster & Risk Evaluation and Management Center, E-DREaM” from the
Featured Areas Research Center program within the framework of the Higher Education
Sprout Project by the Ministry of Education (MOE) in Taiwan. The Taiwan Earthquake
Research Center (TEC) contribution number for this article is 00145. We thank Eric L.
Geist for discussing the stochastic slip model. The authors very much appreciate the
constructive comments received on the manuscript by two anonymous reviewers, as well
as the thoughtful insights of the editor, Piero Lionello.

Edited by: Piero Lionello

Reviewed by: two anonymous referees

References

Back to toptop
Aki, K.: Scaling law of seismic spectrum, J. Geophys. Res., 72, 1217–1231, https://doi.org/10.1029/JZ072i004p01217, 1967.

Amante, C. and Eakins, B. W.: ETOPO1 1 arc-minute global relief model: procedures, data sources and analysis, US Department of Commerce, National Oceanic and Atmospheric Administration, National Environmental Satellite, Data, and Information Service, National Geophysical Data Center, Marine Geology and Geophysics Division Colorado, https://doi.org/10.7289/V5C8276M, 2009.

Andrews, D. J.: A stochastic fault model: 1. Static case, J. Geophys. Res.-Sol. Ea., 85, 3867–3877, https://doi.org/10.1029/JB085iB07p03867, 1980.

Burroughs, S. M. and Tebbens, S. F.: Power-law scaling and probabilistic forecasting of tsunami runup heights, Pure Appl. Geophys., 162, 331–342, https://doi.org/10.1007/s00024-004-2603-5, 2005.

Chen, P.-F., Newman, A. V., Wu, T.-R., and Lin, C.-C.: Earthquake Probabilities and Energy Characteristics of Seismicity Offshore Southwest Taiwan, Terr. Atmos. Ocean. Sci., 19, 697–703, https://doi.org/10.3319/TAO.2008.19.6.697(PT), 2008.

Cheng, S.-N., Shaw, C.-F., and Yeh, Y. T.: Reconstructing the 1867 Keelung Earthquake and Tsunami Based on Historical Documents, Terr. Atmos. Ocean. Sci., 27, 431–449, https://doi.org/10.3319/TAO.2016.03.18.01(TEM), 2016.

Cornell, C. A.: Engineering seismic risk analysis, Bull. Seism. Soc. Am., 58, 1583–1606, 1968.

Davis, T. L. and Namson, J. S.: A Balanced Cross-Section of the 1994 Northridge Earthquake, Southern California, Nature, 372, 167–169, https://doi.org/10.1038/372167a0, 1994.

Dean, R. G. and Dalrymple, R. A.: Water wave mechanics for engineers and scientists, World Scientific Publishing Co Inc, Singapore, https://doi.org/10.1142/9789812385512_0004, 1991.

Fowler, C. M. R.: The Solid Earth: An Introduction to Global Geophysics, Cambridge University Press, Cambridge, 728pp., 2004.

Fritz, H. M., Petroff, C. M., Catalán, P. A., Cienfuegos, R., Winckler, P., Kalligeris, N., Weiss, R., Barrientos, S. E., Meneses, G., Valderas-Bermejo, C., Ebeling, C., Papadopoulos, A., Contreras, M., Almar, R., Dominguez, J. C., and Synolakis, C. E.: Field Survey of the 27 February 2010 Chile Tsunami, Pure Appl. Geophys., 168, 1989–2010, https://doi.org/10.1007/s00024-011-0283-5, 2011.

Geist, E. L.: Complex earthquake rupture and local tsunamis, J. Geophys. Res.-Sol. Ea., 107, ESE 2-1-ESE 2-15, https://doi.org/10.1029/2000JB000139, 2002.

Geist, E. L. and Dmowska, R.: Local Tsunamis and Distributed Slip at the Source, in: Seismogenic and Tsunamigenic Processes in Shallow Subduction Zones, edited by: Sauber, J. and Dmowska, R., Birkhäuser Basel, Basel, https://doi.org/10.1007/978-3-0348-8679-6_6, 1999.

Geist, E. L. and Parsons, T.: Probabilistic Analysis of Tsunami Hazards*, Nat. Hazards, 37, 277–314, https://doi.org/10.1007/s11069-005-4646-z, 2006.

Geist, E. L. and Parsons, T.: Assessment of source probabilities for potential tsunamis affecting the U.S. Atlantic coast, Mar. Geol., 264, 98–108, https://doi.org/10.1016/j.margeo.2008.08.005, 2009.

Goda, K. and Song, J.: Uncertainty modeling and visualization for tsunami hazard and risk mapping: a case study for the 2011 Tohoku earthquake, Stoch. Environ. Res. Risk Assess., 30, 2271–2285, https://doi.org/10.1007/s00477-015-1146-x, 2016.

Goda, K., Yasuda, T., Mori, N., and Mai, P. M.: Variability of tsunami inundation footprints considering stochastic scenarios based on a single rupture model: Application to the 2011 Tohoku earthquake, J. Geophys. Res.-Oceans, 120, 4552–4575, https://doi.org/10.1002/2014JC010626, 2015.

Gutenberg, B. and Richter, C. F.: Frequency of earthquakes in California, Bull. Seism. Soc. Am., 34, 185–188, 1944.

Hanks, T. C. and Kanamori, H.: A moment magnitude scale, J. Geophys. Res.-Sol. Ea., 84, 2348–2350, https://doi.org/10.1029/JB084iB05p02348, 1979.

Herrero, A. and Bernard, P.: A Kinematic Self-Similar Rupture Process for Earthquakes, Bull. Seism. Soc. Am., 84, 1216–1228, 1994.

Horikawa, K. and Shuto, N.: Tsunami disasters and protection measures in Japan, Tsunamis-Their Science and Engineering, Terra Scientific Publishing Company, Tokyo, 9–22, 1983.

Houston, J. R., Carver, R. D., and Markle, D. G.: Tsunami-Wave Elevation Frequency of Occurrence for the Hawaiian Islands, Army Engineer Waterways Experiment Station, Vicksburg, MS, 66 pp., 1977.

Hsu, Y.-J., Yu, S.-B., Simons, M., Kuo, L.-C., and Chen, H.-Y.: Interseismic crustal deformation in the Taiwan plate boundary zone revealed by GPS observations, seismicity, and earthquake focal mechanisms, Tectonophysics, 479, 4–18, https://doi.org/10.1016/j.tecto.2008.11.016, 2009.

Hsu, Y.-J., Ando, M., Yu, S.-B., and Simons, M.: The potential for a great earthquake along the southernmost Ryukyu subduction zone, Geophys. Res. Lett., Geophys. Res. Lett., 39, L14302, https://doi.org/10.1029/2012GL052764, 2012.

Ide, S., Baltay, A., and Beroza, G. C.: Shallow dynamic overshoot and
energetic deep rupture in the 2011 *M*_{w} 9.0 Tohoku-Oki earthquake,
Science, 332, 1426–1429, https://doi.org/10.1126/science.1207020, 2011.

Kanamori, H. and Anderson, D. L.: Theoretical basis of some empirical relations in seismology, Bull. Seism. Soc. Am., 65, 1073–1095, 1975.

Lavallée, D. and Archuleta, R. J.: Stochastic modeling of slip spatial complexities for the 1979 Imperial Valley, California, earthquake, Geophys. Res. Lett., 30, 1245, https://doi.org/10.1029/2002GL015839, 2003.

Lavallée, D., Liu, P., and Archuleta, R. J.: Stochastic model of heterogeneity in earthquake slip spatial distributions, Geophys. J. Int., 165, 622–640, https://doi.org/10.1111/j.1365-246X.2006.02943.x, 2006.

Lay, T. and Wallace, T. C.: Modern global seismology, Academic Press, San Diego, California, 1995.

Lay, T., Kanamori, H., Ammon, C. J., Nettles, M., Ward, S. N., Aster, R. C., Beck, S. L., Bilek, S. L., Brudzinski, M. R., Butler, R., DeShon, H. R., Ekström, G., Satake, K., and Sipkin, S.: The Great Sumatra-Andaman Earthquake of 26 December 2004, Science, 308, 1127–1133, https://doi.org/10.1126/science.1112250, 2005.

Lay, T., Ammon, C. J., Kanamori, H., Koper, K. D., Sufri, O., and Hutko, A.
R.: Teleseismic inversion for rupture process of the 27 February 2010 Chile
(*M*_{w} 8.8) earthquake, Geophys. Res. Lett., 37, L13301,
https://doi.org/10.1029/2010GL043379, 2010.

Lay, T., Ammon, C. J., Kanamori, H., Xue, L., and Kim, M. J.: Possible large
near-trench slip during the 2011 *M*_{w} 9.0 off the Pacific coast of Tohoku
Earthquake, Earth Planets Space, 63, 687–692, https://doi.org/10.5047/eps.2011.05.033,
2011.

Liu, P. L. F., Cho, Y. S., Yoon, S. B., and Seo, S. N.: Numerical Simulations of the 1960 Chilean Tsunami Propagation and Inundation at Hilo, Hawaii, in: Tsunami: Progress in Prediction, Disaster Prevention and Warning, edited by: Tsuchiya, Y. and Shuto, N., Springer Netherlands, Dordrecht, 1995.

Liu, P. L. F., Cho, Y. S., Briggs, M. J., Kanoglu, U., and Synolakis, C. E.: Runup of solitary waves on a circular Island, J. Fluid Mech., 302, 259–285, https://doi.org/10.1017/S0022112095004095, 1995.

Ma, K.-F. and Lee, M.-F.: Simulation of historical tsunamis in the Taiwan region, Terr. Atmos. Ocean. Sci., 8, 13–30, 1997.

Ma, K.-F., Satake, K., and Kanamori, H.: The origin of the tsunami excited by the 1989 Loma Prieta Earthquake – Faulting or slumping?, Geophys. Res. Lett., 18, 637–640, https://doi.org/10.1029/91GL00818, 1991.

Mimura, N., Yasuhara, K., Kawagoe, S., Yokoki, H., and Kazama, S.: Damage from the Great East Japan Earthquake and Tsunami – A quick report, Mitig. Adapt. Strat. Gl., 16, 803–818, https://doi.org/10.1007/s11027-011-9297-7, 2011.

Nakamura, M.: Fault model of the 1771 Yaeyama earthquake along the Ryukyu Trench estimated from the devastating tsunami, Geophys. Res. Lett., 36, L19307, https://doi.org/10.1029/2009GL039730, 2009.

Okada, Y.: Surface deformation due to shear and tensile faults in a half-space, Bull. Seism. Soc. Am., 75, 1135–1154, 1985.

Okal, E. A.: Mode-wave equivalence and other asymptotic problems in tsunami theory, Phys. Earth Planet. Inter., 30, 1–11, https://doi.org/10.1016/0031-9201(82)90123-6, 1982.

Piombo, A., Tallarico, A., and Dragoni, M.: Displacement, strain and stress fields due to shear and tensile dislocations in a viscoelastic half-space, Geophys. J. Int., 170, 1399–1417, https://doi.org/10.1111/j.1365-246X.2007.03283.x, 2007.

Ruiz, J. A., Fuentes, M., Riquelme, S., Campos, J., and Cisternas, A.: Numerical simulation of tsunami runup in northern Chile based on non-uniform k −2 slip distributions, Nat. Hazards, 79, 1177–1198, https://doi.org/10.1007/s11069-015-1901-9, 2015.

Satake, K., Nishimura, Y., Putra, P. S., Gusman, A. R., Sunendar, H., Fujii, Y., Tanioka, Y., Latief, H., and Yulianto, E.: Tsunami Source of the 2010 Mentawai, Indonesia Earthquake Inferred from Tsunami Field Survey and Waveform Modeling, Pure Appl. Geophys., 170, 1567–1582, https://doi.org/10.1007/s00024-012-0536-y, 2013.

Sella, G. F., Dixon, T. H., and Mao, A.: REVEL: A model for Recent plate velocities from space geodesy, J. Geophys. Res.-Sol. Ea., 107, ETG 11-11-ETG 11-30, https://doi.org/10.1029/2000JB000033, 2002.

Senior Seismic Hazard Analysis Committee (SSHAC): Recommendations for probabilistic seismic hazard analysis: guidance on uncertainty and use of experts, US Nuclear Regulatory Commission Washington, DC, 1997.

Seno, T., Stein, S., and Gripp, A. E.: A model for the motion of the Philippine Sea Plate consistent with NUVEL-1 and geological data, J. Geophys. Res.-Sol. Ea., 98, 17941–17948, https://doi.org/10.1029/93JB00782, 1993.

Shao, G., Li, X., Ji, C., and Maeda, T.: Focal mechanism and slip history of
the 2011 *M*_{w} 9.1 off the Pacific coast of Tohoku Earthquake, constrained
with teleseismic body and surface waves, Earth Planets Space, 63, 559–564,
https://doi.org/10.5047/eps.2011.06.028, 2011.

Soloviev, S. L.: Recurrence of tsunamis in the Pacific, in: Tsunamis in the Pacific Ocean, edited by: Adams, W. M., East-West Center Press, 149–163, 1969.

Theunissen, T., Font, Y., Lallemand, S., and Liang, W.-T.: The largest instrumentally recorded earthquake in Taiwan: revised location and magnitude, and tectonic significance of the 1920 event, Geophys. J. Int., 183, 1119–1133, https://doi.org/10.1111/j.1365-246X.2010.04813.x, 2010.

Tsai, C.-C. P.: Slip, Stress Drop and Ground Motion of Earthquakes: A View from the Perspective of Fractional Brownian Motion, Pure Appl. Geophys., 149, 689–706, https://doi.org/10.1007/s000240050047, 1997.

Tsai, Y.-B.: A study of disastrous earthquakes in Taiwan, 1683–1895, Bull. Inst. Earth Sci., Acad. Sin., 5, 1–44, 1985.

Wang, X.: User manual for COMCOT version 1.7 (first draft), Cornel University, 65, 2009.

Wang, X. M. and Liu, P. L. F.: A numerical investigation of Boumerdes-Zemmouri (Algeria) earthquake and tsunami, CMES Comput. Model. Eng. Sci., 10, 171–183, https://doi.org/10.3970/cmes.2005.010.171, 2005.

Wang, X. and Liu, P. L. F.: An analysis of 2004 Sumatra earthquake fault plane mechanisms and Indian Ocean tsunami, J. Hydraul. Res., 44, 147–154, https://doi.org/10.1080/00221686.2006.9521671, 2006.

Wang, X. and Liu, P. L. F.: Numerical simulations of the 2004 Indian Ocean tsunamis – coastal effects, J. Earthquake and Tsunami, 01, 273–297, https://doi.org/10.1142/s179343110700016x, 2007.

Wang, Y.-J., Chan, C.-H., Lee, Y.-T., Ma, K.-F., Shyu, J., Rau, R.-J., and Cheng, C.-T.: Probabilistic seismic hazard assessment for Taiwan, Terr. Atmos. Ocean. Sci., 27, 325–340, https://doi.org/10.3319/TAO.2016.05.03.01(TEM), 2016.

Wu, T.-R., Chen, P.-F., Tsai, W.-T., and Chen, G.-Y.: Numerical Study on Tsunamis Excited by 2006 Pingtung Earthquake Doublet, Terr. Atmos. Ocean. Sci., 19, 705–715, https://doi.org/10.3319/TAO.2008.19.6.705(PT), 2008.

Wu, Y.-H., Chen, C.-C., Turcotte, D. L., and Rundle, J. B.: Quantifying the seismicity on Taiwan, Geophys. J. Int., 194, 465–469, https://doi.org/10.1093/gji/ggt101, 2013.

Yu, N.-T., Yen, J.-Y., Chen, W.-S., Yen, I. C., and Liu, J.-H.: Geological records of western Pacific tsunamis in northern Taiwan: AD 1867 and earlier event deposits, Mar. Geol., 372, 1–16, https://doi.org/10.1016/j.margeo.2015.11.010, 2016.

Yu, S.-B., Chen, H.-Y., and Kuo, L.-C.: Velocity field of GPS stations in the Taiwan area, Tectonophysics, 274, 41–59, https://doi.org/10.1016/S0040-1951(96)00297-1, 1997.

Yue, H. and Lay, T.: Inversion of high-rate (1 sps) GPS data for rupture
process of the 11 March 2011 Tohoku earthquake (*M*_{w} 9.1), Geophys. Res.
Lett., 38, L00G09, https://doi.org/10.1029/2011GL048700, 2011.

Short summary

The maximum possible earthquake magnitude is *M*_{w} 8.15 with an average slip of 8.25 m in the southernmost portion of the Ryukyu Trench. One hundred slip distributions of the seismic rupture surface were generated by a stochastic slip model. The simulated results demonstrate that the complexity of the rupture plane has a significant influence on the near field for local tsunamis. The propagation of tsunami waves and the peak wave heights largely vary in response to the slip distribution.

The maximum possible earthquake magnitude is *M*_{w} 8.15 with an average slip of 8.25 m in the...

Natural Hazards and Earth System Sciences

An interactive open-access journal of the European Geosciences Union