the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
A comparison of a twodimensional depthaveraged flow model and a threedimensional RANS model for predicting tsunami inundation and fluid forces
Xinsheng Qin
Michael Motley
Randall LeVeque
Frank Gonzalez
Kaspar Mueller
The numerical modeling of tsunami inundation that incorporates the built environment of coastal communities is challenging for both 2D and 3D depthintegrated models, not only in modeling the flow but also in predicting forces on coastal structures. For depthintegrated 2D models, inundation and flooding in this region can be very complex with variation in the vertical direction caused by wave breaking on shore and interactions with the built environment, and the model may not be able to produce enough detail. For 3D models, a very fine mesh is required to properly capture the physics, dramatically increasing the computational cost and rendering impractical the modeling of some problems. In this paper, comparisons are made between GeoClaw, a depthintegrated 2D model based on the nonlinear shallowwater equations (NSWEs), and OpenFOAM, a 3D model based on Reynoldsaveraged Navier–Stokes (RANS) equation for tsunami inundation modeling. The two models were first validated against existing experimental data of a bore impinging onto a single square column. Then they were used to simulate tsunami inundation of a physical model of Seaside, Oregon. The resulting flow parameters from the models are compared and discussed, and these results are used to extrapolate tsunamiinduced force predictions. It was found that the 2D model did not accurately capture the important details of the flow near initial impact due to the transiency and large vertical variation of the flow. Tuning the drag coefficient of the 2D model worked well to predict tsunami forces on structures in simple cases, but this approach was not always reliable in complicated cases. The 3D model was able to capture transient characteristic of the flow, but at a much higher computational cost; it was found this cost can be alleviated by subdividing the region into reasonably sized subdomains without loss of accuracy in critical regions.
1.1 Scope and motivation of the study
For many years, researchers have been working on different numerical models that can predict tsunami behavior. Tsunami prediction generally requires modeling at a wide range of spatial scales, including (from large to small scale) offshore wave propagation, beach runup, inland inundation, and impact on individual structures.
Due to the large differences in scale for the different processes, most tsunami models solve twodimensional depthintegrated equations, e.g., the nonlinear shallowwater equations (NSWEs) or some form of Boussinesq wave equations, to predict tsunami behavior, using computational grids that vary several orders of magnitude in spatial resolution, from several kilometers far from the shoreline to 10 m inland. The NSWEs are often used in the nearshore and inundation zone, since they can handle nonlinearities that arise in very shallow water and can be adapted to deal robustly with wetting and drying. However, it is not clear whether these equations are adequate to properly model fully threedimensional turbulent flow, particularly at the scale necessary to determine tsunami impact and corresponding tsunamiinduced forces on individual structures.
It would be preferable to solve the threedimensional Navier–Stokes equations with a proper turbulence closure. However, this is still extremely expensive computationally relative to twodimensional models, and only practical for detailed simulations over small spatial regions.
Although such modeling is challenging, the latest version of ASCE 716 (American Society of Civil Engineers) in the United States has a chapter for tsunami loads and effect for coastal structures. The provision requires sitespecific inundation modeling and analysis be performed for all vertical evacuation structures. One of such examples is the design of the first vertical evacuation structure in the United States (Ash, 2015), the sitespecific inundation analysis of which was conducted by González et al. (2013).
Given the challenges in tsunami inundation modeling, and the necessity of doing so for coastal structures design and tsunami hazard mitigation, the goal of this paper is to (1) study the characteristics of two different types of models when applied to tsunami inundation modeling (with an explicitly represented constructed environment); (2) show some guidance for modeling tsunamis or other flooding events in similar constructed environments; and, by showing the pros and cons of both models, (3) provide some insight to tsunami hazard researchers, coastal structure designers, and relevant government agencies.
1.2 Previous work
Before introducing the two numerical models used in the current study, a brief review of previous research involving different types of models is given below.
The twodimensional depthintegrated equations are the most widely used tsunami models for their simplicity and computational efficiency. Popinet (2012) simulated the 2011 Tōhoku tsunami by solving the 2D NSWE with dynamically adapted spatial resolution that varied from 250 m in flooded areas nearshore up to 250 km offshore. The model accurately predicted longdistance wave and coarsescale flooding; the initial surface elevation was determined from a source model based on seismic inversion (as opposed to inversion of DART buoys and tidal gauge time series). This also showed that an accurate and consistent model of tsunami wave propagation can sometimes be constructed using only seismic wave inversion. Wei et al. (2013) used the Method of Splitting Tsunamis (MOST) model to simulate the same tsunami event. The MOST model solves the shallowwater equations in spherical coordinates with numerical dispersion. Their results demonstrated that it may be possible to forecast nearfield tsunami inundation in real time. Hu et al. (2000) presented an NSWE model that can simulate storm waves propagating in the coastal surf zone and overtopping a sea wall. They found that waves overtopping a vertical wall may be approximately modeled by representing the wall as a steep slope and that the overtopping rate is sensitive to the bottom friction and the minimum friction depth. The twodimensional NSWE model of wave runup and overtopping by Hubbard and Dodd (2002) features an adaptive mesh refinement (AMR) algorithm. Their model can accurately reproduce 1D and 2D wave transformation, runup, and overtopping in physical experiments. Their modeling of seawall overtopping by offnormal incident waves showed that there can be more flooding in such a situation than at normal incidence. Lynett (2007) simulated longwave runup obstructed by an obstacle and concluded that the obstacle can help reduce runup and maximum overland velocity if the wave is highly nonlinear (with a ratio of wave height to shelf water depth ≥0.5). The sensitivity study also showed that in cases of breaking waves the Boussinesq model was more accurate than the nonlinear shallowwater equations in terms of wave runup (maximum differences up to 10 %). For nonbreaking long waves, differences between the two were negligible. Shi et al. (2012) developed a highorder adaptive timestepping total variation diminishing (TVD) solver for a fully nonlinear Boussinesq model and validated it against a series of laboratory experiments for wave shoaling and breaking and a suite of benchmark tests for wave runup. The results showed that the model was able to accurately model wave shoaling, wave breaking, and waveinduced nearshore circulation. With a Boussinesq model, Lynett et al. (2010) simulated overtopping of levees of the Mississippi River–Gulf Outlet (MRGO) during Hurricane Katrina at four characteristic transects along the 20 km long stretch of the levees. The predicted overtopping rates agreed well with the observed data.
As computing power increases, it becomes possible to model the tsunami runup process by solving threedimensional Navier–Stokes equations with a proper turbulence closure. Choi et al. (2007) solved threedimensional Reynoldsaveraged Navier–Stokes (RANS) equations to simulate wave runup on an conical island and compared different turbulence closure models including kϵ, RNG (renormalization group methods; Yakhot et al., 1992), and LES (largeeddy simulation). Their results showed that LES and RNG kϵ are similar and more accurate than kϵ is worse than those two. Williams and Fuhrman (2016) solved incompressible RANS equations with a transitional variant of the standard twoequation kω turbulence closure to study boundary layer flow induced by tsunamiscale waves. Their results indicated that the boundary layer generated by a tsunami is both currentlike due to the long duration and wavelike due to its unsteadiness. The study also indicated that an existing expression for maximum bed shear stress under wind wave scale can be reasonably extrapolated to full tsunami scale. Mayer and Madsen (2000) investigated wave breaking in the surf zone by solving the RANS equations with a kω turbulence model. They found that the volumeoffluid method could be used successfully to simulate wave breaking and that, although some instabilities occurred in applying the RANS equations, they can be eliminated by an ad hoc modification of the turbulence model. Other relevant studies include Biscarini (2010), Montagna et al. (2011), and Larsen et al. (2017), which model landslide generated tsunamis and tsunamiinduced scours around coastal structures.
The prediction of tsunami impact on individual structures is also important because it provides guidance on designing coastal structures in tsunami inundation zones. The twodimensional depthintegrated model may not work properly for these scenarios since the problems are more threedimensional with large variation in the vertical direction and with transient and turbulent flow impacting the structure. In these cases, a threedimensional model that solves the Navier–Stokes equation may give much better results. Researchers at the University of Washington (UW) modeled a series of “dam break” experiments by solving the 3D RANS equations for boretype impact of a wave on a series of 1∕20scale model girder bridges to assess the 3D effects on bridge skew (Motley et al., 2015; Wong, 2015).
The scale of modeling tsunami inundation inland with an explicitly represented constructed environment lies between that of modeling the largescale tsunami wave propagation offshore and the smallscale tsunami impact on individual structures. This process is actually even more challenging to model since, for twodimensional depthintegrated models, inclusion of the constructed environment increases the complexity of the topography, and the flow begins to have more variation in the vertical direction; while for the threedimensional model that solves the Navier–Stokes equations, a fine mesh needs to be generated around each individual structure, which dramatically increases the number of cells in the computational domain.
Some researchers have tried to model this process with twodimensional models. Ozer Sozdinler et al. (2015) used the numerical code NAMI DANCE to investigate hydrodynamic parameters in tsunami inundation zones with idealized structures – three rows of 20 blocks representing threestory concrete buildings. The code solved the NSWE using a finitedifference technique in a staggered leapfrog scheme. The effect of wave period, wave shape, protection structures, building layout, and Manning's coefficient are discussed. Some major conclusions included that the coastal protection structures like seawalls and breakwaters have very limited effect if the waves are able to overtop them and that it is preferable to use different Manning's coefficients for the sea, land, and buildings if more accurate values of hydrodynamic parameters are needed, albeit at the expense of more computational time. Similar conclusions on the Manning's coefficient were presented by Park et al. (2013). They simulated tsunami inundation in part of Seaside, Oregon, and compared flow parameters with their physical experiment. The comparison showed that the flow parameters were sensitive to the friction coefficient, especially for the momentum flux, which is proportional to tsunami loads on structures. For instance, decreasing the friction coefficient by a factor of 10 increased the predicted momentum flux by 208 %. Muhari et al. (2011) compared three different tsunami inundation models to evaluate tsunami impact on coastal communities: (1) the Constant Roughness Model (CRM) which uses a constant friction coefficient, does not include the constructed environment, and assumes that all buildings are not able to withstand the tsunami; (2) the Topographic Model (TM), which includes the constructed environment by incorporating building shape and height information into the topography; and (3) The Equivalent Roughness Model (ERM), which represents the building by using a different equivalent friction coefficient at the site of a building on the original topography (with only terrain information but not building height). Both the TM model and the ERM model gave more reliable prediction than the CRM model did, which confirmed the importance of taking the constructed environment into consideration.
However, few researchers have tried to use a threedimensional model for inundation in a complex built environment. Shin et al. (2012) applied a 3D LES model with twophase flow to simulate inland tsunami inundation in a coastal city with hundreds of buildings and compared the prediction with experimental measurements. However, a fairly coarse mesh was used on land, and each building had only three to five mesh cells along its edge in the alongshore or crossshore direction, so that the resulting agreement in flooding depth can only be considered qualitative. Qin et al. (2018) used 3D RANS equations to predict tsunami inundation process and loads on individual buildings in part of Seaside, Oregon, and demonstrated that the whole part can be modeled using subsections with proper width without loss of accuracy in areas of interest.
In this paper, the results from a twodimensional NSWE model and a 3D Navier–Stokes model are presented for the test case of flow through a scale model of a portion of Seaside, Oregon. Two opensource models are used: the 2D GeoClaw software from Clawpack (Clawpack Development Team, 2015), which is widely used for modeling tsunamis (both global propagation and local inundation), and the 3D OpenFOAM software (The OpenFOAM Foundation, 2014). The two models are first compared and validated against an experiment in which a simple bore impinges on a single column and then compared for the Seaside model.
2.1 Twodimensional model
The nonlinear shallowwater equations can be written as
where $u(x,y,t)$ and $v(x,y,t)$ are the depthaveraged velocities in the two horizontal directions, h is the water depth, g is gravitational acceleration, B(x,y) is elevation of the seabed, and $D=D(h,u,v)$ is the drag coefficient. The drag coefficient D could have many forms; in this study it is represented by
where M is the Manning's coefficient and is set to 0.025 s m${}^{\mathrm{1}/\mathrm{3}}$ for all twodimensional simulations in this study. The subscripts in these equations represent firstorder partial derivatives.
The GeoClaw model (LeVeque et al., 2011; Berger et al., 2011) features AMR and is released as a submodule of the Clawpack software (Clawpack Development Team, 2015), an opensource package for solving hyperbolic systems of partial differential equations (PDEs) of one, two, and three dimensions, through finitevolume implementation of highresolution Godunovtype “wavepropagation algorithms”. Cell averages of the solution variables q are computed over the volume of each cell and updated with waves propagating into the cell from all surrounding cell edges. The wave at each edge is computed by solving a “Riemann problem” with initial piecewise constant data determined by cell averages on each side of the edge. This method is especially good at solving problems with discontinuous solutions like shock waves, which usually arise in the solution of nonlinear hyperbolic equations (e.g., bores in the case of NSWEs).
Specifically, GeoClaw uses a variant of the fwave formulation of the wavepropagation algorithms that allow incorporation of the topography source terms on the righthand side of Eqs. (2) and (3) into the Riemann problem directly. The augmented Riemann solver in GeoClaw combines the desirable qualities of the Roe solver (Roe, 1981), HLLEtype (Harten, Lax, van Leer and Einfeldt) solvers (Einfeldt, 1988; Einfeldt et al., 1991), and the fwave approach (Bale et al., 2003). The Roe solver provides an exact solution for the singleshock Riemann problem. It is also depthpositive semidefinite like the HLLE solvers, has a natural entropy fix by providing more than two waves, and yields a better approximation for problems with large rarefactions. A large class of steady states is also preserved, even for nonstationary steady states with nonzero fluid velocity. In addition, it is able to handle the presence of dry states in the “Riemann problem”, in which one state is wet (h>0) while another is dry (h=0), or both states are dry. It also works robustly in situations where the topography changes abruptly from one cell to another by an arbitrarily large value. For more details on the augmented Riemann solver in GeoClaw, see George (2008).
A typical characteristic of tsunami inundation models, especially those that incorporate the built environment, is that the spatial scale of regions of interest may vary from kilometers to meters. For regions several kilometers offshore, grid cells can be as large as thousands of meters on a side, while for regions near the shoreline or in a built environment onshore, grid cells must be refined to several meters or less, since the size of a building may be only several meters and an adequate number of grid cells are required to achieve acceptable accuracy. In GeoClaw, a patchbased AMR technique can efficiently handle these situations (LeVeque et al., 2011; Berger and Leveque, 1998).
2.2 Threedimensional model
For the threedimensional model, version 2.3.1 of the opensource computational fluid dynamics (CFD) package OpenFOAM was used (The OpenFOAM Foundation, 2014). The package comes with different solvers for different types of flow. For tsunami inundation, in which there are two immiscible fluids (air and water) with a free interface, the interFoam solver can be chosen which uses the PISO algorithm to solve the RANS equations with a volumeoffluid (VOF) approach to model the free surface. For details on these numerical methods, readers can refer to Hirt and Nichols (1981) and Versteeg and Malalasekera (2011). The VOF approach defines a scalar field α_{water}, which represents fractional volume of water in each cell. A cell full of water (ρ=1000 kg m^{−3}, $\mathit{\nu}=\mathrm{1.0}\times {\mathrm{10}}^{\mathrm{6}}$ m^{2} s^{−1}) has α_{water}=1.0, while a cell full of air (ρ=1.22 kg m^{−3}, $\mathit{\nu}=\mathrm{1.48}\times {\mathrm{10}}^{\mathrm{5}}$ m^{2} s^{−1}) has α_{water}=0.0. Here ρ is the mass density of the fluid and ν is the kinematic viscosity. A cell with α_{water} between 0 and 1 contains the interface. A special transport equation is solved to advance the α_{water} field. To close the RANS equations, Menter's kω shear stress transport (SST) model (Menter and Esch, 2001) was applied.
There are many other turbulence closure models, among which the kϵ model is also very popular. It is suitable for fully turbulent and nonseparated flows and has the shortcoming of numerical stiffness in the viscous sublayer, which can result in stability issues (Menter, 1993). It was also applied to model the inundation process in this study but became unstable during the simulation. The kω SST is generally more stable and behaves better in modeling partially separated flows, which is the case in the current study (flow becomes separated after passing around the built environment).
With the assumption of an incompressible fluid, the RANS equations are listed below:
where ${\stackrel{\mathrm{\u203e}}{u}}_{i}$ is the mean velocity in the i direction, ${{u}_{i}}^{\prime}$ is the fluctuating component of velocity in the i direction, and $\stackrel{\mathrm{\u203e}}{p}$ is the mean pressure. If u_{i} is the velocity component in the i direction, then ${u}_{i}={\stackrel{\mathrm{\u203e}}{u}}_{i}+{{u}_{i}}^{\prime}$. The Reynolds stress term in Eq. (6) is
where k is the turbulence kinetic energy and ν_{t} is the turbulence eddy viscosity. The equations above need to be closed with some closure model. Here Menter's kω SST model (Menter and Esch, 2001) was applied:
where ν is the kinematic viscosity of fluid and $\stackrel{\mathrm{\u0303}}{G}$ is defined as $\stackrel{\mathrm{\u0303}}{G}=min\left\{G,{c}_{\mathrm{1}}{\mathit{\beta}}^{\ast}k\mathit{\omega}\right\}$, where G is the production term and defined as
S is the invariant measure of the strain rate, defined by
and S_{ij} is the strain rate tensor, defined by ${S}_{ij}=\frac{\mathrm{1}}{\mathrm{2}}\left(\mathrm{\nabla}\mathbf{U}+{\mathbf{U}}^{T}\right)$. F_{1} is a blending function defined by
where ${\mathrm{CD}}_{k\mathit{\omega}}^{\ast}$ is defined by
and CD_{kω} is defined by
After solving Eqs. (8) and (9), ν_{t} can be calculated by
where F_{2} is a second blending function, defined as
All other constants are computed using a blend from the corresponding constants associated with the kϵ and kω models via blending functions like $\mathit{\varphi}={\mathit{\varphi}}_{\mathrm{1}}{F}_{\mathrm{1}}+{\mathit{\varphi}}_{\mathrm{2}}\left(\mathrm{1}{F}_{\mathrm{1}}\right)$. Values for these constants are α_{k1}=0.85013, α_{k2}=1.0, α_{ω1}=0.5, α_{ω2}=0.85616, β_{1}=0.075, β_{2}=0.0828, γ_{1}=0.5532, γ_{2}=0.4403, ${\mathit{\beta}}^{\ast}=\mathrm{0.09}$, a_{1}=0.31, and c_{1}=10.0 (Menter et al., 2003).
A force vector, F, on a structure is computed by summing forces from pressure, F_{p}, and from viscous stress, F_{v}.
F_{p} and F_{v} are calculated respectively by
where i is the index of cell faces on the building on which forces need to be evaluated; p_{i} is the total pressure on face i; A_{i} is area of face i; ${\left({\mathit{\alpha}}_{\mathrm{water}}\right)}_{i}$ is volume fraction of water in the adjacent cell of face i; n_{i} is the unit normal vector of face i pointing into the computational domain; and τ_{i} is the viscous stress tensor at face i, which can be expressed by ${\mathit{\tau}}_{i}=\left\{\mathit{\rho}\left(\mathit{\nu}+{\mathit{\nu}}_{\mathrm{t}}\right)\left[\mathrm{\nabla}\mathbf{U}+\mathrm{\nabla}{\mathbf{U}}^{T}\right]\right\}$ on face i.
An initial comparison of the two numerical models was conducted by modeling the interaction between a bore and a freestanding coastal structure, with experimental results from Árnason (2005). The experiment was performed at the Charles W. Harris Hydraulics Laboratory at UW, Seattle. In the experiment, a square column was placed in a 16.6 m long, 0.6 m wide, and 0.45 m deep wave tank, and aligned in parallel to the tank side walls (Fig. 1).
A thin gate separated water in the tank into two parts with different depths: 0.02 m deep on the square column side and 0.25 m deep on the other side. When the gate was lifted to the top of the tank in 0.2 s by a 6.4 cm diameter pneumatic piston, a bore formed and propagated toward the square column downstream. The square column with a 12×12 cm squareshaped cross section was placed 5.2 m downstream from the gate. To measure hydrodynamic forces, the column was supported from above and connected with a force sensor.
Both the threedimensional and twodimensional models were developed at model scale to simulate the physical experiment. The threedimensional OpenFOAM model incorporated the column into the computational domain by simply cutting off a block of mesh of the same shape from the computational domain. The mesh was coarse far from the column (1 cm by 1 cm by 0.5 cm in the x, y, and z directions, where z is the vertical direction) and was refined gradually to 0.125 cm by 0.125 cm by 0.0625 cm in the x, y, and z directions near the column surface. These distances are evaluated to be 70, 70, and 35, respectively, in terms of dimensionless wall distance defined as ${y}^{+}=\frac{y{u}^{\ast}}{\mathit{\nu}}$, where y is the distance, u^{∗} is the friction velocity, and ν is the kinematic fluid viscosity. The mesh was finer in the z directions to better capture the water surface. Forces on the column were obtained by integrating pressure and shear forces from fluid on the surface of the column.
In the twodimensional GeoClaw model, the column was incorporated into the computational domain through the topography term B(x,y) on the righthand side of Eqs. (2) and (3). Values for B(x,y) are set to a very large constant value, h_{c}, in the region of the column and to 0 elsewhere. This prevents water from overtopping the area, thus simulating a column. Setting h_{c} to a very large value also made all four side walls of the square column be more “vertical” in the model since they are represented by steep slopes arising from B=0 (outside the column) to B=h_{c} (inside the column). The coarsest level grid had a resolution of 0.02 m by 0.02 m and covered most of the computational domain; the finest mesh near the column was 0.25 cm by 0.25 cm, which corresponds to a dimensionless wall distance of ${y}^{+}=\mathrm{140}$.
First, a case without the column was modeled. Figure 2 shows predictions of water level history, measured at 5.2 m downstream from the gate (i.e., at x=11.1, the center of the column; see Fig. 1 for location of the gauge) by the two numerical models and the experiment. In general, both 2D and 3D models accurately predict the arrival time of the bore, which is t=3.2 s.
The OpenFOAM model matches the measurement better than GeoClaw with a sharp (but not vertical) slope at the front, a gradually rising surface to the peak near t=8 s, then a downward slope, followed by interactions with the reflected wave from the back wall that creates the second jump in water level at around t=14 s.
OpenFOAM includes water viscosity, which diffuses sharp discontinuities. In contrast, solving the nonlinear shallowwater equations with an initial discontinuity yields a shock wave (discontinuity) propagating to the right as a vertical bore front followed by a region with constant water depth; as a consequence, GeoClaw slightly overestimates the initial height of the bore front, underestimates the height at t=8 s, and presents the reflected wave as a second sharp discontinuity at t=13.1 s.
At the same location, streamwise (the alongchannel direction) components of the velocity at different depths were also predicted. Figure 3 shows time histories of streamwise velocity at nine different distances from the bottom. Note that, since the twodimensional model is depthaveraged, its predicted velocity is constant with depth. The prediction from the twodimensional model matches the measurements very well near the water surface, except for the spike at the front, which is better captured by the threedimensional model. The threedimensional model underestimating flow velocity near the bottom might be due to our nearwall treatment being imperfect. The velocities in the upper region are also hard to predict well because of air entrained in the water near the free surface as well as the fact that the velocimeter may not be immersed in water at times, causing the measurement to oscillate dramatically.
Figure 4 shows a comparison of total forces on the square column from the experiment, the threedimensional model, and the twodimensional model. The force predicted by the threedimensional model was obtained by integrating the pressure and viscous fluid forces on the surface of the column (see Eq. 17). The threedimensional model predicts the force very well in terms of magnitude and is able to capture even the small spike near t=4 s. In the twodimensional model, no hydrodynamic pressure field is available for force prediction. To predict forces from the twodimensional model, data from the previous case without the column were used instead. The water level, h, and streamwise velocity, u, were first sampled at the center of the footprint of the column that was removed from the domain, to compute the momentum flux, M=hu^{2}. As recommended by FEMA P646 (Applied Technology Council, 2012), the hydrodynamic forces on such a structure can be computed as
where C_{d} is the drag coefficient and may be conservatively chosen to be 2.0 as recommended by FEMA P646 (Applied Technology Council, 2012), F_{d} is the streamwise component of the fluid forces, ρ is the density of the fluids, h is the water depth, u is the fluid velocity at the location of the structure, and b is the breadth of the structure in the plane normal to the direction of flow. Note that the hu^{2} term in parentheses is the momentum flux, M.
Note that in the experiment or threedimensional model the water level on the upstream side of the column is different from that on the downstream side of the column. This causes a difference in hydrostatic pressure and thus a hydrostatic force on the column. For this reason, it may be more appropriate to refer to this value as the coefficient of resistance instead of solely as a drag coefficient. Using a drag coefficient of 2.0 overestimates the force by 13 % in general. This is as expected since it is said to be “conservative” according to FEMA P646 (Applied Technology Council, 2012). Figure 4 also shows that, if a drag coefficient of 1.76 is used instead, the force prediction from the twodimensional model matches the measurement more closely.
4.1 The physical experiment
A 1 : 50 scale physical model of part of Seaside, Oregon, adjacent to the Cascadia Subduction Zone (CSZ), was constructed in the Tsunami Wave Basin at the O.H. Hinsdale Wave Research Laboratory at Oregon State University, and a series of experiments were conducted to measure flow velocities and water levels at 31 locations within the modelscale community. For full details on the experiment, one can refer to Park et al. (2013).
The rectangular basin for the experiment is 48.8 m long, 26.5 m wide, and 2.1 m deep. Figure 5 shows the top and side view of the basin. The still water depth at the wave maker is 0.97 m and decreases as it approaches the shoreline. A 0.04 m height (modelscale) seawall was also constructed between all idealized buildings and the shoreline and was parallel to the wave maker. Figures 6 and 7 show the locations of the 31 gauges where water level and flow velocity were measured at a sampling frequency of 50 Hz. The gauges are grouped into four groups – A, B, C, and D (from bottom to top) – and marked by different symbols. Buildings in blue are large commercial buildings like hotels and hospitals. All red buildings are of the same size and represent small commercial buildings. Buildings in yellow are residential structures and are also all the same size.
In the experiment, the pistontype wave maker was designed to generate an initial wave with a wave height of approximately 0.2 m (model scale) at the lower horizontal section of the basin; this is equivalent to 10 m at full scale, which corresponds to a 500year CSZ tsunami for this region (Tsunami Pilot Study Working Group, 2006). Note that this is not a solitary wave but a long singlepeak wave. Experimental measurement of the wave maker speed was fit with a Gaussian function of the form
with β=0.25, t_{0}=14.75, and amplitude A=0.51. Equation (21) was used as input to generate numerical wave in current simulation. The experiment was repeated many times with identical initial conditions. Data from multiple trials were averaged to obtain the results presented here to smooth out stochastic features of the experiment, more details of which were presented in Park et al. (2013).
4.2 Setup of numerical models
4.2.1 OpenFOAM model
In the threedimensional OpenFOAM model, a numerical wave basin was developed to simulate the experiments. It was built at the model scale instead of full scale to exclude scaling effects. This facilitated the comparison between the numerical model and the physical experiment.
To generate the required waves, a numerical wave generator was previously developed in OpenFOAM (Motley et al., 2014), and it was validated against available data from a pair of experiments. Two steps are taken by the numerical wave generator to simulate the wavegenerating procedure of a pistontype wave maker. First, a short subsection of the wave basin adjacent to the wave maker is modeled. This step is conducted with the wave maker as the reference frame, eliminating the need for a moving mesh, and fluid is forced to enter the domain at the wave maker's speed from the end of the domain that is opposite to the wave maker, to simulate the movement of the wave maker. A timevarying acceleration vector field is also embedded in the solver to compensate for the noninertial frame. The second step is to map all field data in this domain (the generated wave) to a full model of the basin with the mapFields utility in OpenFOAM after the wave maker stops moving. Further simulations can then start from here.
One disadvantage of the threedimensional model is that it requires heavy computational resources. Even with four dualeightcore 2 GHz Intel Xeon e52650 machines (64 total processors), it was not possible to model the entire basin. Instead, the entire domain was divided into four different subsections of equal width to predict flow parameters at different groups of gauges (see Fig. 7). For clarity, only the onshore domain is shown in the figure; however, the numerical domain spans the entire 48.8 m from the wave maker to the back wall of the basin. For each simulation, approximately 60 million cells were used, and the solver was run in parallel with the 64 processor cores mentioned above for ∼10 days (including wave generation), which is equivalent to a total CPU time of ∼640 days.
The boundary conditions for each boundary in the numerical wave basin are
listed in Table 1. The term all walls and floor in the
table includes the bottom, side walls, two end walls, and surfaces of internal
buildings. Another term, atmosphere, refers to the upper boundary of
the computational domain. A zeroGradient
boundary condition
specifies zero gradient on the boundary. A fixedValue
boundary
condition sets the value of a quantity to a constant specified value on the
boundary. The velocity field on a wall is set to 0. An inletOutlet
boundary condition is identical to the zeroGradient
boundary
condition if the flux is out of domain but is switched to apply a
fixedValue
boundary condition if the flux is into the domain. The
pressureInletOutletVelocity
condition at the top of the domain is
essentially identical to a zeroGradient
boundary condition in our
current model. On all walls and floor, p_{rgh} is defined such
that there is zero flux, using the fixedFluxPressure
boundary
condition, while the atmosphere was defined with a uniform reference
pressure p_{0} using the totalPressure
boundary condition:
Here p_{rgh} is pressure subtracted by static pressure ρgh, where ρ is the water density, g is the gravitational acceleration, and h is relative depth under the initial free surface. The turbulence quantities near solid walls are obtained with wall functions that model them as functions of distance from the boundary.
Centers of the first layer of cells near the wall are chosen as positions in
the loglaw region of the boundary layer where the wall functions are
applied. A kqRWallFunction
boundary condition can be expressed as
$\frac{\partial k}{\partial n}=\mathrm{0}$ for k on a wall, where n is a unit
normal vector to the wall. An omegaWallFunction
boundary condition
provides a wall function for the turbulence specific dissipation, ω,
with default model coefficients E=9.8, κ=0.41, and C_{μ}=0.09.
It is computed with
where ω_{vis} is the value of ω in the viscous region, and
ω_{log} is the value of ω in the logarithmic region
(Menter and Esch, 2001). The nutUSpaldingWallFunction
boundary condition
for ν_{t} is used for smooth walls. It computes a continuous ν_{t}
profile to the wall based on Spalding's law (Spalding, 1961), which
is essentially a unified law of the wall which works for the viscous
sublayer, buffer layer, and the logarithmic region in a boundary layer.
The initial condition for α_{water} is set to 1 for cells where there is water at the beginning and to 0 for the rest. The initial values of U and p_{rgh} were zero since the flow is initially at rest. Although the fluid is at rest at the beginning, a small value of the turbulent kinetic energy k must be “seeded” in the domain, because the production term in the governing equation of the turbulent kinetic energy k is zero and thus will produce no turbulence if initially k is 0.
In the 3D model, the turbulent kinetic energy k is estimated from the crossshore velocity u_{1}, with a factor of 1.25 to take into account the fact that turbulence is threedimensional (Scott et al., 2005): $k\approx \frac{\mathrm{1.25}}{\mathrm{2}}{u}_{\mathrm{1}}^{\prime \mathrm{2}}$. The velocity fluctuation ${u}_{\mathrm{1}}^{\prime}$ is computed from $I=\frac{{u}^{\prime}}{U}$, where I is the turbulence intensity, ${u}^{\prime}=\sqrt{\frac{\mathrm{1}}{\mathrm{3}}({u}_{\mathrm{1}}^{\prime \mathrm{2}}+{u}_{\mathrm{2}}^{\prime \mathrm{2}}+{u}_{\mathrm{3}}^{\prime \mathrm{2}})}$, and U can be chosen as wave celerity in this case. This approach is the same as Svendsen (1987) and Lin and Liu (1998). Several choices of initial turbulence intensity were tested, and an turbulence intensity of 1 % is chosen. For the specific dissipation rate, ω, $\mathit{\omega}=\frac{\sqrt{k}}{l}$ is used, where l is the turbulent length scale and is set to 7 % of the hydraulic diameter of the channellike computational domain, according to Pope (2000).
Two computational meshes with refinement focused on different regions were used during the simulation to mitigate computational demand. The first mesh was used in the first phase, from the beginning of the simulation to the time when the wave almost started to break. In this mesh, all buildings were removed from the domain, leaving only the flat bottom of the wave basin, which reduced the number of cells needed onshore significantly, allowing for use of a much finer mesh offshore. The mesh size was approximately 0.08 m × 0.08 m × 0.01 m (length × width × height) near the wave maker and was gradually reduced to 0.08 m × 0.08 m × 0.004 m onshore due to changes in topography. Note that the mesh cells onshore seem to have large aspect ratio but there is no water onshore at all during this time period. The second mesh was used until the end of simulation. The buildings were added into the domain, and very fine mesh was generated around the the onshore bathymetry. The mesh size was 0.3 m × 0.015 m × 0.035 m near the wave maker and refined to 0.0075 m × 0.0075 m × 0.0025 m near the flat bottom of the onshore segment of the basin and at the edges and corners of the buildings. It should be noted that this was only the size of a structured background mesh, which was further refined by a factor of 2 and deformed by a mesh tool, snappyHexMesh, from OpenFOAM near the buildings and seawalls to make the mesh accurately represent the complex and irregular geometry of the boundaries. Simulation results from the end of the first phase were mapped to the second phase, and the simulation continued. This strategy is similar to the dynamic AMR feature in the 2D GeoClaw model. Here, however, statically refined meshes were used instead of dynamically refined grids as used in the 2D GeoClaw model. The average Courant number across the entire computational domain during these simulations is approximately 0.01. While this is considerably low for a typical analysis, this is due to the fact that grid sizes vary by several orders of magnitude.
4.2.2 GeoClaw model
With GeoClaw, it is possible to model the entire basin. Thus, the computational domain is a 48.8 m by 26.5 m rectangle. The geometry of the basin bottom and built environment is described by topography files of different resolution, which specify B(x,y) on the righthand side of Eqs. (2) and (3). A typical wall time for one simulation is approximately 6 h with a single core in an Intel® Core™ i74790 CPU processor, which means the CPU time is also 6 h (0.25 day). Note that the computational resources required by the GeoClaw model are only $\mathrm{0.25}\xf7\mathrm{640}\approx \frac{\mathrm{1}}{\mathrm{2500}}$ of what is required by the threedimensional OpenFOAM model in this study.
To generate tsunami waves in GeoClaw, userdefined timevarying boundary conditions can be specified at the inlet of the computational domain, based on data for the wave maker speed s(t) in the physical experiment. The data from the physical experiment can be fit quite well with Eq. (21). However, the way we imposed velocity boundary conditions at a fixed location rather than having a moving boundary, we found better agreement with the observed wave at several offshore wave gauges by setting A=0.6 in Eq. (21), which was therefore used for all simulations.
The AMR feature of GeoClaw was used, with a mesh size for the baselevel grid of 0.5 m (corresponding to 25 m in full scale) in both crossshore direction and alongshore direction. The term crossshore is used to refer to the direction that the wave propagates from the wave maker to the structures onshore, while the direction perpendicular to the crossshore direction is referred to as the alongshore direction. The mesh is refined in the nearshore region up to four levels, with specified refinement ratios: 4 from level 1 to 2, 5 from level 2 to 3, and 2 from level 3 to 4. The finest mesh in the domain with this setup for AMR is 0.0125 m by 0.0125 m (corresponding to 0.625 m in full scale) and eventually covers the entire onshore region. The desired Courant number is set to 0.9 to guarantee the stability of the explicit numerical scheme.
One thing to be noted is that, for both numerical models described above, all coastal structures, including different types of buildings and the seawall, are assumed to be undamaged and thus fixed and rigid during the inundation.
4.3 Comparison of flow parameters
The predicted freesurface elevation, crossshore velocity, and corresponding momentum flux from the two numerical models will be compared and discussed in this section. All experimental data in this study were provided by the NTHMP Mapping and Modeling Benchmarking Workshop: Tsunami Currents (University of Southern California, 2015), and descriptions of the physical experiments to gather the data are provided by Park et al. (2013) and Rueben et al. (2011).
Gauges were positioned as shown in Figs. 5–7. Ultrasonic surface wave gauges were used to measure the free surface. The bore front propagation speed was obtained by analysis of imagery gathered by two highresolution video cameras located above the wave basin (Rueben et al., 2011). Fluid velocity measurements were acquired by acoustic Doppler velocimeter (ADV) only after peaks; air entrainment in the bore at and shortly after the initial impact rendered the ADV measurements inconsistent in repeated trials (Park et al., 2013). Park et al. (2013) then assumed that the propagation speed and fluid velocity at the bore front are equal and fit a secondorder polynomial to that value and ensembleaveraged ADV measurements in this region.
Time histories of the freesurface elevation, crossshore velocity, and corresponding momentum flux at selected gauges are shown in Figs. 8–11. After the peak (initial impact), there appears to be a significant drop in discrepancies between modeled and measured water level and fluid velocity; therefore, the discussion that follows will separately compare the results before and after the peak.
4.3.1 Onshore time series near initial impact
Water level amplitude by OpenFOAM and arrival time by both OpenFOAM and GeoClaw agree fairly well with measurements at many of the gauges in groups A, B, and C, but GeoClaw underestimates the amplitude at many gauges. These differences reflect the challenge of modeling a turbulent and rapidly varying bore front. An additional factor is that the gauges in groups A, B, and C are placed along straight lines, representing roads within the community, whereas those in group D are set behind buildings. As a consequence, flow around group A, B, and C gauges is dominated by flow in the crossshore direction, while flow around group D gauges is more complex and challenging to model.
Fluid velocity experimental values derived by optical means are significantly lower than the modeled OpenFOAM and GeoClaw velocity in many of the 16 cases presented in Figs. 8–11. This is because the optical measurement of the bore front is not necessarily representative of flow velocity. Here the animation of GeoClaw numerical results was analyzed to obtain estimates of 1.3 m s^{−1} for peak velocity: Fig. 12 showed modeled velocity distributions in the bore at two consecutive time steps in the GeoClaw simulation at gauge A4, illustrating that the modeled maximum fluid velocity occurs at some point behind the bore front.
Momentum flux modeled by OpenFOAM and GeoClaw do not agree well with experimental estimates, due to the discrepancies in fluid velocity estimates, discussed above. This is critical, since momentum flux is often used to compute the tsunami forces on structure, as discussed in detail in Sect. 5.
In summary, predictions near the initial impact are challenging for both models, but the threedimensional OpenFOAM model performs better than the twodimensional GeoClaw model because it models turbulence and the variation of velocity with depth.
4.3.2 Onshore time series in postimpact region
Water level agreement among both models and the experimental data is significantly improved after initial impact. Note that some gauges are quite far from the shoreline (for example, gauges A6, B8, and C8), where the inundation depth is very shallow compared to the peak value near the shoreline (less than 20 % of the peak value). Even at these locations, however, both numerical models provide reasonable predictions. It is also of interest that, as noted above, GeoClaw predicts a lower bore front propagation speed than OpenFOAM; as a result, arrival of the OpenFOAM bore front agrees well with experiment, but the GeoClaw bore front is significantly delayed at gauges farther inland, such as B8 and C8 (see Figs. 9 and 10).
Fluid velocity measurements by the ADV are more stable after 30 s, and both OpenFOAM and GeoClaw velocity time series agree much better with the experimental data at gauges in groups A, B, and C. Agreement does degrade significantly in group D, especially in the case of GeoClaw; this is no doubt due to the more complicated fluid flow in the group D environment, behind buildings, compared to the relatively simpler crossshore flow in the street environments of groups A, B, and C (Fig. 7).
Momentum flux from both numerical models are in better agreement with the measurements at most gauges, since water level and velocity agreements are better than in the t<30 s time period.
Figure 13 compares snapshots of the simulation near line A from the two models at three different times. The threedimensional model provides substantial detail about the complex flow among buildings, including the strong channeling effect along line A, aligned with the street, and among the buildings on both sides of the street. These channeling effects can alter the forces exerted on both sides of that street, so that any differences between OpenFOAM and GeoClaw in modeling such effects may result in different prediction of forces on the buildings.
Some representative buildings along line A were selected for preliminary analysis of fluid forces on the coastal infrastructure, as shown in Fig. 14. Building I is one of the two large structures adjacent to gauge A1. It directly faces the shoreline, with dimensions of 0.29 m by 0.78 m by 0.246 m (length in crossshore direction by length in 10 alongshore directions by height). Building II has dimensions of 0.39 m by 0.39 m by 0.091 m. Buildings III and IV, representing small houses within the community, are identical but placed in different directions, which have a length, width, and height of 0.17 m, 0.26, and 0.154 m.
In terms of force measurements, the singlecolumn case presented in Sect. 3 was the only dataset available with experimental measurements of wave impact forces on similar structures. Through validation against that data, it was shown that, provided the water height and fluid velocity are properly modeled, the fluidinduced forces could also be properly predicted. This could be generally extrapolated and applied to the Seaside problem, where the only available measured data included flow parameters (water depth and velocity).
Figure 15 shows predicted forces in the crossshore direction from the two models on selected buildings. Note that these forces are normalized by the width of the western (left) wall of the buildings. Since no pressure field exists in the twodimensional GeoClaw model, the same approach as was used in Sect. 3 is applied here to compute forces on these selected buildings for the GeoClaw model (C_{d} chosen as 2.0 as well). In this case, note that not all the buildings are removed to get the momentum flux for a specific building. Instead, only the building at the center of which the momentum flux is to be predicted is removed, with all other constructed environment unchanged. This minimizes the influence of removing that building on the flow overall.
Peak values of forces predicted by the GeoClaw model on all buildings are only approximately half of those predicted by the OpenFOAM model, except for building III. This indicates that this approach for predicting tsunami load can be off by as large as 100 % in such a complex scenario where multiple objects are present, although it is recommended by FEMA P646 (Applied Technology Council, 2012) as an empirical method when only velocity and surface elevation map is available and is prevalent in tsunami inundation problems.
In this paper, two different types of numerical models of tsunami inundation were developed and compared. They were first validated by comparing water level, velocity profile, and forces on a single column impacted by a bore from a dam break. Then the two models were used to predict freesurface elevation, velocity, and momentum flux of a tsunami inundation on a modelscale constructed environment. The predicted flow parameters agree well with experimental measurements in the postimpact region at most gauges. During initial impact, however, the twodimensional GeoClaw model has difficulty in capturing transient characteristic of the flow. The threedimensional OpenFOAM model can solve this challenge better, albeit at the expense of many more computational resources being required. This is because the variation in the vertical direction is “eliminated” by the integration in the twodimensional model, while all threedimensional characteristics of the flow as well as turbulence are modeled by the threedimensional model. Several primary conclusions can be drawn from this work:

The threedimensional RANS model can predict flow parameters and forces on structures by modeling a subsection of 1∕3 of the width of the entire basin, while the twodimensional NSWE model can model the entire basin at one time, with much less computational resources. Both models agree well with experimental measurements at most locations considered after the initial impact. The RANS model, however, can provide more details on the flow, especially near the initial impact region.

The fluid dynamics in the bore front are transient and turbulent. Thus near the initial impact, prediction of flow parameters and forces is challenging but also the most critical since the flow parameters and forces have maximum value near this point. The threedimensional RANS model solves this challenge better than the twodimensional NSWE model but needs many more computational resources.

Using the approach recommended by FEMA P646 to predict fluid forces on structures from the twodimensional model works well in the simple case of flow around a column but becomes less reliable in a complex constructed environment. Although choosing a drag coefficient of 2.0 is considered conservative, the 2D model with this value was still seen to significantly underestimate fluid forces (in some cases giving only half of the prediction from the 3D model as discussed in Sect. 5) because the 2D model underestimates peak velocities in this complex flow. For this reason, it is recommended that a 3D model should be used to determine the tsunami loads on structures when possible, which eliminates the necessity of choosing a large safety factor when only flow velocity is available from the 2D model as done in Ash (2015).
This research compares different characteristics of a twodimensional model and a threedimensional model of tsunami inundation with the constructed environment. Challenges in prediction of flow parameters and forces are revealed, and the capabilities of the two numerical models in solving this type of problem are analyzed. A tradeoff needs to be made between the two models due to their different levels of accuracy and required computational resources. The comparisons in the current study can provide a reference when choosing between a twodimensional model and threedimensional model to predict required information in tsunami inundation.
The setups for both models described in the paper have been uploaded to https://doi.org/10.5281/zenodo.1419317 (Qin, 2018).
Xinsheng Qin was the primary modeler and author of this work. All coauthors provided assistance in the modeling and subsequent interpretation of results.
The authors declare that they have no conflict of interest.
The authors would like to thank the National Science Foundation for their
financial support through grants EAR1331412 and CMMI1536198. This work was
facilitated through the use of advanced computational, storage, and
networking infrastructure provided by the Hyak supercomputer system.
Edited by: Maria Ana Baptista
Reviewed by: three anonymous referees
Applied Technology Council: Guidelines for Design of Structures for Vertical Evacuation from Tsunamis. Second Edition (FEMA P646), FEMA P646 Publication, https://doi.org/10.1061/40978(313)7, 2012. a, b, c, d
Árnason, H.: Interactions between an incident bore and a freestanding coastal structure, PhD thesis, University of Washington, 2005. a
Ash, C.: Design of a Tsunami Vertical Evacuation Refuge Structure in Westport, Washington, in: Structures Congress, Portland, Oregon, USA, 23–25 April 2015, 1530–1537, 2015. a, b
Bale, D. S., LeVeque, R. J., Mitran, S., and Rossmanith, J. A.: A wave propagation method for conservation laws and balance laws with spatially varying flux functions, SIAM J. Sci. Comput., 24, 955–978, 2003. a
Berger, M. J. and Leveque, R. J.: Adaptive Mesh Refinement Using WavePropagation Algorithms for Hyperbolic Systems, SIAM J. Numer. Anal., 35, 2298–2316, https://doi.org/10.1137/S0036142997315974, 1998. a
Berger, M. J., George, D. L., LeVeque, R. J., and Mandli, K. T.: The GeoClaw software for depthaveraged flows with adaptive refinement, Adv. Water Resour., 34, 1195–1206, 2011. a
Biscarini, C.: Computational fluid dynamics modelling of landslide generated water waves, Landslides, 7, 117–124, 2010. a
Choi, B. H., Kim, D. C., Pelinovsky, E., and Woo, S. B.: Threedimensional simulation of tsunami runup around conical island, Coast. Eng., 54, 618–629, https://doi.org/10.1016/j.coastaleng.2007.02.001, 2007. a
Clawpack Development Team: Clawpack software, version 5.3.1, available at: http://www.clawpack.org (last access: June 2018), 2015. a, b
Einfeldt, B.: On Godunovtype methods for gas dynamics, SIAM J. Numer. Anal., 25, 294–318, 1988. a
Einfeldt, B., Munz, C.D., Roe, P. L., and Sjögreen, B.: On Godunovtype methods near low densities, J. Comput. Phys., 92, 273–295, 1991. a
George, D. L.: Augmented Riemann solvers for the shallow water equations over variable topography with steady states and inundation, J. Comput. Phys., 227, 3089–3113, 2008. a
González, F., LeVeque, R., and Adams, L.: Tsunami Hazard Assessment of the Ocosta School Site in Westport, WA, Technical Report, Washington State Emergency Management Division, 2013. a
Hirt, C. W. and Nichols, B. D.: Volume of fluid (VOF) method for the dynamics of free boundaries, J. Comput. Phys., 39, 201–225, 1981. a
Hu, K., Mingham, C., and Causon, D.: Numerical simulation of wave overtopping of coastal structures using the nonlinear shallow water equations, Coast. Eng., 41, 433–465, https://doi.org/10.1016/S03783839(00)000405, 2000. a
Hubbard, M. E. and Dodd, N.: A 2D numerical model of wave runup and overtopping, Coast. Eng., 47, 1–26, https://doi.org/10.1016/S03783839(02)000947, 2002. a
Larsen, B. E., Fuhrman, D. R., Baykal, C., and Sumer, B. M.: Tsunamiinduced scour around monopile foundations, Coast. Eng., 129, 36–49, 2017. a
LeVeque, R. J., George, D. L., and Berger, M. J.: Tsunami modelling with adaptively refined finite volume methods, Acta Numer., 20, 211–289, 2011. a, b
Lin, P. and Liu, P. L. F.: A numerical study of breaking waves in the surf zone, J. Fluid Mech., 359, 239–264, https://doi.org/10.1017/S002211209700846X, 1998. a
Lynett, P. J.: Effect of a Shallow Water Obstruction on Long Wave Runup and Overland Flow Velocity, J. Waterw. Port C., 133, 455–462, https://doi.org/10.1061/(ASCE)0733950X(2007)133:6(455), 2007. a
Lynett, P. J., Melby, J. A., and Kim, D. H.: An application of Boussinesq modeling to Hurricane wave overtopping and inundation, Ocean Eng., 37, 135–153, https://doi.org/10.1016/j.oceaneng.2009.08.021, 2010. a
Mayer, S. and Madsen, P. A.: Simulation of Breaking Waves in the Surf Zone using a NavierStokes Solver, Proceeding to Coast. Eng. Conference, I, Sydney, Australia, 16–21 July 2000, 928–941, 2000. a
Menter, F.: Zonal two equation kturbulence models for aerodynamic flows, AIAA paper, in: 23rd fluid dynamics, plasmadynamics, and lasers conference, Orlando, Florida, USA, July 1993, p. 2906, 1993. a
Menter, F. and Esch, T.: Elements of Industrial Heat Transfer Predictions, 16th Brazilian Congress of Mechanical Engineering (COBEM) Uberlândia, Minas Gerais, Brazil, 26–30 November 2011, 117–127, 2001. a, b, c
Menter, F. R., Kuntz, M., and Langtry, R.: Ten Years of Industrial Experience with the SST Turbulence Model, Turbulence Heat and Mass Transfer, 4, 625–632, 2003. a
Montagna, F., Bellotti, G., and Di Risio, M.: 3D numerical modeling of landslidegenerated tsunamis around a conical island, Nat. Hazards, 58, 591–608, 2011. a
Motley, M., Lemoine, G., and Livermore, S.: Three Dimensional Loading Effects of Tsunamis on Bridge Superstructures, ASCE Structures Congress 2014, Boston, Massachusetts, USA, 3–5 April 2014, 1348–1358, 2014. a
Motley, M. R., Wong, H. K., Qin, X., Winter, A. O., and Eberhard, M. O.: TsunamiInduced Forces on Skewed Bridges, J. Waterw. Port C., 142, 04015025, https://doi.org/10.1061/(ASCE)WW.19435460.0000328, 2015. a, b
Muhari, A., Imamura, F., Koshimura, S., and Post, J.: Examination of three practical runup models for assessing tsunami impact on highly populated areas, Nat. Hazards Earth Syst. Sci., 11, 3107–3123, https://doi.org/10.5194/nhess1131072011, 2011. a
Ozer Sozdinler, C., Yalciner, A. C., and Zaytsev, A.: Investigation of Tsunami Hydrodynamic Parameters in Inundation Zones with Different Structural Layouts, Pure Appl. Geophys., 172, 931–952, https://doi.org/10.1007/s000240140947z, 2015. a
Park, H., Cox, D. T., Lynett, P. J., Wiebe, D. M., and Shin, S.: Tsunami inundation modeling in constructed environments: A physical and numerical comparison of freesurface elevation, velocity, and momentum flux, Coast. Eng., 79, 9–21, https://doi.org/10.1016/j.coastaleng.2013.04.002, 2013. a, b, c, d, e, f
Pope, S. B.: Turbulent flows, Cambridge University Press, https://doi.org/10.1017/CBO9780511840531 2000. a
Popinet, S.: Adaptive modelling of longdistance wave propagation and finescale flooding during the Tohoku tsunami, Nat. Hazards Earth Syst. Sci., 12, 1213–1227, https://doi.org/10.5194/nhess1212132012, 2012. a
Qin, X.: Model setups for the paper: A comparison of a twodimensional depth averaged flow model and a threedimensional RANS model for predicting tsunami inundation and fluid forces, https://doi.org/10.5281/zenodo.1419317, 2018.
Qin, X., Motley, M. R., and Marafi, N.: Threedimensional modeling of tsunami forces on coastal communities, Coast. Eng., 140, 43–59, 2018. a
Roe, P. L.: Approximate Riemann solvers, parameter vectors, and difference schemes, J. Comput. Phys., 43, 357–372, 1981. a
Rueben, M., Holman, R., Cox, D., Shin, S., Killian, J., and Stanley, J.: Optical measurements of tsunami inundation through an urban waterfront modeled in a largescale laboratory basin, Coast. Eng., 58, 229–238, https://doi.org/10.1016/j.coastaleng.2010.10.005, 2011. a, b
Scott, C. P., Cox, D. T., Maddux, T. B., and Long, J. W.: Largescale laboratory observations of turbulence on a fixed barred beach, Meas. Sci. Technol., 16, 1903, 2005. a
Shi, F., Kirby, J. T., Harris, J. C., Geiman, J. D., and Grilli, S. T.: A highorder adaptive timestepping TVD solver for Boussinesq modeling of breaking waves and coastal inundation, Ocean Model., 43–44, 36–51, https://doi.org/10.1016/j.ocemod.2011.12.004, 2012. a
Shin, S., Lee, K.H., Park, H., Cox, D. T., and Kim, K.: Influence of a Infrastructure on Tsunami Inundation in a Coastal City, Coast. Eng., 1, https://doi.org/10.9753/icce.v33.currents.8, 2012. a
Spalding, D.: A single formula for the “law of the wall”, J. Appl. Mech., 28, 455–458, 1961. a
Svendsen, I. A.: Analysis of surf zone turbulence, J. Geophys. Res.Oceans, 92, 5115–5124, https://doi.org/10.1029/JC092iC05p05115, 1987. a
The OpenFOAM Foundation: OpenFOAM v2.3.1, available at: https://openfoam.org/release/231/ (last access: September 2018), 2014. a, b
Tsunami Pilot Study Working Group: Seaside, Oregon Tsunami Pilot Study – modernization of FEMA flood hazard maps, Technical Report, NOAA OAR Special Report, NOAA/OAR/PMEL, Seattle, WA, 94 pp. + 7 appendices, available at: http://archives.pdx.edu/ds/psu/8333, 2006. a
University of Southern California: NTHMP Mapping and Modeling Benchmarking Wrokshop: Tsunami Currents, Portland, Oregon, 9–10 February 2015, available at: http://coastal.usc.edu/currents_workshop/ (last access: September 2018), 2015. a
Versteeg, H. K. and Malalasekera, W. W.: An introduction to computational fluid dynamics: the finite volume method, Pearson Education Ltd, 2011. a
Wei, Y., Chamberlin, C., Titov, V. V., Tang, L., and Bernard, E. N.: Modeling of the 2011 Japan Tsunami: Lessons for NearField Forecast, Pure Appl. Geophys., 170, 1309–1331, https://doi.org/10.1007/s000240120519z, 2013. a
Williams, I. A. and Fuhrman, D. R.: Numerical simulation of tsunamiscale wave boundary layers, Coast. Eng., 110, 17–31, https://doi.org/10.1016/j.coastaleng.2015.12.002, 2016. a
Wong, H. K.: ThreeDimensional Effects of Tsunami Impact on Bridges, Master thesis, University of Washington, available at: https://digital.lib.washington.edu/researchworks/handle/1773/33666 (last access: September 2018), 2015. a
Yakhot, V., Orszag, S., Thangam, S., Gatski, T., and Speziale, C.: Development of turbulence models for shear flows by a double expansion technique, Phys. Fluids AFluid, 4, 1510–1520, 1992. a