FOR TSUNAMI FLOOD RISK ASSESSMENT IN URBAN AREAS

Several recent tsunami events worldwide have increased the general public ́s awareness of the inherent risk many coastal communities face. Accordingly, academics and professionals have sought to acquire a better understanding of the underlying physical processes involved, and have developed increasingly sophisticated numerical tools to simulate them. One important characteristic of the tsunami flood tragedies that must be acknowledged is that they occur in urban areas, where the built environment is likely to determine the flow, and hence, the flood risk. Despite this fact, relatively little literature covers the particular challenges of modelling the hydrodynamic interactions of the tsunami flood and the buildings and infrastructure. In this work we propose the indirect inclusion of the built environment's effect by actually excluding every building from the numerical domain, enclosing them in an internal reflective boundary condition. The objective is two-fold: it allows a reasonable inclusion of the flow funneling, and it reduces the number of cells in the numerical domain, allowing the use of a higher grid resolution in the expected wet domain, and to achieve shorter simulation runtimes. The methodology was validated against measurements of flow depth, velocity, and momentum on a scale model of the city of Seaside, Oregon. The model and methodology were subsequently used to produce flood depth and velocity maps for 11 coastal cities in Chile, caused by the likely maximum mega-thrust earthquake for each of them.


INTRODUCTION
Within the development of a model to determine the Probable Maximum Loss (PML) for the Chilean Insurance Association (AACH), a tsunami risk assessment was carried out and inundation maps were developed for 11 coastal cities of Chile (PRDW, 2014): Arica, Iquique, Mejillones, Antofagasta, La Serena, Coquimbo, Quintero, Viña del Mar, Valparaíso, San Antonio, and Talcahuano.These areas are of interest given their high population density and the presence of public and private infrastructure near the coast.The tsunami threat over these areas was evaluated in terms of their inundation depth and induced current velocities.
A numerical modelling methodology was defined in two separated physical scales; for the regional scale the propagation over oceanic waters of the free surface's initial deformation induced by the earthquake was modelled, while for the local scale, propagation over waters with a depth shallower than 100 m and inundation over originally dry land was modelled.This difference is justified by the physical processes fundamentally different between both scales, and the usage of specific tools for the modelling process.
Additionally, buildings were schematized as reflective boundaries in the inundation model, which allowed the reproduction of the funneling effect through the streets and a proper characterization of flooding velocities.The latter parameter is highly relevant in various engineering applications, such as hydrodynamic loads calculation and scour estimation, among others.At the same time, reduction of the effective calculation area of the model domain resulted in an increased resolution of active areas, improving precision of the obtained results.

METHODOLOGY
For each analyzed city a team specialized in geophysics determined a maximum credible event, as presented in Table 1 (PRDW, 2014).For the cities of Arica and Iquique, the fault was divided into 5 sub-faults in order to reproduce the curved shape of the fault in that area.
Cortical initial deformation produced by maximum credible events was calculated according to (Okada, 1985) methodology.This deformation was used as an initial condition on the water surface elevation.All 5 initial conditions are presented from Fig. 5 to Fig. 10.This perturbation is then propagated radially, like a long lineal wave of low amplitude, without suffering big changes while it propagates over the open deep ocean.Interaction with the continental shelf slope is the main phenomenon at this stage, since it determines the generation of edge waves and it could control the tsunami behavior at a certain location.A second stage is produced when the tsunami wave approaches the shore.With shallower depths, non-linear interactions between celerity, wave height and local depth become relevant.Also, interaction with the shoreline and possible land inundation is produced, which is a numerical challenge by itself.
Both processes present different physical characteristics and temporal/spatial scales, which is the reason why they were modelled separately in this study.Results from the regional scale or propagation model were used as a boundary condition for the local scale or inundation model.

Propagation model
For propagation over the open ocean a hydrodynamic model of the MIKE 21 numerical suite developed by the Danish Hydraulic Institute (DHI, 2016) was used.The model utilizes geographical coordinates to solve the shallow water non-linear equations over an irregular mesh, which allows a usage over extended domains.It must be noted that numerical schemes tend to produce numerical dispersion, which can be relevant in the process of far field tsunami propagation, however they are not relevant for the scale modelled in this study, since all tsunamis are considered of nearfield generation.
For the topo-bathymetric construction, data was obtained from a global database (GEBCO, 2009), nautical charts from the Chilean Navy Hydrographic and Oceanographic Service (SHOA), Topography from the NASA project Shuttle Radar Topography Mission (Jarvis et al., 2008) and local bathymetries and topographies available.
Shallow water equations consider the following equations of conservation of mass and momentum: where η is the instant level of the water free surface, h is the instant depth, u and v are the orthogonal velocity components along the x and y axes and g is the gravity acceleration.The Coriolis effect is considered in parameter f = 2Ω sinϕ, where Ω is the Earth's angular rate of revolution and ϕ the geographic latitude.Bottom friction is considered proportional to the squared magnitude of the flow velocity, characterized by the Manning roughness parameter n.
Finite volumes method is used to discretize and resolve nonlinear shallow water equations.Domain was discretized with bi-dimensional triangular elements, considering equations are integrated over the vertical axis.
Along the coastline a closed boundary condition was used, where the velocity component perpendicular to it is forced to be zero.In other words, the domain is considered fixed, no inundation is modelled, increasing computation velocity.For the sea boundary a Flather (1976) condition was used, one of the most efficient open boundary conditions (DHI, 2016), where velocity components u and v, and the sea level on the boundary are forced to zero.

Inundation model
Wave tsunami interaction with the coastline and originally dry areas is modelled using ANUGA (Nielsen et al., 2005), developed by the Australian National University and Geoscience Australia.The model simulates the propagation of tsunami waves over a flexible irregular mesh.The model is particularly efficient modelling inundation since it does not consider the coastline as an explicit boundary condition but it is naturally computed by the numerical scheme.
To incorporate urban morphology in the model, buildings that could be an obstacle for the flooding were identified.These buildings were excluded from the inundation area and were assigned a reflective boundary condition.Information was obtained from Coastal Usage Drawings provided by (SHOA, n.d.), and complemented with Google Earth®.
The model resolves the same shallow water equations presented in equations ( 1), ( 2) and ( 3), but presented in its conservative form, which is in terms of the instant depth (h) and momentum (hu and hv).For the discretization and resolution of the equations a finite volumes method is used, following an explicit shock capturing numerical scheme, to ensure the flux depth remains positive.This scheme also allows discontinuities of the water free surface, and dissipates energy through breaking, although it carries a computational cost.
In these equations, conserved variables evolution depends on the mass and momentum fluxes through the boundaries of each cell of the domain and on the source terms (S on equation ( 5)), which refer to external forces, in this case associated to gravity force and bottom friction.The terms zx and zy represent slopes on x and y directions.
Boundary condition of the local model along the sea side, correspond to the time series of sea level and momentum obtained from the bigger propagation model, interpolated along the domain's boundary.The reflective condition used along the buildings boundary is assigned with the value of the free surface of the neighboring volume and the same momentum magnitude but with opposite direction.

VALIDATION
The urban inundation model was compared against laboratory tests results produced at Oregon State Uniersity (OSU) (Park et al. 2013).The experiment consisted on a long wave propagating over a bathymetry with two consecutive slopes, representing a 1:50 scale idealization of Seaside in Oregon, as presented on the left panel in Fig. 1.The proceeding included wave breaking before urban inundation is produced.
Water depth and flux velocity were measured on the locations indicated on the right panel of Fig. 1.Using these data momentum was calculated.These locations are placed along streets of the scale model and are representative of maximum velocities expected during an urban inundation.Bathymetry, including buildings height was obtained from the authors so researchers and engineers can prove and validate their models.The tsunami wave was generated with a paddle as wide as the tank itself at a depth of 0.97 m.Water free elevation measured at x =5 m is presented in Fig. 2.This time series (and its corresponding celerity) was used in this work to force the numerical model.In Fig. 3 the numerical model domain used to simulate the experiment with ANUGA is presented.It can be observed how buildings have been explicitly removed from the computation domain.It must be noted that it is assumed that there is no flux above the level of houses and buildings or that its occurrence is unlikely to happen, enough to be discarded.Numerical domain was constructed on a triangular mesh with a resolution of 100 cm 2 over the wet domain and 50 cm 2 over the urban area, equivalent to triangles with equivalent side of approximately 15.2 cm and 10.7 cm respectively.This resolution is enough to include several triangles across each street.Regarding the calibration, only a sensibility analysis of the bottom roughness, using Manning coefficients n = 0.013, 0.011, 0.008, and 0.005 was considered.This type of shock capture explicit models in finite volumes, and particularly ANUGA, present very little numerical parameters, therefore the modelling quality greatly depends on the resolution adopted, boundary conditions and bottom friction.This is a desirable behavior since the operation of the model is more transparent for the modeler.Obtained results were compared against available measurements, showing a high correspondence in most locations analyzed.Bigger differences were presented for points located far from the coastline towards inland, as could be expected.Some of the time series compared are presented in Fig. 4, for points A6 and B5.In general, it is observed that bottom friction, at least at this scale, presents a low incidence in results, except for the arrival time and the amplitude of the first wave.The best agreement is obtained for n = 0.011, which corresponds to the recommended value for well finished concrete, which is the material reported for the tank used for the scale model.In presented time series, and most of those time series compared, there is a tendency to underestimate the momentum flux, even though the velocity magnitude is generally well represented.

RESULTS
In following figures (Fig. 5 to Fig. 10) the 5 initial conditions used are presented, with maximum inundation heights and maximum velocity magnitudes for the 11 cities of interest (PRDW, 2014).
On the upper left corner of each initial condition image, the names of the cities modelled for this condition are presented.The geographic location of each city is also presented with the corresponding number.In general, inundations with a high horizontal intrusion in relatively low areas are observed, such as the northern part of Arica, La Serena and Talcahuano.This phenomenon is also observed in estuaries, as observed in the northern part of Quintero (Campiche creek) and the southern part (landing track).The same effect is observed in Viña del Mar, along the Marga Marga creek.
In Talcahuano it is observed that, differently from occurred in 2010, a significant horizontal intrusion is produced from the west (San Vicente bay).This is explained due to a relevant cortical displacement produced in front of that area.
Regarding maximum velocities, values between 2 m/s to 6 m/s are observed.Exceptionally some higher values are presented, as the Arica situation, especially on the norther part which is the most inundated.

CONCLUSIONS
A methodology of tsunami modelling considering two scales is presented.A regional scale considered for propagation over oceanic waters and a local scale for simulating inundation on originally dry areas.Each scale was modelled with the numerical scheme most suitable to associated physical phenomena.
The relevance of the presented methodology is that it considers the need to incorporate the effect of buildings in the inundation of urban areas, with the objective of simulating the complex flow patterns expected around buildings.Additionally, the exclusion of areas of the urban area analyzed allows higher resolutions without compromising the operability of the model, i.e. within a reasonable computation time.This focus and the model utilized were validated with satisfactory comparisons between data obtained from laboratory tests over a physical model with reduced scale.
Despite the significant improvement that this methodology brings to tsunami modelling of urban areas, it is evident that mathematical modelling of tsunami inundation still presents high uncertainty (seismological, hydrodynamic, mathematical), especially in urban inundation scenarios, where additionally to funneling effects on streets, scouring processes and sediment and debris transport processes can be very relevant.

Figure 2 .
Figure 2. Time series of incident wave at x = 5 m, used as input for the numerical model (Park et al. 2013).

Figure 5 .
Figure 5.Initial condition, inundation height and velocity magnitude for Arica and Iquique.

Figure 6 .
Figure 6.Initial condition, inundation height and velocity magnitude for Mejillones and Antofagasta.

Figure 7 .
Figure 7. Initial condition, inundation height and velocity magnitude for La Serena and Coquimbo.

Figure 8 .
Figure 8.Initial condition, inundation height and velocity magnitude for Quintero, Viña del Mar and Valparaíso.

Figure 9 .
Figure 9. Inundation height and velocity magnitude for San Antonio.