Tsunami interaction with coastlines and elevation predictions

PDF ( 6 MBs ) ( .pdf )

TSUNAMI INTERACTION WITH COASTLINES
AND ELEVATION PREDICTIONS

By

JAMES ROBERT HOUSTON


A Dissertation Presented to the Graduate Council of
The University of Flori'da
In Partial Fulfillment of the Requlrerents for the
Degree of Doctor of Philosophy









UNIVERSITY OF FLORIDA


1978

















ACKNOWLEDGMENTS


The writer wishes to express his sincere gratitude to his

supervisory committee chairman, Professor Knox Millsaps, and co-

chairman, Professor Joseph Hammack, for their support and efforts to

provide him the opportunity to complete his studies at the University

of Florida.

A special debt of gratitude is owed to Dr. R. W. hhalin, Chief of

the Wave Dynamics Division, Hydraulics Laboratory, U. S. Army Waterways

Experiment Station for his encouragement and active support during

the research and preparation periods of this dissertation.

The writer wishes to thank Professor Umit Un)uata for the insight

he provided and for organizing a supervisory committee prior to his

departure from ths University. Appreciation also is extended to

Dr. D. L. Durham for many helpful discussions at the Waterways Experi-

ment Station.

The research upon which this dissertation is based was funded by

the Pacific Ocean Division of the U. S. Corps of Engineers and the

Federal Insurance Administration of the Department of Housing and

Urban Development through the Office of the Chief of Engineers.

















TABLE OF CONTENTS

Page

ACKNOWLEDGMENTS ...... ................. i

KEY TO SYMBOLS . . . .. v

ABSTRACT . ................... .viii

CHAPTER I: INTRODUCTION .. ........ ..... . 1

Tsunamis ... . ..... .. 1
Objective and Scope of Study . ... . ....... 2

CHAPTER II: TSUNAMI INTERACTION WITH FTHE HAWAIIAN ISLANDS 4

Introduction .... ... .. ...... . .... 4
Governing Equations . . . .. . ..... 7
Finite Element Model . . ............ 12
Numerical Grid . .......... ......... 21
Model Verification ................. 25

CHAPTER III: TSUNAMI INTERACTION WITH CONTINENTAL
COASTLINES . .. .. ... .... ...... . 43

Introduction.. .. . .............. . 43
Governing Equations .. .... ... . . . . 45
Finite Difference Model .... . . .. . . 48
Model Verification ...... . . . . . 51

CHAPTER IV: ELEVATION PREDICTIONS . . . . . 58

Hawaiian Island Predictions . . . ..... .. . 58

Introduction .. .. . . . . . 58
Interpolation of Recent Historical Data . . 61
Time Period Analysis . .. ... . . .. 63
Reconstruction of Historical Data . ... .. . 65
Frequency of Occurrence Distribution . . . 69
Results .................. 71











Page

West Coast of the United States Predictions . . .. 75

Introduction . ... .... . .. ... 75
Tsunami Occurrence Probabilities .. . .... 77
Use of Deterministic Numerical Models . . . 80
Effect of the Astronomical Tides . . . . 85
Comparison with Local Observation Predictions . 88

CHAPTER V: CONCLUSIONS ...... . . . . . 95

APPENDIX: SUPERPOSITION OF LARGE tSUNAMIS AND
ASTRONOMICAL 'ITDES ................. 97

REFERENCES .................. .. 99

BIOGRAPHICAL SKETCH .... ............. . 105


















KEY TO SYMBOLS


a Wave amplitude

b Incident wave amplitude
0

b(m) Amplitude of frequency component j

C Chezy's coefficient
th
C Amplitude of m tidal coefficient

D Period of time

f Frequency of occurrence
-->-
F External force vector

F( ) Functional

g Acceleration of gravity

h W'ater depth

S Average runup, meters
avg
H Hankel function of first kind of order m
n

i (-1)1/2

I Tsunami intensity

j Integer

k Wave number

[K] Coefficient matrix

L Wavelength

m Integer

n Unit normal vector

n Manning's coefficient
ma










n( ) Frequency of occurrence

N Interpolation function

P Pressure
-?-
q Velocity vector

{Q} Total load vector

r Spherical coordinate

t Time

T Wave period

u Velocity component in x-direction

v Velocity component in y-direction

w Velocity component in z-direction

x Cartesian coordinate

Ax Small increment in x

y Cartesian coordinate

Ay Small increment in y

z Cartesian coordinate (vertical direction)



a Constant

a Constant coefficient
m

B Constant coefficient
m

S Area of element

cE Coefficient of horizontal eddy viscosity

E2 Coefficient of vertical eddy viscosity

n Free-surface elevation

















P
SM

p



p(u)
2
o

Tb

0


Velocity potential

Velocity potential

Velocity potential

Velocity potential

Angular frequency

Gradient operator


in region A

in region B

of incident wave

of scattered wave


Spherical coordinate

Dynamic viscosity constant

Response to arbitrary tsunami

Mass density

Phase angle

Tidal variance

Bottom stress

Velocity potential

















Abstract of Dissertation Presented to the Graduate Council of
The University of Florida in Partial Fulfillment of the Requirements
For the Degree of Doctor of Philosophy



TSUNAMI INTERACTION WITH COASTLINES
AND ELEVATION PREDICTIONS

By

James Robert Houston

June 1978

Chairman: K. Millsaps
Cochairman: J. L. Hammack
Major Department: Engineering Science

The development is described of two numerical models that accurately

simulate the propagation of tsunamis to nearshore regions and the

interaction of tsunamis with coastlines. One of these models is a

finite element model which uses a telescoping numerical grid to cover

a section of the Pacific Ocean, including all eight islands of the

Hawaiian Islands. The second model is a finite difference scheme that

uses four rectilinear grids to cover most of the west coast of the

continental United States. The finite element model solves linear and

dissipationless long-wave equations. Such equations govern nearshore

propagation in the Hawaiian Islands since the short continental shelf

of the islands limits the time available for nonlinearities and

dissipation to cause significant effects. The finite difference

model solves nonlinear long-wave equations that include bottom stress

terms which are important for the long continental shelf of the west

viii










coast of the United States. Both models are verified by hindcasting

actual historical tsunamis and comparing the numerical model calcu-

lations with tide-gage recordings. A frequency of occurrence analysis

of tsunami elevations at the shoreline in the Hawaiian Islands is

described. This analysis is based upon local historical data with the

finite element model used to interpolate between historical data re-

corded during the period of accurate survey measurements since 1946.

Historical data recorded at Hilo, Hawaii, and dating to 1837 is used

in conjuction with data recorded since 1946 in Hilo and throughout the

Hawaiian Islands to reconstruct elevations at locations in the islands

lacking data prior to 1946. Frequency of occurrence curves are deter-

mined for the entire coastline of the Hawaiian Islands using these

reconstructed elevations. Since most of the west coast of the United

States lacks local data of tsunami activity, a frequency of occurrence

analysis of tsunami elevations is based upon historical data of tsunami

occurrence in tsunamigenic regions in addition to numerical model cal-

culations. A generation and deep-ocean propagation numerical model is

used to propagate tsunamis with varying intensities from locations

throughout the Aleutian-Alaskan and Peru-Chile regions to a water depth

of 500 meters off the west coast. The nearshore finite difference model

propagates these tsunamis from the 500 meter depth to shore. The

frequency of occurrence of combined tsunami and astronomical tide

elevations is determined by an analysis involving the numerical super-

position of tsunamis and local tides.
















CHIAPT',R I

INTRODUCTION

Tsunamis



Of all water waves that occur in nature, or.e of the most destruc-

tive is the tsunami. The term "tsunami," originating from the Japanese

words "tsu" (harbor) and "nami" (wave), is used to describe sea waves

of seismic origin. Tectonic earthquakes, ie., earthquakes that cause

a deformation of the sea bed, appear to be the principal seismic mech-

anism responsible for the generation of tsunamis. Coastal and submarine

landslides and volcanic eruptions also have triggered tsunamis.

Tsunamis are principally generated by undersea earthquakes of the

dip-slip type with magnitudes greater than 6.5 on the Richter scale and

focal depths less tha n 50 kilometers. They are very long period waves

(5 minutes to several hours) of low height (a few feet or less) when

transversing water of oceanic depth. Consequently, they are not dis-

cernible in the deep ccean and go unnoticed by ships. Tsunamis travel

at the shallow-water wave celerity equal to the square root of the

acceleration due to gravity times water depth even in the deepest oceans

because of their very long wavelengths. This speed of propagation can

be in excess of 500 miles per hour in the deep ocean.

When tsunamis approach a coastal region where the water depth

decreases rapidly, wave refraction, shoaling, and bay or harbor reso-

nance may result in significantly increased wave heights. The great











period and wavelength of tsunamis preclude their dissipating energy

as a breaking surf; instead, they are apt to appear as rapidly rising

water levels and only occasionally as bores.

Over 500 tsunamis have been reported within recorded history.

Virtually all of these tsunamis have occurred in the Pacific Basin.

This is because most tsunamis are associated with earthquakes, and most

seismic activity beneath the oceans is concentrated in the narrow

fault zones adjacent to the great oceanic trench systems that are

predominantly confined to the Pacific Ocean.

The loss of life and destruction of property due to tsunamis have

been immense. 'Ihe Great Ioei Tokaido-Nanhaido tsunami of Japan killed

30,000 people in 1707. In 1868, the Great Peru tsunami caused 25,000

deaths and carried the frigate U.S.S. Waterlee 1,300 feet inland.

The Great Meiji Sanriku tsunami of 1896 killed 27,122 persons in Japan

and washed away over 10,000 houses.

In recent times, three tsunamis have caused major destruction in

areas of the United States. The Great Aleutian tsunami of 1946 killed

173 persons in Hawaii, where runup heights as great as 55 feet were

recorded. The 1960 Chilean tsunami killed 330 people in Chile, 61 in

Hawaii, and 199 in distant Japan. The most recent major tsunami to

affect the United States, the 1964 Alaskan tsunami, killed 107 people in

Alaska, 4 in Oregon, and 11 in Crescent City, California, and caused

over 100 million dollars in damage on the west coast of North America.










Objective and Scope of Study

The objective of the present study is to determine tsunami eleva-

tion frequencies of occurrence at the shoreline for the Hawaiian Islands

and the west coast of the continental United States. Such elevation pre-

dictions are required by the Pacific Ocean Division of the U. S. Corps

of Engineers for use in tsunami flood hazard evaluations for floodplain

management and the Federal Insurance Administration of the Department of

Housing and Urban Development for flood insurance rate calculations.

The current trend of general population migration to coastal

regions and the growth of the tourist industry makes it essential to

know possible tsunami inundation levels even for coastal regions that

are currently undeveloped. There are only a few isolated locations in

the United States that have sufficient local historical data of tsunami

activity to allow accurate predictions solely based upon historical

data. This lack of historical data of tsunami activity on most of the

coastline of the United States makes it necessary to use numerical

models in addition to existing historical data to predict tsunami

elevation frequency of occurrence. However, numerical models have not

been successfully applied in the past to the simulation of actual

tsunami nearshore propagation and interaction with coastlines. There

fore, two numerical models are developed in the present study to

accurately simulate tsunami nearshore activity. Calculations of these

deterministic numerical models are combined with probabilistic analyses

based upon historical data to predict frequencies of occurrence of

tsunami elevations for the Hawaiian Islands and the west coast of the

United States.
















CHAPTER II

TSUNAMI INTERACTION WITH THE HAWAIIAN ISLANDS


Introduction



Early studies of tsunami interaction with islands involved analytic

investigations of the scattering of monochromatic plane waves by a

single island of simple geometric shape. Oner and Hall (1949) deter-

mined the theoretical diffraction pattern for long-wave scattering by

a circular cylinder in water of constant depth by adopting results from

acoustic theory. They found qualitative agreement between this dif-

fraction pattern and smoothed runup heights observed for the 1946

tsunami on the island of Kauai in the Hawaiian Islands. Hom-ma (1950)

determined the diffraction pattern for long wave scattering off a circu-

lar cylinder surrounded by a parabolic bathymetry that extended to a

prescribed distance offshore, beyond which the depth was constant.

Hydraulic models have been used to study tsunami interaction with

single islands that are realistically shaped and surrounded by variable

bathymetry. Van Dorn (1970) studied tsunami interaction with Wake

Island using a 1:57000 undistorted scale model. Jordaan and Adams

(1968) studied tsunami interaction with the island of Oahu in the

Hawaiian Islands using a 1:20000 undistorted scale model. They

found poor agreement between historical measure-ents of tsunami runup

and the hydraulic model data. Scale effects (e.g. viscous effects)

and the effects of the arbitrary boundaries that confine the nodel are










problems that are inherent in hydraulic modeling of tsunamis. The small

heights of tsunamis relative to deep ocean depths is also difficult to

simulate. Jordaan and Adams (1968) modeled waves to a vertical scale

of 1:2000; thus, the waves had heights ten times the normal proportion.

Even with this distortion, waves typically had heights in the model of

only 0.3 millimeters.

Vastano and Reid (1966) studied the problem of the response of

a single island to monochromatic plane waves of tsunami period using

a finite difference numerical model. A transformation of coordinates

allowed a mapping of the island shore as a circle in the image plane.

The finite difference solution employed a grid which allowed greater

resolution in the vicinity of the island than in the deep ocean.

Vastano and Bernard (1973) extended the techniques developed by

Vastano and Reid (1966) to multiple-island systems. However, the

transformation of coordinates technique allows high resolution only

in the vicinity of one island of a multiple-island system. Thus, for

the three-island system of Kauai, Oahu, and Niihau in the Hawaiian

Islands, the two islands of Oahu and Niihau had to be represented by

cylinders with vertical walls whose cross sections were truncated

wedges. Kauai was represented by a circular cylinder with the sur-

rounding bathymetry increasing linearly in depth with distance radially

from the island until a constant depth was attained. A single Gaussian-

shaped plane wave was used as input to the model. No comparisons were

made with historical tsunami data for the three islands.

Hwang and Divoky (1975) presented a simulation of the 1964 tsunami

at Hilo, Hawaii, using a finite difference numerical model. The 1964

Alaskan tsunami was generated and propagated across the deep ocean using










a finite difference numerical model employing a coarse grid. The

Hawaiian Islands were too small to be represented by this numerical grid.

A waveform was calculated at a point in water 5000 meters deep that was

approximately 200 kilometers from where Hilo, Hawaii, would have been if

the Hawaiian Islands could have been modeled. This waveform was used as

input to a finite difference grid with fine grid cells that covered an

area approximately 22 kilometers by 67 kilometers in the immediate

vicinity of Hilo. The water depth along the input boundary of the fine

grid varied from approximately 250 meters to 1500 meters. Hwang's

approach neglected the influence of the Hawaiian Islands and of the

varying bathymetry and depth changes between the deep ocean point and the

fine grid input boundary on the tsunami arriving in Hilo. The numerical

model calculations for the 1964 tsunami were in poor agreement with

the Hilo tide gage recording of this tsunami.

A finite difference model employing a grid covering the eight

major islands of the Hawaiian Island chain was used by Bernard and

Vastano (1977) to study'the interaction of a plane Gaussian pulse with

the Hawaiian Islands. The square grid cells were 5.5 kilometers on a

side and close to the minimum feasible (because of computer time and

cost limitations) size for a rectangular cell finite difference grid

covering the major islands of Hawaii. Yet historical data indicate that

significant variations of tsunami elevations occur over distances much

less than 5.5 kilometers. The islands of Hawaii are relatively small

and very poorly represented by a 5.5 kilometer grid. For example,

Oahu has a diameter of only approximately 30 kilometers and the land-

water boundary of the island has characteristic direction changes










that occur over distances much less than 5.5 kilometers. The off-

shore bathymetery of the islands also varies rapidly, with depth changes

of more than 1500 meters frequently occurring in distances of 5.5 kil-

ometers. Furthermore, if a resolution of eight grid cells per tsunami

wavelength is maintained for tsunami periods as low as 15 minutes, a

5.5 kilometer grid cannot be used to propagate waves into depths less

than approximately 150 meters. The processes that cause significant

wave modifications and subsequent rapid variations of elevations along

the coastline (that are known from historical observations to occur

during tsunami activity in the Hawaiian Islands) probably occur within

this region extending from water at a depth of 150 meters to the

shoreline.

In the present study, a finite element numerical model is used to

propagate tsunamis from the deep ocean to the shoreline of the Hawaiian

Islands. Since finite element techniques allow dramatic changes in

element sizes and shapes, this model has great advantages relative to

finite difference numerical models in accurately representing land

shapes, ocean bathymetry, and tsunami waveforms. Element size can be

large in deep water where bathymetric variations are gradual and wave

lengths are long. As a wave enters shallow water its length decreases

and the elements of the grid can be telescoped to smaller sizes with no

accompanying loss of resolution.


Governing Equations


Typical tsunamis have lengths much greater than the water depths

over which they propagate; for such waves, fluid motions are










approximately two-dimensional, i.e., vertical fluid velocities and

accelerations are small in comparison with horizontal velocities and

accelerations, respectively.

Consider a small element of fluid as shown in Figure 1. The rate

of increase of volume of the element is given by


S{n+h) dxdy = dxdy (2.1)

where,

n = a free-surface elevation,

h = the water depth,

t = time,

and (x,y) are Cartesian coordinates.




7--,4





dy x


Figure 1. Element of Fluid


In Equation (2.1), it is assumed that h=h(x,y), i.e. h has no time

dependence. The net rate of volume flux is given in difference form by

{u(n+h)} x + dx dy {u(n+h)}xdy

+ (v(n+h)}y+dy dx {v(n+h) }dx

+ vn+h)}I + 0(d-x dy2
= dxdy { u(n+h) + y v

= dxdy V q(n+h) (2.2)










where u,v are horizontal velocity components, w is the vertical

velocity component (w=0 here), and q = (u,v,w) is a velocity vector

(a Taylor expansion is used in Equation (2.2) and terms of order

dx2dy2 have been dropped). Equating the rate of increase of volume to

the rate of volume flux yields the following continuity equation:



S+ q (Tl+h) = 0 (2.3)

The governing momentum equation for incompressible and constant

viscosity fluid flow is the Navier-Stokes equation in the form

2q V, V2(
p t + q V = pF Vp + q (2.4)
where

p = the fluid pressure,

p = the mass density,

p = the dynamic viscosity constant,

I = an external force vector.

Since vertical fluid velocities and accelerations are small and

assuming that gravity is the only external force and the fluid is in-

viscid so that the last term of Equation (2.4) can be dropped, Equation

(2.4) becomes in component form


uu u + u 1 u2.
xu + v .. (2.5)
wt X Sy p x


__ yv av I arp
S+ u + v Ip (2.6)
3t 3x Uy p S

0 = g (2.7)
p 1(










WIhen integrated, Equation (2.7) yields (with atmospheric pressure

equal to zero)

p = -pg(S-n) (2.8)

i.e., the pressure is hydrostatic, and

Vp = pgVn (2.9)

Equations (2.5) and (2.6) can be recombined to yield

t + (q V)q = -gn (2.10)

Since along most of the path of propagation of a tsunami the

wave amplitude is much less than the water depth, Equation (2.3) can be

approximated as

+ V (hq) = 0 (2.11)


It can also be shown that the second term in Equation (2.10) is

generally less than the other terms by a factor of the order of n/h

Therefore, for n/h small, Equation (2.10) is approximated by
->
aq
= -gVn (2.12)

Taking a time derivative of Equation (2.12) and substituting Equation

(2.11) into the right-hand side of the subsequent equation yields
2->
t =2 gV(V (hq)) (2.13)
dt

Equation (2.13) is known as the linear long wave equation. It is a

governing equation for small amplitude long waves in an inviscid and

incompressible variable depth fluid and is valid when the assumptions

made in its derivation remain valid.

Assuming that the fluid is irrotational q can be expressed

as a gradient of a scalar function t(x,y,t) known as the










velocity potential, i.e,
q= V (2.14)

Assuming all motions to be simple harmonic in time yields

5(x,yt) = Re [(x,y) e t] (2.15)

where Re [ ] indicates the real part of the braced quantity, and a

is the frequency of the motion. Substituting Equations (2.14) and

(2.15) into (2.13) yields
2
7 (h7V) + 0 = 0 (2.16)

Equation (2.16) is the well known generalized Helmholtz equation.

The wave amplitude n can be related to the velocity potential D

by evaluating the Bernoulli equation at the free surface. The Bernoulli

equation is given by

p g7 + V = 0 (2.17)

however; the last term in Equation (2.17) is generally smaller than the

others by a factor of the order of n/h and can be neglected. Also,

taking p=0 at the free surface where =n Bernoulli's equation

becomes to the same order of approximation as Equation (2.16)

n = -a (2.18)
2 = 0
Several assumptions are made in the derivation of Equations (2.16)

and (2.18) and the adequacy of these assumptions for the particular

problem of tsunami interaction with the Hawaiian Islands must be con-

sidered. An important factor in this consideration is the fact that the

continental shelf of the Hawaiian Islands is exceptionally short.

Thus, during most of the time of tsunami propagation, fluid velocities

are much less than wave propagation velocities since wave amplitudes










are much less than water depths and bottom stress forces as well as

nonlinearities have little time to cause significant effects.

Ursell (1953) and Hammack (1972) have shown that nonlinear effects

grow with time. However, for a typical large tsunami, such as the 1964

Alaskan tsunami; Hammack and Segur (1977) have shown that neither non-

linearity nor frequency dispersion have any significant effect on the

lead wave as it propagates across the ocean. In addition, when the lead

wave reaches a continental shelf a certain time must pass before non-

linearity or frequency dispersion become significant. For a typical

large tsunami, Hammack and Segur (1977) show that linear and nondis-

persive theory holds for the lead wave over a distance of approximately
2h
a ( a is a wave amplitude). Taking the average depth of the region

from the edge of the continental shelf (600 foot depth) to shore as 300

feet and the average amplitude as 3 feet (the amplitude becomes large
2h2
only near shore) yields -- = 12 miles. This distance is much greater
a
than the length of the continental shelf of the Hawaiian Islands

(usually less than 1 nile). Since the lead wave is followed by waves

which have been reflected by land areas in the generation region, linear

and nondispersive theory should be adequate, at least for the major

waves of a tsuna i.


Finite Element Model


Figure 2 shows the problem of tsunami interaction with islands

with the process divided into two domains. One domain (A) includes the

islands and the surrounding variable depth region. The second domain

(B) is a constant depth region extending to infinity. This constant

















/ -BOUNDARY B
(AT INFINITY)


/
/
/
/
/


\
K


REGION B


Figure 2. Regions of Computation.










depth region is realistic since the Pacific Ocean has a remarkably

constant depth beyond the immediate region of the Hawaiian Islands.

The boundary value problem expressed in ter s of these domains becomes

the following:
2
V (, (hV) + = 0 in Region A (2.19)
g
2
V (VI) + O- = 0 in Region B (2.20)


SA = B on 3A (2.21)

(h nA ) (h )B on .A (2.22)
A A


-= 0 on 3C = 3C1+3C2.. .3C (2.23)
1/2 a
lim (r)112 (y-ik)%S = 0 (2.24)
r ->

where

j = number of land masses,

r = a spherical coordinate,

k = the wave number,

S = the scattered wave velocity potential,

i = (-1)1/2

Equations (2.21) an4 (2.22) express the continuity of the velocity

potential and its deriva:ive along the boundary separating the two

domains and Equation (2.23) expresses the impermeability of solid bound-

aries. Equation (2.24) is the Sommerfeld radiation condition that

requires the scattered wave to be an outgoing wave at infinity.

A calculus of variations approach can be used to obtain an Euler-

Lagrange formulation of boundary value problems. This variational

formulation is based upon the principle that of all possible displace-










ment configurations a body can assume that satisfy compatibility and

the constraints of cinematic boundary conditions, the configuration

satisfying equilibrium makes the potential energy a minimum. The Euler-

Lagrange formulation of boundary value problems has often been used in

classical dynamics (Hamilton's principle of least action) and solid

mechanics (structural, solid, and rock mechanics). Of course, it is not

possible to obtain analytical solutions for many engineering problems.

The finite element method is a numerical method first introduced by

Turner et al. (1956) (although some of the underlying ideas were dis-

cussed by Courant (1943)) which has found wide application in solid

mechanics for problems based upon an Euler-Lagrange formulation involving

complex material properties or boundary conditions (see Zienkiewicz,

1971 and Desai and Abel, 1972). The finite element method is a discrete

approximation procedure applicable whenever a variational principle can

be formed. This method has a very short history in fluid mechanic

applications despite advantages over conventional numerical schemes such

as finite difference methods. The first application of the finite

element method to fluid mechanic problems involved a study by

Zienkiewicz and Cheung (1966) of seepage through porous media.

Chen and Mei (1974) studied the problem of the forced oscillation

of a small basin protected by a breakwater and containing floating

nuclear power plants by using the finite element method. They dis-

cretized a small constant depth region using approximately 250 triangular

elements that were all about the same size and shape. Since small ampli-

tude long waves were used as a forcing function, Equation (2.20) was the











governing equation (the water depth being constant in their application)

In the present study an approach similar to that of Chen and Mei is

used to solve the problem of tsunami interaction with islands in a

variable depth region.

For the finite element method to be useful for tsunami problems, it

is necessary to take advantage of two important properties of 'this

algorithm. First, the finite element grid which discretizes a domain

can have grid cells (elements) of arbitrary size and shape. Thus, the

grid can telescope from a large cell size in one section of a variable

depth domain to a very fine cell size in another. Secondly, if the

form of the spatial variation of a physical parameter is known a priori,

the finite element method can employ a very small number of elements in

the discretization used to estimate the exact variation. For example, a

single finite element can be used to solve the problem of the seiching

of a lake for a constant depth rectangular lake if the basic form

(unknown coefficients) of the solution is known. The finite element

method degenerates into the well-known Rayleigh-Ritz variational method

in this case. This property of the finite element method will be im-

portant in modeling the infinite region B shown in Figure 2.

The variational principle for the boundary value problem given

by Equations (2.19) through (2.24) requires a certain functional

F(4) to be stationary with respect to arbitrary first variation of .

The first variation of the functional F(4) for this problem is well-

known and can be expressed as









2 ,6
,F (j) f.[ (hVi) + C- ] ; h;- a B



S/[(I A h- @1 DB A C A h 'A (2.25
A B 'O -

The functional can he derived from this expression for its first vari-

ation (Chen and Mei, 1974) and is


2 2 13 5 )
F(,) = fI [h(VD) + h' ) 3n
A A 2 I mA

Sa( B-I) a I
h, -i h- -
MA A 3A A A


+ / ho (2.26)
A D nA
where 0 is the velocity potential of the incident waves (forcing

function).

The finite element method (unlike the fin-te difference method)

discretizes the domain itself rather than the governing equations.

The region under consideration is partitioned into small regions known

as finite elements and the assemblage of all such elements represents

the original region. Instead of solving the problem for the entire body

in one operation, the solutions are formulated for each constituent

element and combined to obtain the solution for the original region.

Thus, for example, the first term of Equation :2.26) is evaluated on the

element level and the results for all elements are summed. This sum is

equivalent to the original integral over the complete region A.

In order to evaluate the first term in Equation (2.26) on the

element level, it is necessary to approximate the field variable 6










within each element by a polynomial N called a trial or interpolation

function, in the coordinates x and y

For each element

o = {N} {C) (2.27)

where

{0} = (2.28)



is a column vector with ., 0j, k0 the unknown values of the field
3 k
variable at the nodal points (vertices) of a triangular element and


{N} = [N.,N ,Nk] (2.29)

is a row vector with

N. = (a. + bix + c.y)/26ri=1,2,3) (2.30)

ai = xYj xjyk (2.31)

bi = -y j (2.32)

ci = xk (2.33)

and A is the area of the element. The coordinates of the nodal point

i are (xi,yi) and the N. and N' are similarly expressed by the

cyclical permutation of i, j, and k .

A linear interpolation function can be used for tsunamis propa-

gating within region A since elements within this region will always be

much smaller than local wavelengths (long waves are adequately repre-

sented as linear variations over such short distances). The convergence

of the solution is assured as the elements are made infinitesimal

(Strang and Fix, 1973). Since each element has its own polynomial











expansion independent of all other elements, the treatment of the entire

domain is systematically handled by surmming the contributions from each

element.

The expression g- 1 in Equation '2.26) is the scattered wave

velocity potential %S Since region B has a constant depth, :S can

be solved analytically. The governing equation for this case is the

Ilelmholtz equation (Equation (2.20)). The boundary condition at in-

finity requires the scattered wave obey the Sonr.erfeld radiation

condition given by Equation (2.24). The ell-kr.nown solution to this

boundary value problem is as follows:

s = Zi H(kr) (a cos I: sin ) (2.34)
s m=0 m m n
m=0
whore

11 are Hankel functions of the first kind of order m,
m
a and P are constant coefficients,
m m
and 8 is a spherical coordinate.

In the finite element formulation, the infinite region B can

be considered to be a single element with the interpolation function

given by Equation (2.34). Thus, within region A the interpolation

functions are linear, but within region B the interpolation function

is an infinite series involving unknown coefficients. The terms in

Equation (2.34) are oscillatory but have a ronotonically decreasing

modulation; hence, the series can be truncated after a finite number of

terms depending upon the desired accuracy (see next section).

The required stationarity of the functional F(P) with respect to

arbitrary variations of P at each nodal point yields a set of linear

algebraic equations when expressed on the element level that can be

solved using matrix methods. For each specified incident wave period,










the velocity potential 0 can be determined at all nodal points. The

linearized dynamic Free surface boundary condition (Equation (2.18)) can

then be used to calculate the surface elevation r .

The set of linear algebraic equations evolving from the element

level can be written in the simple matrix forn

[K] {4} = (Q} (2.35)

where [K] is the large coefficient matrix, {f} is a column vector

representing a combination of nodal unknowns for the velocity potential

and the coefficient unknowns from the single element covering region B,

and {Q} is the total load column vector related to the forcing function.

In the finite element solution of structural mechanics problems, one

obtains a similar equilibrium relation with [K] being a stiffness matrix,

{)} a nodal displacement vector, and {Q) a nodal force vector (Desai

and Abel, 1972). The coefficient matrix [K] is symmetric, banded, and

complex (complex notation was used to represent waveforns, e.g., the

Hankel functions are complex). For the case of tsunami interaction with

the Hawaiian Islands, [K] will be exceptionally large and sparse (many

zero terms). Of course, only half of the banded part of the matrix

alone needs to be manipulated in the solution process. The banded part

of the matrix [K] is solved using standard Gaussian elimination methods.

The matrix solution used by Chen and Mei (1974) was modified in this

study to take advantage of the sparseness of the matrix in order to re-

duce computational time. Also, the elimination method was modified so

that calculations involved only small blocks of terms at a time with the

remainder of the matrix kept in peripheral storage. This matrix parti-

tioning was necessary for the finite element grid shown in the next










section since the banded part of the matrix had '.-er Sno,000 terms.

Partitioning is possible because the syi.r.eric eSfficient matrix is

positive definite; hence, a solution is possible b.; elimination methods

without pivoting. Without pivoting, cli,.inatior. erforned using one

row, affects only the triangle of elements -ithin the band below that

row. Thus, only two triangular submatrices neecde to be stored in

computer memory at a given time.


Numerical Grid



The finite element grid covering the variable depth region A

is shown in Figure 3. Elements of the grid telescope from large sizes

in the deep ocean to triangles with areas as snall as 0.5 square

kilometers in shallow coastal waters. The grid covers an area of

approximately 450,000 square kilometers. In general, finite element

techniques allow elements to be any arbitrary shape (e.g. quadrilateral)

provided that the ratio of the lengths of the shortest and longest

sides of an element is not extreme. Triangular shapes were used

for convenience in the present study.

The geometric shapes of the eight islands co.-prising the Hawaiian

Island chain and the rapid bathynetric variations surrounding the islands

are modeled very precisely by the finite element grid of Figure 3.

Also, the number of node points along the shorelines of the islands

is very dense. Since wave heights are calculated at node points by

the finite element model, the model can adequately represent the rapid

wave-height variations along coastlines that are known to occur

during tsunami activity in the Hawaiian Islands.





































































Figure 3. Finite Eleiment Grid.










The finite element grid was drawn by hand with the element sizes

dictated by the local water depth. The lengths of element sides were

allowed to be no greater than one-eighth of the local wavelength for the

shortest period waves considered. Such a cell size results in a wave

resolution error of no more than one percent ('ei 1978) There are

computer programs (Cole and Rjech, 1976) that generate finite element

grids and number the nodes to minimize bandwidths; however, such pro-

grams only can be applied to fairly homogeneous domains. The rapid

depth variations in region A require a rapidly telescoping grid that is

too complicated for automatic grid generation techniques.

Each node of the finite element grid must be numbered and this

numbering is extremely important since it determines the bandwidth of

the large coefficient matrix discussed in the last section. The band-

width of this matrix is controlled by the largest difference (among all

elements) in the nodal numbers between two nodes of the same element.

The computational ti:e of the numerical model is proportional to the

cube of the bandwidth of the coefficient matrix.

For a simple region covered by homogeneous elements, the bandwidth

is minimized by numbering nodal points in cyclic sequence following a

spiraling pattern :here one spiral is shielded by the following spiral

from a third spiral. The grid pictured in Figure 3, however, is so

large and complicated that considerable judgement and effort is involved

in determining a numbering pattern that minimizes the bandwidth, thus

resulting in reasonable computational times for the numerical model.

The numerical model also requires the spatial coordinates of all

node points and information concerning how the elements are connected











relative to each other (connectiveness). ]he coordinates were deter-

mined using an electronic digitizer and a simple computer program

developed to relate the internal coordinate system of the digitizer

to a prototype coordinate system through coordinate rotation and

scaling. The elements are numbered in an arbitrary fashion and then

connectiveness information assembled for the numerical model by listing

the three nodes of each element.

The nodal numbering, digitizing, and connectiveness assembly

necessarily introduce errors. For example, the connectiveness assembly

involves determining over 17,000 numbers and placing them on computer

cards. Many connectiveness mistakes do not result in error conditions

which terminate the numerical model's computations. Instead, incorrect

answers are produced by the model. To eliminate all errors in the

input data, computer plotting programs were developed to draw the

grid and number nodes and elements based upon the connectiveness and

node coordinate information used as input to the inuerical model. Such

plots make errors readily apparent.

The number of terms retained in the truncation of Equation (2.34)

was determined by trial and error. The final number of terms was such

that addition of further terns produced negligible effects. It was

found that this occurs for terms approximately 10-5 times smaller than

the initial terms.

The computational time requirements of the finite element model for

the grid shown in Figure 3 are very modest, making it economically

feasible to determine the interaction of arbitrary tsunamis with the

Hawaiian Islands (discussed in next section). The reason that the










computationtal time required by the finite element model is reasonably

small is that the grid uses small elements only in areas where they

are necessary. The grid of Figure 3 has approximately 2500 points,

whereas the finite difference grid of Bernard and Vastano (1977) had

26,000 grid points. Even so, some of the elements of Figure 3 are as

much as sixty times smaller than the finite difference grid cells used

by Bernard and Vastano.



Model Verification



The finite element model is verified in the present study by

comparing numerical simulations of the 1960 Chilean and 1964 Alaskan

tsunamis with tide gage recordings of these tsunamis in the Hawaiian

Islands. These two tsunamis are the only major tsunamis for which

reliable information exists concerning characteristics of the ground

uplift which generated the waves (much more information exists for the

Alaskan source).

A deep ocean recording of a tsunami near the Hawaiian Islands

has never been .sad. Therefore, prototype wave records in the deep

ocean of the 1960 Chilean and 1964 Alaskan tsunamis are not available

for use as input to verify the finite element model. However, a finite

difference numerical model which solves the linear long wave equations

employing a spherical coordinate grid has been developed (Hiwang et al.,

1972), and used in several studies (Houston and Garcia, 1974; Houston

and Garcia, 1975; and Houston et al., 1975b) to generate tsunamis and

propagate them across the deep ocean. This deep ocean numerical model










is used in this study to determine deep ocean waveforms of the 1960 and

1964 tsunamis that were needed to verify the finite element model.

The deep-ocean numerical model simulates the uplift deformation of

the ocean water surface caused by the permanent vertical displacement of

the ocean bottom during an earthquake and the subsequent propagation of

the resulting tsunami across the deep ocean. The permanent deformation

(permanent in the sense that the time scale associated with it is much

longer than the period of the tsunami) and not the transient movements

within the time-history of the ground motion is considered to be the

important parameter governing far-field wave characteristics. Tran-

sient movements occur for a period of the order of tens of seconds,

whereas, tsunami wave periods are of the order of tens of minutes.

Experimental investigations of tsunami generation by !lammack (1972) show

that the transient ground movements do not influence far-field char-

acteristics of resulting tsunamis for spatially large ground displace-

ments occurring over a short period of time.

Hwang et al. (1972) verified the deep ocean model by simulating

the 1964 Alaskan tsunami. The permanent ground deformation in the

source region was taken from Plafker (1969) and a 1/40 by 1/6 grid was

used for the simulation. Since this grid was much too coarse to allow

accurate modeling of shallow-water propagation, a comparison was made

between the numerical model calculations and a gage recording of the

1964 tsunami at Wake Island where Van Dorn (1970) had a pressure trans-

ducer in relatively deep water (800 foot depth) some distance offshore.

The amplitude and length of the first wave of this tsunami were shown

to be in good agreement. (Actually, the grid used by Hwang et al.










(1972), did not extend as far as Wake Island. The deep ocean model

propagated the tsunami to a point in the Pacific Ocean and then simple

refraction techniques -were used to estimate the change in the waveform

during propagation from this point to ;ake Island.)

The simulation of the 1960 and 1964 tsunamis performed in this

study uses a 1/30 by 1/30 grid. Prior simulations of the 1964 tsunami

using both a 1/4 by 1/60 grid and a 1/3' by 1/30 grid demonstrated

negligible differences between use of the two grid spacings (Houston et

al., 1975a). The permanent ground deformation used as input is also

taken from Plafker (1969). The numerical grid was large enough to cover

a substantial portion of the Pacific Ocean including Wake Island. The

numerical model calculations at W'ake Island were found to be similar to

those presented by Hwang et al. (1972). Figure 4 shows a comparison of

the numerical model calculations of this study and the Wake Island

recording of the 1964 tsunami. The lack of agreement in the trailing

wave region may be due to local effects (e.g. wave trapping by the Wake

Island seamount) or reflections from the nearby transmission boundaries

of the numerical grid (these boundaries are not perfectly transmitting).

Figures 5 and 6 show time-histories of the 1960 Chilean and 1964

Alaskan tsunamis in deep water near the Hawaiian Islands calculated

using the deep ocean propagation model. The permanent ground defor-

mation in the source region of the Chilean earthquake of 1960 is taken

from Plafker and Savage (1970). Calculations for the Alaskan tsunami

involve a grid with approximately 30,000 cells and the calculations for

the Chilean tsunami a grid with approximately 90,000 cells.

The finite element numerical model is a time-harmonic solution of

































































































0 0 0
133J


-13









0 .-








C
Uo













OHo
o n










4 0
O CG


0 -0







0
O -C
00




o 0


u-^


0



L-






































I (CJ


LC)







C)





1.).





< 2














0 -~ ,- -


I-A 2dA.L00Qd'SNOI-VA--1'





































































































10 N 0
o d

IJ 3dA1010Od 'SNO!IVA313


*n

o
O *'
- 31





HO


OO
z -



an o




So=
3-
I- 0:3










I-



C
0^

















ot
o o








O
:3



('











the boundary value problem. The response of a groupp of islands to an

arbitrary tsunami can easily be determined within the framework of a

linearized theory using the theory of superposition. For example, an

arbitrary tsunami in the deep ocean can be Fourier decomposed as

follows:

b (t) = i b(a) e-i(t+p[l)] d. (2.36)

where

b = incident wave amplitude

b(w) = amplitude of frequency component L

p(m) = phase angle

If n(x,y,u) is the response amplitude at any point (x,y) along the

island coasts due to an incident plane wave of unit amplitude and

frequency w then the response of the islands to the arbitrary tsunami

time history b (t) is given by


S(x,y,t) = Re [ fb(j)n(x,vy,w) e-it-' ) )d,] (2.37)


Therefore, when n(x,y,') is known for all u the island response to

an arbitrary tsunami can be calculated. Of course, it is not feasible

to calculate the integrals of Equations (2.36) and (2.37) over all

frequencies. Instead, the frequency range must be discretized and the

integrals replaced by sums over a frequency range containing most of the

energy of the tsunami.

Equation (2.36) involves a Fourier decomposition of a time series.

This decomposition is accomplished for the time-histories of the 1960

Chilean and 1964 Alaskan tsunamis (Figures S and 6) using a least

squares harmonic fitting procedure. The time-history of the Alaskan










tsunami in deep water was decomposed into 18 components with periods

ranging from 14.5 minutes to the time length of the record (260

minutes). The variance of the residual (difference between the actual

record and a recomposition of the 18 components) was approximately

0.2 percent of the variance of the record. Therefore, virtually all the

energy of the wave record was contained in the 18 components. The

time-history of the Chilean tsunami in deep water was decomposed into

11 components with periods ranging from 15.5 minutes to the time length

of the record (170 minutes). The variance of the residual was less than

0.1 percent of the variance of the record. For both cases, the original

time-history and a time-history constructed from a recomposition of the

components were virtually indistinguishable.

Equations (2.36) and (2.37) take the following form when discre-

tized:
m (2.3
b (t) i b(n ) e- )} (2.38)
boft) = I b(,n)
n=]

(x,y,t) = Re b(mn) n(x,yn) ei nt+p(n}] (2.39)
n=l

where

n = number of components,
th
a = frequency of n component

The b(n ) tern is determined by the least-squares harmonic

fitting procedure and n(x,y,Yn) is determined by the finite element

numerical model for each frequency n and at each location (x,y).

Therefore, the time-history E(x,y,t), which represents the response of

location (x,y) to the deepwater tsunami time-history given by bo(t),










can be calculated for any location along the coastline of the Hawaiian

Islands.

Figures 7, S, and 9 present comparisons between tide gage re-

cordings and numerical model calculations of the 196- Alaskan tsunami at

Kahului, Maui; Honolulu, Oahu; and Hilo, Hawaii, respectively. The

largest waves recorded at each of these sites are shown. Waves

arriving at later times are all much smaller than those shown. When-

ever tide gage limits were encountered, the recordings were linearly

extended (the tide gage locations are shown in Figure 10).

The wave records shown in these figures are in remarkable agree-

ment, especially considering the fact that the ground displacement of

the 1964 earthquake was not precisely known. Since the Hilo breakwater

was not included in the numerical model, the numerical model calcula-

tions for Hilo are probably too large; however, this breakwater was

undoubtedly highly permeable to the 1964 tsunami. The numerical results

appear to be too large in Hilo by some constant factor, since the tide

gage recording and the numerical model calculations have the same form

with the first wave crest and trough having approximately the same

amplitude and being proportionately smaller than the second crest.

Figure 11 shows a comparison between the tide gage recording of the

1960 Chilean tsunami at Honolulu and the numerical model calculations.

Again the largest waves recorded are shown (the Honolulu gage was the

only tide gage in the Hawaiian Islands not destroyed by the 1960

tsunami). The permanent ground motion of the 1960 earthquake is not

known nearly as well as the ground motion for the 1964 earthquake.


~































; // 1 \\
-I 1


I-1 \







GAGE LIMIT GAGE LIIT
'\ i 'I










-7 I I
-S
4(





8.45 9:00 10:00
APPROXIMATE HOURS GREENWICH MEAN TIME
Figure 7. 1964 Tsina'mi from \Aaska Recorded at Kahului, Mui



























































B:45 9:00 10:00
APPROXIMATE HOURS GREENWICH MEAN TIME


Figure 8. 1964 Tsunami from Alaska Rcorded at .Ionolulu, Oabu.


















































U 0

1
Z



l'-

-2


-3


-4
LEGEND
S -- TIDE GAGE RECORD
-5
SI NUMERICAL MODEL
CALCULATIONS
-6



-7

-e
9:00 10.00
APPROXIMATE HOURS GREENWICH MEAN TIME




Figure 9. 1961 Tsunami from Alaska Recorded at iilo, Hawaii.


















































HILO


HAWAII


1:igirc 10. (;;e Lo locations.
































































































S-


i '33NV1SIC


N N NO m
I I











Hence, the agreement indicated in Figure 11 is quite good considering

the uncertainty of the uplift in the source region.

Tide gages, of course, do not record tsunamis perfectly since they

are nonlinear devices with responses which depend upon both the period

and amplitude of the disturbances they measure. However, simple

calculations based upon the paper of Noye (1970) show that the distortion

is small for the tsunamis shown in Figures 7 through 9 and Figure 11.

Figure 12 shows the response amplitude as a function of wave period

at the locations of tide gages in Kahului, Honolulu, and Hilo. The

response versus wave period is greatly dependent upon the precise

location of recording for bays such as Hilo. For example, the Hilo

response (Figure 12) at the tide gage location shown in Figure 10 has

a single large resonant peak at 28.9 minutes. This peak occurs at the

fundamental period of oscillation of Hilo Bay and, consequently, all

locations in the bay have a resonant peak at this incident wave period

(the amplitude of the peak depends on location in the bay). Other

locations in Hilo Bay also have significant peaks at lower periods.

Loomis (1970) calculated a free oscillation fundamental period (node

assumed at bay mouth) of 26.2 minutes for Hilo Bay using a finite

difference grid spacing of approximately one-half a kilometer.

Bernard and Vastano (1977) show a contour nap of a normalized

energy parameter versus wave period and location around the island of

Hawaii. The map shows Hilo displaying an energy peak from 10 to 16

minutes centered at approximately 12.5 minutes; no significant peak at

28.9 minutes is observed. The difference between Bernard and Vastano's

calculations and the results presented in this paper is attributed to




















KAHULUI TIDE GAGE


HONOLULU TIDE GAGE


HILO TIDE GAGE


10 20 30 40 50 60
PERIOD,MIN

Figure 12. Amplification Fnctor Versus WI'ave Period.










the poor resolution of the spatial features of Hilo Bay by their finite

difference grid spacing of 5.5 kilometers. HilD Bay has the shape of a

right triangle with sides that are approximately 9 kilometers and 13.5

kilometers long. The grid used by Bernard and '.astano represented Hilo

Bay as a right triangle with sides that were 5.3 kilometers and 11

kilometers long; thus, the surface area of Hilo Bay was approximately

one half its actual area. This distortion of the size of Hilo Bay and

also, perhaps, the distortion of the bathymetry of the bay by the coarse

grid accounts for the lower response periods for Hilo Bay presented by

Bernard and Vastano.

The importance of accurately representing island shapes and off-

shore bathymetry and resolving tsunami wavelengths can be illustrated by

considering a study by Adams (1975) that used the numerical model

developed by Bernard and Vastano (197") to interpolate between his-

torical measurements of tsunami elevations prod:ceI by the 1946 Aleutian

tsunami. A single Gaussian wave crest was used as input for the

numerical model. If a resolution of four grid cells per crest wave-

length is maintained, this numerical model cannot propagate waves into

water with a depth less than about 500 feet. The main processes in-

volved in the transformation of the form of a tsunami probably occur in

the region over the continental shelf (depths less than 600 feet).

Consequently, the elevations calculated by Ada:s were small and did not

display the great variation along coastlines known to occur during

historical tsunamis in the Hawaiian Islands.

Adams included the effects of propagation over continental shelves

by multiplying the numerical model elevations by a factor (different

for each island) that forced the best agreement with the historical






412



data. However, even with artifact, agreement with historical data was

poor. For example, on the northeast coast of Oahu, the numerical model

elevations multiplied by the factor resulted in predicted elevations

ranging from 14 to 21.4I feet. Loomis (1976) presents 36 historical

elevations for the 19i6 Aleutian tsunami on the northeast coast of Oahu

and only 7 even fall within this range of elevations. The historical

elevations actually varied from 1 to 37 feet on this coast.

















CHAPTER III

TSUNAMI INTERACTION 1'ITH CONTINEN_-AL COASTLINES

Introduction


Experimental investigations of long wave interactions with

continental coastlines have generally been ore-dimensional and have

involved idealized conditions. Flume tests by Savage (1958) and

Saville (1956) measured runup of periodic wavere on constant slope

beaches ad 'all and Watts (1953), Kaplan (1953), and Kishi and

Saeki (1966) measured runup of solitary waves on constant slope

beaches. All of these flume experiments used waves much steeper

than typical tsunanis due to the difficulties involved in the

generation of small amplitude long waves.

Many theoretical studies of one-dimensional propagation of

weakly dispersive long waves have been perfo-red in recent years.

Tsunamis become weakly dispersive (amplitude and frequency dispersive)

long waves at so.e point during their propagation, although (as was

the case for the Hawaiian Islands) they .-ay not be weakly dispersive

for a time period sufficiently long for related effects to develop.

The proper equations for the investigation of weakly dispersive

long wave propagation over a steeply sloping bottom were first

derived by Mei and LeMehaute (1966) and in a slightly different

form by Peregrine (1967) who also gave the first numerical solution

in the case of a solitary wave using an implicit finite difference










scheme. A more complete solution was later given by Madsen and Mel

(1969) who transformed the equations to characteristic form before

they were solved numerically. Svendsen (1974) studied weakly dis-

persive wave propagation over a gently sloping bottom. Hammack

(1972) made comparisons of long wave generation by bottom displace-

ments in the laboratory and numerical solutions of the Korteweg-

DeVries equation using the numerical approach of Peregrine (1966).

Several studies of one-dimensional tsunami propagation over

continental shelves of constant slope have been performed. Keller

and Keller (1964) solved nonlinear long wave equations numerically

for a constant slope bathymetry. Takahasi (1964) studied tsunami

propagation over a simple bathymetry through analytic solutions of

the linear long wave equation. Carrier (1966) assumed a tsunami

was generated by a point source and studied nonlinear propagation

over a constant slope bathymetry. Heitner and Housner (1970) used

a finite element numerical method to solve nonlinear equations for

propagation up constant slope beaches. More recently, Mader (1974)

studied tsunami propagation over linear slopes using a finite

difference solution of the complete Navier-Stokes equations (not

vertically integrated).

Occasional two-dimensional studies have been performed for

tsunami propagation in nearshore regions. Grace (1969) studied

tsunamis affecting the reef runway of the Honolulu International

Airport using a hydraulic model. However, the wavemakers used were

only a small fraction of a tsunami wavelength from the region to be

studied and waves reflected by land were reflected again almost










immediately by the wavemakers. This is a common problem in hydraulic

model tests involving long waves. The wavemakers can be moved a

few wavelengths from the region of interest by reducing the model

scale; however, scale effects may then become very significant.

Wave absorber screens may also be placed in front of the wavemakers;

however, it is very difficult to absorb very long waves and such an

approach is often not very helpful. Similar problems with re-

reflected waves can occur in numerical studies. The calculations

of Hwang and Divoky (1975) were in poor agreement with the tide gage

recording of the 1964 tsunami at Hilo not only for the reasons

presented in Chapter II, but also because their algorithm allowed

waves to leave the numerical grid efficiently only for waves approach-

ing the input boundary normally. Since Hilo Bay is triangularly

shaped, the reflected waves approached the input boundary at various

angles and were partially re-reflected toward Hilo.


Governing Equations


The governing equations for long wave propagation over the

continental shelf of the west coast of the United States are Equations

(2.3) and (2.4). Thus, the long wave assumption is made once again

with vertical velocities and accelerations neglected relative to

horizontal components. The continuity equation is given by the

equation

-+ V q(n + h) = 0 (3.1)

and the momentum equation is

p t + q v = pF Vp +V (3.2)











If the only external force is gravity and the pressure is hydro-

static, then Equation (3.2) becomes


p -+ q Vq = pgVn + Iv~q (3.3)
It

The west coast of the United States has a much broader continental

shelf than that of the Hawaiian Islands; thus, the nonlinear advective

and bottom stress effects may be more significant, and the nonlinear

and viscous terms of Equation (3.3) must be retained.

The dynamic viscosity p in the viscous term of Equation (3.3)

is replaced by an eddy viscosity constant E for the case of

turbulent flow. Using this approach the viscous term can be expressed

as
2t 2 + 2

Ix ly Iz


where the coefficients s and E2 represent coefficients of

horizontal and vertical eddy viscosity, respectively. The term

7 2

Sx 3y
2
is usually small in comparison with the E 0- term (Dronkers,
z2 2
1964) and will be neglected in this study. If C1 s
az
integrated with respect to z from the bottom to the surface

and divided by the distance from bottom to surface n + h we find


c2 au (3.5)
n+-h Ia z=-h s z =

The terms e2-J zh and c2(- ) are components of the
SsZ z=-h n Iz 2 zZ=T










tangential stresses at the sea bottom and surface, respectively.

The shear stress at the sea surface is quite important for storm

surge problems where wind velocities are very large but can be

neglected relative to the bottom stress for tsunamis.

From experiments and theoretical considerations for one-

dimensional flow in a river, it has been shown (Dronkers, 1964), that

the bottom stress Tb can be expressed as follows


T= -pguul (3.6)
b 2


where C is Chezy's coefficient and is given by the expression


C = 8 h1/6 (3.7)
n


where n is Manning's n Thus, Equation (3.3) may be written

->+ 2 2 1/2
p Vq = pgV pgq(uv) (3.8)
C (n+h)


Equation (3.8) contains the nonlinear advective term but neglects

frequency dispersion since vertical accelerations have been neglected.

This equation adequately describes wave propagation for Ursell

numbers much greater than unity (Ursell, 1953). The Ursell number

is proportional to the ratio of nonlinear effects to linear effects
n2
(i.e., frequency dispersion) and is equal to where L is a
1/2 h
typical wavelength. Since L=(gh) 12T for long waves, where T is
2
the wave period, an Ursell number may be chosen as gn(-)2

The numerical model described in the next section uses Equation

(3.8) to propagate waves from depths of 500 meters to shore. These










waves have periods of the order of 1800 seconds and initial have

amplitudes of the order of a foot. 'Thus, the Ursell number, gn( )2

is approximately 30 and increases as the wave enters shallow water

with h decreasing and n increasing.

Mader (1974) found that the one-dimensional form of Equation

(3.8), without the viscous term, produced very similar results to a

one-dimensional Navier-Stokes equation (viscosity neglected) for

propagation of long wavelength tsunamis over linear continental

slopes. The pressure term in Mader's solution was not hydrostatic;

thus, frequency dispersion was permitted.



Finite Difference Model



To solve Equations (3.1) and (3.8) numerically, the differential

equations are discretized and replaced by a system of finite difference

equations using central differences on a space-staggered grid. The

space-staggered scheme describes velocities, water levels, and depths

at different grid points. This scheme, first used by Platzman (1958)

for lake surge problems, has the advantage that in the equation for the

variable operated upon in time, there is a centrally located spatial

derivative for the linear term. Leendertse's (1967) implicit-explicit

multi-operational method is employed in determining the solution

for q and n as functions of time. The stability and convergence

of this method is described in detail by Leendertse (1967). This

method for solving finite-difference equations has been applied in

studies analyzing the tidal hydraulics of harbors and inlets (Raney,











1976). Houston and Garcia (1975) used the method to calculate the

decay of a tsunami as it enters and spreads throughout San Francisco

Bay.

The nearshore finite-difference numerical model uses the time-

history calculated by the generation and deep ocean propagation

model described in Chapter II as input in the following manner. A

tsunami is generated in the Aleutian-Alaskan area or the west coast

of South America and propagated across the deep ocean to a 500-

meter depth off the west coast of the United States. Waveforms

calculated at this depth by the deep ocean numerical model are

recorded all along the west coast. These waveforms are then used

as input to the nearshore numerical model, which propagates the

tsunamis from the 500-meter depths across the continental slope and

shelf to shore.

Figure 13 shows outlines of the four numerical grids used to

cover the west coast. The grids have square grid cells two miles

on a side (southern California and some of San Francisco Bay were

not covered by numerical grids because elevation predictions were

made in previous studies (Houston and Garcia, 1974; and Houston and

Garcia, 1975) using simple one-dimensional analytical solutions for

propagation over the continental shelf). The offshore bathymetry

was modeled from the 500-meter contour to shore. Beyond the 500-

meter contour, the ocean was assumed to have a constant 500-meter

depth. The input boundary of each grid was located approximately

one and one-half wavelengths of a 30-minute wave from the shore.

Therefore, at least three typical waves could arrive at the shore

before waves reflected from the input boundary became a problem.




























OREGON


CALIFORNIA


SCALE IN MILES
50 0 50 100












BAJA
CALIFORNIA


Figure 13. Finite Differcnce Grid L.ocations.











The input boundaries of the grids were oriented approximately

parallel to the shoreline since refraction will bend the wavefronts

to such an orientation before they reach the 500-meter contour

(Wilson and Torum, 1968). Lateral boundaries of the grids were

taken to be impermeable vertical walls.



Model Verification



The finite difference nearshore numerical model was also veri-

fied by numerical simulations of the 1964 Alaskan tsunami. The

deep ocean propagation model discussed in Chapter II was used to

generate the 1964 tsunami and propagate it across the deep ocean.

Figure 14 shows surface elevation contours 3-1/2 hours after the

1964 Alaskan earthquake, as calculated by the deep-ocean propagation

numerical model. This figure illustrates the concentration of energy

on the northern California, Oregon, and Washington coasts due to the

directional radiation of energy from the source region. The upper

section of Figure 15 shows a time history of the 1964 Alaskan tsunami

calculated by the deer-ocean model near Crescent City, California, in

a water depth of 500 meters. This waveform was used as input to the

finite difference nearshore model.

The lower section of Figure 15 shows a comparison between a tide

gage recording of the 1964 tsunami at Crescent City and the nearshore

model simulation of this tsunami. The numerical model calculations

agree well with the tide gage recording; periods, phases, and

amplitudes are accurately reproduced. The main disagreement occurs


































L



o 0


\I ,e '

d P L0
i ~I




-! 6





/'i /
N Li"
~ / / L
/ I 2i ,
U) 1 6 C~i

O D 4


1.!i
Z IS/ N


Liii
I, IL -"
:1 L! -
'-I.-"--':~ Li

J I IL o
o U) U)
U) 0 0 ~ L'S


Un o
Uo t0
















I-

I 2.0


0
- 1.0
O


z 0

- i.o-1I

> -1.0
j 0 20 40 60 BO
l TIME, PROTOTYPE MIN


100 120


PROTOTYPE GAGE RECORD
INFERRED PROTOTYPE GAGE
RECORD (FROM WILSON, 1968)
NUMERICAL SIMULATION


8 9 10
APPROX TIME, HRS, GMT


11 12


Figure 15. 1964 Tsunami from Alaska Recorded at Crescent City, California











for the amplitude of the third wave. With regard to this disagree-

ment, it should be noted that there is also a lack of agreement

on the amplitudes of the actual historical third and fourth wave

crests. Figure 15 shows a reconstruction of the 1964 tsunami at

Crescent City inferred by Wilson and Torum (196S) from the prototype

gage record and later survey measurements by Magoon (1965). The

elevation of the fourth wave (not shown) was estimated by Wilson

and Torum to have been a little less than 14 feet above mean lower

low water (MLLW) (the zero elevation shown in Figure 15 is MLLW

datum); however, Wiegel (1965) estimated that the elevation of the

third wave was approximately 16 feet above MLL' and that the fourth

wave attained the highest elevation of 18 or 19 feet above MLLW at the

tide gage location. The numerical model predicted a fourth wave

elevation approximately the same as that estimated by Wilson and

Torum (14 feet above MLLIW). However, reflections from the input

boundary were probably growing in importance in the numerical model

calculations during the arrival of the fourth wave.

There may be several reasons why the maximum wave elevation

predicted by the numerical model does not agree with the historical

record as well as the elevations of the initial waves. First, the

later waves in the deepwater waveform are probably waves which have

been reflected from land areas in the source region. The fairly

large grid cell spacing used in the deep ocean propagation model

to represent the land-water boundary in this region probably distorts

the reflected waveform. Total reflection of waves in this region

also is not completely realistic. Furthermore, the nearshore grid may










not have fine enough cells to completely model the resonance

phenomenon leading to the large third or fourth wave.

Since 1964, there has been considerable speculation concerning

the reasons that the effects of the 1964 tsunami were so great at

Crescent City. .The finite difference numerical models show that

both the directional radiation of energy from the source region

(see Figure 14) and a local resonance caused the relatively large

elevations at Crescent City. Actually, large elevations at Crescent

City were not unique since directional radiation of energy from the

source also caused large elevations along the Oregon and northern

California coasts. Runup 10 to 15 feet above the high tide level

occurred all along the Oregon coast south of the Columbia River

(runup at Crescent City was approximately 15 feet above the high

tide level). However, the Oregon coast is very sparsely populated

and there were few damage reports. The severity of structural

damage at Crescent City, which has a large logging industry, also

was apparently due to the impact of logs carried by the tsunami

(Wilson and Torum, 1968). The nearshore numerical model indicates

that the resonant effects at Crescent City were fairly local,

extending over 2 to 4 miles of coastline. This behavior contradicts

Wilson and Torum's (1968) speculation that a bowl-shaped section of

the continental shelf with a diameter of approximately 50 miles

experienced a resonant oscillation. Historical data support the

numerical model calculations. At the mouth of the Klamath River,

approximately 15 miles south of Crescent City, elevations observed

during the 1964 tsunami were only 2 to 3 feet above normal high










tide. Wilson and Torum's hypothesis of a shelf oscillation would

predict elevations at the mouth of the Klamath River greater than

those that occurred at Crescent City.

Figure 16 shows a comparison between a tide gage recording of

the 1964 tsunami at Avila Beach, California, and the nearshore

numerical model calculations. The elevations recorded by the Avila

Beach tide gage were larger than those recorded at any tide gage on

the west coast except the Crescent City gage. The historical and

calculated waveforms shown in Figure 16 are in good general agreement.

The tide gage record obviously has higher frequency components which

are not predicted by the numerical model; these components may be local

oscillations of water areas which were too small to be accurately

represented by the numerical grid. Important features such as the wave

amplitudes are in good agreement.

The waveforms recorded for the 1964 tsunami by the Crescent

City and Avila Beach tide gages have elevations larger than those

recorded by any other gage on the west coast. Other tide gage

recordings of this tsunami on the west coast (except at Astoria,

Oregon) were in areas not covered by the grids shown in Figure 13.

A comparison with the tide gage at Astoria, Oregon, was not made

because Astoria is approximately 12 miles away from the coast in

the estuary of the Columbia River, and amplitudes of the 1964

tsunami were small at Astoria.































8-







z 2

>0

Lh
J
-2

-4 GAGE ---- -- LIMIT
I LEGEND
-6 \ PROTOTYPE GAGE RECORD
"J --- NUMERICAL SIMULATION
-8 I I
8.5 9 9.5 10 10.5 11
APPROX. TIME, HRS., GMT.


Fitwre 16. 196-i Tsunami from Alaska Recordcd at Availa Beach,
Cal ifornia.

















CHAPTER IV

ELEVATION PREDICTIONS

Hawaiian Island Predictions


Introduction

The Hawaiian Islands, a chain of eight islands as shown in

Figure 17, have a history of destructive tsunamis generated both in

distant areas and locally. The earliest recording of a severe tsunami

in the Hawaiian Islands was in 1837 when a tsunami from Chile reached

an elevation of 20 feet at Hilo, Hawaii, and killed 46 people in the

Kau District of the big island of Hawaii. Prior to 1837, a number

of severe tsunamis undoubtedly reached the islands but, unfortu-

nately no detailed records were kept. Since 1837, there have been 16

tsunamis that have caused significant damage.

Most of the destructive tsunamis in the Hawaiian Islands have

been generated along the coast of South America, the Aleutian Islands,

and the Kanchatkan Peninsula of the Soviet Union. Approximately one-

fourth of all the tsunamis recorded in the Hawaiian Islands have

originated along the coast of South America, while more than one-

half have originated in the Kuril-Kamchatka-Aleutian region of the

north and northwestern Pacific. Tsunamis generated by local seismic

events have caused large runup in the islands, especially on the

southeast coast of the big island of Hawaii. The 1868 tsunami

produced the largest waves of record in the Hawaiian Islands with













A' *



/99 /5 < \ S
:9 :90








,191
/ ? ^
':9-













11 4




a ,g g
A. j -
.9;, :9 -- :










I i h^ I5

:sZ1 I \
^~ ~ Ar
-r ->= <0
:: 'CI











^h r n\ I
'"











s ;*
-:9- '-I









z s
.sJA.
4 99'









N i '* 5
? I ^



S -J 2 o



0 < .9;
^___________-___________ ___________)










60-foot waves reported on the South Puna coast of the island of

Hawaii. The most recent major tsunami in Hawaiian history occurred

on November 29, 1975, when waves generated by an earthquake with

an epicenter on the South Puna Coast may have reached elevations as

great as 45 feet along the southeast coast. The most destructive

tsunami to ever hit the islands in terms of both loss of life and

property destruction was the Great Aleutian tsunami of 1946, which

killed 173 people and produced waves over 55 feet in elevation. Hilo

incurred S26 million in property damage attributable to this tsunami.

The 1960 Chilean tsunami is the most recent distantly generated

tsunami that produced major effects in the Hawaiian Islands. Sixty-

one lives, all at Hilo, were claimed by the tsunami. Damage through-

out the state was estimated to be $23.5 million. Inspection of the

damage at Hilo revealed much evidence of the tremendous forces

developed by the waves. Twenty-ton boulders had been moved hundreds

of feet, asphaltic concrete pavements were peeled from their sub-

base, and hundreds of automobiles were moved and crushed.

The following sections discuss the methods used to obtain

sufficient data to predict runup elevations for the entire coastline

of the Hawaiian Islands. The finite element numerical model discussed

in Chapter II is used to interpolate between local historical data

for the years of accurate survey measurements from 1946 through 1977.

Such a method of interpolation is needed so that elevations can be

predicted for the entire coastline even for locations having no

actual historical data. A method is discussed herein which allows

reconstruction of historical data for the period of time from 1837 to










1946 when elevations for the period from 1946 through 1977 are known.

Such an approach is necessary in order to have sufficient data for

accurate elevation predictions since it is shown herein that reason-

able predictions cannot be made solely based upon known historical

data for the short time span from 1946 through 1977.

Interpolation of Recent Historical Data

Historical data of tsunami activity in the Hawaiian Islands

are, of course, often limited to certain locations. Information on

tsunami activity in the islands prior to the 1946 tsunami is concentrated

in Ililo (Hawaii) and to a lesser extent in Honolulu (Oahu) and

Kahului and Lahaina (Maui). Even data for tsunamis occurring during

the period of accurate survey measurements from 1946 through the

present are absent or fragmentary (i.e., data exist only for certain

of the events) for most of the coastline of the islands. Therefore,

it is necessary to rely on more than just available historical data

to determine tsunami occurrence frequencies.

The numerical model is used to fill in historical data gaps for

tsunamis from 1946 through 1977 by providing relative responses of

the Hawaiian Islands to tsunamis. For example, although the deepwater

waveform of a tsunami such as that of 1946 is not known, the direction

of approach of this tsunami and its range of wave periods are known.

The average response of the islands to a tsunami similar to the 1946

tsunami is determined by inputting sinusoidal waves of unit amplitude

from the direction of the 1946 tsunami into the numerical model over

a band of wave periods. If historical data exist at one location,

wave elevations at a nearby location for which historical data do

not exist are determined by multiplying the historical data at the










first location by the ratio of the response calculated by the numer-

ical model at the second location and the response calculated at the

first. The numerical model takes into account the major processes

that would cause different relative wave elevations at the two

locations. That is, the model calculates shoaling, refraction,

diffraction, reflection, resonance, shielding of the backside of an

island by the front side, and interactions betwee-n islands. The

historical data account for the unknown absolute heights of the

deepwater tsunamis. The "correction factors" provided by the

historical data are local instead of global (single factor for all

islands) because there are certain local factors that cannot be

determined by the numerical model. For example, coral reefs extend over

parts of some of the islands and are known to protect against tsunami

attack (e.g., the extensive reefs in Kaneohe Bay, Oahu, which are

powerful tsunami dissipators). The effects of the reefs on tsunanis

are implicitly contained in the historical data itself. The historical

data from 1946 through the present were taken from the most recent

and complete compilation of these data (Loomis, 1976).

As an example of the use of the numerical model to fill in

historical data gaps for tsunamis from 1946 through the present,

consider the 1946 tsunami recorded on the island of Lanai only at

Kaumalapua Harbor and Manele Bay. Wave elevations can be calculated

for the 1946 tsunami at other locations on the coastline of Lanai by

multiplying these historical values by the ratios of the responses

calculated by the numerical model at locations of interest and the

responses calculated at Kaumalapua Harbor and Manele Bay. In this










way, elevations for the 1946 tsunami can be determined for all of

Lanai.

Time Period Analysis

It is necessary to use historical data of tsunami occurrence

covering the greatest possible time span to properly determine

frequency of occurrence of tsunami elevations for the Hawaiian

Islands. Tsunami activity has not been uniform in these islands

during recorded history. For example, the two largest and four of

the ten largest tsunamis striking Hilo from 1837 through 1977 occurred

during the 15-year period from 1946 through 1960. Two of the tsunamis

from 1946 through 1960 originated in the Aleutian Islands, one in

Kamchatka and one in Chile. However, six of the ten largest tsunamis

occurred during the 109-year period from 1837 to 1945 with three

originating in Chile, two in Kamchatka and one in Hawaii. Therefore,

both the frequency of occurrence and place of origin of tsunamis

have been remarkably variable. Any study (e.g., Towill Corporation,

1975) basing frequency calculations on a short time span that includes

the period from 1946 through 1960 will predict a significantly more

frequent occurrence of large tsunamis than is warranted by historical

data from 1837 through 1977.

The errors introduced in frequency of occurrence calculations

by consideration only of a short time period that includes the

unrepresentative years from 1946 through 1960 will be greater than

the errors resulting from possible observational inaccuracies of the

19th Century. This is easily demonstrated for Hilo, Hawaii, which

is the location in the Hawaiian Islands having the most complete










data of tsunami activity. Using the data compiled by Cox (1961) for

the ten largest tsunamis in Hilo from 1837 through 1977 (his data

were actually through 1964, but no large tsunamis have occurred in

Hilo since 1964), a least-squares fit of the data was made employing

the logarithmic frequency distribution discussed later in this

chapter. The 100-year tsunami elevation for Hilo calculated using

this method was 27.3 feet. If the elevations recorded during the

nineteenth century (five of the ten largest tsunamis), are increased

by 50 percent (which is much larger than possible observational

error), the 100-year elevation is found to be 30.4 feet. Similar

calculations just based upon data taken during the period of accurate

survey measurements in Hilo from 1946 through 1977 yields a 100-year

elevation of 44.2 feet. Since the largest elevation in Cox's data

for the 141-year period from 1I37 through 1977 was just 28 feet, the

100-year elevation of 44.2 feet is obviously too large. Similar

overestimates will occur at any location in the Hawaiian Islands if

an analysis is based up3n a short period of time that includes the

years 1946 through 1964 since the exceptionally frequent occurrence

of major tsunamis in Hilo from 1946 through 1964 is a property of

the unusual activity of tsunami generation areas and not of special

properties of Hilo.

The unusual tsunamigenic activity in generation regions in

recent years is reflected in the most up-to-date catalog of tsunami

occurrence in the circumpacific area complied by Soloviev and Go

(1969). This catalog shows, for example, that the three greatest

intensity tsunamis generated in the Aleutian-Alaskan region since










1788 occurred in 1946, 1957, and 1964. The catalog also lists the

1960 Chilean tsunami as the greatest intensity tsunami generated in

South America in recorded history (116-year period from 1562 through

1977). Reporting of tsunamis on the west coast of South America is

quite good with,105 tsunamis reported during the 416-year period.

Kanamori (1977) also notes the unusual seismic activity that occurred

during a relatively short time span in recent years. He shows that

the four earthquakes having the greatest seismic moments (true

measure of magnitude for great earthquakes) in this century occurred

during a thirteen year period from 1952 through 1964 (1960 Chilean,

1964 Alaskan, 1957 Aleutian, and 1952 namchatkan earthquakes).

Evidence that tsunamis were unusually active in the Hawaiian

Islands from 1946 through 1964 is also apparent from historical data

at Kahului and Lahaina, Maui. Kahului, Lahaina, luilo, and Honolulu

are the only locations in the islands with historical data extending

to 1837. Four of the five largest tsunami elevations recorded at

Kahului occurred during the years from 1946 through 1964. Three of

the four largest tsunami elevations recorded at Lahaina occurred

during the same period. The historical data for Honolulu is not

very useful because the commonly reported elevations for the 1837

and 1841 tsunamis were, in fact, drops in water level (Pararas-

Carayannis, 1977) and there are no reports for the 1868 and 1877

tsunamis from Chile.

Reconstruction of Historical Data

In order to have data covering a sufficiently long time span to

make reasonable predictions of tsunami elevations for the Hawaiian

Islands, it is necessary to reconstruct historical data for the










years prior to 1946. To do this, it is assumed in this study that

tsunamis generated in a single source region (Kamchatka or Chile but

not the Aleutian-Alaskan region) approach the islands from approximately

the same direction and have energy lying within the same band of

wave periods. The difference in wave elevations at the shoreline in

the Hawaiian Islands produced by tsunamis generated at different

times in the same region is attributed mainly to differences in

deepwater wave amplitudes. For example, the 1841 tsunami from

Kamchatka produced a wave elevation in lilo that was approximately

25 percent greater than that of the 1952 tsunami from Kamchatka. It

is assumed that these two tsunamis had the same relative magnitudes

throughout the islands. Since the finite element model allows a

determination of elevations for the 1952 tsunami along the entire

coastline of the islands, the elevations of the 1S41 tsunami and all

ether Kamchatkan tsunamis occurring during the period of record in

Hilo from 1837 to 1946 can be reconstructed for the entire coastline.

Similarly, Chilean tsunami elevations can be reconstructed knowing

elevations for the 1960 Chilean tsunami. No major tsunamis were

generated in the Aleutian-Alaskan region from 1837 through 1946; hence,

no reconstruction is necessary for tsunamis from this region.

Fortunately, tsunamis from this region since 1946 were representative

of the entire region since these tsunamis were generated in the

western Aleutians (1957), central Aleutians (1946) and eastern

Alaskan area (1964).

The assumption that tsunamis generated in a single source

region (Kamchatka or Chile) approach the Hawaiian Islands from










nearly the same direction is justified in the case of Kamchatka by

the small spatial extent of the known generation area. However, the

tsunamigenic area of Chile is considerably larger than that of

Kamchatka. Still, the direction of approach of different tsunamis

from Chile is nearly the same because the tsunamigenic region in

Chile subtends a relatively small angle with respect to the islands

as a result of the great distance between Chile and the Hawaiian

Islands and the axis of the Chile Trench maintains a constant orienta-

tion relative to the islands (thus directional radiation effects are

not important).

The assumption that tsunamis generated in a single source region

have energy lying within the same band of wave periods is supported

by the elevation patterns produced along coastlines in the Hawaiian

Islands by different historical tsunamis. These patterns are such

that it has been concluded that "tsunamis of diverse geographic

origin are strikingly different, whereas those from nearly the same

origin are remarkably similar" (Eaton, et al., 1961). A recent

study by Wybro (1977) shows that even the distributions of normalized

elevations (i.e. the elevation patterns) produced in the Hawaiian

Islands by different Aleutian-Alaskan tsunamis are nearly the same

yet quite different front the distributions for tsunamis of other

origins. This agreement occurs despite the fact that the relatively

close proximity to Hawaii of the Aleutian-Alaskan Trench and the

varying orientation of the trench axis relative to Hawaii introduces

important directional effects for tsunamis generated in the Aleutian-

Alaskan area. Apparently, these directional effects influence the










magnitude of the elevations in Hawaii but do not greatly alter ele-

vation patterns. Therefore, historical observations support the

approach used in this study which estimates the elevations produced

by tsunamis from Chile or Kamchatka prior to 1946 based upon data

for tsunamis from these tsunamigenic regions recorded during the

years of accurate survey measurements from 1946 to 1977.

In this study elevations are calculated along the coastline

of the Hawaiian Islands for the ten tsunamis which produced the

greatest elevations in Hilo, Hawaii, from 1837 through 1977. The

finite element numerical model fills in all data gaps for those tsunamis

of the ten that occurred from 1946 through 1977 and then the Hilo

data (Cox, 1964) is used to determine relative magnitudes to reconstruct

the elevations for those tsunamis that occurred prior to 1946. For

much of the coast of the Hawaiian Islands, the ten largest tsunamis

since 1837 would be the same as the ten largest tsunamis in Hilo for

the same period. That is, the ten largest would be the tsunamis of

1960, 1946, 1923, 1837, 1877, 1841, 1957, 1952, and the two in 1868.

Of course, the order of the ten largest would vary from location to

location. For example, the largest tsunami on the South Puna coast

of the big island of Hawaii was the locally generated tsunami of

1868.

There are a few locations where historical data show that one

of the ten largest tsunamis was not among the ten largest in Hilo

from 1837 through 1977. For example, the 1896 tsunami generated

near Japan produced small elevations in the Hawaiian Islands except

at Keauhou on the Kona coast of the big island of Hawaii. The










observation at Keauhou was, therefore, included in the frequency of

occurrence analysis for that location. Similarly, all other histor-

ical observations of significant tsunami elevations at isolated

locations were included in the frequency of occurrence analyses for

those locations,. Historical data for tsunamis prior to 1946 were

taken from Pararas-Carayannis (1977). Data for the local tsunami of

November 29, 1975, were taken from the compilation by Loomis (1976).

Since the local tsunami of 1868 had on earthquake epicenter very

near that of the 1975 tsunami and the two apparently produced very

similar elevation patterns (Tilling et al., 1976), historical ele-

vations for the two tsunamis recorded at the same location were used to

determine a relative magnitude of the tsunamis. The many observations

for the 1975 tsunami were then used to reconstruct elevations for the

1868 tsunami.

Frequency of Occurrence Distribution

After the ten largest tsunamis from 1837 through 1977 were

determined at locations all along the coasts of the Hawaiian Islands

using the methods discussed in the previous sections, elevation

versus frequency of occurrence curves were determined at each location

by least-squares fitting of the data by the following expression:



h = -B A log0 f (4.1)

where

h = elevation of maximum wave at the shoreline,

f = frequency per year of occurrence.

Equation (4.1) was used as the frequency of occurrence distribution










since it has been found to agree with historical data at several

locations. Cox (1964) found that the logarithm of tsunami frequency

of occurrence was linearly related to tsunami elevations for the ten

largest tsunamis occurring from 1837 through 1964 at Hilo, Hawaii.

Soloviev (1969) .has shown a similar relationship between tsunami

frequency of generation and intensity for moderate to large tsunamis.

Also, earthquake frequency of occurrence and magnitude have been

similarly related by Gutenberg and Richter (1965); hiegel (1965)

found the same relationship for historical tsunamis at San Francisco

and Crescent City, California and Adams (1970) for tsunamis at Kahuku

Point, Oahu. A recent study by Rascon and Villarreal (1975) revealed

the same relationship for historical tsunamis on the west coast of

America, excluding Mexico.

It is possible that other distributions may agree with the

historical data equally well as the logarithmic distribution. For

example, the Gunbel distribution has been used in the past to study

annual streamflow extremes (Gumbel, 1955). Borgman and Resio (1977)

illustrate the use of this distribution to determine frequency curves

for non-annual events in wave climatology. To investigate the sensi-

tivity of calculations of 100-year elevations on the assumed frequency

distribution, the approach of Borgman and Resio was applied to the

Hilo data of Cox. The Gumbel distribution yields a 100-year elevation

of 28.8 feet for data from 1837 through 1977 and an elevation of 42.5

feet for data from 1946 through 1977. This compares with the elevations

of 27.3 feet and 44.2 feet calculated for the same time periods using

the logarithmic distribution (Equation (4.1)). Clearly, the arguments

used earlier concerning the period of time that must be considered for










a valid analysis are not dependent upon the assumed frequency distri-

bution. Since ample precedent exists for using the logarithmic distri-

bution for analyzing tsunami frequency of occurrence, this distribution

is used in this stud)'.

Results

The methods described in the previous sections allow frequency of

occurrence curves to be constructed for points all along the coast-

line of the Hawaiian Islands. The coefficients A and B of Equation

(4.1) contain all necessary information to determine frequency of

occurrence of tsunami elevations at a location. Plots of the coeffi-

cients A and B versus location along the coasts were constructed for

all of the islands. Smooth curves were drawn through these coeffi-

cients to allow elevation predictions at any location. Figures 18

ard 19 show typical plots of A and B for the southeastern coast of

the big island of Hawaii (Figure 20). For example, at Hoopuloa

(location 11 in Figure 20) A = 6 and B = 4; hence, the 100-year

elevation (F = .01) is only 8 feet. The 100-year elevations for

most locations in the Hawaiian Islands are considerably larger.

For such situations involving human life, one is less interested

in the mean exceeda-:e frequency f (the average frequency per year

of tsunamis of equal or greater elevations) than in the chance of a

given elevation being exceeded in a certain period of time. To

calculate risk, one first notes that tsunamis are usually caused by

earthquakes, and earthquakes are often idealized as a generalized

Poisson process (Newmark and Rosenblueth, 1971 and Der Kiureghian

and Ang, 1977). It was assumed by Wiegel (1965) and Rascon and










72


















m
















C~
rj











CD


0





% W O7 CD-



































SS31INOISN3MG~ V



























--- )- N,
N
m















I N






---4---- -- T ~
ii ,U.
0






















- ----










o o u N c o
rc




























BAY


S93-99


PEE P-EE
P A A


HAWA I I


20 NAPOOPOO


Figure 20. Location Map for !awaii.


PT










Villarreal (1975) that tsunamis also follow such a stochastic process.

The probability that a tsunami with an average frequency of occur-

rence of f is exceeded in D years, assuming that tsunamis follow

a Poisson process, is given by the following equation:


P = 1 e (4.2)


Thus, if an acceptable risk is one chance in one thousand that an

elevation be exceeded during a 100-year period, P = .001 and

D = 100 years. Substituting these values in Equation (4.2) yields

the elevation h which has one chance in one thousand of being

exceeded during the next 100 years.


West Coast of the United States Predictions


Introduction

Unlike the Hawaiian Islands, the west coast of the continental

United States lacks sufficient data to allow tsunami elevation

predictions based upon local historical records of tsunami activity.

Virtually all of the west coast is completely without data of tsunami

occurrence, even for the prominent 1964 tsunami. Only a handful of

locations have historical data for tsunamis other than the 1964 tsunami.

However, the Federal Insurance Administration requires information on

tsunami elevations for the entire west coast of the continental United

States; even for the many locations that have no known historical

data of tsunami activity and for coastal areas that are currently

not developed (since these areas may be developed in the future).

The lack of historical data of tsunami activity on the west coast










necessitates the use of numerical models to predict runup eleva-

tions. The Aleutian-Alaskai area and the west coast of South America

were found by Houston and Garcia (1974) to be the tsunamigenic

regions of concern to the west coast of the United States. Both

regions have sufficient data on the generation of major tsunamis to

allow a statistical investigation of tsunami generation. The genera-

tion and deep-ocean propagation model described in Chapter III can

then be used to generate representative tsunamis and propagate them

across the deep ocean. The nearshore numerical model also described

in Chapter III, can then be used to propagate tsunamis from the deep

ocean over the continental slope and shelf to shore.

In this study, only tsunamis of distant origin are considered

in the analysis. Hammack (1972) has shown that near the generation

region of an impulsively generated tsunami, the waveform is dependent

upon details of the time-dependent movement of the ground during the

earthquake. Little is known about the actual time-dependent ground

motion during earthquakes generating tsunamis and this motion cannot

be predicted in advance. Also, there is not enough historical data

concerning locally generated tsunamis on the west coast to allow

predictions of tsunami occurrence. Thus, reasonable predictions on

the west coast of the properties of locally generated tsunamis or

their likelihood of occurrence are not possible at this time.

The probability is not considered very great that a destructive,

locally generated tsunami will occur on the west coast of the conti-

nental United States. Tsunamis are generally produced by earthquakes

having fault movements that exhibit a pronounced "dip-slip," or










vertical component of motion. "Strike-slip," or horizontal displace-

ment, fault movements are inefficient generators of tsunamis.

Faults on the west coast of the United States characteristically

exhibit "strike-slip" motion since the Pacific block of the earth's

crust is moving horizontally relative to the North American block.

The west coast of the United States does not share the characteristics

(ocean trenches and island arcs) of known-generating areas and, in

fact, has not historically been one. Relatively small locally

generated tsunamis have been known to occur on the west coast, but

there are no reliable reports of major locally generated tsunamis.

Heights of 6 feet in the immediate vicinity of the 1927 Point

Arguello earthquake are the largest authenticated heights produced

by local tsunamis on the west coast). There could be a few locations

on the west coast for which locally generated tsunamis pose a greater

hazard than do distantly generated tsunamis because the elevations

produced by distantly generated tsunamis are small. However, pre-

dictions of elevations produced by locally generated tsunamis are

beyond the scope of this study.

Tsunami Occurrence Probabilities

Historical data of tsunami generation must be the basis for an

analysis that considers the probability of tsunami generation in the

two tsunamigenic regions of the Pacific Ocean of concern to the west

coast of the continental United States--the Aleutian and Peru-Chile

Trench regions. A satisfactory correlation between earthquake

magnitude and tsunami intensity has never been demonstrated. Not

all large earthquakes occurring in the ocean even generate noticeable










tsunamis. Furthermore, earthquake parameters of importance to

tsunami generation, such as focal depth, rise time, and vertical

ground displacement, have only been measured for earthquakes occurring

in recent years. Therefore, data of earthquake occurrence cannot be

used to determine occurrence probabilities of tsunamis. Instead,

historical data of tsunami occurrence in generation regions must be

used to determine these probabilities.

In South America, a wealth of information exists concerning

tsunami generation. Reliable data (grouped in intensity increments

of one-half) exist for tsunamis with intensity greater than or equal

to 0 for a 171-year period and greater than or equal to 2-1/2 for a

416-year period (Soloviev and Go, 1969). The intensity scale used

is a modification (by Soloviev and Go) of the standard Imamura-lida

tsunami intensity. Intensity is defined as

S I= log2(2 Ha) (4.3)

This definition in terms of an average runup (in meters) over a

coast instead of a maximum runup elevation at a single location

(used for the standard Imamura-Iida scale) tends to eliminate any

spurious intensity magnitudes caused by often observed anomalous

responses (due, for example, to local resonances) of single isolated

locations.

Using the most recent and complete catalog of tsunami occurrence.

in the Pacific Ocean (Soloviev and Go, 1969) a relationship between

tsunami intensity and frequency of occurrence was determined for the

tsunamigenic trench running the length of the Peru-Chile coast.

Tsunamis with intensity greater than or equal to 0 were considered.










It was assumed that the logarithm of the tsunami frequency of occur-

rence was linearly related to the tsunami intensity. Earlier in

this chapter, it was shown that such a relationship holds between

both earthquake magnitude and frequency of occurrence and tsunami

intensity and frequency of occurrence.

Letting n(I) equal the probability of a tsunami with an

intensity I being generated during any given year and using

statistics for the entire trench along the Peru-Chile coast, a

least-squares analysis results in the following expression:


n(I) = 0.74e-0.631 (4.4)

In using statistics for the entire trench area along the Peru-Chile

coast, it is assumed that the probability of tsunami occurrence is

uniform along the trench. This is a standard assumption for earth-

quake frequency analysis (Gutenberg and Richter, 1965). The tectonic

justification of this assumption lies in the fact that a single Sialic

block or plate of the earth's crust or lithosphere is dipping into

the Peru-Chile Trench (Wilson, 1959). It can reasonably be expected

that the movement of this single plate is similar along its entire

length.

In the Aleutian Trench region, only large tsunamis occurring

in relatively recent years (since 1788) have been recorded due to the

isolation of the area. Assuming an exponential coefficient of -0.71

for this trench area (determined by Soloviev, 1969) as a mean value

for regions of the Pacific with the most data on tsunamis) and using

only the reliable data for large tsunamis (intensity greater than or

equal to 3.5) from Soloviev and Go (1969), the following relationship










is determined by a least-squares analysis:

n(l) = 0.l e-0.711

Again, the probability of tsunami occurrence is assumed to be uniform

along the trench.

Use of Deterministic Numerical Models

To relate the probability distributions of tsunami intensities to

source characteristics, it is assumed that the ratio of the source

uplift heights producing two tsunamis of different intensity (as

defined in the previous section) is equal to the ratio of the average

runup heights produced on the coasts near these tsunami sources. This

ratio is equal to 2 for two tsunamis with intensities I1 and

12

If H is the wave height in the direction parallel to the major
a
axis of length a of a tsunami source with an elliptical shape

(large tsunamis have historically had elliptically-shaped uplifts)

and Hb is tha wave height in the direction parallel to the minor

axis of length b, then experimental research of tsunami generation

has shown that Hb/Ha is approximately equal to a/b (Hatori,

1963). For a large tsunami, Hb can be larger than H by a factor
b a
of as much as S or 6. Thus, the orientation of the tsunami source

relative to the area where runup is to be determined is very important,

i.e., the runup at a distant site due to the generation of a tsunami

at one location along a trench cannot be considered as being repre-

sentative of all possible placements of the tsunami source in the

entire trench region. Hence, the Aleutian and Peru-Chile Trenches

had to be segmented and runup along the west coast of the United






81



States determined for tsunami sources located at the center of each

of the segments.

The spatial size of a tsunami source was standardized since

there is not an apparent correlation between tsunami intensity and

spatial size of a tsunami source. For example, the 1946 Aleutian

tsunami had an uplift region of very small spatial extent, whereas

the 1957 Aleutian tsunami had an uplift region that covered perhaps

the greatest spatial extent of any known earthquake (Kelleher,

1972). Yet the 1946 tsunami had the greater intensity, producing,

in general, greater runup elevations in the near and distant regions.

The standard source used is discussed in detail by Houston and

Garcia (1974). It is elliptical in plan view with a major axis length

of 600 miles and a minor axis length of 130 miles. The vertical dis-

placement increases linearly from a zero elevation on the ellipse

perimeter to a maximum displacement of 30 feet. The vertical dis-

placement has a parabolic crest (concave downward) parallel to the

direction of the major axis. The standard source represents a large

tsunami with intensity 4 on the modified Imamura-lida scale.

Certainly, tsunamis of low intensity may have smaller spatial extents;

however, large tsunamis pose the greatest threat to a distant area

such as the west coast of the United States. These large tsunamis

can be expected to have similar spatial extents, with any spatial

differences being unimportant in the far-field compared with the effects

of source orientation and vertical uplift. Vertical uplifts are assigned

to different intensities in accordance with the convention discussed

at the beginning of this section.










Figures 21 and 22 show the Aleutian Trench divided into 12 seg-

ments and the Peru-Chile Trench into 3 segments. The segments in the

Aleutian Trench were approximately one-quarter -he length of the major

axis of the standard source, whereas the segments in the Peru-Chile

Trench were approximately the length of the major axis of the standard

source. The standard source was centered in each segment such that

the major axis of the source was parallel to the trench axis. Uplift

regions historically have had such an orientation relative to trench

systems. The Aleutian Trench is segmented much finer than the Peru-

Chile Trench since the Aleutian Trench is oriented relative to the

west coast such that elevations produced on the west coast are very

sensitive to the exact location of a source along the Trench. Uplifts

along the Peru-Chile Trench do not radiate energy directly toward the

west coast regardless of their position along the Trench. The Peru and

Chile sections of the Feru-Chile Trench have constant orientations

relative to the west coast of the United States; therefore, eleva-

tions on the west coast -are not very sensitive to source location

within these sections.

In each of the segments of the Aleutian and Peru-Chile Trenches

tsunamis with intensities from 2 to 5 in steps of one-half intensity

are generated and propagated across the deep ocean using the deep-

ocean numerical model discussed in Chapter II. Tsunamis with inten-

sity less than 2 are too small to produce significant runup on the

west coast. An upper limit of 5 was chosen because the largest

tsunami intensity ever reported was less than 5 (Soloviev and Go,

1969). Gutenberg and Richter (1965) indicate that there is an upper

























[- -






o
O
0- 1

H- z

































1 >
< *4























--
Jc*
r3
44J







0

o-
C3



C










1 \ C,

1 \ 44
\ \ 1Cs;


\ \

\ 1 cs)
^\v U


C -











v-f3


6'

















ENTURA


pAC/FIC 0c




lCUlouE


Figure 22. SegmcnLed Chilean Trench.










limit to the strain that car be supporl od by rock before fracture.

Thus, earthquakes only reach certain: s;axJmum inagn i udes and tsunamis

can be expected to have sii ;iar uppi.; limits to intensity. Perkins

(19721 and McGarr (1976) have denor, .fo ':ted thai. future earthquaklce

cannot have seismic moments fmeasuce of earthquake magnitude for

large earthquakes) much larger than hose of earthquakes that occurred

in recorded history.

f'ne waveforms propagated to c.be icst coast by the deep-ocean

model are used as before as input to Ihe nearshore propagation model.

Each waveform is propagated from a '*: cr depth of 5>00 meters to

shore using the nearslore model and c.o, of the grids shown in FPiure

13. Thus, at each grid location ou the shoreline of the wen.t coast

there is a group of 105 waveforms seven waveforms (for intensities

fro- 2 to 5 in one-half intensity i:nciremerts) for each segment of

the Aleutait'n and Peru-Chile Trenchi.. Each of ihese waveforms has an

associated probabilityi equal to the 'irrbi)ijty that a certain

intensity tsunami will be generated ;ji a parltcular segment of a

trench region.

Effect of the Astrononical Tides

The maxiun "still-water" elevation produced during tsunami

activity is the result of a superposition of tsunamis and tides.

Therefore, -:he statistical effect of the astronomical tides on

total tsun mi runup must be include i in the predictive scheme pre-

sented in this study. Since the -;aveforms calculated by the nearshore

model do not have a simple form (e.g. sinusoidal), the statistical

effect of the astronomical tide or tsunami runup nust be determined






86



through a numerical approach. (The effect of the astronomical tides

was not considered in the section of predictions for the Hawaiian

Islands because the tidal range is quite small for these islands and

the local historical data implicitly contain the effects of the

tides.)

The waveforms calculated by the nearshore numerical model

extend over a period of time of approximately two hours. Three or

four wave crests (the largest waves in the tsunami) arrive during

this time. Smaller waves arriving at later times, however, often

persist for days at a coastal location. An analysis of tide gage

records of the 1960 and 1964 tsunamis on the west coast indicates

that these smaller waves have amplitudes on the average of 40 percent

of the maximum wave amplitude of the tsunami. Therefore, a sinu-

soidal group of these smaller waves were added to each of the calculated

waveforms so that the total waveform extended over a 24 hour period.

These smaller waves are important for locations where tsunamis

are fairly small compared with tidal variations. At such locations

the maximum combined tsunami and astronomical elevation occurs

during the maximum tidal elevation.

A computer program was developed to predict time-histories of

the astronomical tides at all grid locations on the west coast of

the United States. The program was based upon the harmonic analysis

methods used in the past by the Coast and Geodetic Survey for mechan-

ical tide-predicting machines (Schureman, 1948). Tidal constants

available from the Coast and Geodetic Survey were used as input to

the computer program; and a year of tidal elevations was then predicted

for grid locations all along the west coast. The year 1964 was


w~










selected because all the major tidal components had a node factor

equal to approximately 1.00 during this year making it an average

year (node factors can vary by about 10 percent). The node factor

is associated with the revolution of the moon's node and has an

18.6-year cycle. Since a tsunami can arrive at any time during this

18.6-year period (arrival at a low of the node factor is equally

likely as an arrival at a high), the statistical effect of the

varying node factor is small and an average value should be used.

The statistical effect of the varying node factor on the predicted

runup elevations can be shown to be a small fraction of an inch,

using the approach discussed in the next section and Appendix A

(with the variance of the nodal variation equal to approximately

0.1 feet squared).

The year of tidal elevations calculated at each of the nearshore

numerical model grid points along the west coast was then subdivided

into 15-minute segments. The 24-hour waveforms were allowed to

arrive at the beginning of each of these 15-minute segments and

then superposed upon the astronomical tide for the 24-hour period.

The maximum combined tsunami and astronomical tide elevation over

the 24-hour period was determined for tsunamis arriving during each

of these 15-minute starting times during a year. All of the maximum

elevations had an associated probability equal to the probability

that a certain intensity tsunami would be generated in a particular

segment of the two trench regions and arrive during a particular 15-

minute period of a year.










The many maximum elevations with associated probabilities can

be used to determine cumulative probability distributions of combined

tsunami and astronomical tide elevations. The maximum elevations

are ordered and probabilities summed, starting with the largest

elevations, until a desired probability is obtained. The elevation

encountered when the summed probabilities reach a desired value P

is the elevation that is equaled or exceeded with an average fre-

quency of once every 1/P years. Thus, when the summed probabilities

reach the value .01, the elevation associated with the last prob-

ability summed is the 100-year elevation.

The 100-year and 500-year elevations are determined at all

grid points of the nearshore numerical grids using techniques des-

cribed previously in this chapter. Smooth curves are then used to

connect all discrete elevations so that continuous predictions of

elevation versus coastline location can be made for all the west

coast. Figure 23 shows 100-year and 500-year elevations in the

Crescent City, California, area.

Comparison with Local Observation Predictions

Crescent City and San Francisco, California, are the only

locations on the portion of the west coast of the United States

considered in this study that have sufficient historical data of

tsunami activity to allow frequency of occurrence predictions based upon

local historical observations. Wiegel (1965) made such predictions of

tsunami height (trough to crest height) for Crescent City based upon

the period from 1900 to 1965. He predicted 100-year and 500-year

heights at Crescent City of approximately 25.6 feet and 43.2 feet,










respectively. If the crest amplitude is taken to be one-half the total

height, the 100-year and 500-year elevations are 12.8 feet and 21.6

feet, respectively. However, the crest amplitude at Crescent City is

typically greater than one-half the wave height (e.g. the crest ampli-

tude of the largest wave of the 1964 tsunami was approximately 60 per-

cent of the total height). If Wiegel's analysis is applied to historical

crest elevations instead of heights, the 100-year and 500-year elevations

are found to be 15.4 and 26.4 feet, respectively. Furthermore, the

analysis can now be applied to the longer time period from 1900 through

1977. The 100-year and 500-year elevations based upon this longer time

span are 14.5 feet and 25.5 feet respectively.

Figure 23 shows 100-year and 500-year elevations of 13.1 feet

and 24.9 feet, respectively. These values compare very favorably

with the 14.5 feet and 25.5 feet elevations determined from historical

data for the period of time from 1900 through 1977. The elevations

predicted by the analysis based upon the local historical data are

probably somewhat larger than the elevations predicted in this

report because the short time period (relative to 100 or 500 years)

from 1900 through 1977 includes the exceptionally active years of

tsunami generation from 1946 through 1964 (see Chapter III).

The 13.1-foot and 24.9-foot elevations predicted for Crescent

City using the techniques discussed in this study are not totally

comparable to the elevations predicted using Wiegel's analysis based

upon local historical data since the effect of the astronomical tide

has been included in the elevations predicted in this study. However,






































I' _


I [i K


i--iI- F\ ii~

I I





Cci





i t" ____ ____ ___


L


b
0




LA


0


4J

0






CD
0













0
L)
0CD
z
'4f Ci

0
I- LC





CO
o ci

oo
-(

CCC
'ax
0 I





-UU


cW c o Nc 0 cO Co N 0 co o
N N N N .
1A 'NOIiVA313










the statistical effect of the astronomical tide on the total ele-

vation is not significant at Crescent City due to the large amplitude

of tsunamis there. Appendix A shows that the statistical effect of

the astronomical tide for a location where tsunamis are large is to
2
increase the predicted elevation by an amount equal to 2

where c- is the tidal variance and equals E C C is
m m m
m=l
equal to the mth tidal constituent, and a is given by the following

expression (since it was shown earlier in this chapter that elevations

are linearly related to the logarithm of probabilities):

P(Z) = Ae-aZ

with

Z = the elevation above local mean sea level, and

P(Z) = the cumulative probability distribution for the elevation

at a given site being equal to or exceeding Z due only to the

maximum wave of the tsunami. Based upon the tidal constituents

predicted by the Coast and Geodetic Survey, o = 7.1 for Crescent

City. Using Wiegel's data for the period of time from 1900 through

1977 gives a = .145 Therefore, the astronomical tide contributes

approximately 0.5 feet to the elevations predicted for Crescent

City.

Frequency of occurrence calculations just inside San Francisco

Bay were made in an earlier study (Houston and Garcia, 1975) using a

simple one-dimensional analytical solution for nearshore propagation.

However, portions of San Francisco Bay are included in the nearshore

numerical model (see Figure 13) used in this study in order that the

effect of the Bay on elevations outside the Bay can be properly

PAGE 1

TSUNAMI INTERACTION WITH COASTLINES AND ELEVATION PREDICTIONS By JAMES ROBERT HOUSTON A Dissertation Presented to the Graduate Council of The University of Florida In Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy UNIVERSITY OF FLORIDA 197S

PAGE 2

ACKNOWLEDGMENTS The writer wishes to express his sincere gratitude to his supervisory committee chairman, Professor Knox Millsaps, and cochairman, Professor Joseph Hammack, for their support and efforts to provide him the opportunity to complete his studies at the University of Florida. A special debt of gratitude is owed to Dr. R. IV. Khalin, Chief of the Wave Dynamics Division, Hydraulics Laboratory, U. S. Army Waterways Experiment Station for his encouragement and active support during the research and preparation periods of this dissertation. The writer wishes to thank Professor limit Unluata for the insight he provided and for organizing a supervisory committee prior to his departure from the University. Appreciation also is extended to Dr. D. L. Durham for many helpful discussions at the Waterways Experiment Station. The research upon which this dissertation is based was funded by the Pacific Ocean Division of the U. S. Corps of Engineers and the Federal Insurance Administration of the Department of Housing and Urban Development through the Office of the Chief of Engineers.

PAGE 3

TABLE OF CONTENTS ACKNOWLEDGMENTS KEY TO SYMBOLS ABSTRACT Page VI 11 CHAPTER I: INTRODUCTION 1 Tsunamis \ Objective and Scope of Study 2 CHAPTER II: TSUNAMI INTERACTION KITH THE HAWAIIAN ISLANDS . 4 Introduction 4 Governing Equations . 7 Finite Element Model 12 Numerical Grid 21 Model Verification 25 CHAPTER III: TSUNAMI INTERACTION WITH CONTINENTAL COASTLINES 43 Introduction . . 43 Governing Equations 45 Finite Difference Model 43 Model Verification 51 CHAPTER IV: ELEVATION PREDICTIONS 58 Hawaiian Island Predictions 58 Introduction 58 Interpolation of Recent Historical Data 61 Time Period Analysis 63 Reconstruction of Historical Data 65 Frequency of Occurrence Distribution 69 Results 71

PAGE 4

Page West Coast of the United States Predictions 75 Introduction 75 Tsunami Occurrence Probabilities 77 Use of Deterministic Numerical Models 80 Effect of the Astronomical Tides 85 Comparison with Local Observation Predictions . . 88 CHAPTER V: CONCLUSIONS 95 APPENDIX: SUPERPOSITION OF LARGE TSUNAMIS AND ASTRONOMICAL TIDES 97 REFERENCES 99 BIOGRAPHICAL SKETCH 105

PAGE 5

KEY TO SYMROLS a Wave amplitude b Incident wave amplitude o t b(co) Amplitude of frequency component u C Chezy's coefficient C Amplitude of m tidal coefficient m r D Period of time f Frequency of occurrence F External force vector F( ) Functional g Acceleration of gravity h Water depth H Average runun, meters avg Hankel function of first kind of order m m i (-1J I Tsunami intensity j Integer k Wave number [K] Coefficient matrix L Wavelength m Integer n Unit normal vector n Manning's coefficient

PAGE 6

n( ) Frequency of occurrence N Interpolation function P Pressure -> q Velocity vector {Q} Total load vector r Spherical coordinate t Time T Wave period u Velocity component in x-direction v Velocity component in y-direction w Velocity component in z-direction x Cartesian coordinate Ax Small increment in x y Cartesian coordinate Ay Small increment in y z Cartesian coordinate (vertical direction) a Constant a Constant coefficient m 3 Constant coefficient m A Area of element £ Coefficient of horizontal eddy viscosity e„ Coefficient of vertical eddy viscosity n Free-surface elevation

PAGE 7

6 Spherical coordinate \\ Dynamic viscosity constant £ Response to arbitrary tsunami p Mass density p (u>) Phase angle 2 a Tidal variance t. Bottom stress b $ Velocity potential <£> . A Velocity potential m region A $ Velocity potential in region B B $ Velocity potential of incident wave $ Velocity potential of scattered wave to Angular frequency V Gradient operator

PAGE 8

Abstract of Dissertation Presented to the Graduate Council of The University of Florida in Partial Fulfillment of the Requirements For the Degree of Doctor of Philosophy TSUNAMI INTERACTION WITH COASTLINES AND ELEVATION PREDICTIONS By James Robert Houston June 197S Chairman: K. Millsaps Cochairman: J. L. Hammack Major Department: Engineering Science The development is described of two numerical models that accurately simulate the propagation of tsunamis to nearshore regions and the interaction of tsunamis with coastlines. One of these models is a finite element model which uses a telescoping numerical grid to cover a section of the Pacific Ocean, including all eight islands of the Hawaiian Islands. The second model is a finite difference scheme that uses four rectilinear grids to cover most of the west coast of the continental United States. The finite element model solves linear and dissipationless long-wave equations. Such equations govern nearshore propagation in the Hawaiian Islands since the short continental shelf of the islands limits the time available for nonlinearities and dissipation to cause significant effects. The finite difference model solves nonlinear long-wave equations that include bottom stress terms which are important for the long continental shelf of the west

PAGE 9

coast of the United States. Both models are verified by hindcasting actual historical tsunamis and comparing the numerical model calculations with tide-gage recordings. A frequency of occurrence analysis of tsunami elevations at the shoreline in the Hawaiian Islands is described. This analysis is based upon local historical data with the finite element model used to interpolate between historical data recorded during the period of accurate survey measurements since 1946. Historical data recorded at Hilo, Hawaii, and dating to 1337 is used in conduction with data recorded since 1946 in Hilo and throughout the Hawaiian Islands to reconstruct elevations at locations in the islands lacking data prior to 1946. Frequency of occurrence curves are determined for the entire coastline of the Hawaiian Islands using these reconstructed elevations. Since most of the west coast of the United States lacks local data of tsunami activity, a frequency of occurrence analysis of tsunami elevations is based upon historical data of tsunami occurrence in tsunamigenic regions in addition to numerical model calculations. A generation and deep-ocean propagation numerical model is used to propagate tsunamis with varying intensities from locations throughout the Aleutian-Alaskan and Peru-Chile regions to a water depth of 500 meters off the west coast. The nearshore finite difference model propagates these tsunamis from the 500 meter depth to shore. The frequency of occurrence of combined tsunami and astronomical tide elevations is determined by an analysis involving the numerical superposition of tsunamis and local tides.

PAGE 10

CHAPTER I INTRODUCTION Tsuna mis Of all water waves that occur in nature, ore of the most destructive is the tsunami. The term "tsunami," originating from the Japanese words "tsu" (harbor) and "nami" (wave) , is used to describe sea waves of seismic origin. Tectonic earthquakes, i.e., earthquakes that cause a deformation of the sea bed, appear to be the principal seismic mechanism responsible for the generation of tsunamis. Coastal and submarine landslides and volcanic eruptions also have triggered tsunamis. Tsunamis are principally generated by undersea earthquakes of the dip-slip type with magnitudes greater than 6.5 en the Richter scale and focal depths less than 59 kilometers. They are very long period waves (5 minutes to several hoars) of low height (a few feet or less) when transversing water cf oceanic depth, Consequently, they are not discernible in the deep ocean and go unnoticed by ships. Tsunamis travel at the shallow-water wave celerity equal to the square root of the acceleration due to gravity times water depth even in the deepest ocean; because of their very long wavelengths. This speed of propagation can be in excess of 500 miles per hour in the deep ocean. When tsunamis approach a coastal region where the water depth decreases rapidly, wave refraction, shoaling, and bay or harbor resonance may result in significantly increased wave heights. The great

PAGE 11

period and wavelength of tsunamis preclude their dissipating energy as a breaking surf; instead, they are apt to appear as rapidly rising water levels and only occasionally as bores. Over 500 tsunamis have been reported within recorded history. Virtually all of these tsunamis have occurred in the Pacific Basin. This is because most tsunamis are associated with earthquakes, and most seismic activity beneath the oceans is concentrated in the narrow fault zones adjacent to the great oceanic trench systems that are predominantly confined to the Pacific Ocean. The loss of life and destruction of property due to tsunamis have been immense. The Great Iloei Tokaido-Nanhaido tsunami of Japan killed 30,000 people in 1707. In 186S, the Great Peru tsunami caused 25,000 deaths and carried the frigate U.S.S. Waterlee 1,500 feet inland. The Great Meiji Sanriku tsunami of 1896 killed 27,122 persons in Japan and washed away over 10,000 houses. In recent times, three tsunamis have caused major destruction in areas of the United States. The Great Aleutian tsunami of 1946 killed 173 persons in Hawaii, where runup heights as great as 55 feet were recorded. The 1960 Chilean tsunami killed 330 people in Chile, 61 in Hawaii, and 199 in distant Japan. The most recent major tsunami to affect the United States, the 1964 Alaskan tsunami, killed 10 7 people in Alaska, 4 in Oregon, and 11 in Crescent City, California, and caused over 100 million dollars in damage on the west coast of North America.

PAGE 12

Objectiv e and Scope of Stud}' The objective of the present study is to determine tsunami elevation frequencies of occurrence at the shoreline for the Hawaiian Islands and the west coast of the continental United States, Such elevation predictions are required by the Pacific Ocean Division of the U. S. Corps of Engineers for use in tsunami flood hazard evaluations for floodplain management and the Federal Insurance Administration of the Department of Housing and Urban Development for flood insurance rate calculations. The current trend of general population migration to coastal regions and the growth of the tourist industry makes it essential to know possible tsunami inundation levels even for coastal regions that are currently undeveloped. There arc only a few isolated locations in the United States that have sufficient local historical data of tsunami activity to allow accurate predictions solely based upon historical data. This lack of historical data of tsunami activity on most of the coastline of the United States makes it necessary to use numerical models in addition to existing historical data to predict tsunami elevation frequency of occurrence. However, numerical models have not been successfully applied in the past to the simulation of actual tsunami nearshore propagation and interaction with coastlines. There fore, two numerical models are developed in the present study to accurately simulate tsunami nearshore activity. Calculations of these deterministic numerical models are combined with probabilistic analyses based upon historical data to predict frequencies of occurrence of tsunami elevations for the Hawaiian Islands and the west coast of the United States.

PAGE 13

CHAPTER II TSUNAMI INTERACTION KITH THE HAWAIIAN ISLANDS Int rodu ction Early studies of tsunami interaction with islands involved analytic investigations of the scattering of monochromatic plane waves by a single island of simple geometric shape. Oner and Hall (1949) determined the theoretical diffraction pattern for long-wave scattering by a circular cylinder in water of constant depth by adopting results from acoustic theory. They found qualitative agreement between this diffraction pattern and smoothed runup heights observed for the 1946 tsunami on the island of Kauai in the Hawaiian Islands. Hom-ma (1950) determined the diffraction pattern for long wave scattering off a circular cylinder surrounded by a parabolic bathymetry that extended to a prescribed distance offshore, beyond which the depth was constant. Hydraulic models have been used to study tsunami interaction with single islands that are realistically shaped and surrounded by variable bathymetry. Van Porn (1970) studied tsunami interaction with Viake Island using a 1:57000 undistorted scale model. Jordaan and Adams (1968) studied tsunami interaction with the island of Oahu in the Hawaiian Islands using a 1:20000 undistorted scale model. They found poor agreement between historical measurements of tsunami runup and the hydraulic model data. Scale effects (e.g. viscous effects) and the effects of the arbitrary boundaries that confine the model are 4

PAGE 14

problems that are inherent in hydraulic modeling of tsunamis. The small heights of tsunamis relative to deep ocean depths is also difficult to simulate. Jordaan and Adams [1968) modeled waves to a vertical scale of 1:2000; thus, the waves had heights ten times the normal proportion. Even with this distortion, waves typically had heights in the model of only 0.3 millimeters. Vastano and Reid (1966) studied the problem of the response of a single island to monochromatic plane waves of tsunami period using a finite difference numerical model. A transformation of coordinates allowed a mapping of the island shore as a circle in the image plane. The finite difference solution employed a grid which allowed greater resolution in the vicinity of the island than in the deep ocean. Vastano and Bernard (1973) extended the techniques developed by Vastano and Reid (1966) to multipleisland systems. However, the transformation of coordinates technique allows high resolution only in the vicinity of one island of a multiple-island system. Thus, for the three-island system of Kauai, Oahu, and Niihau in the Hawaiian Islands, the two islands of Oahu and Niihau had to be represented by cylinders with vertical walls whose cross sections were truncated wedges. Kauai was represented by a circular cylinder with the surrounding bathymetry increasing linearly in depth with distance radially from the island until a constant depth was attained. A single Gaussianshaped plane wave was used as input to the model. No comparisons were made with historical tsunami data for the three islands. Hwang and Divoky (1975) presented a simulation of the 1964 tsunami at Hilo, Hawaii, using a finite difference numerical model. The 1964 Alaskan tsunami was generated and propagated across the deep ocean using

PAGE 15

a finite difference numerical model employing a coarse grid. The Hawaiian Islands were too small to be represented by this numerical grid. A waveform was calculated at a point in water 5000 meters deep that was approximately 200 kilometers from where Hilo, Hawaii, would have been if the Hawaiian Islands could have been modeled. This waveform was used as input to a finite difference grid with fine grid cells that covered an area approximately 22 kilometers by b7 kilometers in the immediate vicinity of Hilo. The water depth along the input boundary of the fine grid varied from approximately 250 meters to 1500 meters. Hwang's approach neglected the influence of the Hawaiian Islands and of the varying bathymetry and depth changes between the deep ocean point and the fine grid input boundary on the tsunami arriving in Hilo. The numerical model calculations for the 1964 tsunami were in poor agreement with the Hilo tide gage recording of this tsunami. A finite difference model employing a grid covering the eight major islands of the Hawaiian Island chain was used by Bernard and Vastano (1977) to study 'the interaction of a plane Gaussian pulse with the Hawaiian Islands. The square grid cells were 5.5 kilometers on a side and close to the minimum feasible (because of computer time and cost limitations) size for a rectangular cell finite difference grid covering the major islands of Hawaii. Yet historical data indicate that significant variations of tsunami elevations occur over distances much less than 5.5 kilometers. The islands of Hawaii are relatively small and very poorly represented by a 5.5 kilometer grid. For example, Oahu has a diameter of only approximately 50 kilometers and the landwater boundary of the island has characteristic direction changes

PAGE 16

that occur over distances much less than 5.5 kilometers. The offshore bathymetery of the islands also varies rapidly, with depth changes of more than 1500 meters frequently occurring in distances of 5.5 kilometers. Furthermore, if a resolution of eight grid cells per tsunami wavelength is maintained for tsunami periods as low as 15 minutes, a 5.5 kilometer grid cannot be used to propagate waves into depths less than approximately 150 meters. The processes that cause significant wave modifications and subsequent rapid variations of elevations along the coastline (that are known from historical observations to occur during tsunami activity in the Hawaiian Islands) probably occur within this region extending from water at a depth of 150 meters to the shoreline. In the present study, a finite element numerical model is used to propagate tsunamis from the deep ocean to the shoreline of the Hawaiian Islands. Since finite element techniques allow dramatic changes in element sizes and shapes, this model has great advantages relative to finite difference numerical models in accurately representing land shapes, ocean bathymetry, and tsunami waveforms. Element size can be large in deep water where bathymetric variations are gradual and wave lengths are long. As a wave enters shallow water its length decreases and the elements of the grid can be telescoped to smaller sizes with no accompanying less of resolution. Governing Equations Typical tsunamis have lengths much greater than the water depths over which they propagate; for such waves, fluid motions are

PAGE 17

approximately two-dimensional, i.e., vertical fluid velocities and accelerations are small in comparison with horizontal velocities and accelerations, respectively. Consider a small element of fluid as shown in Figure 1. The rate of increase of volume of the element is given by — {(n+h) dxdy} = jdxdy where, n = a free-surface elevation, h = the water depth, t = time, and (x,y) are Cartesian coordinates. (2.1) r,4h Figure 1. Element of Fluid In Equation (2.1), it is assumed that h=h(x,y), i.e. h has no time dependence. The net rate of volume flux is given in difference form by {u(n+h)> A , dy {u(n+h)} dy + {v(n+h)> , dx {v(n+h)> dx y+dy )' = dxdy {-7— u(n+h) + 3y oX = dxdy V • q(n+h) v(n+h)l + 0(dx 2 dy ) (2.2)

PAGE 18

where u,v are horizontal velocity components, w is the vertical velocity component (w=0 here), and q = (u,v,w) is a velocity vector (a Taylor expansion is used in liquation (2.2) and terms of order 2 2 dx dy have been dropped). Equating the rate of increase of volume to the rate of volume flux yields the following continuity equation: f£ + V • q(mh) =0 (2.3) at The governing momentum equation for incompressible and constant viscosity fluid flow is the Navier-Stokes equation in the form p ||+ q • Vq = pF Vp + yV 2 q (2.4) where p = the fluid pressure, p = the mass density, u the dynamic viscosity constant, F = an external force vector. Since vertical fluid velocities and accelerations are small and assuming that gravity is the only external force and the fluid is inviscid so that the last term of Equation (2.4) can be dropped, Equation (2.4) becomes in component form 3u 3u 3u 1 3p ,. .. — + u — + v — = -^(2.5) 3t 3x 3y p 9x 3v 3v 3v 1 3p ,-,,-, -— + u t+ v — = -^ (2.6) 3t 3x 3y P 3y

PAGE 19

10 When integrated, Equation (2.7) yields (with atmospheric pressure equal to zero) p = -pg(S-Ti) (2.8) i.e., the pressure is hydrostatic, and Vp = pgVn (2-9) Equations (2.5) and (2.6) can be recombined to yield ||+ (q ' V)q = -gVn "(2.10) Since along most of the path of propagation of a tsunami the wave amplitude is much less than the. water depth, Equation (2.3) can be approximated as |£+ V • (hq) = (2.11) It can also be shown that the second term in Equation (2.10) is generally less than the other terms by a factor of the order of n/h . Therefore, for n/h snail, Equation (2.10) is approximated by 8q 9t Taking a time derivative of Equation (2.12) and substituting Equation g 7n (2.12) (2.11) into the right-hand side of the subsequent equation yields i-| = gV(V • (hq)) (2.13) dt Equation (2.13) is known as the linear long wave equation. It is a governing equation for small amplitude long waves in an inviscid and incompressible variable depth fluid and is valid when the assumptions made in its derivation remain valid. Assuming that the fluid is irrotational q can be expressed as a gradient of a scalar function (x,y,t) , known as the

PAGE 20

11 velocity potential, i.e, q = VC (2.14) Assuming all motions to be simple harmonic in time yields Hx,y,t) = Re [OC-x.yJe -1 ^] (2.15) where Re [ ] indicates the real part of the braced quantity, and w is the frequency of the motion. Substituting Equations (2.1.4) and (2.15) into (2.13) yields 2 7 • (hVO) + ~ i (2.16) Equation (2.16) is the well known generalized Helmholtz equation. The wave amplitude r, can be related to the velocity potential $ by evaluating the Bernoulli equation at the free surface. The Bernoulli equation is given by 3 ; P 1 1 1 2 _ = _ g 7 + vol = , (2.17) at p 2 ' ' v ' however; the last term in Equation (2.17) is generally smaller than the others by a factor of the order of n/h and can be neglected. Also, taking p=0 at the free surface where I-n , Bernoulli's equation becomes to the same order of approximation as Equation (2.16) Several assumptions are made in the derivation of Equations (2.16) and (2.18) and the adequacy of these assumptions for the particular problem of tsunami interaction with the Hawaiian Islands must be considered. An important factor in this consideration is the fact that the continental shelf of the Hawaiian Islands is exceptionally short. Thus, during most of the time of tsunami propagation, fluid velocities are much less than wave propagation velocities since wave amplitudes

PAGE 21

12 arc much loss than water depths and bottom stress forces as well as nonlinearities have little time to cause significant effects. Ursell (1953) and Hammack (1972) have shown that nonlinear effects grow with time. However, for a typical large tsunami, such as the 1964 Alaskan tsunami} Hammack and Segur (1977) have shown that neither nonlinearity nor frequency dispersion have any significant effect on the lead wave as it propagates across the ocean. In addition, when the lead wave reaches a continental shelf a certain time must pass before nonlinearity or frequency dispersion become significant. For a typical large tsunami, Hammack and Segur (1977) show that linear and nondispersive theory holds for the lead wave over a distance of approximately 2 ±L ^— (a is a wave amplitude). Taking the average depth of the region from the edge of the continental shelf (600 foot depth) to shore as 300 feet and the average amplitude as 3 feet (the amplitude becomes large ?h 2 only near shorej yields = 12 miles. This distance is much greater than the length of the continental shelf of the Hawaiian Islands (usually less than 1 mile). Since the lead wave is followed by waves which have been reflected by land areas in the generation region, linear and nondispersive theory should be adequate, at least for the major waves of a tsunami. Finite Element Model Figure 2 shows the problem of tsunami interaction with islands with the process divided into two domains. One domain (A) includes the islands and the surrounding variable depth region. The second domain (B) is a constant depth region extending to infinity. This constant

PAGE 22

15 Figure 2. Regions of Computation

PAGE 23

14 depth region is realistic since the Pacific Ocean has a remarkably constant depth beyond the immediate region of the Hawaiian Islands. The boundary value problem expressed in terms of these domains becomes the following: 2 V ., (hr:) + — 5 = in Region A (2.19) 2 V • (Vo) + !L co where j = number of land masses, r = a spherical coordinate, k = the wave number, § = the scattered wave velocity potential, i = (-D • Equations (2.21) and (2.22) express the continuity of the velocity potential and its derivative along the boundary separating the two domains and Equation (2.23) expresses the impermeability of solid boundaries. Equation (2.24) is the Sommerfeld radiation condition that requires the scattered wave to be an outgoing wave at infinity. A calculus of variations approach can be used to obtain an EulerLagrange formulation of boundary value problems. This variational formulation is based upon the principle that of all possible displace-

PAGE 24

13 merit configurations a body can assume that satisfy compatibility and the constraints of kinematic boundary conditions, the configuration satisfying equilibrium makes the potential energy a minimum. The EulerLagrange formulation of boundary value problems has often been used in classical dynamics (Hamilton's principle of least action) and solid mechanics (structural, solid, and rock mechanics). Of course, it is not possible to obtain analytical solutions for many engineering problems. The finite element method is a numerical method first introduced by Turner et al. (1956) (although some of the underlying ideas were discussed by Courant (1943)) which has found wide application in solid mechanics for problems based upon an EulerLagrange formulation involving complex material properties or boundary conditions (see Zienkiewicz, 1971 and Desai and Abel, 1972). The finite element method is a discrete approximation procedure applicable whenever a variational principle can be formed. This method has a very short history in fluid mechanic applications despite advantages over conventional numerical schemes such as finite difference methods. The first application of the finite element method to fluid mechanic problems involved a study by Zienkiewicz and Cheung (1966) of seepage through porous media. Chen and Mei (1974) studied the problem of the forced oscillation of a small basin protected by a breakwater and containing floating nuclear power plants by using the finite element method. They discretized a small constant depth region using approximately 250 triangular elements that were all about the same size and shape. Since small amplitude long waves were used as a forcing function, Equation (2.20) was the

PAGE 25

16 governing equation (the water depth being constant in their application) . In the present study an approach similar to that of Chen and Mei is used to solve the problem of tsunami interaction with islands in a variable depth region. For the finite element method to be useful for tsunami problems, it is necessary to take advantage of two important properties of 'this algorithm. First, the finite element grid which discretizes a domain can have grid cells (elements) of arbitrary size and shape. Thus, the grid can telescope from a large cell size in one section of a variable depth domain to a very fine cell size in another. Secondly, if the form of the spatial variation of a physical parameter is known a priori, the finite element method can employ a very small number of elements in the discretization used to estimate the exact variation. For example, a single finite element can be used to solve the problem of the seiching of a lake for a constant depth rectangular lake if the basic form (unknown coefficients) of the solution is known. The finite element method degenerates into the well-known Rayleigh-Ritz variational method in this case. This property of the finite element method will be important in modeling the infinite region B shown in Figure 2. The variational principle for the boundary value problem given by Equations (2.19) through (2.24) requires a certain functional F(-}>) to be stationary with respect to arbitrary first variation of . The first variation of the functional F(4') for this problem is wellknown and can be expressed as

PAGE 26

17 2 3<54> ;f co) = //[v • (iwo) £*] 6? + ; h(i B -; A J H A ° "-A A The functional can be derived from this expression for its first variation (Chen and Mei, 1974) and is ? 3($ -$ T ) Ft.) // i[M«) z 1* 2 i * f * «VV thHA fe -A A 3(VV 3i 1 SA A 3n A 3A A 3n A 3($ R -$ T ) + / ho gg T , (2.26) 3A l Sn A where $ is the velocity potential of the incident waves (forcing function) . The finite element method (unlike the. finite difference method) discretizes the domain itself rather than the governing equations. The region under consideration is partitioned into small regions known as finite elements and the assemblage of all such elements represents the original region. Instead of solving the problem for the entire body in one operation, the solutions are formulated for each constituent element and combined to obtain the solution for the original region. Thus, for example, the first term of Equation [2.26) is evaluated on the element level and the results for all elements are summed. This sum is equivalent to the original integral over the complete region A. In order to evaluate the first term in Equation (2.26) on the element level, it is necessary to approximate the field variable

PAGE 27

within each element by a polynomial N , called a trial or interpolation function, in the coordinates x and y . For each element $ = {N} {e(2.27) where {£} I

PAGE 28

19 expansion independent of all other elements, the treatment of the entire domain is systematically handled by summing the contributions from each element. The expression $ R -* T i n Equation (2. 26 j is the scattered wave velocity potential $„ . Since region B has a constant depth, £> can be solved analytically. The governing equation for this case is the Helmholtz equation (Equation (2.20)). The boundary condition at infinity requires the scattered wave obey the Sommerfeld radiation condition given by Equation (2.24). The well-known solution to this boundary value problem is as follows: $ I H (kr) (a cos in? + s sin m9) (2.34) s „ m m m m=0 where H are Hankel functions of the first kind of order m, m a and 6 are constant coefficients, m in and 9 is a spherical coordinate. In the finite element formulation, the infinite region B can be considered to be a single element with the interpolation function given by Equation (2.54). Thus, within region A the interpolation functions are linear, but within region B the interpolation function is an infinite series involving unknown coefficients. The terms in Equation (2.54) are oscillatory but have a monotonically decreasing modulation; hence, the series can be truncated after a finite number of terms depending upon the desired accuracy (see next section) . The required stationarity of the functional F(e) with respect to arbitrary variations of $ at each nodal point yields a set of linear algebraic equations when expressed on the element level that can be solved using matrix methods. For each specified incident wave period,

PAGE 29

the velocity potential can be determined at all nodal points. The linearized dynamic free surface boundary condition (Equation (2.18)) can then be used to calculate the surface elevation n . The set of linear algebraic equations evolving from the element level can be written in the simple matrix forn [K] W = {Q} (2.35) where [K] is the large coefficient matrix, {0} is a column vector representing a combination of nodal unknowns for the velocity potential and the coefficient unknowns from the single element covering region B, and {Q} is the total load column vector related to the forcing function. in the finite element solution of structural mechanics problems, one obtains a similar equilibrium relation with [K] being a stiffness matrix, {^} a nodal displacement vector, and {Q} a nodal force vector (Desai and Abel, 1972). The coefficient matrix [K] is symmetric, banded, and complex (complex notation was used to represent waveforms, e.g., the Hankel functions are complex) . For the case of tsunami interaction with the Hawaiian Islands, [K] will be exceptionally large and sparse (many zero terms) . Of course, only half of the banded part of the matrix alone needs to be manipulated in the solution process. The banded part of the matrix [K] is solved using standard Gaussian elimination methods. The matrix solution used by Chen and Mei (1974) was modified in this study to take advantage of the sparseness of the matrix in order to reduce computational time. Also, the elimination method was modified so that calculations involved only small blocks of terms at a time with the remainder of the matrix kept in pe7 % ipheral storage. This matrix partitioning was necessary for the finite element grid shown in the next

PAGE 30

21 section since the banded part of the matrix had over S00,000 terms. Partitioning is possible because the symmetric coefficient matrix is positive definite; hence, a solution is possible by elimination methods without pivoting. Without pivoting, elimination performed using one row, affects only the triangle of elements within the band below that row. Thus, only two triangular submatrices needed to be stored in computer memory at a given time. Numerical G rid The finite element grid covering the variable depth region A is shown in Figure 3. Elements of the grid telescope from large sizes in the deep ocean to triangles with areas as small as 0.5 square kilometers in shallow coastal waters. The grid covers an area of approximately 450,000 square kilometers. In general, finite element techniques allow elements to be any arbitrary shape (e.g. quadrilateral) provided that the ratio of the lengths of the shortest and longest sides of an element is not extreme. Triangular shapes were used for convenience in the present study. The geometric shapes of the eight islands comprising the Hawaiian Island chain and the rapid bathymetric variations surrounding the. islands are modeled very precisely by the finite element grid of Figure 3. Also, the number of node points along the shorelines of the islands is very dense. Since wave heights are calculated at node points by the finite element model, the model can adequately represent the rapid wave-height variations along coastlines that are known to occur during tsunami activity in the Hawaiian Islands.

PAGE 31

Figure 3. Finite I; lenient Grid.

PAGE 32

The finite element grid was drawn by hand with the element sizes dictated by the local water depth. The lengths of element sides were allowed to be no greater than one-eighth of the local wavelength for the shortest period waves considered. Such a cell size results in a wave resolution error of no more than one percent f.Mei 197S) . There are computer programs (Cole and Riech, 1976) that generate finite element grids and number the nodes to minimize bandwidths; however, such programs only can be applied to fairly homogeneous domains. The rapid depth variations in region A require a rapidly telescoping grid that is too complicated for automatic grid generation techniques. Each node of the finite element grid must be numbered and this numbering is extremely important since it determines the bandwidth of the large coefficient matrix discussed in the last section. The bandwidth of this matrix is controlled by the largest difference (among all elements) in the nodal numbers between two nodes of the same element. The computational time of the numerical model is proportional to the cube of the bandwidth of the coefficient matrix. For a simple region covered by homogeneous elements, the bandwidth is minimized by numbering nodal points in cyclic sequence following a spiraling pattern '..here one spiral is shielded by the following spiral from a third spiral. The grid pictured in Figure 3, however, is so large and complicated that considerable judgement and effort is involved in determining a numbering pattern that minimizes the bandwidth, thus resulting in reasonable computational times for the numerical model. The numerical model also requires the spatial coordinates of all node points and information concerning how the elements are connected

PAGE 33

relative to each other (connectiveness) . The coordinates were determined using an electronic digitizer and a simple computer program developed to relate the internal coordinate system of the digitizer to a prototype coordinate system through coordinate rotation and scaling. The elements are numbered in an arbitrary fashion and then connectiveness information assembled for the numerical model by listing the three nodes of each element. The nodal numbering, digitizing, and connectiveness assembly necessarily introduce errors. For example, the connectiveness assembly involves determining over 17,000 numbers and placing them on computer cards. Many connectiveness mistakes do not result in error conditions which terminate the numerical model's computations. Instead, incorrect answers are produced by the model. To eliminate all errors in the input data, computer plotting programs were developed to draw the grid and number nodes and elements based upon the connectiveness and node coordinate information used as input to the numerical model. Such plots make errors readily apparent. The number of terms retained in the truncation of Equation (2.34) was determined by trial and error. The final number of terms was such that addition of further terms produced negligible effects. It was found that this occurs for terms approximately 10 " times smaller than the initial terms. The computational time requirements of the finite element model for the grid shown in Figure 3 are very modest, making it economically feasible to determine the interaction of arbitrary tsunamis with the Hawaiian Islands (discussed in next section) . The reason that the

PAGE 34

computationtal time required by the finite element model is reasonably small is that the grid uses small elements only in areas where they are necessary. The grid of Figure 3 has approximately 2500 points, whereas the finite difference grid of Bernard and Vastano (1977) had 26,000 grid points. Even so, some of the elements of Figure 3 are as much as sixty times smaller than the finite difference grid cells used by Bernard and Vastano. Model Verification The finite element model is verified in the present study by comparing numerical simulations of the 1960 Chilean and 1964 Alaskan tsunamis with tide gage recordings of these tsunamis in the Hawaiian Islands. These two tsunamis are the only major tsunamis for which reliable information exists concerning characteristics of the ground uplift which generated the waves (much more information exists for the Alaskan source) . A deep ocean recording of a tsunami near the Hawaiian Islands has never been made. Therefore, prototype wave records in the deep ocean cf the 1960 Chilean and 1964 Alaskan tsunamis are not available for use as input to verify the finite element model. However, a finite difference numerical model which solves the linear long wave equations employing a spherical coordinate grid has been developed (Hwang et al., 1972), and used in several studies (Houston and Garcia, 1974; Houston and Garcia, 1975; and Houston et al . , 1975b) to generate tsunamis and propagate them across the deep ocean. This deep ocean numerical model

PAGE 35

26 is used in this study to determine deep ocean waveforms of the 1960 and 1964 tsunamis that were needed to verify the finite element model. The deep-ocean numerical model simulates the uplift deformation of the ocean water surface caused by the permanent vertical displacement of the ocean bottom during an earthquake and the subsequent propagation of the resulting tsunami across the deep ocean. The permanent deformation (permanent in the sense that the time scale associated with it is much longer than the period of the tsunami) and not the transient movements within the time-history of the ground motion is considered to be the important parameter governing farfield wave characteristics. Transient movements occur for a period of the order of tens of seconds, whereas, tsunami wave periods are of the order of tens of minutes. Experimental investigations of tsunami generation by Ilammack (1972) show that the transient ground movements do not influence farfield characteristics of resulting tsunamis for spatially large ground displacements occurring over a short period of time. Hwang et al , (1972) verified the deep ocean model by simulating the 1964 Alaskan tsunami. The permanent ground deformation in the source region was taken from Plafker (1969) and a 1/4° by 1/6° grid was used for the simulation. Since this grid was much too coarse to allow accurate modeling of shallow-water propagation, a comparison was made between the numerical model calculations and a gage recording of the 1964 tsunami at Wake Island where Van Dorn (1970) had a pressure transducer in relatively deep water (800 foot depth) some distance offshore. The amplitude and length of the first wave of this tsunami were shown to be in good agreement. (Actually, the grid used by Hwang et al.

PAGE 36

27 (1972), did not extend as far as Wake Island. The deep ocean model propagated the tsunami to a point in the Pacific Ocean and then simple refraction techniques were used to estimate the change in the waveform during propagation from this point to Wake Island.) The simulation of the 1960 and 1964 tsunamis performed in this study uses a 1/3° by 1/3° grid. Prior simulations of the 1964 tsunami using both a 1/4° by 1/6° grid and a 1/5 D by 1/3° grid demonstrated negligible differences between use of the two grid spacings (Houston et al . , 1975a). The permanent ground deformation used as input is also taken from Plafker (1969). The numerical grid was large enough to cover a substantial portion of the Pacific Ocean including Wake Island. The numerical model calculations at Wake Island were found to be similar to those presented by Hwang et al. (1972). Figure 4 shows a comparison of the numerical model calculations of this study and the Wake Island recording of the 1964 tsunami. The lack of agreement in the trailing wave region may be due to local effects (e.g. wave trapping by the Wake Island seamount) or reflections from the nearby transmission boundaries of the numerical grid (these boundaries are not perfectly transmitting) . Figures 5 and 6 show time-histories of the 1960 Chilean and 1964 Alaskan tsunamis in deep water near the Hawaiian Islands calculated using the deep ocean propagation model. The permanent ground deformation in the source region of the Chilean earthquake of 1960 is taken from Plafker and Savage (1970). Calculations for the Alaskan tsunami involve a grid with approximately 30,000 cells and the calculations for the Chilean tsunami a grid with approximately 90,000 cells. The finite element numerical model is a time-harmonic solution of

PAGE 37

28 J,

PAGE 38

1 I I I I C^T ~ — ^ ^— « „. I C 5 o z O 2 c •.— -* c r. — JLJ 3d^lOXOad 'SNOIXVA313

PAGE 39

— _ o U O li BdAlOlOHd sno;xva3~i^

PAGE 40

the boundary value problem. The response of a group of islands to an arbitrary tsunami can easily be determined within the framework of a linearized theory using the theory of superposition. For example, an arbitrary tsunami in the deep ocean can be Fourier decomposed as follows: b (tO = f b(w) e" l(a)t+pfu) ^ du (2.36) where b = incident wave amplitude o r b(w) = amplitude of frequency component tu p (u) ~ phase angle If n(x,y,(j) is the response amplitude at any point (x,y) along the island coasts due to an incident plane wave of unit amplitude and frequency cj , then the response of the islands to the arbitrary tsunami time history b ft) is given by o f(x,y,t) = Re [ /b(w)n(x,y,oO e" 1 ^^ ^ *du] (2.37) Therefore, when n(x,y,co) is known for all u , the island response to an arbitrary tsunami can be calculated. Of course, it is not feasible to calculate the integrals of Equations (2.56) and (2.37) over all frequencies. Instead, the frequency range must be discretized and the integrals replaced by sums over a frequency range containing most of the energy of the tsunami. Equation (2.36) involves a Fourier decomposition of a time series. This decomposition is accomplished for the time-histories of the 1960 Chilean and 1964 Alaskan tsunamis (Figures 5 and 6) using a least squares harmonic fitting procedure. The time-history of the Alaskan

PAGE 41

tsunami in deep water was decomposed into IS components with periods ranging from 14.5 minutes to the time length of the record (260 minutes) . The variance of the residual (difference between the actual record and a recomposition of the 18 components) was approximately 0.2 percent of the variance of the record. Therefore, virtually all the energy of the wave record was contained in the 18 components. The time-history of the Chilean tsunami in deep water was decomposed into 11 components with periods ranging from 15.5 minutes to the time length of the record (170 minutes). The variance of the residual was less than 0.1 percent of the variance of the record. For both cases, the original time-history and a time-history constructed from a recomposition of the components were virtual ly indistinguishable. Equations (2.56) and (2.57) take the following form when discretized : b o (t) = I b(u) n ) e i{ V +p (V } ( 2 5g ) C(x,y,t) = Re[ E b(u ) n(x,y,ai ) e" i{u n t+pCw n 5} ] (2.59) n=l n where m = number of components, th to = frequency or n component The M" ) term is determined by the leastsquares harmonic fitting procedure and n(x,y, u ) is determined by the finite element numerical model for each frequency to and at each location (x,y). Therefore, the time-history 5(x,y,t), which represents the response of location (x,y) to the deepwater tsunami time-history given by b (t) ,

PAGE 42

,s5 can be calculated for any location along the coastline of the Hawaiian Islands. Figures 7, S, and 9 present comparisons between tide gage recordings and numerical model calculations of the 1964 Alaskan tsunami at Kahului, Maui; Honolulu, Oahu; and Hilo, Hawaii, respectively. The largest waves recorded at each of these sites are shown. Waves arriving at later times are all much smaller than those shown. Whenever tide gage limits were encountered, the recordings were linearly extended (the tide gage locations are shown in Figure 10) . The wave records shown in these figures are in remarkable agreement, especially considering the fact that the ground displacement of the 1964 earthquake was not precisely known. Since the Hilo breakwater was not included in the numerical model, the numerical model calculations for Hilo are probably too large; however, this breakwater was undoubtedly highly permeable to the 1964 tsunami. The numerical results appear to be too large in Hilo by some constant factor, since the tide gage recording and the numerical model calculations have the same form with the first wave crest and trough having approximately the same amplitude and being proportionately smaller than the second crest. Figure 11 shows a comparison between the tide gage recording of the 1960 Chilean tsunami at Honolulu and the numerical model calculations. Again the largest waves recorded are shown (the Honolulu gage was the only tide gage in the Hawaiian Islands not destroyed by the 1960 tsunami) . The permanent ground motion of the 1960 earthquake is not known nearly as well as the ground motion for the 1964 earthquake.

PAGE 43

LEGEND NUMERICAL MODEL CALCULATIONS ESTIMATED (BEYOND GAGE LIMIT) APPROXIMATE HOURS GREENWICH MEAN TIME Fig ur c 7. 1964 Tsunami from Alaska Recorded at Kahului, Maui

PAGE 44

:45 9:00 10:00 APPROXIMATE HOURS GREENWICH MEAN TIME : igure 8. 1964 Tsunami from Alaska Recorded at Honolulu, Oahu.

PAGE 45

TIDE GAGE RECORD NUMERICAL MODEL CALCULATIONS 9:00 10-00 APPROXIMATE HOURS GREENWICH MEAN TIME Figure 9. 1961 Tsunami from Alaska Recorded at fiiio, Hawaii.

PAGE 46

Figure 10. Gage Location:

PAGE 48

39 Hence, the agreement indicated in Figure 11 is quite good considering the uncertainty of the uplift in the source region. Tide gages, of course, do not record tsunamis perfectly since they are nonlinear devices with responses which depend upon both the period and amplitude of the disturbances they measure. However, simple calculations based upon the paper of Noye (1970) show that the distortion is small for the tsunamis shown in Figures 7 through 9 and Figure 11. Figure 12 shows the response amplitude as a function of wave period at the locations of tide gages in Kahului, Honolulu, and Hilo. The response versus wave period is greatly dependent upon the precise location of recording for bays such as Hilo. For example, the Hilo response (Figure 12) at the tide gage location shown in Figure 10 has a single large resonant peak at 2S.9 minutes. This peak occurs at the fundamental period of oscillation of Hilo Bay and, consequently, all locations in the bay have a resonant peak at this incident wave period (the amplitude of the peak depends on location in the bay). Other locations in Hilo Bay also have significant peaks at lower periods. Loomis (1970) calculated a free oscillation fundamental period (node assumed at bay mouth) of 26.2 minutes for Hilo Bay using a finite difference grid spacing of approximately one-half a kilometer. Bernard and Vastano (1977) show a contour map of a normalized energy parameter versus wave period and location around the island of Hawaii. The map shows Hilo displaying an energy peak from 10 to 16 minutes centered at approximately 12.5 minutes; no significant peak at 28.9 minutes is observed. The difference between Bernard and Vastano 's calculations and the results presented in this paper is attributed to

PAGE 49

KAHULUI TIDE GAGE 15 CE

PAGE 50

41 the poor resolution of the spatial features of Hilo Pay by their finite difference grid spacing of 5.5 kilometers. Hi Id Bay has the shape of a right triangle with sides that are approximately 9 kilometers and 13.5 kilometers long. The grid used by Bernard and Vastano represented Hilo Bay as a right triangle with sides that were 5.5 kilometers and 11 kilometers long; thus, the surface area of Hilo Bay was approximately one half its actual area. This distortion of the sir.e of Hilo Bay and also, perhaps, the distortion of the bathymetry of the bay by the coarse grid accounts for the lower response periods for Hi lo Bay presented by Bernard and Vastano. The importance of accurately representing island shapes and offshore bathymetry and resolving tsunami wavelengths can be illustrated by considering a study by Adams (1975) that used the numerical model developed by Bernard and Vastano (19""") to interpolate between historical measurements of tsunami elevations produced by the 1946 Aleutian tsunami. A single Gaussian wave crest was used as input for the numerical model. If a resolution of four grid cells per crest wavelength is maintained, this numerical model cannot propagate waves into water with a depth less than about 500 feet. The main processes involved in the transformation of the form of a tsunami probably occur in the region over the continental shelf (depths less than 600 feet). Consequently, the elevations calculated by Adams were small and did not display the great variation along coastlines known to occur during historical tsunamis in the Hawaiian Islands. Adams included the effects of propagation over continental shelves by multiplying the numerical model elevations by a factor (different for each island) that forced the best agreement with the historical

PAGE 51

42 data. However, even with artiface, agreement with historical data was poor. For example, on the northeast coast of Oahu, the numerical model elevations multiplied by the factor resulted in predicted elevations ranging from 14 to 21.4 feet. Loomis (1976) presents 36 historical elevations for the 1946 Aleutian tsunami on the northeast coast of Oahu and only 7 even fall within this range of elevations. The historical elevations actually varied from 1 to 37 feet on this coast.

PAGE 52

CHAPTER III TSUNAMI INTERACTION WITH CONTINENTAL COASTLINES In troduc tion Experimental investigations of long wave interactions with continental coastlines have generally been one-dimensional and have involved idealized conditions. Flume tests by Savage [1358) and Saville (1956) measured runup of periodic waves on constant slope beaches and Mail and Watts (1953), Kaplan (1935), and Kishi and Saeki (1966) measured runup of solitary waves :>n constant slope beaches. All of these flume experiments used waves much steeper than typical tsunamis due to the difficulties involved in the generation of small amplitude long waves. Many theoretical studies of one-dimensional propagation of weakly dispersive long waves have been performed in recent, years. Tsunamis become weakly dispersive (amplitude and frequency dispersive) long waves at some point during their propagation, although (as was the case for the Hawaiian Islands) they may not be weakly dispersive for a time period sufficiently long for related effects to develop. The proper equations for the investigation of weakly dispersive long wave propagation over a steeply sloping bottom were first derived by Mei and LeMehaute (1966) and in a slightly different form by Peregrine (1967) who also gave the first numerical solution in the case of a solitary wave using an implicit finite difference 43

PAGE 53

44 scheme. A more complete solution was later given by Madsen and Mei (1969) who transformed the equations to characteristic form before they were solved numerically. Svendsen (1974) studied weakly dispersive wave propagation over a gently sloping bottom. Hammack (1972) made comparisons of long wave generation by bottom displacements in the laboratory and numerical solutions of the KortewegDeVries equation using the numerical approach of Peregrine (1966). Several studies of one-dimensional tsunami propagation over continental shelves of constant slope have been performed. Keller and Keller (1964) solved nonlinear long wave equations numerically for a constant slope bathymetry. Takahasi (1964) studied tsunami propagation over a simple bathymetry through analytic solutions of the linear long wave equation. Carrier (1966) assumed a tsunami was generated by a point source and studied nonlinear propagation over a constant slope bathymetry. Heitner and Housner (1970) used a finite element numerical method to solve nonlinear equations for propagation up constant' slope beaches. More recently, Mader (1974) studied tsunami propagation over linear slopes using a finite difference solution of the complete Navier-Stokes equations (not vertically integrated). Occasional two-dimensional studies have been performed for tsunami propagation in ncarshore regions. Grace (1969) studied tsunamis affecting the reef runway of the Honolulu International Airport using a hydraulic model. However, the wavemakers used were only a small fraction of a tsunami wavelength from the region to be studied and waves reflected by land were reflected again almost

PAGE 54

45 immediately by the wavemakers. This is a common problem in hydraulic model tests involving long waves. The wavemakers can be moved a few wavelengths from the region of interest by reducing the model scale; however, scale effects may then become very significant. Wave absorber screens may also be placed in front of the wavemakers; however, it is very difficult to absorb very long waves and such an approach is often not very helpful. Similar problems with rereflected waves can occur in numerical studies. The calculations of Hwang and Divoky (1975) were in poor agreement with the tide gage recording of the 1964 tsunami at Hilo not only for the reasons presented in Chapter II, but also because their algorithm allowed waves to leave the numerical grid efficiently only for waves approaching the input boundary normally. Since Hilo Bay is triangularly shaped, the reflected waves approached the input boundary at various angles and were partially re-reflected toward Hilo. Governing Equations The governing equations for long wave propagation over the continental shelf of the west coast of the United States are Equations (2.5) and (2.4). Thus, the long wave assumption is made once again with vertical velocities and accelerations neglected relative to horizontal components. The continuity equation is given by the equation ~ + V • qCn + h) = (3.1) and the momentum equation is 8q -> -v 2->yf + q • Vq = pF Vp + pV q (3.2)

PAGE 55

46 If the only external force is gravity and the pressure is hydrostatic, then liquation (3.2) becomes p S+ q • Vc i = ps Vr i + ^ iV " c i ( 3 -^) c)t The west coast of the United States has a much broader continental shelf than that of the Hawaiian Islands; thus, the nonlinear advective and bottom stress effects may be more significant, and the nonlinear and viscous terms of Equation (3.3) must be retained. The dynamic viscosity u in the viscous term of Equation (3.3) is replaced by an eddy viscosity constant e for the case of turbulent flow. Using this approach the viscous term can be expressed as *, c^4 ^4) + b 2 M C3.4) 1 3x" 3y ^z where the coefficients e. and e_ represent coefficients of horizontal and vertical eddy viscosity, respectively. The term 2 2 3x 3y 8 2 is usuallv small in comparison with the e — § term (Dronkers, 1964) and will be neglected in this study. If e, — | is 9z integrated with respect to z from the bottom to the surface and divided by the distance from bottom to surface n + h , we find [(ItK-.k (It) . J V-V n+h L 9z z=-h 3z z = n The terms £:_(-—) _ . and E 2 ( '^z=n &re cora P onents of the

PAGE 56

47 tangential stresses at the sea bottom and surface, respectively. The shear stress at the sea surface is quite important for storm surge problems where wind velocities are very large but can be neglected relative to the bottom stress for tsunamis. From experiments and theoretical considerations for onedimensional flow in a river, it has been shown (Dronkers, 1964), that the bottom stress t, can be expressed as follows T =l£SHM (3.6) b c 2 where C is Chezy's coefficient and is given by the expression 1,486 h l/6 (37) n m where n is Manning's n . Thus, Equation (3.3) may be written m -*-, 2 2,1/2 3q •*" ^ n p gq(u +v ) f ~ R1 p —x + q • Vq = pgVr. ~^ — ' ' U6 J CT(n+h) Equation (3.8) contains the nonlinear advective term but neglects frequency dispersion since vertical accelerations have been neglected, This equation adequately describes wave propagation for Ursell numbers much greater than unity (Ursell, 1953). The Ursell number is proportional to the ratio of nonlinear effects to linear effects nL 2 (i.e., frequency dispersion) and is equal to — =— , where L is a 1/2 h typical wavelength. Since L=(gh) T for long waves, where T is T 2 the wave period, an Ursell number may be chosen as gn(^) The numerical model described in the next section uses Equation (3.8) to propagate waves from depths of 500 meters to shore. These

PAGE 57

48 waves have periods of the order of 1800 seconds and initial wave T 2 amplitudes of the order of a foot. Thus, the Ursell number, gnGr) is aproximately 50 and increases as the wave enters shallow water with h decreasing and n increasing. Mader (1974) found that the one-dimensional form of Equation (3.8), without the viscous term, produced very similar results to a one-dimensional N'avier-Stokes equation (viscosity neglected) for propagation of long wavelength tsunamis over linear continental slopes. The pressure term in Mader' s solution was not hydrostatic; thus, frequency dispersion was permitted. Finite Difference Model To solve Equations (5.1) and (5.8) numerically, the differential equations are ciscretized and replaced by a system of finite difference equations using central differences on a space-staggered grid. The space-staggered scheme describes velocities, water levels, and depths at different grid points. This scheme, first used by Platzman (195S) for lake surge problems, has the advantage that in the equation for the variable operated upon in time, there is a centrally located spatial derivative for the linear term. Leendertse's (1967) implicit-explicit multi-operational method is employed in determining the solution for q and n as functions of time. The stability and convergence of this method is described in detail by Leendert.se (1967) . This method for solving finite-difference equations has been applied in studies analyzing the tidal hydraulics of harbors and inlets (Raney,

PAGE 58

49 1976). Houston and Garcia (1975) used the method to calculate the decay of a tsunami as it enters and spreads throughout San Francisco Bay. The nearshore finite-difference numerical model uses the timehistory calculated by the generation and deep ocean propagation model described in Chapter II as input in the following manner. A tsunami is generated in the Aleutian-Alaskan area or the west coast of South America and propagated across the deep ocean to a 500meter depth off the west coast of the United States. Waveforms calculated at this depth by the deep ocean numerical model are recorded all along the west coast. These waveforms are then used as input to the nearshore numerical model, which propagates the tsunamis from the 500-meter depths across the continental slope and shelf to shore. Figure 13 shows outlines of the four numerical grids used to cover the west coast. The grids have square grid cells two miles on a side (southern California and some of San Francisco Bay were not covered by numerical grids because elevation predictions were made in previous studies (Houston and Garcia, 1974; and Houston and Garcia, 1975) using simple one-dimensional analytical solutions for propagation over the continental shelf). The offshore bathymetry was modeled from the 500-meter contour to shore. Beyond the 500meter contour, the ocean was assumed to have a constant 500-meter depth. The input boundary of each grid was located approximately one and one-half wavelengths of a 30-minutc wave from the shore. Therefore, at least three typical waves could arrive at the shore before waves reflected from the input boundary became a problem.

PAGE 59

'QGruc. WASHINGTON ' x "baja CALIFORNIA Figure 15. Finite Difference Grid Locations,

PAGE 60

51 The input boundaries of the grids were oriented approximately parallel to the shoreline since refraction will bend the wavefronts to such an orientation before they reach the 500-meter contour (Wilson and Torum, 1968). Lateral boundaries of the grids were taken to be impermeable vertical walls. Model Verification The finite difference nearshcre numerical model was also verified by numerical simulations of the 1964 Alaskan tsunami. The deep ocean propagation model discussed in Chapter II was used to generate the 1964 tsunami and propagate it across the deep ocean. Figure 14 shows surface elevation contours 3-1/2 hours after the 1964 Alaskan earthquake, as calculated by the deep-ocean propagation numerical model. This figure illustrates the concentration of energy on the northern California, Oregon, and Washington coasts due to the directional radiation of energy from the source region. The upper section of Figure 15 shows a time history of the 1964 Alaskan tsunami calculated by the deep-ocean model near Crescent City, California, in a water depth of 500 r.eters. This waveform was used as input to the finite difference nearshore model. The lower section of Figure 15 shows a comparison between a tide gage recording of the 1964 tsunami at Crescent City and the nearshore model simulation of this tsunami. The numerical model calculations agree well with the tide gage recording; periods, phases, and amplitudes are accurately reproduced. The main disagreement occurs

PAGE 61

52

PAGE 62

40 60 80 TIME, PROTOTYPE MIN LEGEND PROTOTYPE GAGE RECORD INFERRED PROTOTYPE GAGE RECORD (FROM WILSON, 196B) NUMERICAL SIMULATION APPROX TIME, HRS, GMT Figure 15. 1964 Tsunami from Alaska Recorded at Crescent City, California

PAGE 63

54 for the amplitude of the third wave. With regard to this disagreement, it should be noted that there is also a lack of agreement on the amplitudes of the actual historical third and fourth wave crests. Figure 15 shows a reconsti action of the 1964 tsunami at Crescent City inferred by Wilson and Torum (1963) from the prototype gage record and later survey measurements by Magoon (1965). The elevation of the fourth wave (not shown) was estimated by Wilson and Torum to have been a little less than 14 feet above mean lower low water (MLLW) (the zero elevation shown in Figure 15 is MLLW datum); however, Wiegel (1965) estimated that the elevation of the third wave was approximately 16 feet above MLLW and that the fourth wave attained the highest elevation of 18 or 19 feet above MLLW at the tide gage location. The numerical mode] predicted a fourth wave elevation approximately the same as that estimated by Wilson and Torum (14 feet above MLLW) . However, reflections from the input boundary were probably growing in importance in the numerical model calculations during the arrival of the fourth wave. There may be several reasons why the maximum wave elevation predicted by the numerical model does not agree with the historical record as well as the elevations of the initial waves. First, the later waves in the deepwater waveform are probably waves which have been reflected from land areas in the source region. The fairly large grid cell spacing used in the deep ocean propagation model to represent the land-water boundary in this region probably distorts the reflected waveform. Total reflection of waves in this region also is not completely realistic. Furthermore, the nearshore grid may

PAGE 64

55 not have fine enough cells to completely model the resonance phenomenon leading to the large third or fourth wave. Since 1964, there has been considerable speculation concerning the reasons that the effects of the 1964 tsunami were so great at Crescent City. .The finite difference numerical models show that both the directional radiation of energy from the source region (see Figure 14) and a local resonance caused the relatively large elevations at Crescent City. Actually, large elevations at Crescent City were not unique since directional radiation of energy from the source also caused large elevations along the Oregon and northern California coasts. Runup 10 to 15 feet above the high tide level occurred all along the Oregon coast south of the Columbia River (runup at Crescent City was approximately 15 feet above the high tide level). However, the Oregon coast is very sparsely populated and there were few damage reports. The severity of structural damage at Crescent City, which has a large logging industry, also was apparently due to the impact of logs carried by the tsunami (Wilson and Torum, 1968) . The nearshore numerical model indicates that the resonant effects at Crescent City were fairly local, extending over 2 to 4 miles of coastline. This behavior contradicts Wilson and Torum' s (1968) speculation that a bowl-shaped section of the continental shelf with a diameter of approximately 50 miles experienced a resonant oscillation. Historical data support the numerical model calculations. At the mouth of the Klamath River, approximately 15 miles south of Crescent City, elevations observed during the 1964 tsunami were only 2 to 3 feet above normal high

PAGE 65

56 tide. Wilson and Torum's hypothesis of a shelf oscillation would predict elevations at the mouth of the Klamath River greater than those that occurred at Crescent City. Figure 16 shows a comparison between a tide gage recording of the 1964 tsunami at Avila Beach, California, and the nearshore numerical model calculations. The elevations recorded by the Avila Beach tide gage were larger than those recorded at any tide gage on the west coast except the Crescent City gage. The historical and calculated waveforms shown in Figure 16 are in good general agreement. The tide gage record obviously has higher frequency components which are not predicted by the numerical model; these components may be local oscillations of water areas which were too small to be accurately represented by the numerical grid. Important features such as the wave amplitudes are in good agreement. The waveforms recorded for the 1964 tsunami by the Crescent City and Avila Beach tide gages have elevations larger than those recorded by any other gage on the west coast. Other tide gage recordings of this tsunami on the west coast (except at Astoria, Oregon) were in areas not covered by the grids shown in Figure 13. A comparison with the tide gage at Astoria, Oregon, was not made because Astoria is approximately 12 miles away from the coast in the estuary of the Columbia River, and amplitudes of the 1964 tsunami were small at Astoria.

PAGE 66

GAGE \ .' PROTOTYPE GAGE RECORD \S NUMERICAL SIMULATION 8.5 9 9.5 10 10.5 APPROX. TIME, HRS., GMT. Pi "ure J 6. 'sunami from Alaska Recorded at Avaiia Beach. .alii ornia .

PAGE 67

CHAPTER IV ELEVATION PREDICTIONS Hawaiian Island Predictions Introduction The Hawaiian Islands, a chain of eight islands as shown in Figure 17, have a history of destructive tsunamis generated both in distant areas and locally. The earliest recording of a severe tsunami in the Hawaiian Islands was in 1S37 when a tsunami from Chile reached an elevation of 20 feet at Hilo, Hawaii, and killed 46 people in the Kau District of the big island of Hawaii. Prior to 1S37, a number of severe tsunamis undoubtedly reached the islands but, unfortunately no detailed records were kept. Since 1S37, there have been 16 tsunamis that have caused significant damage. Most of the destructive tsunamis in the Hawaiian Islands have been generated along the coast of South America, the Aleutian Islands, and the Kamchatkan Peninsula of the Soviet Union. Approximately onefourth of all the tsunamis recorded in the Hawaiian Islands have originated along the coast of South America, while more than onehalf have originated in the Kuril-Kamchatka-Aleutian region of the north and northwestern Pacific. Tsunamis generated by local seismic events have caused large runup in the islands, especially on the southeast coast of the big island of Hawaii. The 1868 tsunami produced the largest waves of record in the Hawaiian Islands with 58

PAGE 68

59

PAGE 69

60 60-foot waves reported on the South Puna coast of the island of Hawaii. The most recent major tsunami in Hawaiian history occurred on November 29, 1975, when waves generated by an earthquake with an epicenter on the South Puna Coast may have reached elevations as great as 45 feet along the southeast coast. The most destructivetsunami to ever hit the islands in terms of both loss of life and property destruction was the Great Aleutian tsunami of 1946, which killed 175 people and produced waves over 55 feet in elevation. Hilo incurred $26 million in property damage attributable to this tsunami. The 1960 Chilean tsunami is the most recent distantly generated tsunami that produced major effects in the Hawaiian Islands. Sixtyone lives, all at Hilo, were claimed by the tsunami. Damage throughout the state was estimated to be $25.5 million. Inspection of the damage at Hilo revealed much evidence of the tremendous forces developed by the waves. Twenty-ton boulders had beer, moved hundreds of feet, asphaltic concrete pavements were peeled from their subbase, and hundreds of automobiles were moved and crushed. The following sections discuss the methods used to obtain sufficient data to predict runup elevations for the entire coastline of the Hawaiian Islands. The finite element numerical model discussed in Chapter II is used to interpolate between local historical data for the years of accurate survey measurements from 1946 through 1977. Such a method of interpolation is needed so that elevations can be predicted for the entire coastline even for locations having no actual historical data. A method is discussed herein which allows reconstruction of historical data for the period of time from 1857 to

PAGE 70

61 1946 when elevations for the period from 1946 through 1977 are known. Such an approach is necessary in order to have sufficient data for accurate elevation predictions since it is shown herein that reasonable predictions cannot be made solely based upon known historical data for the short time span from 1946 through 1977. Interp ol ation of Recent Historical Pat a Historical data of tsunami activity in the Hawaiian Islands are, of course, often limited to certain locations. Information on tsunami activity in the islands prior to the 1946 tsunami is concentrated in Hilo (Hawaii) and to a lesser extent in Honolulu (Oahu) and Kahului and Lahaina (Maui). Even data for tsunamis occurring during the period of accurate survey measurements from 194 6 through the present are absent or fragmentary (i.e., data exist only for certain of the events) for most of the coastline of the islands. Therefore, it is necessary to rely on more than just available historical data to determine tsunami occurrence frequencies. The numerical model is used to fill in historical data gaps for tsunamis from 1946 through 1977 by providing relative responses of the Hawaiian Islands to tsunamis. For example, although the deepwater waveform of a tsunami such as that of 1946 is not known, the direction of approach of this tsunami and its range of wave periods are known. The average response of the islands to a tsunami similar to the 1946 tsunami is determined by inputting sinusoidal waves of unit amplitude from the direction of the 1946 tsunami into the numerical model over a band of wave periods. If historical data exist at one location, wave elevations at a nearby location for which historical data do not exist are determined by multiplying the historical data at the

PAGE 71

b2 first location by the ratio of the response calculated by the numerical model at the second location and the response calculated at the first. The numerical model takes into account the major processes that would cause different relative wave elevations at the two locations. That is, the model calculates shoaling, refraction, diffraction, reflection, resonance, shielding of the backside of an island by the front side, and interactions between islands. The historical data account for the unknown absolute heights of the deepwater tsunamis. The "correction factors" provided by the historical data are local instead of global (single factor for all islands) because there are certain local factors that cannot be determined by the numerical model. For example, coral reefs extend over parts of some of the islands and are known, to protect against tsunami attack (e.g., the extensive reefs in Kaneohe Bay, Oahu, which are powerful tsunami dissipators) . The effects of the reefs on tsunamis are implicitly contained in the historical data itself. The historical data from 1946 through the present were taken from the most recent and complete compilation of these data (Loomis, 1976) . As an example of the use of the numerical model to fill in historical data gaps for tsunamis from 194 6 through the present, consider the 1946 tsunami recorded on the island of Lanai only at Kaumalapua Harbor and Manele Bay. Wave elevations can be calculated for the 1946 tsunami at other locations on the coastline of Lanai by multiplying these historical values by the ratios of the responses calculated by the numerical model at locations of interest and the responses calculated at Kaumalapua Harbor and Manele Bay. In this

PAGE 72

way, elevations for the 1946 tsunami can be determined for all of Lanai . Time Period Analysis It is necessary to use historical data of tsunami occurrence covering the greatest possible time span to properly determine frequency of occurrence of tsunami elevations for the Hawaiian Islands. Tsunami activity has not been uniform in these islands during recorded history. For example, the two largest and four of the ten largest tsunamis striking Hilo from 1857 through 1977 occurred during the 15-year period from 1946 through 1960. Two of the tsunamis from 1946 through 1960 originated in the Aleutian Islands, one in Kamchatka and one in Chile. However, six of the ten largest tsunamis occurred during the 109-year period from 1857 to 1945 with three originating in Chile, two in Kamchatka and one in Hawaii. Therefore, both the frequency of occurrence and place of origin of tsunamis have been remarkably variable. Any study (e.g., Towill Corporation, 1975) basing frequency calculations on a short time span that includes the period from 1946 through 1960 will predict a significantly more frequent occurrence of large tsunamis than is warranted by historical data from 1857 through 1977. The errors introduced in frequency of occurrence calculations by consideration only of a short time period that includes the unrepresentative years from 1946 through 1960 will be greater than the errors resulting from possible observational inaccuracies of the 19th Century. This is easily demonstrated for Hilo, Hawaii, which is the location in the Hawaiian Islands having the most complete

PAGE 73

61 data of tsunami activity. Using the data compiled by Cox (1964) for the ten largest tsunamis in Hilo from 1857 through 197" (his data were actually through 1964, but no large tsunamis have occurred in Hilo since 1964), a least-squares fit of the data was made employing the logarithmic frequency distribution discussed later in this chapter. The 100-year tsunami elevation for Hilo calculated using this method was 27.3 feet. If the elevations recorded during the nineteenth centry (five of the ten largest tsunamis), are increased by 50 percent (which is much larger than possible observational error), the 100-year elevation is found to be 50.4 feet. Similar calculations just based upon data taken during the period of accurate survey measurements in Hilo from 194 6 through 1977 yields a 100-year elevation of 44.2 feet. Since the largest elevation in Cox's data for the 141-year period from 1S57 through 1977 was just 28 feet, the 100-year elevation of 44.2 feet is obviously too large. Similar overestimates will occur at any location in the Hawaiian Islands if an analysis is based upon a short period of time that includes the years 1946 through 1964 since the exceptionally frequent occurrence of major tsunamis in Hilo from 1946 through 1964 is a property of the unusual activity of tsunami generation areas and not of special properties of Hilo. The unusual tsunamigenic activity in generation regions in recent years is reflected in the most up-to-date catalog of tsunami occurrence in the circumpacific area complied by Soloviev and Go (1969). This catalog shows, for example, that the three greatest intensity tsunamis generated in the Aleutian-Alaskan region since

PAGE 74

65 1788 occurred in 1946, 1957, and 1964. The catalog also lists the 1960 Chilean tsunami as the greatest intensity tsunami generated in South America in recorded history [416-year period from 1562 through 1977) . Reporting of tsunamis on the west coast of South America is quite good with, 105 tsunamis reported during the 416-year period. Kanamori (1977] also notes the unusual seismic activity that occurred during a relatively short time span in recent years. He shows that the four earthquakes having the greatest seismic moments (true measure of magnitude for great earthquakes) in this century occurred during a thirteen year period from 1952 through 1964 (1960 Chilean, 1964 Alaskan, 1957 Aleutian, and 1952 Kamchatkan earthquakes). Evidence that tsunamis were unusually active in the Hawaiian Islands from 1946 through 1964 is also apparent from historical data at Kahului and Lahaina, Maui. Kahului, Lahaina, iiilo, and Honolulu are the only locations in the islands with historical data extending to 1837. Four of the five largest tsunami elevations recorded at Kahului occurred during the years from 1946 through 1964. Three of the four largest tsunami elevations recorded at Lahaina occurred during the same period. The historical data for Honolulu is not very useful because the commonly reported elevations for the 1837 and 1841 tsunamis were, in fact, drops in water level (PararasCarayannis, 1977) and there are no reports for the 1868 and 1877 tsunamis from Chile. Reconstruction of Historic al Data In order to have data covering a sufficiently long time span to make reasonable predictions of tsunami elevations for the Hawaiian Islands, it is necessary to reconstruct historical data for the

PAGE 75

66 years prior to 1946. To do this, it is assumed in this study that tsunamis generated in a single source region (Kamchatka or Chile but not the Aleutian-Alaskan region) approach the islands from approximately the same direction and have energy lying within the same band of wave periods. The difference in wave elevations at the shoreline in the Hawaiian Islands produced by tsunamis generated at different times in the same region is attributed mainly to differences in deepwater wave amplitudes. For example, the 1841 tsunami from Kamchatka produced a wave elevation in Hilo that was approximately 25 percent greater than that of the 1952 tsunami from Kamchatka. It is assumed that these two tsunamis had the same relative magnitudes throughout the islands. Since the finite element model allows a determination of elevations for the 1952 tsunami along the entire coastline of the islands, the elevations of the 1S41 tsunami and all ether Kamchatkan tsunamis occurring during the period of record in Hilo from 1857 to 1946 can be reconstructed for the entire coastline. Similarly, Chilean tsunami elevations can be reconstructed knowing elevations for the 1960 Chilean tsunami. No major tsunamis were generated in the Aleutian-Alaskan region from 1857 through 1946; hence, no reconstruction is necessary for tsunamis from this region. Fortunately, tsunamis from this region since 1946 were representative of the entire region since these tsunamis were generated in the western Aleutians (1957), central Aleutians (1946) and eastern Alaskan area (1964) . The assumption that tsunamis generated in a single source region (Kamchatka or Chile) approach, the Hawaiian Islands from

PAGE 76

67 nearly the same direction is justified in the case of Kamchatka by the small spatial extent of the known generation area. However, the tsunamigenic area of Chile is considerably larger than that of Kamchatka. Still, the direction oP approach of different tsunamis from Chile is nearly the same because the tsunamigenic region in Chile subtends a relatively small angle with respect to the islands as a result of the great distance between Chile and the Hawaiian Islands and the axis of the Chile Trench maintains a constant orientation relative to the islands (thus directional radiation effects are not important) . The assumption that tsunamis generated in a single source region have energy lying within the same band of wave periods is supported by the elevation patterns produced along coastlines in the Hawaiian Islands by different historical tsunamis. These patterns are such that it has been concluded that "tsunamis of diverse geographic origin are strikingly different, whereas those from nearly the same origin are remarkably similar" (Eaton, et al., 1961). A recent study by IVybro (1977) shows that even the distributions of normalized elevations (i.e. the elevation patterns) produced in the Hawaiian Islands by different Aleutian-Alaskan tsunamis are nearly the same yet quite different from the distributions for tsunamis of otherorigins. This agreement occurs dispite the fact that the relatively close proximity to Hawaii of the Aleutian-Alaskan Trench and the varying orientation of the trench axis relative to Hawaii introduces important directional effects for tsunamis generated in the AleutianAlaskan area. Apparently, these directional effects influence the

PAGE 77

magnitude of the elevations in Hawaii but do not greatly alter elevation patterns. Therefore, historical observations support the approach used in this study which estimates the elevations produced by tsunamis from Chile or Kamchatka prior to 1946 based upon data for tsunamis from these tsunamigenic regions recorded during the years of accurate survey measurements from 1946 to 1977. In this study elevations are calculated along the coastline of the Hawaiian Islands for the ten tsunamis which produced the greatest elevations in Hilo, Hawaii, from 1S57 through 1977. The finite element numerical model fills in all data gaps for those tsunamis of the ten that occurred from 1946 through 1977 and then the Hilo data (Cox, 1964) is used to determine relative magnitudes to reconstruct the elevations for those tsunamis that occurred prior to 1946. For much of the coast of the Hawaiian Islands, the ten largest tsunamis since 1837 would be the same as the ten largest tsunamis in Hilo for the same period. That is, the ten largest wouJd be the tsunamis of 1960, 1946, 1925, lS57,'lS77, 1841, 1957, 1952, and the two in 1868. Of course, the order of the ten largest would vary from location to location. For example, the largest tsunami on the South Puna coast of the big island of Hawaii was the locally generated tsunami of 1868. There are a few locations where historical data show that one of the ten largest tsunamis was not among the ten largest in Hilo from 1837 through 1977. For example, the 1896 tsunami generated near Japan produced small elevations in the Hawaiian Islands except at Keauhou on the Kona coast of the bis island of Hawaii. The

PAGE 78

69 observation at Keauhou was, therefore, included in the frequency of occurrence analysis for that location. Similarly, all other historical observations of significant tsunami elevations at isolated locations were included in the frequency of occurrence analyses for those locations,. Historical data for tsunamis prior to 1946 were taken from Pararas-Carayannis (1977) . Data for the local tsunami of .November 29, 1975, were taken from the compilation by Loomis (1976). Since the local tsunami of 1868 had an earthquake epicenter very near that of the 1975 tsunami and the two apparently produced very similar elevation patterns (Tilling et al . , 1976), historical elevations for the two tsunamis recorded at the same location were used to determine a relative magnitude of the tsunamis. The many observations for the 1975 tsunami were then used to reconstruct elevations for the 186S tsunami. Frequency of Occurre nce Distribution After the ten largest tsunamis from 1837 through 1977 were determined at locations all along the coasts of the Hawaiian Islands using the methods discussed in the previous sections, elevation versus frequency of occurrence curves were determined at each location by least-squares fitting of the data by the following expression: h = -B A log f (4.1) where h = elevation of maximum wave at the shoreline, f = frequency per year of occurrence. Equation (4.1) was used as the frequency of occurrence distribution

PAGE 79

70 since it has been found to agree with historical data at several locations. Cox (1964) found that the logarithm of tsunami frequency of occurrence was linearly related to tsunami elevations for the ten largest tsunamis occurring from 1837 through 1964 at Hilo, Hawaii. Soloviev (1969] has shown a similar relationship between tsunami frequency of generation and intensity for moderate to large tsunamis. Also, earthquake frequency of occurrence and magnitude have been similarly related by Gutenberg and Richter (1965) ; IViegel (1965) found the same relationship for historical tsunamis at San Francisco and Crescent City, California and Adams (1970) for tsunamis at Kahuku Point, Oahu. A recent study by Rascon and Villarreal (1975) revealed the same relationship for historical tsunamis on the west coast of America, excluding Mexico. It is possible that other distributions may agree with the historical data equally well as the logarithmic distribution. For example, the Gumbel distribution has been used in the past to study annual streamflow extremes (Gumbel, 1955). Borgman and Resio (1977) illustrate the use of this distribution to determine frequency curves for non-annual events in wave climatology. To investigate the sensitivity of calculations of 100-year elevations on the assumed frequency distribution, the approach of Borgman and Resio was applied to the Ililo data of Cox. The Gumbel distribution yields a 100-year elevation of 28.8 feet for data from 1837 through 1977 and an elevation of 42.5 feet for data from 1946 through 1977. This compares with the elevations of 27.3 feet and 44.2 feet calculated for the same time periods using the logarithmic distribution (Equation (4.1)). Clearly, the arguments used earlier concerning the period of time that must be considered for

PAGE 80

"1 a valid analysis are not dependent upon the assumed frequency distribution. Since ample precedent exists for using the logarithmic distribution for analysing tsunami frequency of occurrence, this distribution is used in this study. Results The methods described in the previous sections allow frequency of occurrence curves to be constructed for points all along the coastline of the Hawaiian Islands. The coefficients A and B of Equation (4.1) contain all necessary information to determine frequency of occurrence of tsunami elevations at a location. Plots of the coefficients A and B versus location along the coasts were constructed for all of the islands. Smooth curves were drawn through these coefficients to allow elevation predictions at any location. Figures 18 and 19 show typical plots of A and B for the southeastern coast of the big island of Hawaii (Figure 20). For example, at Hoopuloa (location 11 in Figure 20) A = 6 and B = 4; hence, the 109-year elevation (F = .01) is only 8 feet. The 100-year elevations for most locations in the Hawaiian Islands are considerably larger. For such situations involving human life, one is less interested in the mean exceedar.ee frequency f (the average frequency per year of tsunamis of equal or greater elevations) than in the chance of a given elevation being exceeded in a certain period of time. To calculate risk, one first notes that tsunamis are usually caused by earthquakes, and earthquakes are often idealized as a generalized Poisson process (Ncwmark and Rosenblueth, 1971 and Der Kiureghian and Ang, 1977). It was assumed by Wiegel (1965) and Rascon and

PAGE 82

73

PAGE 83

74 Figure 20. Location Map for Hawaii.

PAGE 84

75 Villarreal (1975) that tsunamis also follow such a stochastic process. The probability that a tsunami with an average frequency of occurrence of f is exceeded in D years, assuming that tsunamis follow a Poisson process, is given by the following equation: P = 1 e" fD (4.2) Thus, if an acceptable risk is one chance in one thousand that an elevation be exceeded during a 100-year period, P = .001 and D = 100 years. Substituting these values in Equation (4.2) yields the elevation h which has one chance in one thousand of being exceeded during the next 100 years. West Coast of the United States Predictions Introductio n Unlike the Hawaiian Islands, the west coast of the continental United States lacks sufficient data to allow tsunami elevation predictions based upon local historical records of tsunami activity. Virtually all of the west coast is completely without data of tsunami occurrence, even for the prominent 1964 tsunami. Only a handful of locations have historical data for tsunamis other than the 1964 tsunami However, the Federal Insurance Administration requires information on tsunami elevations for the entire west coast of the continental United States; even for the many locations that have no known historical data of tsunami activity and for coastal areas that are currently not developed (since these areas may be developed in the future) . The lack of historical data of tsunami activity on the west coast

PAGE 85

76 necessitates the use of numerical models to predict runup elevations. The Aleutian-Alaskan area and the west coast of South .America were found by Houston and Garcia (1974) to be the tsunamigenic regions of concern to the west coast of the United States. Both regions have sufficient data on the generation of major tsunamis to allow a statistical investigation of tsunami generation. The generation and deep-ocean propagation model described in Chapter III can then be used to generate representative tsunamis and propagate them across the deep ocean. The nearshore numerical model also described in Chapter III, can then be used to propagate tsunamis from the deep ocean over the continental slope and shelf to shore. In this study, only tsunamis of distant origin are considered in the analysis. Hammack (1972) has shown that near the generation region of an impulsively generated tsunami, the waveform is dependent upon details of the time-dependent movement of the ground during the earthquake. Little is known about the actual timedependent ground motion during earthquakes generating tsunamis and this motion cannot be predicted in advance. Also, there is not enough historical data concerning locally generated tsunamis on the west coast to allow predictions of tsunami occurrence. Thus, reasonable predictions on the west coast of the properties of locally generated tsunamis or their likelihood of occurrence are not possible at this time. The probability is not considered very great that a destructive, locally generated tsunami will occur on the west coast of the continental United States. Tsunamis are generally produced by earthquakes having fault movements that exhibit a pronounced "dip-slip," or

PAGE 86

77 vertical component of motion. "Strike-slip/' or horizontal displacement, fault movements are inefficient generators of tsunamis. Faults on the west coast of the United States characteristically exhibit "strike-slip" motion since the Pacific block of the earth's crust is moving, horizontally relative to the North American block. The west coast of the United States does not share the characteristics (ocean trenches and island arcs) of known-generating areas and, in fact, has not historically been one. Relatively small locally generated tsunamis have been known to occur on the west coast, but there are no reliable reports of major locally generated tsunamis. Heights of 6 feet in the immediate vicinity of the 1927 Point Arguello earthquake are the largest authenticated heights produced by local tsunamis on the west coast) . There could be a few locations on the west coast for which locally generated tsunamis pose a greater hazard than do distantly generated tsunamis because the elevations produced by distantly generated tsunamis are small. However, predictions of elevations produced by locally generated tsunamis are beyond the scope of this study. Tsu nami Occurrence Probabilities Historical data of tsunami generation must be the basis for an analysis that considers the probability of tsunami generation in the two tsunamigenic regions of the Pacific Ocean of concern to the west coast of the continental United States--the Aleutian and Peru-Chile Trench regions. A satisfactory correlation between earthquake magnitude and tsunami intensity has never been demonstrated. Not all large earthquakes occurring in the ocean even generate noticeable

PAGE 87

78 tsunamis. Furthermore, earthquake parameters of importance to tsunami generation, such as focal depth, rise time, and vertical ground displacement, have only been measured for earthquakes occurring in recent years. Therefore, data of earthquake occurrence cannot be used to determine occurrence probabilities of tsunamis. Instead, historical data of tsunami occurrence in generation regions must be used to determine these probabilities. In South America, a wealth of information exists concerning tsunami generation. Reliable data (grouped in intensity increments of one-half) exist for tsunamis with intensity greater than or equal to for a 171 -year period and greater than or equal to 2-1/2 for a 416-year period (Soloviev and Go, 1969). The intensity scale used is a modification (by Soloviev and Go) of the standard Imamura-Iida tsunami intensity. Intensity is defined as I = log (/2 H ) (4.3) This definition in terns of an average runup (in meters) over a coast instead of a maximum runup elevation at a single location (used for the standard Imamura-Iida scale) tends to eliminate any spurious intensity magnitudes caused by often observed anomalous responses (due, for example, to local resonances) of single isolated locations. Using the most recent and complete catalog of tsunami occurrence, in the Pacific Ocean (Soloviev and Go, 1969) a relationship between tsunami intensity and frequency of occurrence was determined for the tsunamigenic trench running the length of the Peru-Chile coast. Tsunamis with intensity greater than or equal to were considered.

PAGE 88

79 It was assumed that the logarithm of the tsunami frequency of occurrence was linearly related to the tsunami intensity. Earlier in this chapter, it was shown that such a relationship holds between both earthquake magnitude and frequency of occurrence and tsunami intensity and frequency of occurrence. Letting n(I) equal the probability of a tsunami with an intensity 1 being generated during any given year and using statistics for the entire trench along the Peru-Chile coast, a least-squares analysis results in the following expression: n(I) = 0.74e~°6iI (4.4) In using statistics for the entire trench area along the Peru-Chile coast, it is assumed that the probability of tsunami occurrence is uniform along the trench. This is a standard assumption for earthquake frequency analysis (Gutenberg and Richter, 1965) . The tectonic justification of this assumption lies in the fact that a single Sialic block or plate of the earth's crust or lithosphere is dipping into the Peru-Chile Trench (Wilson, 1959) . It can reasonably be expected that the movement of this single plate is similar along its entire length. In the Aleutian Trench region, only large tsunamis occurring in relatively recent years (since 17SS) have been recorded due to the isolation of the area. Assuming an exponential coefficient of -0.71 for this trench area (determined by Soloviev, 1969) as a mean value for regions of the Pacific with the most data on tsunamis) and using only the reliable data for large tsunamis (intensity greater than or equal to 3.5) from Soloviev and Go (1969), the following relationship

PAGE 89

80 is determined by a least-squares analysis: nl n it-0.711 n ( I j = . 1 1 3 e Again, the probability of tsunami occurrence is assumed to be uniform along the trench. Use of Deterministic Numerical M odels To relate the probability distributions of tsunami intensities to source characteristics, it is assumed that the ratio of the source uplift heights producing two tsunamis of different intensity (as defined in the previous section) is equal to the ratio of the average runup heights produced on the coasts near these tsunami sources. This |l r I 2 l ratio is equal to 2 for two tsunamis with intensities I and l 2 If H is the wave height in the direction parallel to the maior a axis of length a cf a tsunami source with an elliptical shape (large tsunamis have historically had elliptically-shaped uplifts) and H. is the wave height in the direction parallel to the minor b axis of length b, then experimental research of tsunami generation has shown that H, /H is approximately equal to a/b (Hatori, b a 1963). For a large tsunami. H, can be larger than H by a factor J a k o g of as much as 5 or 6. Thus, the orientation of the tsunami source relative to the area where runup is to be determined is very important, i.e., the runup at a distant site due to the generation of a tsunami at one location along a trench cannot be considered as being representative of all possible placements of the tsunami source in the entire trench region. Hence, the Aleutian and Peru-Chile Trenches had to be segmented and runup along the west coast of the United

PAGE 90

81 States determined for tsunami sources located at the center of each of the segments. The spatial size of a tsunami source was standardized since there is not an apparent correlation between tsunami intensity and spatial size of a tsunami source. for example, the 1946 Aleutian tsunami had an uplift region of very small spatial extent, whereas the 1957 Aleutian tsunami had an uplift region that covered perhaps the greatest spatial extent of any known earthquake (Kclleher, 1972) . Yet the 1946 tsunami had the greater intensity, producing, in general, greater runup elevations in the near and distant regions. The standard source used is discussed in detail by Houston and Garcia (1974). It is elliptical in plan view with a major axis length of 600 miles and a minor axis length of 130 miles. The vertical displacement increases linearly from a zero elevation on the ellipse perimeter to a maximum displacement of 30 feet. The vertical displacement has a parabolic crest (concave downward) parallel to the direction of the major axis. The standard source represents a large tsunami with intensity 4 on the modified Imamura-Iida scale. Certainly, tsunamis of low intensity may have smaller spatial extents; however, large tsunamis pose the greatest threat to a distant area such as the west coast of the United States. These large tsunamis can be expected to have similar spatial extents, with any spatial differences being unimportant in the farfield compared with the effects of source orientation and vertical uplift. Vertical uplifts are assigned to different intensities in accordance with the convention discussed at the beginning of this section.

PAGE 91

82 Figures 21 and 22 show the Aleutian Trench divided into 12 segments and the Peru-Chile Trench into 3 segments. The segments in the Aleutian Trench were approximately one-quarter the length of the major axis of the standard source, whereas the segments in the Peru-Chile Trench were approximately the length of the major axis of the standard source. The standard source was centered in each segment such that the major axis of the source was parallel to the trench axis. Uplift regions historically have had such an orientation relative to trench systems. The Aleutian Trench is segmented much finer than the PeruChile Trench since the Aleutian Trench is oriented relative to the west coast such that elevations produced on the west coast are very sensitive to the exact location of a source along the Trench. Uplifts along the Peru-Chile Trench do not radiate energy directly toward the west coast regardless of their position along the Trench. The Peru and Chile sections of the Peru-Chile Trench have constant orientations relative to the west coast of the United States; therefore, elevations on the west coast -are not very sensitive to source location within these sections. In each of the segments of the Aleutian and Peru-Chile Trenches tsunamis with intensities from 2 to 5 in steps of one-half intensity are generated and propagated across the deep ocean using the deepocean numerical model discussed in Chapter II. Tsunamis with intensity less than 2 are too small to produce significant runup on the west coast. An upper limit of 5 was chosen because the largest tsunami intensity ever reported was less than 5 (Soloviev and Go, 1969). Gutenberg and Richter (1965) indicate that there is an upper

PAGE 92

rJ kj o o k

PAGE 93

SOUT/y P A i F I C C e <4 4/ S4 COLUMBIA BUEM.ENTURA ECUADOR iGUAYAauti PERU Figure 22. Segmented Chilean Trench (CHILE

PAGE 94

85 limit to the strain that can be supported by rock before fracture. Thus, earthquakes only reach certain maximum magnitudes and tsunamis can be expected to have similar upp . r limits to intensity. Perkins (1972) and McGarr (1976) have demon ^tro Led that future earthquakes cannot have seismic moments (measure of earthquake magnitude for large earthquakes) much larger than, those of earthquakes that occurred in recorded history. The waveforms propagated to the west coast by the deep-ocean model are used as before as input to the nearshore propagation model. Each waveform is propagated from a water depth of 500 meters to shore using the nearshore model and. one of the grids shown in Figure 13. Thus, at each grid location on the shoreline of the west coast, there is a group of 105 waveforms seven waveforms (for intensities from 2 to 5 in onehalf intensity increments) for each segment, of the Aleutian and Peru-Chile Trenches. Each of these waveforms has an associated probability equal to the probability that a certain intensity tsunami will be generated in a particular segment of a trench region. Effect of the Astronomic al Tide s The maximum "still-water" elevation produced during tsunami activity is the result of a superposition of tsunamis and tides. Therefore, she statistical affect of the astronomical tides on total tsunami runup must be include.? in the predictive scheme presented in this study. Since the waveforms calculated by the nearshore model go not have a simple form (e.g. sinusoidal), the statistical effect of the astronomical tide on tsunami runup must be determined

PAGE 95

86 through a numerical approach. (The effect of the astronomical tides was not considered in the section of predictions for the Hawaiian Islands because the tidal range is quite snail for these islands and the local historical data implicitly contain the effects of the tides. ) The waveforms calculated by the nearshore numerical model extend over a period of time of approximately two hours. Three or four wave crests (the largest waves in the tsunami) arrive during this time. Smaller waves arriving at later times, however, often persist for days at a coastal location. An analysis of tide gage records of the 1960 and 1964 tsunamis on the west coast indicates that these smaller waves have amplitudes on the average of 40 percent of the maximum wave amplitude of the tsunami. Therefore, a sinusoidal group of these smaller waves were added to each of the calculated waveforms so that the total waveform extended over a 24 hour period. These smaller waves are important for locations where tsunamis are fairly small compared with tidal variations. At such locations the maximum combined tsunami and astronomical elevation occurs during the maximum tidal elevation. A computer program was developed to predict time-histories of the astronomical tides at all grid locations on the west coast of the United States. The program was based upon the harmonic analysis methods used in the past by the Coast and Geodetic Survey for mechanical tide-predicting machines (Schureman, 1948). Tidal constants available from the Coast and Geodetic Survey were used as input to the computer program; and a year of tidal elevations was then predicted for grid locations all along the west coast. The year 1964 was

PAGE 96

87 selected because all the major tidal components had a node factor equal to approximately 1.00 during this year making it an average year (node factors can vary by about 10 percent). The node factor is associated with the revolution of the moon's node and has an 18.6-year cycle. Since a tsunami can arrive at any time during this 18.6-year period (arrival at a low of the node factor is equally likely as an arrival at a high], the statistical effect of the varying node factor is small and an average value should be used. The statistical effect of the varying node factor on the predicted runup elevations can be shown to be a small fraction of an inch, using the approach discussed in the next section and Appendix A (with the variance of the nodal variation equal to approximately 0.1 feet squared). The year of tidal elevations calculated at each of the nearshore numerical model grid points along the west coast was then subdivided into 15-minute segments. The 24-hour waveforms were allowed to arrive at the beginning of each of these 15-minute segments and then superposed upon the astronomical tide for the 24-hour period. The maximum combined tsunami and astronomical tide elevation over the 24-hour period was determined for tsunamis arriving during each of these 15-minute starting times during a year. All of the maximum elevations had an associated probability equal to the probability that a certain intensity tsunami would be generated in a particular segment of the two trench regions and arrive during a particular 15minute period of a year.

PAGE 97

The many maximum elevations with associated probabilities can be used to determine cumulative probability distributions of combined tsunami and astronomical tide elevations. The maximum elevations are ordered and probabilities summed, starting with the largest elevations, until a desired probability is obtained. The elevation encountered when the summed probabilities reach a desired value P is the elevation that is equaled or exceeded with an average frequency of once every 1/P years. Thus, when the summed probabilities reach the value .01, the elevation associated with the last probability summed is the 100-year elevation. The 100-year and 500-year elevations are determined at all grid points of the nearshore numerical grids using techniques described previously in this chapter. Smooth curves are then used to connect all discrete elevations so that continuous predictions of elevation versus coastline location can be made for all the west coast. Figure 25 shows 100-year and 500-year elevations in the Crescent City, California, area. Comparison with Local O bservation Predictions Crescent City and San Francisco, California, are the only locations on the portion of the west coast of the United States considered in this study that have sufficient historical data of tsunami activity to allow frequency of occurrence predictions based upon local historical observations. Wiegel (1965) made such predictions of tsunami height (trough to crest height) for Crescent City based upon the period from 1900 to 1965. He predicted 100-year and 500-year heights at Crescent City of approximately 25.6 feet and 43.2 feet,

PAGE 98

89 respectively. If the crest amplitude is taken to be one-half the total height, the 100-year and 500-year elevations are 12.8 feet and 21.6 feet, respectively. However, the crest amplitude at Crescent City is typically greater than one-half the wave height (e.g. the crest amplitude of the largest wave of the 1964 tsunami was approximately 60 percent of the total height). If Wiegel's analysis is applied to historical crest elevations instead of heights, the 100-year and 500-year elevations are found to be 15.4 and 26.4 feet, respectively. Furthermore, the analysis can now be applied to the longer time period from 1900 through 1977. The 100-year and 500-year elevations based upon this longer time span are 14.5 feet and 25.5 feet respectively. Figure 23 shows 100-year and 500-year elevations of 13.1 feet and 24.9 feet, respectively. These values compare very favorably with the 14.5 feet and 25.5 feet elevations determined from historical data for the period of time from 1909 through 1977. The elevations predicted by the analysis based upon the local historical data are probably somewhat larger than the elevations predicted in this report because the short time period (relative to 100 or 500 years) from 1900 through 1977 includes the exceptionally active years of tsunami generation from 1946 through 1964 (see Chapter III) . The 13.1-foot and 24.9-foot elevations predicted for Crescent City using the techniques discussed in this study are not totally comparable to the elevations predicted using Wiegel's analysis based upon local historical data since the effect of the astronomical tide has been included in the elevations predicted in this study. However,

PAGE 99

90 z UJ H h < _J 1J'N011VA313

PAGE 100

91 the statistical effect of the astronomical tide on the total elevation is not significant at Crescent City due to the large amplitude of tsunamis there. Appendix A shows that the statistical effect of the astronomical tide for a location where tsunamis are large is to 2 ceo increase the predicted elevation by an amount equal to — x— , 2 °° 2 where a is the tidal variance and equals E C , C is 1 m m m=l equal to the m tidal constituent, and a is given by the following expression [since it was shown earlier in this chapter that elevations are linearly related to the logarithm of probabilities): P(Z) Ae~ a ~ with Z = the elevation above local mean sea level, and P(Z) = the cumulative probability distribution for the elevation at a given site being equal to or exceeding Z due only to the maximum wave of the tsunami. Based upon the tidal constituents 2 predicted by the Coast and Geodetic Survey, o =7.1 fox Crescent City. Using Wiegel's data for the period of time from 1900 through 1977 gives a = .145 . Therefore, the astronomical tide contributes approximately 0.5 feet to the elevations predicted for Crescent City. Frequency of occurrence calculations just inside San Francisco Bay were made in an earlier study (Houston and Garcia, 1975) using a simple onedimensional analytical solution for nearshore propagation. However, portions of San Francisco Bay are included in the nearshore numerical model (see Figure 13) used in this study in order that the effect of the Bay on elevations outside the Bay can be properly

PAGE 101

92 simulated. Thus, approximate elevation calculations can be made in this study at the location of the Presidio tide gage within San Francisco Bay. Figure 24 shows the 100-year and 500-year elevation predictions at the Presidio to be 6.1 feet and 10.0 feet, respectively. Local historical data of tsunami activity at the Presidio was compiled by Wiegel (1965). Using a logarithmic distribution and Wiegel's wave height data for San Francisco over the period 1900 through 1977 yields 100-year and 500-year heights (trough to crest) of 8.0 feet and 14.2 feet respectively. If the crest amplitude is one-half the height, this yields 100-year and 500-year elevations of 4.0 feet and 7.1 feet respectively. The same analysis based upon data of crest amplitudes yields 100-year and 500-year elevations of 4.4 feet and 7.8 feet, respectively. It is difficult to compare the elevations predicted using local historical data with the elevations predicted in this study for San Francisco, since the effect of the astronomical tides (included in the calculations of this study) cannot easily be estimated, as was the case for Crescent City. Tsunamis recorded at the Presidio in San Francisco are known to persist at fairly substantial levels for extended periods of time apparently due to an oscillation phenomenon. For example, during the 1964 tsunami there were six waves of approximately equal or greater amplitude than the initial wave (in addition to smaller waves) during the first 8 hours of the tsunami. Oscillations approximately 40 percent of the amplitude of the initial wave persisted for at least another 24 hours. The tidal range at San Francisco is 5.7 feet with mean higher high water (MHHW) 2.7 feet

PAGE 102

U'NO!JLVA313

PAGE 103

94 above mean sea level (MSL) and mean high water (MHW) 1.0 feet above MSL. At least one of the main waves of a tsunami at San Francisco probably superposes upon either mj iiv' (with, the smaller oscillations then adding to MHHW) or MJIIIW. Thus, the effect of the tides is to contribute 1 to. 2. 7 feet to the total elevation. These elevations when added to the elevations predicted from local historical data yield 100-year and 500-year elevation:; of 5.4 feet to 7.1 feet and 8.8 feet to 10.5 feet, respectively. These elevations compare favorably with the 6.1 feet and 10.0 feet elevations determined in this study.

PAGE 104

CHAPTER V CONCLUSIONS The major objective of the present study has been to determine tsunami elevation frequencies of occurrence at the shoreline for the Hawaiian Islands and the west coast of the continental United States. Numerical models were developed to simulate nearshore tsunami propagation and interaction with coastlines. A finite element numerical model was used to interpolate between modern historical measurements of tsunami elevations in the Hawaiian Islands. Finite difference numerical models were used to propagate tsunamis from the AleutianAlaskan and Peru-Chile regions to the west coast of the United States so that frequency of occurrence predictions could be made. The good agreement between tide gage recordings and numerical simulations of the 1960 Chilean and 1964 Alaskan tsunamis in the Hawaiian Islands indicates that linear, nondispersive, and dissipationfree equations govern nearshore propagation of the main leading waves of long period tsunamis for coastal regions having very short continental shelves. A short continental shelf limits the time available for the development of effects governed by the terms neglected in the equation of motion. Long wave equations that include nonlinear and bottom stress effects were shown to govern propagation of the main leading waves of long period tsunamis over the shallow water region of the west coast of the United States. 95

PAGE 105

9u Good agreement between the numerical model simulations and tide gage recordings of historical tsunamis confirms the experimental flume tests by Hammack (1972) that showed that the far-field waveform generated by a large impulsive ground displacement does not depend upon the transient motions accompanying the permanent ground motion. The good agreement also shows that the accepted permanent ground deformations of the 1960 Chilean and 1964 Alaskan tsunamis, which were partly measured and partly estimated by Plafker and Savage (1970) and Plafker (1969) respectively, are fairly accurate. Finally, this study shows that numerical models can be used in conjunction with limited available historical information to make reasonable predictions of tsunami elevations on coastlines. For an area having local historical data of tsunami occurrence, such as the Hawaiian Islands, numerical models can be used to predict elevations at particular locations not having data. For an area lacking almost any local historical data, such as the west coast of the United States, elevation predictions can be made based upon data of tsunami occurrence in tsunamigenic regions and numerical model propagation of tsunamis from generation regions to the area of interest.

PAGE 106

APPENDIX SUPERPOSITION OF LARGE TSUNAMIS AND ASTRONOMICAL TIDES If the largest wave of a tsunami produces the largest total elevation when the tsunami is superposed upon the astronomical tide (regardless of where the largest wave falls in the tidal cycle), the statistical effect of the astronomical tide on the total elevation can be estimated analytically. Tsunamis satisfy this condition when their maximum waves have amplitudes larger than the sum of the maximum tidal elevation and the elevation of one of the smaller waves which follow the major waves and persist over periods of days. Let Z = the elevation at any time above local mean sea level. P„(Z) = the cumulative probability distribution for the elevation at a given site being equal to or exceeding Z due only to the maximum wave of the tsunami. P (Z) = the probability of the elevation at the same location being equal to or exceeding Z clue only to the astronomical tide. P(Z) = the cumulative probability distribution for the elevation at a given site being equal to or exceeding Z due to the maximum wave of the tsunami and the astronomical tide. According to Chandrasekhar (1943), P(Z) can be calculated from the equation P(Z) = / f g (A)P s (Z-A)dX (Al) 97

PAGE 107

98 where dP R (Z) 'fr ' dz and f fi^ is the P robabilit >' density for the astronomical tide. Pg(Z) can be approximated by a Gaussian distribution. Therefore, c cn> 1 -Z /2a f B ^) -— e (A2) /2ir c where the variance is 2 °° 2 a = E C and C equals the m tidal constituent. P S (Z) can be approximated by an exponential function (see Chapter II) as follows: P S C) = Ae~ uZ (A3) Substituting Equations (A2) and (A3) into (Al) and performing the integration yields P(Z) = Ae 2 a(Z =-) 7 1 or p(z) = Ae" a " J where z 1 = Z °°-2 Thus, the net effect of the astronomical tide is to produce a P(Z) 2 identical with P„(Z) except for a shift of Z by an amount ^~ o ' 2 2 The variance a can be determined for locations in the United States from tidal harmonic coefficients available from the Coast and Geodetic Survey. The exponential constant a in Equation (A3) is obtained from a least-squares exponential fit of tsunami elevation data at a location.

PAGE 108

REFERENCES Adams, W. M., 1970: Tsunami effects and risk at Kahuku Point, Oahu, Hawaii. E ngineering Geology Case Histories, Geo. Soc. of Amer., No. 8, 63-70. ' """ " "' Adams, W. M. , 1975: Conditional expected tsunami inundation for Hawaii. J. of Water ways, Harbors, and Coa stal Eng. Div., Amer. Soc. Civ. Eng. , 101, No. WW4, 319-529. ~ ~ Bernard, E. N. and A. C. Vastano, 1977; Numerical computation of tsunami response for island systems. J. of P lrys. Oceanog r., 7, 389-595. Borgman, L. E. and D. T. Resio, 1977: Extremal prediction of wave climatology. Ports 77, 4th Ann. Symp. of the Waterways, Port , Coastal, and Ocean Div., Amer. Soc. Civ. Eng., 1, 394-412. Carrier, G. F., 1966: Gravity waves on water of variable depth, J. of F luid Me c h., 24, Pt. 4, 641-659. Chandrasekhar, S., 1945: Stochastic problems in physics and astronomy. Rev, of Mod. Phys. 15, 1-89. Chen, H. S. and C. C. Mei, 1974: Oscillations and wave forces in an offshore harbor. Rept. No. 190, Mass. Inst, of Tech., Cambridge, Mass., 215 pp. Cole, R. A. and P. J. Riech, 1976: Computer program for drawing finite element grids. Tech. Rept. N-76-2, U. S. Army Eng. Waterways Exp. Sta., Vicksburg, Miss. 110 pp. Courant, R., 19^5: Variational methods for the solution of problems of equilibrium and vibration. Bull, of the Amer. Math. Sec, 49, 1-25. Cox, D. C, 1964: Tsunami height-frequency relationship of Hilo. Nov. 1964, Hawaii Inst, of Geophys., Univ. of Hawaii, Honolulu, Hawaii, 7 pp. Der Kiureghian, A. and A. H-S. Ang, 1977: A fault-rupture model for seismic risk analysis. Bull, of the Seis. Soc. of A mer., 64, No. 4, 1173-1194. " ' "" " " " Desai, C. S. and J. F. Abel, 1972: Introd u ction t o the Finite Element Method; A Numer ical Meth od for En gineering Analysis, Van Nostrand Reinhold" Co., 477 pp. " " 99

PAGE 109

100 Dronkers, J. J. 1964: Tidal Computations in Rivers and C oastal Waters, North-Holland Publishing Co., ' Amsterd am, 51 S pp . ' Eaton, J. P., D. H. Richter, and W. V. Aults, 1961: The tsunami of May 23, 1960, on the island of Hawaii. Bull, of the Seis. Soc. of Airier . , 51, No. 2, 135-157. Grace, R. A., 1969: Reef runway hydraulic model study: Honolulu International Airport, three-dimensional tsunami study. Dept. of Transp., St. of -Hawaii, Honolulu, 23 pp. Gumbel, E. J., 1955: The calculated risk in flood control, Add. Sci. Res., 5, Sec. A, 273-280. ~ ~ Gutenberg, B. and C. F. Richter, 1965: Seismicity of the Earth and Associated Phenomena, 2nd ed., Hafner Publishing Co., 310 pp~ Hall, J. V. and G. M. Watts, 1953: Laboratory investigation of the vertical rise of solitary waves on impermeable slopes; Tech. mem. 35, U. S. Corps of Eng., Beach Erosion Board, 75 pp. Hammack, J. L., 1972: Tsunamis; a model of their generation and propagation. Rept. No. KH-R-28, Calif. Inst, of Tech. Pasadena, Calif., 261 pp. Hammack, J. L., and H. Segur, 1977: Modeling criteria for long water waves. Univ. of Fla., Gainesville, Fla. (to be published). Hatori, T., 1965: Directivity of tsunamis, B ull, of the Earthq uake Res. Inst., Univ., of Tokoyo, 4T, 61-81. Heitner, K. L. and Housner, G. W., 1970: Numerical model for tsunami runup. J. of the Waterways, Harbors, and Coastal F.ng. Div., Amer. Sec. Civ . Eng. 96, No. WW3, 701-719. " Hom-ma, S., 1950: On the behavior of seismic sea waves around a circuJar island. Geophys. Mag., 21, 199-208. Houston, J. R., R. P. Cax-ver, and D. G. Markle, 1977: Tsunami-wave elevation frequency of occurrence for the Hawaiian Islands. Tech. Rept. H-77-16, U. S. Army Eng. Waterways Exp. Sta., Vicksburg, Miss . , 109 pp. Houston, J. R. and A. W. Garcia, 1974: Type 16 flood insurance study: tsunami predictions for Pacific coastal communities . Tech. Rept. H-74-3, U. S. Army Eng. Waterways Experiment Sta., Vicksburg, Miss., 128 pp. Houston, J. R. and. A. W. Garcia, 1975: Type 16 flood insurance study: tsunami predictions for Monterey and San Francisco Bays and Puget Sound. Tech. Rept. H-75-17, U. S. Army Eng. Waterways Exp. Sta., Vicksburg, Miss. , 265 pp.

PAGE 110

101 Houston, J. R., R. K. Whalin, H. L. Butler, D. Raney, 1975 a: Probable maximum tsunami runup for distant seismic events NORCO-NP-1, PSAR, Fugro Inc., Long Beach, Calif., 253 pp. Houston, J. R., R. IV. Khalin, A. W. Garcia, H. L. Butler, 1975 b: Effect of source orientation and location in the Aleutian Trench on tsunami amplitude along the Pacific coast of the continental United States. Res. Rept . H-75-4, U. S. Army Eng . Waterways Exp. Sta., Vicksburg, Miss., 48 pp. Hwang, L. S., H. L. Butler, and D. J. Divoky, 1975: Tsunami model: generation and open-sea characteristics. Bull, of the Seis. Soc. of A mer., 62, No. 6, ^9-1596. Hwang, L. S. and D. Divoky, 1975: Numerical investigations of tsunami behavior. Tetra Tech, Inc., Pasadena, Calif., 37 pp. Jordaan, J. M. and W. M. Adams, 1968: Tsunami height, Oahu, Hawaii: Model and nature. Proc. of the Eleventh Conf. on C oastal Eng., Inst, of Civ. Eng. and the .Amer. Soc. of Civil Eng. 7 1555-1574. Kanamori, H., 1977: The energy release in great earthquakes. J. of Geoph ys. Res., 82, No. 20, 2981-2987. Kaplan, K. , 1955: Generalized laboratory study of tsunami runup. Tech. Memo. No. 60, U. S. Army Corps of Eng., Beach Erosion Board, 30 pp. Kelleher. J., J. Savino, H. Rowlett, and W. McCann, 1974: Why and where great thrust earthquakes occur along island arcs. J. of Geophys. Res., 79, No. 32, 4SS9-4899. ~ Keller, J. B. and H. B. Keller, 1964: Water wave run-up on a beach. Service Bureau Corp., New York, 34 pp. Kishi, J. and H. Saeki, 1966: The shoaling, breaking, and runup of solitary waves on impermeable rough slopes. Proceed, of the Tenth Conf. on Coastal Eng., .Amer. Soc. of Civ. EngT~7 1, 322-348. Leendertse, J. J., 196": Aspects of a computational model for longperiod water-wave propagation. Erne. Rm-5294-PR, Rand Corp., Santa Monica, Calif., 165 pp. Loomis H. G. 1970: A method of setting up the eigenvalue problem for the linear, shallow water wave equation for irregular bodies of water with variable water depth and application to bays and harbors in Hawaii. NOAA-JTRE-16, Hawaii Inst, of Geophys., Univ. of Hawaii, Honolulu, Hawaii, 77 pp. Loomis, H. B. 1976: The tsunami of November 29, 1975 in Hawaii, HIG-75-21, Hawaii Inst, of Geophys., Univ. of Hawaii, Honolulu, Hawaii, 95 pp.

PAGE 111

io: McGarr, A., 1976: Upper limit to earthquake size. Nature, 262, No. 5567, 378-379. Mader, C. L., 1974: Numerical simulation of tsunamis. J_. o f Phys. , 4, 74-82. Madsen, 0. S. and C. C. Mei, 1969: The transformation of a solitary wave over an uneven bottom. J. of Fl uid Mech. , 59, 781-791. Magoon, 0. T., 1965: Structural damage by tsunamis. Coastal Eng. Santa Barbara Spec, Con f . , Amer. Soc. of Civ. Eng. 356S. Mei, C. C, 1978: Numerical methods in water-wave diffraction and radiation. An n. Rev, of Fluid Mech ., 10. Mei, C. C. and B. LeMehaute, 1966: Note on the equations of Long waves over an uneven bottom. J. of Geophys. Res., 71, 2. New-mark, N. M. and E. Rosenblueth, 1971: Fundamentals of Earthquake Engineering, Prentice-Hall, 528 pp. Noye, B. J., 1970: The frequency response of a tide-well: Proc. of the Third Australasian Conf. on Hyd. and Fluid Mech., Inst, of Fng. , Sydney, Australia, No. 2630. Omer, G. C. and H. H. Hall, 1949: The scattering of tsunamis by a cylinderical island. Bull. Seis. Soc. Ame r., 59 , 257-260. Parraras-Carayannis, G., 1977: Catalog of tsunamis in Hawaii. Rept. SE-4, World Data Center A for Solid Earth Geophysics, Boulder, Colo., 77 pp. Peregrine, D. H., 1967: .Long waves on a beach. J_. o f Fluid Mechanics, 27, pt. 4, 815-827. Perkins, D., 1972: The search for maximum magnitude. Nat. Oceanogr. and Atmos. Admin. Earthquake infor. Bull., July, 1972, 18-23. Piafker, G. , 1955: Tectonics of the March 27, 1964, Alaskan earthquake. U. S. Geologica l Survey Professional Pa p er, 543-1 , 11-174. Piafker, G. and J. C. Savage, 1970: Mechanism of the Chilean earthquake of May 21 and 22, 1960. Bull, of Geo. Soc. of .Amer., 81 No. 4, 1001-1030. Platzman, G. W., 1958: A numerical computation of the surge of 26 June 1954 on Lake Michigan. Geophysica, 6, 409-43S. Raney, D. C., 1976: Numerical analysis of tidal circulation for Long Beach Harbor-existing conditions and alternate plans for Pier J completion, Misc. Paper H-76-4, U. S. Army Eng. Waterways Exp. Sta., Vicksburg, Miss., 138 pp.

PAGE 112

103 Rascon, 0. A. and A. G. Villarreal, 19 7 5: 0a a stochastic model to estimate tsunami risk. J. o f Hyd. Re s., 15, No. 4, 583-405. Savage, R. P., 1958: Wave run-up on roughened and permeable slopes. J. of Wate rways and Harbors Piv., Amcr. Soc. of Civ. Eng., 84, No. WW3, 1640-1641. Saville, T. , 1956: Wave run-up on shoreline structures. J. of Waterway s and Harbors Div., Amer. Soc. of Civ. Eng., 82, No. WW2, 925-950. Schureman P., 1948: Manual of harmonic analysis and prediction of tides. Spec. Public S6, U. S. Coast and Geod. Sur., Washington, D. C. 317 pp. Soloviev, S. L., 1969: Recurrence of tsunamis in the Pacific. Proc. of the inter. Sym p. on Tsunamis and T suna mi Res. , 35-63. Soloviev, S. L. and No. Go, 1969: Ca talog of Tsunamis in the Pacific (Main Data), in Russian, Union of Soviet Socialist Republics, Moscow. Strang, G. and G. J. Fix, 1975: An Analysis of the Finite Element Method, Prentice-Hall, 528 pp. Svendsen, I. A., 1974: Cnoidal waves over a gently sloping bottom. Inst, of Hydrodyn. and Hydr. Eng., Univ. of Denmark, Lyngby, Denmark, 181 pp Takahasi, R. 1964: Coastal effects upon tsunamis and storm surges. Bull, of the Earthquake Res. Inst., Univ. of Tokoyo, 42, 175-180. Tilling, R. I., R. V. Koyanagi, R. W. Lipman, J. P. Lockwood, J. G. Moore, and D. A. Swanson, 1976: Earthquake and related catastrophic events; island of Hawaii, November 29, 1975: a preliminary report. Geological Survey Circular, No. 740, U. S. Geological Survey, Arlington, Virginia, 33 pp. Towill Corp., 1975: Tsunami studies for Oahu, Hawaii. Towill Corp., Honolulu, Hawaii, 53 pp. Turner, M. J., R. W. Clough, H. C. Martin, and L. P. Topp, 1956; Stiffness and deflection analysis of complex structures, J. of the Aeronautical Sci . , 23, No. 9, 805-823, 854. Ursell, F., 1955: The long-wave paradox in the theory of gravity waves. Proceedings, Cambridge Phi lsophical Society, 49, 685-694. Van Dorn, W. G., 1970: Tsunami response at Wake Island: A model study. J. Mar. Res., 28, 356-344.

PAGE 113

104 Vastano, A. C. and E. N. Bernard, 1973: Transient long wave response for a multiple island system. J. Phys . Oceanogr. , 5, 406-418. Vastano, A. C. and R. 0. Rcid, 1966: A numerical study of the tsunami response of an island. Ref. 66-26T, Dept. of Oceanography, Texas A§M University, 141 pp. Wiegel, R. L., 1965: Protection of Crescent City, California from tsunami waves. Redevelopment Agency of the City of Crescent City, 114 pp. Wilson, B. V,'., 1959: Earthquake occurrence and effects in ocean area. CR 69.027, Naval Civ. Eng. Lab., Port Hueneme, Calif., 118 pp. Wilson, B. W. and A. Torum, 1968: The tsunami of the Alaskan earthquake, 1964: engineering evaluation, Tech. Mem. 25, U. S. Army Coastal Eng. Research Center, Washington, D. C, 401 pp. Wybro. P. C, 1977: M. S. Plan B, Univ. of Hawaii, Honolulu, Hawaii, (to be published) . Zienkiewicz, 0. C. and Y. K. Cheung, 1966: Finite elements in the solution of field problems. The Enginee r, 507-510. Zienkiewicz, 0. C, 1971: The finite element method in engineering science, 2d. ed., McGraw Hill, 521 pp.

PAGE 114

BIOGRAPHICAL SKETCH James Robert Houston was born in Berkeley, California, on March 31, 1947. He was awarded a Bachelor of Arts degree in Physics at the University of California at Berkeley in June 1969 and a Masters of Science degree in Physics at the University of Chicago in June 1970. His doctoral studies at the University of Chicago were interrupted in November 1970 when he was inducted into the U. S. Army. For a year and a half he was involved in the study of explosively generated water waves in the Nuclear Weapons Effects Laboratory of the U. S. Corps of Engineers Waterways Experiment Station in Vicksburg, Mississippi. Following an honorable discharge from the U. S. Army in May 1972, he was engaged as a Research Physicist by the Wave Dynamics Division, Hydraulics Laboratory, Waterways Experiment Station. In September 1974, he was awarded a Master of Science degree in Coastal Engineering at the University of Florida. He is the author of several reports and papers dealing with numerical investigations of long wave phenomena. 105

PAGE 115

I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy % •l -J^A&i^fr Knox Mi 11 saps, Chairman^' Professor and Chairman of Engineering Sciences I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Joseph L. Hammack, Cochairman Associate Professor of Engineering Sciences I certify that I have read this study and that, in my opinion it conforms tc acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Sheppard Associate Professor of Coastal and Oceanographic Engineering and Engineering Sciences

PAGE 116

I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. m. A. Eisenberg Professor of Engineering and Engineering Sciences I certify that 1 have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. w *i ?t&r>4L0 Michael E. Thomas Professor and Chairman of Industrial and Systems Engineering This dissertation was submitted to the Graduate Faculty of the College of Engineering and to the Graduate Council, and was accepted as partial fulfillment of the requirements for the degree of Doctor of Philosophy. June 1978 U^i cu. V 7 Dean, College of Engineering Dean, Graduate School

PAGE 117

UNIVERSITY OF FLORIDA 3 1262 08553 2686