QUANTIFYING NUMERICAL MODEL ACCURACY AND VARIABILITY

On March 11, 2011 the Tohoku tsunami event caused the death of thousands of people and generated billions of dollars in damages. Following this and previous tsunami events, there has been an effort to improve tsunami risk mitigation for coastal communities. Currently, numerical models are being used in a state-of-the-art methodology to estimate tsunami risk from near and far field sources. Model predictions are essential for the development of Tsunami Hazard Assessments (THA). By better understanding model bias and uncertainties and if possible minimizing them, a more accurate and reliable THA will result. In this study we compare runup height, inundation lines and flow velocity field measurements between GeoClaw and the Method Of Splitting Tsunami (MOST) predictions in the Sendai plain. Runup elevation and average inundation distance was in general over predicted by the models. However, both models agree relatively well with each other when predicting maximum sea surface elevation and maximum flow velocities. Model predictions show that the flow velocity increases as the tsunami wave front reaches the shoreline and makes its way inland for a couple of kilometers, contrary to what it is generally assumed. The tsunami models used in this study show much more variability when predicting flow velocity than predicting runup elevations and inundation lines. The results provided in this study will help understand the uncertainties in model predictions and locate possible sources of errors within a model.


INTRODUCTION AND BACKGROUND
On March 11, 2011 a   =9.0 earthquake generated a tsunami 130 km off the coast of Sendai, Japan (Mori et al. 2011).This event is one of the worst in Japan history, killing more than 15,000 people and causing more than $200 billion dollars in damage.Available data shows that in some areas runup elevations reached 40 m and flow velocities reached more than 14 m/s (Mori et al. 2011, Koshimura andHayashi 2012).This event raised the safety concerns of many coastal communities.Tsunami Hazard Assessment (THA) has been redefined to improve hazard mitigation and risk assessment in coastal areas susceptible to tsunamis.This event presents a great opportunity for researchers to model, study and understand the nearshore and onshore hydrodynamics of tsunamis.
Available field survey data and video footage analysis measurements are used to compare model runup and flow velocity predictions.This study includes detail comparison between observations and numerical simulations in Sendai focusing more on the Sendai plain area.At this location the tsunami wave traveled more than 5 km inland with a maximum measured runup of approximately 9.4 m with the average being 2.5 m from Mean Sea Level (MSL) (Mori et al. 2011).The tsunami front velocities measured by Koshimura and Hayashi 2012 at different locations on the Sendai plain ranged from 2-8 m/s.Strong currents are also responsible for the damage and loss of life during tsunami events.It has been shown that currents greater than 3-6 knots can cause minor/moderate damage at ports and harbors (Lynett et al. 2014 andLynett 2012).
Numerical model predictions are important for the development of probabilistic and deterministic methods which are used for the development of Tsunami Hazard Assessments (THA).Generally, methodologies for Probabilistic Tsunami Hazard Assessment (PTHA) use numerical models to predict various products (runup, inundation and flow velocity) which are used to perform a statistical analysis to compute recurrence interval rates at a given location (Geist and Parsons 2006, Gonzalez et al. 2009, Thio et al. 2012and Gonzalez et al. 2013).In a deterministic approach worstcase scenarios are modeled by considering physical limitations of the natural occurrences (Gonzalez et al. 2007).So, better understanding of model bias, variability and uncertainties, will lead to more accurate and reliable THA and PTHA development.This will lead to improvement in risk assessment and hazard mitigation in coastal areas susceptible to tsunamis and better evacuation, maritime, land-use and construction planning can be achieved.

FIELD MEASUREMENTS AND TSUNAMI MODELING
The field survey data published in Mori et al. (2011) and the flow velocity measurements from Koshimura and Hayashi (2012) are used to compare accuracy and reliability of numerical model predictions.In Senadai the maximum measured runup was 9.4 m with the average being 2.5 m and standard deviation of 1.7 m (The 2011 Tohoku Earthquake Tsunami Joint Survey Group, 2011).Only 10 percent of the runup measurements were greater than 5 m.For this study we will be focusing more on the Sendai plain (particularly from 38.13° N to 38.28° N).Along this section of the Sendai plain 50 runup measurements are used from which an inundation line can be obtained.Runup averages increase from 38. 13N to 38.18N and from 38.2º N to 38.23º N. The wave front at the Sendai plain reached more than 5 km inland from the shoreline with the average being 4.2 km.
The grids used in this study are from the M7000 digital contoured bathymetric data and the GSI 10-m digital elevation models (http://fgd.gsi.go.jp/download/).Five nested grids were used in the numerical models.The propagation grid or grid A was the coarsest grid consisting of 3 minutes.Four additional nested grids (1 min, 20 arc sec, 4 arc sec and 1 arc sec (where 1 arc sec ≈ 30m)) were used covering the area of interest.All the grids are referenced to Mean Sea Level (MSL) vertical datum and to the World Geodetic System of 1984 (WGS 84) horizontal datum.
The models used in this study were the Method of Splitting Tsunami (MOST) and GeoClaw.MOST was developed as part of the Early Detection and Forecast of Tsunami (EDFT) project and introduced by Titov andSynolakis (1995 and1998) (Titov and Gonzalez 1997).This model is currently used by NOAA for propagation and inundation forecasting (Titov 2009).GeoClaw, released by LeVeque (1997LeVeque ( , 2002)), is an open source tsunami model approved by the United States National Tsunami Hazard Mitigation Program (NTHMP).Both models have been validated extensively and are commonly used for tsunami predictions.For detail description of equations and methods used in MOST and GeoClaw refer to Titov andSynolakis (1995 and1998) and LeVeque 2002 respectively.
The initial condition used in both models was an initial sea-surface deformation based on Yokota et al. (2011).This source model was determined by carrying out a quadruple joint inversion of the strong motion, teleseismic, geodetic and tsunami datasets.The resulting model has a maximum coeseismic slip of approximately 35 m and a seismic moment of 4.2*10^22 Nm, which yields Mw = 9.0.

RESULTS AND DISCUSSION
The 30m resolution grid was used for the inter-model comparison analysis.The Manning coefficient, n, used in both models for this study was n=0.025 constant throughout the grid.For consistency purposes 7 runup measurements from the field data were removed from the analysis because they seem inconsistent with the other data points, see Table 1.These measurements were very close to the shoreline therefore misleading the "possible" inundation line.A total of 43 measurements were used from the Sendai plain in this analysis.The predicted average distance the wave front traveled inland was 4721m and 5251m from MOST and GeoClaw respectively, compared to 4526m calculated from the field data.Table 2 presents some statistics from the field data and model predictions.The average absolute difference between the field data and the models are also presented in Table 2. GeoClaw overestimated the inundation distance by almost 50% more than MOST in the Sendai plain.GeoClaw algorithm might be more susceptible to the Manning coefficient.Figures 2 shows the maximum free surface elevation during the Tohoku tsunami event for MOST (left) and GeoClaw (right).Both models predict that higher free surface elevations occur at the central part of the Sendai plain, around 38.2° N 140.975°E, with maximum wave amplitudes ranging from 8-12m.Both models agree relatively well when predicting sea surface elevation near the shoreline, with MOST making slightly higher predictions.Furthermore, Figure 3 shows the probability density function of the runup height measurements from the interpolated field data and the models.The models distributions have a similar shape but it can be clearly noted that they overestimate the true value of runup along the coast.Tsunami flow velocities have to be carefully analyzed and considered to better study the tsunami hazard at a particular location.Synolakis ( 2004) stated that currents are more destructive than wave height amplitudes during tsunami events.Figures 4 present the maximum flow velocities predicted by MOST (left) and GeoClaw (right), in the Sendai plain.Both models agree on their predictions and locations of high flow velocities.At 30 meter resolution MOST resolves possible eddies better than GeoClaw formed at 38.26° N 141.05°E.
Figures 5 show a comparison between the two models probability density functions of maximum shoreline flow velocity (top) and 1 meter depth maximum flow velocity (bottom).Since there is no available data at these locations it is very difficult to measure the accuracy of the models but instead we present a comparison between the models predictions.Most of the models flow velocity predictions at the shore line fell between 5-9 m/s and with averages of 7.40 and 7.34 m/s for GeoClaw and MOST respectively.The shapes of the shoreline flow velocity pdfs tend to agree well with very small differences in their average and standard deviation.On the contrary, there is not much agreement between the models when estimating flow velocity at the 1 meter depth.Figure 5 (bottom) shows that MOST tends to calculate slightly higher flow velocities than GeoClaw at the 1 meter inundation depth contour line.The average of the GeoClaw distribution is 1.60 m/s while the average for MOST is 1.63 m/s.Also it is important to point out that the standard deviations are very high relative to their average.These discrepancies between the models might be due to the closeness of the 1 meter contour line to dry land.
It is reasonable to expect a steady decrease in velocity, if using the assumption of a simple beach, as the wave front makes its way inland; however the simulation results provided in Figure 6 (top panel) show otherwise.Between 400m-1200m inland, maximum velocities (again the mean of the maximum velocities across the studied Sendia Plain) appear to be constant.After analyzing numerical output, we observe that small bathymetry/topography features can cause large changes in predicted flow velocities, producing secondary peaks.Analysis of animations from three different transects revealed that small bathymetry/topography features caused big changes in flow velocities.Further inland, results from MOST show a steady decline in velocity from 1200m-2800m with a negative slope of about -0.002 for all grid resolutions.Peak flow velocities fluctuate from 6-17 m/s and standard deviations from 0.9-1.7 m/s.

CONCLUSIONS
This study presents a comparison between field data and model predictions during the 2011 Tohoku tsunami.Runup elevation and average inundation distance was in general over predicted by the models.However, both models agree relatively well with each other when predicting maximum sea surface elevation and maximum flow velocities, with MOST making slightly higher predictions.Low variability was found between the models when predicting runup elevation, inundation distance and shoreline flow velocities.
It is important to note that the two models predict similar maximum velocity at the shoreline (Figure 9), which is typically the most important place to do so.Also, there was a noticeable difference between the models maximum flow velocity predictions at the 1-meter depth contour relative to their standard deviations.These results indicate high variability between the models when estimating flow velocities.This might indicate possible uncertainties in the methods employed by the models.Finally, it was expected that the maximum flow velocity would decrease as the wave front reaches the shoreline and makes its way inland but the contrary was found.The analysis of animations from three different transects revealed that small bathymetry/topography features caused big changes in flow velocities, causing no decrease until few kilometers inland.Also, a valley is located in between the 400 to 1200m inland contour which cause the flow to remain constant.Inaccuracies in bathymetry, topography, and bottom friction can cause big changes in the flow which can lead to errors in the calculations.Some of the uncertainties in flow predictions can be acknowledged if higher resolution grids are used.

Figure 1 (
Figure1(left) shows the field data runup measurements and the predicted runup by both models.Average runup from the 43 field data measurements is 1.91 m with a standard deviation of 0.67 m while the average runup calculated by MOST and GeoClaw is 3.09 and 3.29 respectively.The runup standard deviation for both models was 0.16 m and 0.33 m for MOST and GeoClaw respectively.Most of the models runup results lay approximately 2 to 3 standard deviations away from the mean of the field data.Figure1(right) shows the inundation line predicted by both the models and the field data runup height measurements at the Sendai plain (38.12°N to 38.28° N).Both models provide a reasonably accurate prediction of the inundation line.Figures 2 shows the maximum free surface elevation during the Tohoku tsunami event for MOST (left) and GeoClaw (right).Both models predict that higher free surface elevations occur at the central part of the Sendai plain, around 38.2° N 140.975°E, with maximum wave amplitudes ranging from 8-12m.Both models agree relatively well when predicting sea surface elevation near the shoreline, with MOST making slightly higher predictions.Furthermore, Figure3shows the probability density function of the runup height measurements from the interpolated field data and the models.The models distributions have a similar shape but it can be clearly noted that they overestimate the true value of runup along the coast.

Figure 1 :
Figure 1: Comparison of runup height measurements and inundation line between data, MOST and GeoClaw during the 2011 Tohoku event in the Sendai plain.

Figure 2 :
Figure 2: Maximum tsunami amplitudes (m) predicted by MOST (left panel) and GeoClaw (right panel) in the Sendai plain.

Figure 3 :
Figure 3: Comparison of the runup heights probability density functions between the interpolated field data, MOST and GeoClaw models.

Figure 4 :
Figure 4: Maximum flow velocities predicted by MOST (left panel) and GeoClaw (right panel).

Figure 5 :
Figure 5: (top panel) Comparison between GeoClaw and MOST distributions of maximum shoreline flow velocities and (bottom panel) 1 meter depth maximum flow velocities at the Sendai plain.The distributions use a velocity interval spacing of 0.05 m/s.

Figure 6 :
Figure 6: Inland maximum flow velocities across shore in the Sendai plain, (top panel) comparison of the average flow velocities between GeoClaw and MOST; (bottom left panel) comparison of peak flow velocities; (bottom right panel) comparison of standard deviations.