Open access peer-reviewed chapter

# Tsunami Propagation from the Open Sea to the Coast

Written By

Hiroshi Kanayama and Hiroshi Dan

Submitted: 09 December 2015 Reviewed: 20 April 2016 Published: 12 October 2016

DOI: 10.5772/63814

From the Edited Volume

## Tsunami

Chapter metrics overview

View Full Metrics

## Abstract

The tsunami generated by the Great East Japan Earthquake (GEJE) caused serious damage to the coastal areas of the Tohoku district. Numerical simulations are used to predict damage caused by tsunamis. Shallow-water equations are generally used in numerical simulations of tsunami propagation from the open sea to the coast. This research focuses on viscous shallow-water equations and attempts to generate a computational method using finite-element techniques based on the previous investigations of Kanayama and Ohtsuka (1978). First, the viscous shallow-water equation system is derived from the Navier-Stokes equations, based on the assumption of hydrostatic pressure in the direction of gravity. The derived equations have the horizontal viscosity term hereditary from the original Navier-Stokes equations. Next, a numerical finite element scheme is shown. Finally, tsunami simulations of Tohoku-Oki are shown using the above-mentioned approach. Our main concern in this chapter is how to set the boundary condition on the open boundary.

### Keywords

• Tsunami
• the Great East Japan Earthquake (GEJE)
• Viscous shallow-water equations
• Navier-Stokes equations
• Finite-element method

## 1. Introduction

The coastal areas of Tohoku district suffered serious damage from the tsunami caused by the 2011 off the Pacific Coast of Tohoku Earthquake that occurred on March 11, 2011 [1, 2]. Numerical simulations are used to develop disaster-prevention measures to deal with such tsunami disasters. They are also used to predict potential future tsunami disasters, to design disaster-prevention facilities such as coastal breakwaters and levees and to predict tsunami attacks immediately after an earthquake occurs [3, 4]. The Central Disaster Prevention Council [5] prepares a basic disaster prevention plan and participates in the determination of important disaster-prevention matters. If the Tonankai-Nankai earthquake occurs, it may cause considerable damage. Consequently, the council has been performing numerical calculations to predict wave height and arrival time when tsunami reaches the coast.

In this chapter, to let it be self-contained, the viscous shallow-water equations are again derived from the Navier-Stokes equations in which the hydrostatic pressure in the direction of gravity is assumed. In the numerical analysis of tsunamis, the viscosity term is often omitted or simply added [6, 7]. In this study, however, a computational method that does not omit the viscosity term is adopted, that is, because it enables more rigorous analysis to be performed and we intend to include the viscosity term in future stress analysis of tsunamis. The approximation scheme [8] is given below and simulation results [9] for a Tohoku-Oki model are again presented for self-containing. Our main concern in this chapter is how to set the boundary condition on the open boundary, which is completely changed to a new one.

## 2. Viscous shallow-water equations

The tsunami generated by the GEJE caused serious damage to the coastal areas of the Tohoku district. Numerical simulations are used to predict damage caused by tsunamis. Shallow-water equations are generally used in numerical simulations of tsunami propagation from the open sea to the coast. This subsection focuses on viscous shallow-water equations and attempts to generate a computational method using finite-element techniques based on the previous investigations of Kanayama and Ohtsuka [8].

In the numerical analysis of tsunami, a viscosity term is often omitted or simply added [6, 7]. In this subsection, however, a computational method that does not omit the viscosity term is adopted, that is, because it enables more rigorous analysis to be performed, and the authors intend to include the viscosity term in future stress analysis of tsunami. Recently, the authors have found an interesting book [10] that deals with sound mathematical topics related to the viscous shallow-water equations.

### 2.1. Derivation of the viscous shallow-water equations

As shown in Figure 1, Kanayama and Dan [9] considered the shallow-water long-wave flow in which the wavelength is sufficiently long relative to the water depth [11]. First, the viscous shallow-water equations are again derived to let this chapter be self-contained. Detailed derivation comes from the study of Kanayama and Ushijima [11]. Orthogonal coordinates [m] are used, where x1 and x2 represent directions in the horizontal plane and x3 represents the vertical direction, and the time [s] is represented by t. Assuming the hydrostatic pressure in the x3 direction, the following Navier-Stokes equations are used with the external forces, namely the Coriolis forces and the gravity are assumed to act in the x1 and x2 directions and the x3 direction, respectively.

j=13ujxj=0,E1
u1t+j=13uju1xj=1ρpx1+1ρj=13τ1jxj+fu2,E2
u2t+j=13uju2xj=1ρpx2+1ρj=13τ2jxjfu1,E3
0=1ρpx3g.E4

u(x1, x2, x3, t) represents the fluid velocity component [m/s] in the xi (i=1−3) direction, p(x1, x2, x3, t) denotes the pressure [N/m2], ρ is the density [kg/m3], τij is the stress component [N/m2] in the xi direction acting on the xj plane, f is the Coriolis coefficient [1/s], and g is the acceleration [m/s2] due to the gravity. In addition, the ordinate [m] of the water surface (water level) is represented by ζ and the ordinate [m] of the bottom surface is represented by h. Note that h usually has negative input values depending on xi (i = 1–2).

Eqs. (1)(3) [9] are integrated with respect to x3 from the bottom h to the water surface ζ. The velocity component in the normal direction is assumed to be zero at the water and bottom surfaces. By doing this, the viscous shallow-water equations expressed by the water level and the averaged velocity can be derived. If the layer thickness is H(x1, x2, t) and the average velocity in the xi (i = 1, 2) direction is Ui (x1, x2, t), then H and Ui are given by the following equations:

H=ζh,Ui=1Hhζuidx3.E5

Note that H is assumed to be positive in this chapter, whose guarantee is a difficult mathematical problem [9].

The fixed density of the layer is represented by ρ and the stress component in the xi direction acting on the xj plane is represented by τ˜ij. The horizontal viscosity constant [Ns/m2] is represented by μH, the wind effect coefficient is represented by θ, the air density is represented by ρa, the wind velocity component speed in the xi direction is represented by Wi, and the Chezy coefficient [m1/2/s] is represented by C. Then, the following viscous shallow-water equations are derived [9]:

ζt+j=12xj(HUj)=0,E6
H(Uit+j=12UjUixj)=gHζxi+1ρj=12xj(Hτ˜ij)+1ρθρaWi((W1)2+(W2)2)12gC2Ui((U1)2+(U2)2)12+(1)i+1fHUi+1,E7

where

τ˜ij=μH(Uixj+Ujxi),U3=U1. E8

Note that τ˜ij /2 is often used for the definition of τ˜ij. Here, however, we follow the definition in [11] for convenience.

Finite-element approximation is performed for the viscous shallow-water Eqs. (6) and (7) with given initial conditions and boundary conditions.

The computational domain is the two-dimensional (2D) polygonal region Ω surrounded by the boundaries Γc and Γo. In this region, the orthogonal coordinates x = (x1, x2) are used. The boundary conditions on Γc and Γo and the initial conditions at t = 0 are as follows [9]:

Boundary conditions, Eqs. (9) and (10):

Ui(x,t)=0onΓcE9
ζ(x,t)=ζΓo(x,t)j=12τ˜ijnj=0 }onΓo.E10

Here, ni is the component in the xi direction of the unit outward normal vector on the boundary. On Γo, ζΓo(x,t) is specified. Later, we improve Eq. (10).

Initial conditions:

Ui(x,0)=Ui0(x),ζ(x,0)=ζ0(x).E11

Here, Ui0(x) and ζ0(x) are the initial values of Ui and ζ, respectively.

### 2.2. Finite-element approximations for the viscous shallow-water equations

The finite-element approximation [8] is as follows. First, Eqs. (6) and (7) are multiplied by test functions and then integrated over the computational domain Ω. Subsequently, the approximation is carried out for terms containing a derivative with respect to time by using an explicit method. For terms containing spatial derivatives, the finite-element approximation is carried out by using the piecewise linear basis function ϕ^k. For terms without spatial derivatives, the approximation is carried out by using the corresponding step function ϕ¯k. The step function ϕ¯k is 1 in the barycentric region around the node k and is 0 at other places.

(ζ¯n+1ζ¯nΔt,ϕ¯k)+(j=12(H^nU^jn)xj,ϕ^k)=0.E12
(H¯n+1(U¯in+1U¯inΔt),ϕ¯k)+(H^n+1j=12U^jnU^inxj,ϕ^k)=(gH^n+1ζ^n+1xi,ϕ^k)j=12(1ρH^n+1μH(U^inxj+U^jnxi),ϕ^kxj)+(1ρθρaWi((W1)2+(W2)2)12,ϕ¯k)(gC2U¯in((U¯1n)2+(U¯2n)2)12,ϕ¯k)+((1)i+1fH¯n+1U¯i+1n,ϕ¯k).E13

In the above, ζkn+1 is determined from Eq. (12). From the obtained value and Ui,kn, the value Ui,kn+1 is determined by Eq. (13). Here, ζkn and Ui,kn are approximate values of ζ(x,t) and U(x,t) respectively, at the node k after n time steps, and Δt represents the size of time steps. In addition, the following symbols are also used:

(Uixj,ν)=ΩUixjνdx,E14
U^in=kUi,knϕ^k,U¯in=kUi,knϕ¯k.E15

It is well known that integration by parts is used for the viscosity term in (13).

Because of the presence of the basic boundary conditions, Eq. (12) holds for all nodes except for the nodes on Γo and Eq. (13) holds for all the nodes except for the nodes on Γc. However, on the boundary Γo, an approximation method is adopted in which the advective term of Eq. (13) uses Tabata’s upwind approximation [9, 12]. A mathematical justification for the approximation scheme for linearized equations related to the above scheme was given by Kanayama and Ushijima [13, 14].

In general, tsunami is excited in the following two ways. The first one is to consider the tsunami excitation as the initial condition of the water surface, for which we do not have sufficient input information in such an artificial tsunami of Hakata Bay [9]. The second one is to consider it as the boundary condition of the water surface as in the next subsection. In the setting of Hakata Bay, a computational domain is not so wide that the above approach may be the only way. It is also noted that 50 [m] at the open boundary for the later Tohoku-Oki case may be too high. In the computation, the tsunami arrived at Oshika Peninsula after 20 min, and the highest wave height reached 15 [m]. These numerical results should be checked more carefully with data on the open boundary.

## 3. Computational examples

### 3.1. Tohoku-Oki models

Now, numerical results [9] for a Tohoku-Oki model are again shown for self-containing. When not specifically defined, the same physical values as for Hakata Bay were used. Regarding the mesh data, the total number of nodes is 56,562, the total number of elements is 113,013, and the ocean floor ordinate is set to become deep gradually from h = ‒10 to h = ‒1,000 [m]. We set the initial conditions at ζ0 (x) = 0 and Ui0 (x) = 0, and the boundary conditions ζΓo(x,t)=50 [m] (0 < t < 60 [s]) at the tsunami-generating area and ζΓo(x,t)=0 [m] at other areas (see Figure 2). It is noted that 50 [m] at Γo may be too high. In our computation, the tsunami arrived at Oshika Peninsula in Figure 2 after 20 min and the highest wave height reached 15 [m]. These numerical results should more carefully be checked with boundary data on Γo.

Figure 2a shows the ordinates of the water surface (water level) after 18 [min]. Since the computational domain is not wide, it looks that there is a reflection from the northern boundary in Figure 2a. This artificial reflection can be removed by suitable boundary conditions on Γo. Details are later mentioned. Figure 2b shows the water level change at the two points A and B in Figure 2. The wave height at the point B is higher than the point A after about 1500 [s]. When the tsunami reaches coastal areas, the water depth becomes shallow and the wave height becomes high.

Next, new numerical results for the second Tohoku-Oki model are shown. When not specifically defined, the same physical values as for Hakata Bay were used. Regarding the mesh data, the total number of nodes is 563,100, the total number of elements is 1,123,178, and Figure 3 shows the ocean floor ordinate, which is constructed from 30-arc sec interval grid of JTOPO30, provided by the Marine Information Research Center. We set the initial conditions at ζ0 (x) like Figure 4 based on the data of Fujii et al. [15] as a source of tsunami. Figure 5a–f shows the ordinates of the water surface (water level) every 8 [min] from after 8 [min] to after 48 [min]. In our previous paper [9], since the computational domain is not wide, it looks that there is a reflection from the boundary Γo (see Figure 2a). This artificial reflection can be removed by changing boundary conditions from (10) to (16) as follows:

Un(x,t)=cg| h |Hζ(x,t)Ut(x,t)=0 }onΓo,E16

where Un and Ut denote the normal component of velocity and the tangential component of velocity, respectively, and c is a constant. We used (c = 0.9) in this chapter.

Figure 6 shows the water level change at the three points Miyako, Soma and Choshi in Figure 4. Solid curves show observed data [15] and dashed curves show numerical results. Though the arrival time of the first wave is almost same, the wave height is small compared with observed data at the points of Miyako and Soma. Results may be improved by changing initial and boundary conditions.

## 4. Concluding remarks

In this study, the viscous shallow-water equations have again been derived from the Navier-Stokes equations in which the hydrostatic pressure in the direction of gravity is assumed. Tsunami propagation is then simulated by a finite-element computation. In the study of Kanayama and Dan [9] using the Hakata Bay model, in which the tsunami-generating area was taken to be the epicenter of the West off Fukuoka Prefecture Earthquake (2005), tsunami propagation from the open sea to the coastal area might be produced. This is an example in which the conventional open boundary condition Eq. (10) is still valid. In this chapter, however, in the second Tohoku-Oki model, the open boundary condition is changed from Eq. (10) to Eq. (16) to overcome artificial reflection. In addition, we have reconstructed a system in which the tsunami arrival time and the wave height at the time of tsunami attack can be obtained from numerical results. Since our analysis takes the viscosity term into account, this study can also be considered to be a preliminary study for future planned stress analyses of tsunamis (see Figure 7). We still want to check whether the hydrostatic pressure is always the dominant force in the case of tsunami attack compared with other forces including the viscous force.

## Acknowledgments

The first author thanks the JKA RING!RING! Project for the financial support during 2015FY.

Parts of this chapter are reproduced from [9], © 2013 The Author(s); originally published under CC BY license.

## References

1. 1. National Institute for Land and Infrastructure Management, Ministry of Land, Infrastructure, Transport and Tourism, Japan and Building Research Institute, Incorporated Administrative Agency, Japan: Quick report of the field survey and research on the 2011 off the Pacific Coast of Tohoku Earthquake. 2011; Technical Note 636: Building Research Data 132 (in Japanese)
2. 2. Yoshimura S et al.: Understanding tsunami-causing damage mechanism of machines and structures based on a discipline of mechanics. In: Report of JSME Research Committee on the Great East Japan Earthquake Disaster. 2014. Chapter 4.
3. 3. Goto C, Sato K: Development of tsunami numerical simulation system for Sanriku coast in Japan. Report of the Port and Harbour Research Institute. 1993;32(2):3–44 (in Japanese)
4. 4. Takahashi T: Application of numerical simulation to tsunami disaster prevention. Journal of Japan Society of Computational Fluid Dynamics. 2004;12(2):23–32 (in Japanese)
5. 5. The Central Disaster Prevention Council. Special Investigation Council concerning Tounankai and Nankai Earthquakes (in Japanese) http://www.bousai.go.jp/kaigirep/chuobou/senmon/tounankai_nankaijishin/index_nankai.html [Accessed:03-06-2016]
6. 6. Ulutas E: The 2011 off the Pacific Coast of Tohoku-Oki Earthquake and tsunami: influence of the source characteristics on the maximum tsunami heights. In: Proceedings of the International Symposium on Engineering Lessons Learned from the 2011 Great East Japan Earthquake. 2012. pp. 602–611.
7. 7. Imamura F, Yalciner AC, Ozyurt G. Tsunami modelling manual. 2006. http://www.tsunami.civil.tohoku.ac.jp/hokusai3/J/projects/manual-ver-3.1.pdf [Accessed: 27-03-2012]
8. 8. Kanayama H, Ohtsuka K: Finite element analysis on the tidal current and COD distribution in Mikawa Bay. Coastal Engineering in Japan. 1978;21:157–171.
9. 9. Kanayama H, Dan H: A tsunami simulation of Hakata Bay using the viscous shallow-water equations. Japan Journal of Industrial and Applied Mathematics. 2013;30(3):605–624.
10. 10. Bresch D. Shallow-water equations and related topics. In: Dafermos C, Pokorny M, editors. Handbook of Differential Equations: Evolutionary Equations. 5: 2009. pp. 1–104. Elsevier doi:10.1016/S1874-5717(08)00208-9
11. 11. Kanayama H, Ushijima T: On the viscous shallow-water equations I – derivation and conservation laws. Memoirs of Numerical Mathematics. 1981/1982;8/9:39–64.
12. 12. Tabata T: A finite element approximation corresponding to the upwind finite differencing. Memoirs of Numerical Mathematics. 1977;4:47–63.
13. 13. Kanayama H, Ushijima T: On the viscous shallow-water equations II – a linearized system. Bulletin of University of Electro-Communications. 1988;1(2):347–355.
14. 14. Kanayama H, Ushijima T: On the viscous shallow-water equations III – a finite element scheme. Bulletin of University of Electro-Communications. 1989;2(1):47–62.
15. 15. Fujii Y, Satake K, Sakai S, Shinohara M, Kanazawa T: Tsunami source of the 2011 off the pacific coast of Tohoku earthquake. Earth Planets Space. 2011;63:815–820.

Written By

Hiroshi Kanayama and Hiroshi Dan

Submitted: 09 December 2015 Reviewed: 20 April 2016 Published: 12 October 2016