## Abstract

The new mathematical model for tsunami evolution by Tobias & Stiassnie (Tobias & Stiassnie 2011 *J. Geophys. Res. Oceans* **116**, C06026) is used to derive a synthetic tsunami database for the southern part of the Eastern Mediterranean coast. Information about coastal tsunami amplitudes, half-periods, currents and inundation levels is presented.

## 1. Tsunami sources and coastal target sites

The main tsunami sources that threaten the Israeli coast are earthquakes at the Crete–Hellenic arc and Cypriot arc ruptures. Tsunami parameters at 17 coastal sites from Beirut in the north to Alexandria in the south (including three sites in Lebanon, seven in Israel, two in Gaza and five in Egypt) are considered. The tsunami sources and the coastal locations are shown in figure 1.

## 2. The idealized model for tsunami study

The main features of the idealized model for tsunami evolution are given below. The complete derivation and its theoretical consequences can be found in Tobias & Stiassnie [1].

The initial disturbance of the free surface, *η*_{d}, at time *t*=0, is given by
2.1where *η*_{0} and *a* are the initial maximal amplitude and lateral extent of the disturbance, respectively, and *r* is the radial distance from the source of the initial disturbance (located above the earthquake epicentre).

For an ocean of constant depth *H*, the evolution of the flow field (2.1) in time is well known, and is given in Leblond & Mysak [2], §50 by
2.2In (2.2), *J*_{0} is the zero-order Bessel function of the first kind, and the frequency *Ω* is related to the wavenumber *k* through the linear dispersion relationship:
2.3It is assumed that the open ocean extends from the source of the tsunami at *r*=0 to the edge of the continental shelf at *r*=*L*.

The depth *h* of the continental shelf is assumed to vary linearly according to
2.4The evolution of the tsunami on the continental shelf is modelled by the linear shallow-water equation (see Dean & Dalrymple [3]):
2.5In (2.5), the water depth *h* is given by (2.4), *η*_{s} is the free-surface elevation on the continental shelf, and *A* is a friction coefficient.

The solution of (2.5) is required to be finite at the shoreline, *r*=*L*+*l*, for all time, and to be compatible with the deep-sea solution (2.2) at the edge of the continental-shelf:
2.6aIn (2.6a), *H*(*t*) is the Heaviside step function, and *κ*_{t} is the transmission coefficient of long waves over an abrupt depth change, from a depth of *H* to a depth of *h*_{0}, given by
2.6b(see Dean & Dalrymple [3]).

The solution for the continental-shelf region is obtained as follows:

— A Fourier transform from the time domain to the frequency domain is applied to equations (2.2) and (2.5). The Fourier transform leads to a Bessel equation in

*r*, which has to be solved in [*L*,*L*+*l*].— One of the Bessel functions (

*Y*_{0}) is omitted owing to the requirement of boundness at*r*=*L*+*l*. The coefficient in front of the other Bessel function (*J*_{0}) is calculated using the transformed boundary condition, (2.2) and (2.6a).— The inverse Fourier transform, from the frequency domain to the time domain, is applied. The integral over frequency is calculated using the Cauchy residue theorem, to obtain for

*t*≥0: 2.7awhere 2.7band 2.7cwith 2.7dwhere*J*_{1}is the first-order Bessel function of the first kind, and*α*_{j}are the positive zeros of*J*_{0}, i.e.*J*_{0}(*α*_{j})=0,*α*_{j}>0.

The integration over *k* in either (2.2) or (2.7a) is carried out numerically using Matlab. The routine used, QUADGK, is based on the adaptive Gauss–Kronrod quadrature method and is recommended for calculation of oscillatory integrals.

Numerical experimentation revealed that setting the upper bound of the sum in (2.7a) to 100 provides a stable result.

The energy loss over the continental-shelf region has two main components: (i) loss due to bottom friction and (ii) loss due to energy leaking back into the deep sea.

It turns out that the second component is significantly larger than the first. Using the reflection coefficient, given in Dean & Dalrymple [3],
2.8and comparing it to the wave damping of progressive waves, due to the effect of an ‘artificial’ bottom friction, over a distance of 2*l*, and simplifying, gives
2.9The initial amplitude *η*_{0} (*m*) and lateral extent *a* (km) are calculated using equations for slip displacement *D* (m), and rupture area *A* (km^{2}) given by Wells & Coppersmith [4]:
2.10aand
2.10bwhere *M* is the moment magnitude of the earthquake.

The following assumptions give the relationship between the earthquake parameters given by Wells & Coppersmith [4] and the tsunami parameters required for the model calculation:
2.11aassuming an average dip angle of about 45^{°}, and
2.11bwhich is consistent with the lower bound of equation (2.3) in Levin & Nosov [5].

The idealized model provides an analytical solution, based on five input parameters: (i) *M*—the moment magnitude of the earthquake; (ii) *L*—the distance of the target site from the epicentre; (iii) *H*—the average depth of the deep sea; (iv) *h*_{0}—the depth of the continental shelf at its edge; and (v) *l*—the width of the continental shelf.

Various numerical models using different approximations of the full hydrodynamic problem have been used to tackle the tsunami phenomenon. For details, the reader is referred to the recent book of Levin & Nosov [5].

Two recent papers, Madsen *et al*. [6] and Constantin & Henry [7], which have questioned the relevance of the solitary wave paradigm for tsunamis, have encouraged us to investigate the quality of the results produced by linearized tsunami models (see also Hammack & Segur [8]).

A few rather encouraging comparisons with field measurements, during the 2004 Sumatra–Andaman tsunami, are reported in Tobias & Stiassnie [1]. In this paper, we concentrate on the Eastern Mediterranean, and compare our findings with results of more elaborate numerical models.

A rather reassuring quality check, in which the results of the idealized tsunami model are compared with those of Salamon *et al*. [9], Thio [10] and Hamouda [11], is provided in figure 2.

## 3. Tsunami calculations for the Israeli coast

The idealized model for tsunami study by Tobias & Stiassnie [1] is used to calculate the tsunami amplitudes *η* and tsunami half-periods *T* for a wide range of earthquake magnitudes *M*. These results, for *M*=6.5,6.6,…,8.2 in Cyprus, and *M*=7.2,7.3,…,8.9 in Crete, are presented in tables 1 and 2.

Note that the tsunami amplitude is defined as the extreme deviation from the mean sea surface (either elevation or depression); and that the tsunami half-period is defined as the time difference between the arrival of the above-mentioned extreme, either crest or trough, and the arrival of the following trough or crest, respectively.

The tsunami speed at the coast, *u*, can be estimated from
3.1An estimate for inundation level above the mean sea level, *H*, can be calculated from
3.2

## 4. Long-term tsunami estimates

In order to assess the tsunami amplitudes for different return periods, one needs some information regarding the return periods of the magnitudes of the earthquakes. From table 4 in Thio [10], the return periods for *M*=7 and 8 at Cyprus are 1000 and 10 000 years, respectively. Based on the upper left entry of fig. 6 in Tsapanos [12], the return periods for *M*=7.9 and 8.3 at Crete are also 1000 and 10 000 years, respectively.

Assuming a first-type Gumble distribution for both the earthquake magnitudes and the coastal tsunami amplitudes results in the following equation for the tsunami amplitude *η*_{i} (m), at the *i*=1,2,…,17 location, as a function of the return period *t*_{R} (year):
4.1(see Gumbel [13]).

The values of *μ*_{i} (m) and *β*_{i} (m) are given in table 3.

A bar-diagram presentation of the coastal tsunami amplitudes for the return periods *t*_{R}=500, 1000, 2000, 5000, 10 000 and 20 000 years is given in figure 3.

From figure 3, one can see that the most sensitive locations to tsunami attacks, sorted by country, are Sour in Lebanon, Haifa in Israel and Port Said in Egypt, for which the obtained tsunami amplitudes *η* with a return period of 5000 years are 1.5, 2.2 and 4.5 m, respectively. The tsunami half-periods *T*, the tsunami speed, as well as the inundation levels *H*, for these cases are given in table 4.

## 5. Concluding remarks

We have decided to use the term ‘synthetic tsunamis’ in order to stress the fact that our calculations are based on a host of simplifying assumptions. The most severe of them all seem to be the data of the earthquake magnitudes with return periods of 1000 and 10 000 years.

## Acknowledgements

This research was supported by the Israel Science Foundation (Grant 63/09), and by the Israel Ports Company. The research is part of a MSc thesis submitted by Joshua Tobias to the Graduate School of the Technion–Israel Institute of Technology.

## Footnotes

One contribution of 13 to a Theme Issue ‘Nonlinear water waves’.

- This journal is © 2012 The Royal Society