Rational Evaluation Methods of Topographical Change and Building Destruction in the Inundation Area by a Huge Tsunami

Abstract

In the case of huge tsunamis, such as the 2004 Great Indian Ocean Tsunami and 2011 Great East Japan Tsunami, the damage caused by ground scour is serious. Therefore, it is important to improve prediction models for the topographical change of huge tsunamis. For general models that predict topographical change, the flow velocity distribution of a flood region is calculated by a numerical model based on a nonlinear long wave theory, and the distribution of bed-load rates is calculated using this velocity distribution and an equation for evaluating bed-load rates. This bed-load rate equation usually has a coefficient that can be decided using verification simulations. For the purpose, Ribberink’s formula has high reproducibility within an oscillating flow and was chosen by the authors. Ribberink’s formula needs a bed-load transport coefficient that requires sufficient verification simulations, as it consumes plenty of time and money to decide its value. Therefore, the authors generated diagrams that can define the suitable bed-load coefficient simply using the data acquired from hydraulic experiments on a movable bed. Subsequently, for the verification purpose of the model, the authors performed reproduced simulations of topography changes caused by the 2011 Great East Japan Tsunami at some coasts in Northern Japan using suitable coefficients acquired from the generated diagrams. The results of the simulations were in an acceptable range. The authors presented the preliminary generated diagrams of the same methodology but with insubstantial experimental data at the time at the International Society of Offshore and Polar Engineers (ISOPE), (2018 and 2019). However, in this paper, an adequate amount of data was added to the developed diagrams based on many hydraulic experiments to further raise their reliability and their application extent. Furthermore, by reproducing the tsunami simulation on the Sendai Natori coast of Japan, the authors determined that the impact of total bed-load transport was much bigger than that of suspension loads. Besides, the simulation outputs revealed that the mitigation effect of the cemented sand and gravel (CSG) banks and artificial refuge hills reduced tsunami damage on Japan’s Hamamatsu coast. Since a lot of buildings and structures in the inundation area can be destroyed by tsunamis, building destruction design was presented in this paper through an economy and simplified state. Using the proposed tsunami simulation model, we acquired the inundation depth at any specific time and location within the inundated area. Because the inundation breadth due to a huge tsunami can extend kilometers toward the inland area, the evaluation of building destruction is an important measure to consider. Therefore, the authors in this paper presented useful threshold diagrams to evaluate building destruction with an easy and cost-efficient state. The threshold diagrams of “width of a pillar” for buildings or “width of concrete block walls” not breaking to each inundation height were developed using the data of damages due to the 2011 Great East Japan Tsunami.

1. Introduction

Research for predicting coastal scour and erosion by waves or flow has been performed all over the world. In particular, there has been a series of research regarding the prediction of the coastal scour and erosion within a wide area due to tsunamis. Takahashi et al. (1993a, 1993b) [1,2] proposed a formula to estimate bed-load transport to predict scour due to tsunamis. They assumed that the bed-load transport was proportional to power n of the Shields number. However, they also stated that it was not possible to neglect the suspended load transport. Kobayashi et al. (1996) [3] proposed a formula to estimate bed-load transport with the power 1.5 of Shields number and found a good agreement with experimental results. Fujii et al. (1998) [4] and Takahashi et al. (2000) [5] considered suspended load entrainment and deposition that improved the prediction accuracy of the bottom topography change due to tsunamis. Nishihata et al. (2006) [6] considered both bed-load and suspended load when developing a numerical model to estimate sediment deposit flux in Sri Lanka’s Kirinda fishing harbor. Kihara and Matsuyama (2007) [7] developed a three-dimensional (3D) model to compute bottom topography changes caused by tsunamis. Nakamura and Mizutani (2008) [8] proposed a formula to estimate bed-load and suspended-load transport that accounts for the fluctuation of underground stress. Thus, improvement in the prediction model’s accuracy has been advanced. However, the development of a model by which many coastal engineers can easily predict practical accuracy is also desired. Therefore, Ca, Yamamoto, and Charusrojthanadech (2010) [9] presented a two-dimensional numerical simulation model that calculated inundation velocity, inundation depth, and topographical change on an inundation area. For hydraulic calculations, a continuity equation and two-dimensional nonlinear shallow water equations were used. To calculate topographical changes, this model used Ribberink’s equation. This equation requires a bed-load transport coefficient that is usually decided after using verification simulations. To decrease the time and cost needed to decide the value of Ribberink’s bed-load transport coefficient, Ahmadi, Yamamoto, and Hayakawa (2018 and 2019) [10,11] performed many hydraulic experiments and developed useful diagrams to obtain the bed-load coefficient by using the inverse analysis of the experiment results. However, at the time, the amount of experiment data was inadequate. Therefore, in this research, more experimental data was added into the developed diagrams to decide if the bed-load transport coefficient further raised reliability. By reproducing the tsunami simulation on Japan’s Sendai Natori coast, the suspension load showed a smaller effect on topography change as the total bed-load. Moreover, as a mitigation effect to evaluate tsunami countermeasures, we executed prediction simulations for tsunami inundation and topographical change on Japan’s Hamamatsu coast. Since the local government constructed coastal banks of 13 m in height made of cemented sand and gravel (CSG) and some evacuation soil mounds, we examined the effect that CSG banks and artificial refuge hills had on reducing tsunami damage.

As a final result, this paper presented a rational method for predicting building destruction caused by tsunamis. Previously, Yamamoto et al. (2006) [12] showed that the stress analysis using the gate-type Rahmen model could accurately predict whether each building was broken by a tsunami or not. Our proposed hydro-morphodynamics model made it possible to calculate inundation depth and inundation velocity at any desired point and time in a proposed area. Yamamoto, Nunthawath, and Nariyoshi (2011) [13], as well as Charusrojthanadech, Yamamoto, and Ca (2011, 2012) [14,15], showed that the combination of tsunami simulations and diagrams made by stress analysis using the gate-type Rahmen mode easily predicted whether buildings were broken by a tsunami. The formal evaluation method of building destruction due to tsunamis are usually costly and time-consuming, while methods using developed diagrams are economical and simple. The threshold diagrams in their papers were developed using damage surveys from Japan’s 2004 Indian Ocean Tsunami and 1993 Hokkaido-Nansei-Oki Tsunami. Using these diagrams, we evaluated the limit width of a wall or pillar to be safe against a probable tsunami inundation depth. However, since the material quality of Japanese buildings broken by the 2011 Great East Japan Tsunami differed from the quality of buildings in Thailand and Sri Lanka broken by the 2004 Indian Ocean Tsunami, we developed new threshold diagrams to evaluate building destruction using destruction data of the 2011 Great East Japan Tsunami. As a result, we evaluated the limit width of a wall or pillars to be safe against tsunami inundation depth. The evaluation method using threshold diagrams was found to be easy, economical, and useful in countries like Japan.

2. Rational Method for Predicting Topographical Change by Tsunami

2.1. Existing Numerical Simulation Model

To predict topographical change caused by tsunamis, we used Ca et al.’s (2010) [9] numerical simulation model. The model was developed using the following detailed formulas.

2.1.1. Numerical Model for Fluid Motion

The numerical model used for tsunami simulations in this paper was based on a continuity equation of fluid (Equation (1)) and two-dimensional nonlinear long wave equations (Equations (2) and (3)). These governing equations were solved by finite difference methods using the Crank–Nicholson scheme.

$$\frac{\partial f_y{q}_x}{\partial x}+\frac{\partial f_x{q}_y}{\partial y}+\frac{\partial S\eta }{\partial t}=0$$

(1)

$$\frac{\partial q_x}{\partial t}+\frac{1}{S}\frac{\partial }{\partial x}\left(\frac{S{q}_x{}^2}{d}\right)+\frac{1}{S}\frac{\partial }{\partial y}\left(\frac{S{q}_x{q}_y}{d}\right)+gd\frac{\partial \eta }{\partial x}-\frac{1}{S}\frac{\partial }{\partial x}\left\{d{\nu }_{t}S\frac{\partial }{\partial x}\left(\frac{q_x}{d}\right)\right\}-\frac{1}{S}\frac{\partial }{\partial y}\left\{d{\nu }_{t}S\frac{\partial }{\partial y}\left(\frac{q_x}{d}\right)\right\}+\frac{f_c}{d^2}Q{q}_x=0$$

(2)

$$\frac{\partial q_y}{\partial t}+\frac{1}{S}\frac{\partial }{\partial x}\left(\frac{S{q}_y{q}_x}{d}\right)+\frac{1}{S}\frac{\partial }{\partial y}\left(\frac{S{q}_y{}^2}{d}\right)+gd\frac{\partial \eta }{\partial y}-\frac{1}{S}\frac{\partial }{\partial x}\left\{d{\nu }_{t}S\frac{\partial }{\partial x}\left(\frac{q_y}{d}\right)\right\}-\frac{1}{S}\frac{\partial }{\partial y}\left\{d{\nu }_{t}S\frac{\partial }{\partial y}\left(\frac{q_y}{d}\right)\right\}+\frac{f_c}{d^2}Q{q}_y=0$$

(3)

where $q_x$ and $q_y$ are the horizontal fluid fluxes in the x and y directions, respectively. $\eta$ is the water surface elevation.$f_{x }\mathrm{and} f_y$ are the x and y direction ratios of the wet portion in a calculation mesh. $S$ is the area ratio of the wet portion in a calculation mesh. $d$ is the water depth (from the static water surface+$\eta$). g is the gravitational acceleration. ${\nu }_r$ is the eddy viscosity coefficient. $f_c$ is the ground surface friction coefficient and $Q$($=\sqrt{q_x^2+q_y^2}$) is the compound value of $q_x$ and $q_y$.

To calculate the eddy viscosity coefficient and the ground surface friction coefficient, the following equations were used:

$${\nu }_t=\epsilon {\left[{\left(\frac{\partial U}{\partial y}\right)}^2+{\left(\frac{\partial V}{\partial x}\right)}^2\right]}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}·d^2$$

(4)

$$f_c=\frac{g{n}^2}{d^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}}, n^2=n_0{}^2+0.020\frac{B_r}{100-B_r}{d}^{\raisebox{1ex}{$4$}\!\left/ \!\raisebox{-1ex}{$3$}\right.}, n_0{}^2=\frac{n_1{}^2{A}_1+n_2{}^2{A}_2+n_3{}^2{A}_3}{A_1+A_2+A_3}$$

(5)

Here, $U$ and $V$ are the flow velocity in $y$ and $x$ directions. $\epsilon$ equals (0.1). $n$ is the Manning’s roughness coefficient. B_r is the building ratio (= the ratio of the area of all vertical objects like houses and to the mesh area). $n_0$ is the weighted average roughness coefficient of areas like farms, roads, and waste and wetlands, $A_{1, 2, and3}$, with relative roughness coefficients of $n_{1, 2, and3}$.

2.1.2. Numerical Model for Topographical Change

Topographical change based on sediment transport can be expressed using the continuity equation shown in Equation (5):

$$\frac{\partial \zeta }{\partial t}=-\frac{1}{1-{\epsilon }_s}\left(\frac{\partial q_{bx}}{\partial x}+\frac{\partial q_{by}}{\partial y}-C_s+C_{ut}\right)$$

(6)

where $\zeta$ is the ground surface elevation, $q_{bx}$ and $q_{by}$ are the bed-load rate per unit width in x and y directions, respectively, $C_s$ is the deposition rate of the suspended load, $C{ }_{ut}$ is the entrainment rate of the suspended load from the bed, and ${\epsilon }_s$ is the porosity of the sediment.

Modeling of $q_x$ and $q_y$ on Bed-load Transport

To evaluate the bed-load rate, Ribberink’s formula (1998) [16], shown in Equation (6), was used. Yokoyama et al. (2002) [17] performed many flow and wave scouring calculations using indoor and outdoor data on sand and gravel, with diameter ranges of 0.2–10 mm. They found that accurate results could be obtained using this formula.

$$q_{bi}=\{\begin{array}{ll}{C}_b{\left[\left|{\theta }_s\left(t\right)\right|-{\theta }_{sc}\right]}^{1.65}\frac{{\theta }_s\left(t\right)}{\left|{\theta }_s\left(t\right)\right|}\sqrt{\Delta g{(D_{50})}^3}& \left({\theta }_s\left(t\right)\ge {\theta }_{sc}\right)\\ 0& \left({\theta }_s\left(t\right)<{\theta }_{sc}\right)\end{array}$$

(7)

where $q_{bi}$ is the bed-load transport rate per unit width in the i direction, $C_b$ is the bed-load coefficient determined by verification simulations, ${\theta }_s\left(t\right)$ is the Shields parameter in the i direction, ${\theta }_{sc}$ is the critical Shields number, Δ is the relative density of the sand, g is the gravitational acceleration, and D₅₀ is the median diameter of the sediment.

Modeling of $C_s \mathrm{and} C_{ut}$ on Suspended Load Transport

During tsunamis, it is necessary to consider the influence of suspended load transport. The deposition rate of the suspended load C_s and the entrainment rate from the bed C_ut was evaluated using Equation (7) based on the vertical distribution of suspended load concentration:

$$C_s=w_{s}C\left(\frac{w_s}{2}\right),\hspace{1em}{C}_{ut}=-{\nu }_t\frac{\partial C(z)}{\partial z}|{}_{z=z_a},$$

(8)

where $w_s$ is the settling velocity of suspended particles. $C(z)$ is the suspended load concentration and ${\nu }_t$ is the eddy viscosity.

According to Soulsby (1997) [18], when assuming a sheet flow condition for a whole area, the vertical distribution of sediment concentration is expressed as follows:

$$C\left(z\right)=C_a{\left(\frac{z}{z_a}\right)}^b,\hspace{1em}b=-\frac{w_s}{\kappa u_{*}}$$

(9)

where κ is Karman’s constant (=0.4), $u_{*}$ is the friction velocity, $C_a$ and $z_a$ are estimated using the equation of Zyserman and Fredsoe (1994) [19].

2.2. Rational Evaluation Method of Bed-Load Transport Rate

To calculate the bed-load transport using Ribberink’s equation, we needed to perform many verification simulations to obtain a suitable bed-load transport coefficient. However, in this paper, after we performed many hydraulic experiments using inverse analysis, we developed useful diagrams to obtain the bed-load coefficients based on median grain size, uniformity coefficients, and dry soil density for a proposed area. The development of these diagrams saved time and money.

2.2.1. Hydraulic Model Experiment Method

Based on Iwanuma city’s (2015) [20] tsunami protection plan in the Miyagi Prefecture, the model consisted of a typical beach profile, coastal dike, and a general land-side ground. They were set with a scale of 1/20 in the water channel, 0.8 m in height, 19.0 m in length, and 0.5 m in width (as shown in Figure 1). Some tsunami cases were reproduced by collecting water in a tank connected to the left end of the channel and then discharging the water into the channel. The experimental apparatus and measurement methods are illustrated below.

Four cases were tested with water heights of 45 cm, 60 cm, and 65 cm, which produced maximum flow velocities of 0.9 m/s, 1.05 m/s, and 1.1 m/s on the crown of the coastal dike, respectively. They are tabulated in

Table 1. Parameters of scouring experiments on the dike’s land side.

Case	1			2		3		4
Max. Velocity on the Crown (m/s)	0.9 * (0.45)			1.1–1.25 * (0.65)		1.05–1.15 * (0.60)		1.0 * (0.45)
Type of the sediment	Clay	Sand	Soil	Sand	Soil	Sand·Gravel	Soil	Gravel	Soil **
Median grain size (mm)	0.005	0.2	0.2	0.2	0.2	0.1, 0.2, 10	0.2	1.5, 3.4, 5.0, 20, 30	0.2
Uniformity coefficient	--	1.56	20.1	1.56	20.1	1.5–2	11	1.5–3	20.1
Dry density (g/cm3)	Around 1.5							Around 1.5	1.6–2.0

*: There was a difference in the flow velocity due to the effect of slightly different embankment heights depending on the implementation time. **: The concentrations of suspension load rates are around 10% and less than it.

. The number of all experiments was 25. The sediment of 0.2 mm in median grain size and 20.1 in uniformity coefficient was selected as the basic soil because it was common in Japan. In cases 3 and 4, grain sizes up to 10 mm and 30 mm, respectively, were used so that experimental results could be used in areas with larger grain sizes. Here, natural ones (a little round) were used in cases where median grain sizes were smaller than 20 mm. On the other hand, round gravel and crushed gravel were used in the case that the median grain size was 20 mm or 30 mm. Flow velocity on the crown of the dike was measured with an electromagnetic flow velocity meter of disk type and flow velocity in the scour area on the land side of the dike, which were measured using a propeller meter. Then, since the sidewalls of the water channel are made of transparent acrylic boards stiffened with angle irons, the topographical change of the ground was observed from the side by using a video camera, as shown in Figure 2. Moreover, the experiments were performed by a unit that was two times the same case.

The flow velocity and the scour depth of the model experiments were converted to values in the field using Froude’s similarity law. Moreover, “a field friction coefficient/a model friction coefficient = the 3/4th–1st power of the model scale” was obtained from existing experimental data, and this condition was substituted to Shields’ similarity law, i.e., “a model grain size/a field grain size = the 1/4th–0th power of the model scale”. Therefore, we considered the model grain size in our experiments as the same as the field grain size. Indeed, the maximum velocities near the coastal dike in Iwanuma City were 4–7 m/s. On the other hand, the maximum velocities of the model were $÷\sqrt{1/20}=$ 4–5 m/s. The maximum scours depths near the coastal dike were 3–6 m. On the other hand, the maximum scour depths of the model were $÷1/20=$ 3–5 m.

The simulations used to create Figure 3, Figure 4 and Figure 5 are explained by the following three steps:

1)

We reproduced the topographical feature on the right of the section (A) of Figure 1a and inputted time series data of water level from the section (A) so that time-series data of the simulated flow velocity on the crown part of the coastal dike was in agreement with the measured data (refer to Figure 1b).

2)

We setup the first approximate value of the bed-load rate coefficient $C_b$, performed topographical change simulation, and calculated the maximum scour depth, thus forming the average value of scouring width by the dike’s landside. If these two calculated values were mostly in agreement with the measured values, we considered the value of assumed $C_b$ to be a true value.

3)

When these two calculated values were not in agreement with the measured values, we changed the value of $C_b$ and repeated the topographical change simulation until these two calculated values agreed.

2.2.2. Bed-load Coefficient by Inverse Analysis

The bed-load coefficients sufficiently reproduced the scouring depth of the executed hydraulic experiments and were obtained using the numerical model of Ca et al. (2010) [9]. Furthermore, the influence of median grain sizes, uniformity coefficients, and dry densities on bed-load coefficients were also examined. As shown in Figure 3, when median grain size became large (5 mm), the bed-load coefficient “$C_{b0}$” increased from zero to around 60; as the median grain size became larger than 5 mm for the smooth surface grains (from the paper of ISOPE2019 [11]), the coefficient “$C_{b0}$” slightly decreased (shown with the dashed line). Moreover, when the rough surface grains appeared as crushed gravel, the coefficient “$C_{b0}$” rapidly decreased (from the new experiments for the current paper). Meanwhile, Figure 4 shows that as the uniformity coefficient became larger, the reduction coefficient “$C_1$” decreased the bed-load coefficient, changing from 1.0 to 0.8 and corresponding from 1 to 20 in the uniformity coefficient. Moreover, Figure 5 shows that as dry density became larger, the reduction coefficient “$C_2$” decreased the bed-load coefficient, changing from 1.0 to 0.25, thus corresponding from around 1.5 g/cm³ in the dry density to 2.0 g/cm³.

The bed-load coefficient “$C_b$” for Ribberink’s formula can be calculated by the following equation:

$$C_b=C_{bo}\times C_1\times C_2.$$

(10)

2.3. Topography Change Simulation on the Sendai-Natori Coast: A Comparison of Total Bed-Load Transport vs. Suspension Load Only

The target here was to show the impact of suspension load on topographical changes in comparison to the bed-load. For this purpose, we selected for study Japan’s Sendai-Natori coast in the Miyagi Prefecture. The proposed coast is located between Sendai City in the north part and Natori City in the south part where Natori River flows west to east into the Pacific Ocean. It is noteworthy that the proposed area experienced an enormous tsunami after a 2011 earthquake with a magnitude of 9.0. The average incident tsunami height was set to 10 m on the offing boundary line.

In the numerical simulation model, the mesh size was set to 25 m and the bed-load coefficient was set to 18, as described bellows:

• In the land area, since the median grain size was around 0.4 mm, the uniformity coefficient was around 20, and the dry density was around 1.55 g/cm³: $C_b=C_{b0}\times C_1\times C_2$ = 22.5 × 0.8 × 1.0 = 18.
• In the beach area, since the median grain size was around 0.3 mm, the uniformity coefficient was around 10, and the dry density was around 1.55 g/cm³: $C_b=C_{b0}\times C_1\times C_2$ = 20.0 × 0.9 × 1.0 = 18.
• In the sea area, since we could not get sufficient information, $C_b$ = 18 was assumed.

The incident tsunami data on the offing boundary line of this simulation was set with reference to Figure 3.2.3 of Technical Note No.1231 of PARI (2011) [21]. However, since there was missing data for the Miyagi Prefecture tsunami, the average data between the Fukushima Prefecture and Iwate Prefecture were used. Moreover, the maximum tsunami height on the offing boundary line was estimated to be 10 m and the time of the incident phase and whole period of the first wave became 30 min and 48 min, respectively. Figure 6 shows the topography change results of the tsunami simulation considering only the impact of suspension load 48 min after calculation. The scouring depths, depicted with the pink color chart, shows a max depth of less than 0.5 m, while the deposition rate was even less than that. Ahmadi, Yamamoto, and Hayakawa (2019) [11] executed the same simulations considering the impact of both the suspension load and bed-load transport. As shown in Figure 7, the max scouring depth was found to be 4–6 m and the max deposition depth 4 m. Furthermore, to confirm the reproduction accuracy of the model, Figure 8 refers to the topographical change map from the measured data after the 2011 Great East Japan Tsunami and it is originally developed using digital elevation model (DEM). For further clarification, Figure 9 presents random spot elevation-change points along line-A on Figure 7 and Figure 8 for comparison purposes, with the horizontal axis being the distance of the relative point from the coastline and the y-axis being the change in elevation of the point. Considering this comparison, we evaluated the topography change reproduction simulation results to be acceptable.

2.4. Tsunami and Topographical Change Reproduction Simulation in Hamamatsu Coast of Japan

The Hamamatsu coast is located in the southern part of Hamamatsu City in the Shizuoka Prefecture. A tsunami simulation with conditions similar to the 2011 Great East Japan tsunami was executed on the proposed area using a numerical simulation model. Furthermore, the effects of the mitigation of the tsunami countermeasure structure development, which increased the dike height and construction of soil mounds, were also investigated. Then, the suspension load vs. total bed-load effects on topography change was presented and discussed.

The tide level data for the 2011 Great East Japan tsunami, recorded by the Japanese government, is shown in Figure 10a. Because the proposed coast for our simulation (Figure 6, Figure 7, Figure 8 and Figure 9) was the Miyagi Prefecture’s northern part, we used the fourth data (orange curve) from the top of Figure 10a. However, since there was no tidal data after a peak at Miyagi Prefecture’s northern part, we supplied a portion of tidal data after the peak with average values using the data of the Iwate southern part (red curve) and the Fukushima Prefecture part (green curve). At the time, in order to eliminate reflective waves and harassing waves, we eliminated waves where cycles were smaller than the cycle of the original tsunami waves from the epicenter. The results are shown in Figure 10b.

By using the tide waveform shown in Figure 10b and inputting 8 m as the maximum tide from the offshore boundary (water depth 100–110 m), the maximum tsunami heights on the shoreline became 10 m–15 m.

The bed-load coefficient was set to 18 because the sediment characteristics were similar to that of the Sendai-Natori coast. Figure 11 shows the topography of the proposed area and the calculation range. As depicted, the Tenryu River located on the east side of Hamamatsu City flows north to south toward the Pacific Ocean. Figure 12 shows the building ratio in a calculation mesh. The yellow color (the building ratio is 1% in a mesh area) refers to the soil and sand area, whereas the green (5% in a mesh area) refers to the wood and forest area. Further, the light brown (30% in a mesh area) shows a house with a garden, whereas the brown area (70% in a mesh area) shows a building area. The calculation mesh size was set to 12 m.

According to the simulation results, 16 min after beginning calculations, the tsunami height on the shoreline reached its maximum levels for a wide range. Afterward, for about 10 min, the maximum inundation depth was 5 m in a wide range in the inland area.

Figure 13, Figure 14 and Figure 15 show the inundation depth, tsunami velocity, and topographical change results, respectively, 33 min after the start of calculations when the seawater intrusion mode was over. According to these figures, the inundation area extended 3–6 km from the coastline. The forward velocity was calm after 33 min from the start of calculations and the maximum scouring depth was about 6 m.

In the proposed area, as countermeasures to a great tsunami, the crest height of the coastal embankment increased from TP +7 m to TP +13 m after the construction of CSG embankments and evacuation soil mounds; their locations are pointed out with red arrows in Figure 11. Here, including these developments in the simulation criteria, the output results significantly changed. The new CSG embankment did not wash out due to the increased crest height. The maximum inundation depth was halved over a wide range (about 2.5 m). Figure 16, Figure 17 and Figure 18 show the calculation results 33 min after the start, demonstrating that the inundation width from the coastline decreased to about half (1.5–3 km from the shoreline). The maximum inundation depth also decreased to about 2 m except around the Tenryu River and Magome River, where there was no CSG embankment. The tsunami velocity decreased by half later in the calculation time. Furthermore, the 12 soil mounds constructed with a circular shape of 100 m diameter (top) and 150 m diameter (bases) with a height of 8 m were expected to function as evacuees because the scouring depth around them was less than 1 m in depth. However, since they were made of clay and tended to weather easily and lose its shape, plantation cover was thought to potentially overcome this issue.

The next objective was to compare the suspension load and total sediment transport load topography changes. Figure 19 and Figure 20 show that the maximum scouring depth 28 and 33 min from the start of calculations must be around 0.5 m to 1 m, respectively, which was less than that of the total bed-load transport (i.e., 4–6 m in depth, as shown in Figure 17). It is noteworthy that the simulation input criteria for this comparison were similar as the ones defined above.

3. Model of Building Destruction Phenomenon by a Huge Tsunami

Many researchers have already published fragility curves for structures due to tsunami loads. However, since their methods are mostly based on statistics and probability models—such as fragility curves by Nanayakkara and Dias (2016) [23], as well as Suppasri et al. (2012, 2013, 2015) [24,25,26]—their fragility curves can be used to easily evaluate destruction probabilities of various building groups. Further, they are strongly influenced by the characteristics of the collected data.

On the other hand, the fragility curves presented in our paper can be used for a stability check of each building by using its structural parameters, such as pillar width, pillar interval, number of stories, and so on. Moreover, since the fragility curves in our paper were made via stress calculations to main members, our method can consider the influence of the strength for the materials of various main members.

3.1. Field Survey of Building Destruction by the 2011 Great East Japan Tsunami

In April 2011, our team executed an investigation of buildings destroyed by a huge tsunami caused by an earthquake off the Pacific coast of Tohoku where wooden buildings suffered heavy damages when inundation depths exceeded 1.5 m and reinforced concrete buildings broke when inundation depths exceeded 5 m. Using the investigation data presented later in this paper, our team developed diagrams to predict a building’s destruction grade and compared it against the real case site data. Therefore, the proposed diagrams explained the actual building damage. Moreover, the inundation depths taken from the numerical simulation model were used alongside the proposed diagrams to predict the actual building destruction with sufficient accuracy.

The Great East Japan Earthquake tsunami on March 11, 2011 caused enormous inundation damage along the Pacific coast from the Aomori to Chiba prefectures, and many damage investigation results have been already reported. In particular, the Institute of Industrial Science at the University of Tokyo (2011) studied details of the tsunami’s force evaluation method, the destructive tsunami pressure of reinforced concrete buildings and concrete structures, the impact force of drifting objects, and so on, to review the evacuation building design method. However, there was no universal failure limit investigation report to consider the differences in dimensions of a structure’s main components. Therefore, because the Nankai tsunami is expected to occur in the near future, and further research on disaster prevention and mitigation is desired, the authors focused on structural damage caused by the 2011 tsunami. Investigations and comparisons with the compiled structural damage data show that the structural damage limit calculation diagrams of Yamamoto et al. (2011) [13] and two-dimensional run-up numerical simulation model can be used to predict structural damage grade.

From 9–15 April 2011, a building destruction investigation including the destruction status of reinforced concrete buildings, wooden houses, concrete block walls, and steel-framed buildings due to the tsunami inundation was executed. The investigation was accomplished in the southern part of Miyako City, Iwate Prefecture, and the northern part of Hitachi City, Ibaraki Prefecture (excluding the area around Ishinomaki and the off-limits area due to the Fukushima Daiichi nuclear accident).

Table 2. The damage survey for reinforced concrete buildings during the 2011 Great East Japan Earthquake tsunami.

Prefecture	Iwate Prefecture
City Name	Miyako City			Yamada-Cho					Kamaishi City
Latitude	39°38′28″	39°38′25″	39°38′22″	39°28′17″	39°28′17″	39°28′17″	39°28′33″	39°28′33″	39°18′37″	39°15′59″
Longitude	141°52′46″	141°57′56″	141°58′10″	141°57′15″	141°57′14″	141°57′14″	141°57′17″	141°57′17″	141°53′19″	141°53′11″
Inundation depth (m)	4.8	5.0	6.0	6.0	5.5	5.5	6.5	6.5	5.0	6.5
Floor number	2	2	2	3–4	3	2	4	2	2	4
Column width (m)	0.44	0.63	0.85	0.50	0.55			0.49	0.53	0.65
Column height (m)	3.3	3.2	3.5	3.2	4.5			4.0	3.5	3.5
Column spacing (m)	6.0	6.0	5.5	6.5	6.0			4.0	4.5	5.0
Column break	None	None	None	None	None			None	None	None
Wall thickness (m)					0.3	0.25	0.24	0.10	0.17	0.30
Wall height (m)					4.0			3.5	3.2
Wall break	None	None	None	None	None	None	None	Broken	Broken	None
Prefecture	Iwate Prefecture
City name	Ofunato City					Rikuzentakata City
Latitude	39°4′1″	39°4′1″	39°4′1″	39°3′58″	39°4′7″	39°0′36″	39°0′32″	39°0′52″	39°0′51″	39°0′51″
Longitude	141°43′18″	141°43′18″	141°43′14″	141°43′15″	141°43′19″	141°38′38″	141°38′42″	141°38′22″	141°38′02″	141°38′02″
Inundation depth (m)	6.5	6.5	5.0	6.5	8.0	10.0	9.0	8.5	8.5	8.5
Floor number	3	2	2	2	3 – 4	5	2	3	2	2
Column width (m)	0.65	0.45	0.75	0.40	0.75		0.75	0.61	0.60	0.55
Column height (m)	3.5	3.1	3.4	3.2	4.0		3.5	3.5	3.7	3.5
Column spacing (m)	4.8	5.2	5.8	3.5	15.0		6.6	9.4	3.8	5.0
Column break	None	None	None	None	None		None	None	None	None
Wall thickness (m)	0.23	0.21	0.20	0.20	0.20	0.27	0.20	0.17	0.18	0.22
Wall height (m)	3.4								3.5
Wall break	Partial destruction	None	None	None	None	None	None	None	Broken	None
Prefecture	Miyagi Prefecture						Fukushima Prefecture	Ibaraki Prefecture
City name	Kesennuma City			Sendai City			Soma City	Kitaibaraki City		Hitachi City
Latitude	38°54′57″	38°54′58″	38°54′44″	38°16′28″	38°13′23″	38°13′20″	37°49′38″	36°49′49″	36°47′44″	36°34′41″
Longitude	141°34′55″	141°34′58″	141°34′51″	141°0′20″	140°59′05″	140°58′51″	140°58′23″	140°47′32″	140°45′15″	140°39′28″
Inundation depth (m)	8.0	7.0	7.2	5.0	7.5	7.5	6.0	2.7	2.0	3.5
Floor number	3	2	2	2	2	4	2	2 – 3	2	2
Column width (m)	0.40	0.76	0.52	0.82	0.77	0.60	1.20	0.44	0.27	0.90
Column height (m)	2.7	4.0	3.5	3.7	3.0	3.4	3.5	3.4	3.1	4.0
Column spacing (m)	3.5	5.0	3.8	6.5	7.4	4.5	20.0	3.6	3.7	6.0
Column break	None	None	None	None	None	None	None	None	None	None
Wall thickness (m)	No wall				0.21	0.20	Without Walls	0.20		0.23
Wall height (m)								3.2		3.5
Wall break		None	None	None	None	None		None	None	None
Note	The column spacing is a value parallel to the tsunami penetration direction.

,

Table 3. Wooden house damage survey results for the 2011 Great East Japan Earthquake tsunami.

Prefecture	Iwate Prefecture				Miyagi Prefecture		Fukushima Prefecture	Ibaraki Prefecture
City name	Miyako City				Kesen-numa	Sendai City	Iwaki City	Kitaibaraki City			Hitachi City
Latitude	39°38′28″	39°38′28″	39°38′28″	39°38′28″	38°54′18″	38°14′11″	36°51′39″	36°49′49″	36°47′44″	36°47′29″	36°34′41″
Longitude	141°52′45″	141°52′45″	141°52′45″	141°52′46″	141°34′58″	140°57′39″	140°47′25″	140°47′32″	140°45′15″	140°44′59″	140°39′27″
Inundation depth (m)	2.9	2.2	2.9	3.7	2.7	0.95	2.0	2.0	2.0	1.3	1.2
Floor number	2	2	2	2	2	2	2	2	1	2	2
Pillar width (m)	0.34	0.12	0.12	0.15	0.12	0.14	0.20	0.20	0.15	0.12	0.12
Pillar height (m)	3.0	2.7	3.0	3.0	2.7	2.7	2.7	3.1	2.7	2.7	2.7
Pillar interval (m)	6.2	1.8	4.8	2.5	1.8	1.9	3.6	3.6	3.7	3.7	3.6
Pillar break	None	None	Broken	Broken	Broken	None	None	None	None	None	None
Wall structure	Board + metal		Board + metal				Board + metal	Earth + bamboo + board	Earth + bamboo + board		Board + metal
Wall thickness (m)	0.15		0.12				0.10	0.15	0.10		0.10
Wall height (m)	2.9		2.9				2.6	3.0	2.5		2.6
Wall break	Broken	Broken	Broken	Broken	Broken	None	Broken	Partial destruction	Broken	None	Partial destruction
Note	The column spacing is parallel to the tsunami penetration direction.

and

Table 4. Survey on damage to concrete block walls during the 2011 Great East Japan Earthquake tsunami.

City Name	Miyako City	Sendai City				Kitaibaraki City
Latitude	39°36′54″	38°16′41″	38°16′41″	38°14′11″	38°13′46″	36°49′49″	36°47′44″	36°47′29″	36°47′29″
Longitude	141°57′27″	140°59′23″	140°59′23″	140°57′39″	140°58′19″	140°47′32″	140°45′15″	140°44′59″	140°44′59″
Inundation depth (m)	0.5	1.0	1.0	0.9	2.0	2.0	2.0	1.3	1.3
Wall structure	Black Mass Concrete	Perforated rebar Concrete block	Black Mass Concrete	Reinforced concrete Block	Reinforced concrete Block	Reinforced concrete Block	Reinforced concrete Block	Reinforced concrete Block	Black Mass Concrete
Break of the wall	None	None	Broken	None	Broken	Broken	Broken	Broken	Broken

show the details of the investigations including the latitude, longitude, inundation depth, column and wall materials, basic dimensions, and destruction status from the actual devastated area. Figure 21, Figure 22, Figure 23, Figure 24, Figure 25, Figure 26 and Figure 27 show the typical example photographs of heavy destruction taken during the aforementioned surveying.

Considering the above tables and photographs, we considered the breakdown characteristics of various buildings after the tsunami. The results are as follows:

(a)

In the case of a reinforced concrete building, windows and doors were torn when immersed in water at more than half of their surface area, but the walls only broke when the window occupancy ratio was 30% or less, wall thickness was 23 cm or less, and inundation depth was 5 m or more. Although there were no column failure cases for old buildings with insufficient seismic design, as is the case for Onagawa-Cho, the foundation was broken and the entire building was overturned.

(b)

In the case of a wooden house, when the inundation depth became about 1.2 m or more, windows and doors were easily damaged, and the 10 cm thick walls started to break. Pillars were more likely to collapse when the inundation depth exceeded 2.5 m.

(c)

Concrete block walls (standard thickness 10 cm) were more likely to fall if the inundation depth was 1 m or more when there was no reinforcement and 1.3 m or more when there was a reinforcing bar.

(d)

For steel-framed buildings, the data could not be listed due to space limitations, but there were no cases where the main steel frame was broken except for the case where large drifting objects such as ships and automobiles collided with them. However, since the wall was panel-shaped and thin, the wall would likely break when inundation depth was about 3 m or more.

Evaluating tsunami-related building destruction was performed with high accuracy using the conventional building design method, but the cost was high (US$1000 or more per house). The cost was due to the need for experts to conduct such complicated examinations. It is extremely difficult to formulate a disaster prevention and mitigation plan assuming damage from the design method. On the other hand, there is a simple method for examining the relationship between inundation depth (and inundation flow velocity) and the degree of destruction for each type of structure based on actual damage data. The effect on the degree is not known. Therefore, Yamamoto et al. (2011) [13] proposed a simple evaluation method, which is described below.

A typical building is designed to support loads with a rigid frame structure consisting of two columns and one beam. Furthermore, according to a hearing survey on the destruction of buildings caused by the 1993 Hokkaido Nansei-Oki Earthquake tsunami and the 2004 Indian Ocean Tsunami, when a tsunami large enough to destroy buildings occurs, windows and doors are instantaneously damaged, and water damages the inside and outside of buildings. When the static water pressure is offset, walls act openly without a break, and when the column is destroyed, it breaks from the base and becomes “destroyed”. If the column is not destroyed, it is in a “half-broken” state. Since the destruction was the same during the 2011 tsunami, the pillars were integral with the base and were stronger than the walls. Under the condition that walls have windows and doors that break instantaneously, we assumed the presence of a gate-type rigid frame where the tsunami force entered through the top load on the roof, as shown in Figure 28. Thus, the relationship between the inundation depth and unbreakable critical column width can be obtained from the stress at the column base where the tsunami force acts directly.

Tsunami force calculation equations include Asakura et al.’s (2000) [27], which was used near the coastline with a large Froude number, as well as Iizuka and Matsutomi’s (2000) [28] formula in the inland area with reduced momentum. The Institute of Industrial Science at the University of Tokyo (2011) proposed a considerably safe tsunami force calculation standard for evacuation building design. Furthermore, since the dimensions of a building’s main members are standardized by the representative dimensions of humans, combinations of basic dimensions such as column spacing and column height are limited. This chapter introduces threshold diagrams for calculating the column width and wall thickness, which break the limits for a given inundation height and a limited combination of frequently used for column spacing and column height.

As shown in Figure 28, we used a gate-type Rahmen model (i.e., a height of three meters and an interval length between two pillars of five meters as typical values) for the building structure of the on-off shore direction and assumed that the incident tsunami pressure acting on the building’s seaside wall could be obtained using Equation (9) (Iizuka and Matsutomi (2000)) [28] and that stress acted on the base of the seaside pillar.

$$F=\frac{C_f}{2}\rho Bh{v}^2=\frac{C_f}{2}\rho Bh{\left(F_r\sqrt{gh}\right)}^2$$

(11)

where $F$ is the tsunami force, $C_f$ is the hydrodynamic coefficient (= 2), $\rho$ is the seawater density, $B$ is the tsunami force action width, $h$ is the inundation depth in front of the building, $v$ is the inundation velocity, $F_r$ is the Froude number (= 1.1), and $g$ is the gravitational acceleration.

3.2. Confirmation of the Proposed Conventional Method with Disaster Data

3.2.1. Threshold Width of Columns in Reinforced Concrete Buildings

From the frequency characteristics in the actual disaster data, the column spacing of the first floor ramen part is an integral multiple of half (≒ 1 yard) + rounded value. The column height of the first-floor ramen section was 3.35 m and the median value of the high frequency and cross-sectional structure of the columns was a rectangular cross-section of double reinforcement bars in which reinforcing rods were arranged on the seaside and landside because of the high frequency. In the case of two and five stories, the bending/shear stress at the base of the column was calculated by changing the tsunami’s inundation depth and by comparing it with member strength. Thus, we determined the relationship between the inundation depth and the critical width that would not break, as shown in Figure 29. Here, to obtain a result on the safe side, the upper part of the second floor was considered as the upper load and was modeled with a slightly excessive concrete mass of 6.4 m × column spacing × (0.25 m × number of floors). The reinforcement ratio (= the steel ratio) was 0.05 close to the lower limit of the ratios and the clearance thickness (= the cover thickness) was 5.0 cm from the outer concrete surface so that it minimized the effective column width in the collected data. The effect of a band rebar (= a lateral tie, a tie hoop) was ignored. The compressive bending strength of concrete was considered 20 N/mm², which was close to the lower limit. The shear strength was considered 1/10 of that and the tensile bending strength of reinforcing bars was considered 300 N/mm², which was also close to the lower limit.

The non-destructive data of the two- and three-story buildings in Table 2 are plotted in Figure 29a, and the four- and five-story buildings show in Table 2 are plotted in Figure 29b with white circles (unit: meter). The applicability of this calculation diagram was high because the curve corresponded to the column spacing of each white circle and was approximately at the lower right. Under this condition, the type of failure was mainly the tensile bending failure of the reinforcing bar. In the case of seismic design of six stories or more, the possibility that the tsunami force exceeded the seismic force was low, and the necessity of the tsunami stability study became low as well.

3.2.2. Threshold Width of the Pillar in Wooden Buildings

Yamamoto et al. (2011) [13] showed the applicability of the gate-type ramen for the wooden pillars of wooden houses, using the data from the 1993 Hokkaido Nansei-Oki earthquake tsunami, as shown in Figure 30a,b. The height of the column on the first floor was set at 3.3 m, which was higher for safety, and the cross-section of the column was rectangular. The top load (i.e., the roof weight) of a wooden building consisted of the pillar interval × the pillar interval × 0.1 m of a wood portion (density 500 kg/m³) and the pillar interval × the pillar interval × 3.5 cm of a roof tile portion (density of 2000 kg/m³). In the case of a two-story building, a portion other than the roof on the second floor was modeled as a wooden lump with the pillar interval × the pillar interval × 0.4m. Regarding the bending strength and shear strength of wooden pillars, 20 N/mm² and 2.4 N/mm² of intermediate materials were used, respectively, referring to Notification 1452 of the Ministry of Construction (2000) [29].

Figure 30a,b also plots the damage data shown in Table 2, adding the pillar interval (unit: m). The open circles are non-destructive, the gray circles are partial destruction cases, the black circles are all destruction cases, and the curves corresponding to the pillar interval of each white circle are on the lower right, whereas the curves for each gray and black circle are on the upper left. It can be said that the strength of the members was appropriately evaluated.

3.2.3. Threshold Inundation Depth for Concrete Block Walls

For a wall stacked with blocks (standard thickness 10 cm × length 40 cm), the inundation depth of the unbreakable limit for each concrete strength were examined when only the tsunami force acted on the cantilever. In that case, the seaside horizontal hydrostatic pressure could not be ignored and the results are shown Figure 31a,b. Figure 31a shows the case where there was no reinforcement and Figure 31b shows the case where one reinforcing bar with a diameter of 1.3 cm was inserted for each block. Under these conditions, the wall broke by the tensile bending stress of the reinforcing bar and a value 15 times the compressive bending strength of concrete was used for this strength. These figures show the damage data of the Hokkaido Nansei-Oki earthquake tsunami. The open circles are non-destructed, the closed circles are all destructed cases, and the open circles are below the curve, whereas the black circles are above, indicating that the applicability of this calculation diagram is high.

3.3. Verification Examples of Building Destruction due to 2011 Great East-Japan Tsunami

Using the tsunami simulation model and the developed threshold diagrams for building destruction design, some verification examples are presented in

Table 5. Verification examples of building destruction due to the 2011 Great East Japan tsunami.

* RC: reinforced concrete.

. Figure 32, Figure 33 and Figure 34 show the calculated maximum inundation depths in the three proposed areas: Rikuzentakata City, Kesennuma City, and Sendai City. Referring to verification example results, the design destruction state matched the real destruction state of buildings except for the first one. Although the reinforced concrete building located in Rikuzentakata City was stated as “broken” in the design, it was not broken in the real situation and was still on the design’s safe side. Furthermore, when comparing the calculated inundation water depths and measured depths, we observed that the results were acceptable.

4. Conclusions

To acquire the bed-load coefficient of Ribberink’s bed-load transport equation, we developed useful diagrams. These diagrams, illustrated in Figure 3, Figure 4 and Figure 5, present a wide range of grain size, dry density, and uniformity coefficients. Figure 3 shows that the median grain size (i.e., natural soil, sand, and gravel) increased to around 5 mm and the bed-load coefficient increased from zero to around 60. Further, when the median grain size became larger than 5 mm for smooth surface grains, the bed-load coefficient remained almost constant; however, for the rough surface grains, such as crusher-run stone, the bed-load coefficient decreased to around 5 mm. Meanwhile, Figure 4 shows that as the uniformity coefficient became larger, the reduction coefficient for decreasing the bed-load coefficient changed from 1.0 to 0.8, which corresponded from 1 to 20 in the uniformity coefficient. Moreover, Figure 5 shows that as dry density became larger, the reduction coefficient for decreasing the bed-load coefficient changed from 1.0 to 0.25, which corresponded with 1.5 g/cm³ in the dry density to 2.0 g/cm³.

From the Sendai-Natori coast tsunami reproduction simulation results, we observed that the impact of suspension load versus total bed-load was not large. The results of suspension only showed a max depth of scouring less than 0.5 m while that of the total bed-load was around 4–6 m. Furthermore, when we researched to the deposition depths, the deposition depth with suspension loads was negligible, especially when compared to the total bed-load, which was around 5–6 m.

The results of the tsunami simulation on the Hamamatsu coast, where evacuation soil mounds were constructed and the crest height of the coastal dikes increased from TP +6 m to TP +13 m, showed that the evacuation soil mounds served their purpose; however, because they were directly affected by weathering, plantations could overcome this issue. Moreover, the dike improvements significantly decreased inundation depth and height, while inundation water intruded into the area around Magpome River and Tenryu River. Still, however, there were no riverbank improvements.

Lastly, using the building destruction survey data from the 2011 Great East Japan tsunami, useful diagrams were developed and presented to significantly decrease the cost and time needed to evaluate building destruction analysis due to a tsunami. These diagrams are shown in Figure 29, Figure 30 and Figure 31 and confirmed that the tsunami damage prediction method based on the method of calculating the threshold limit of columns and walls was useful.