SECOND EDITION
NUMERICAL MODELING of WATER WAVES
SECOND EDITION
NUMERICAL MODELING of WATER WAVES Charles L.Mader
...

Author:
Charles L. Mader

This content was uploaded by our users and we assume good faith they have the permission to share this book. If you own the copyright to this book and it is wrongfully on our website, we offer a simple DMCA procedure to remove your content from our site. Start by pressing the button below!

SECOND EDITION

NUMERICAL MODELING of WATER WAVES

SECOND EDITION

NUMERICAL MODELING of WATER WAVES Charles L.Mader

CRC PRESS Boca Raton London New York Washington, D.C.

This edition published in the Taylor & Francis e-Library, 2005. To purchase your own copy of this or any of Taylor & Francis or Routledge’s collection of thousands of eBooks please go to www.eBookstore.tandf.co.uk. Library of Congress Cataloging-in-Publication Data Mader, Charles L. Numerical modeling of water waves/Charles L.Mader.—2nd ed. p. cm. Includes bibliographical references and index. ISBN 0-8493-2311-8 (alk. paper) 1. Water waves—Mathematical model. 2. hydrodynamics. I. Title. QA927.M29 2004 532′.593–dc22 2004047812 This book contains information obtained from authentic and highly regarded sources. Reprinted material is quoted with permission, and sources are indicated. A wide variety of references are listed. Reasonable efforts have been made to publish reliable data and information, but the author and the publisher cannot assume responsibility for the validity of all materials or for the consequences of their use. Neither this book nor any part may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, microfilming, and recording, or by any information storage or retrieval system, without prior permission in writing from the publisher. The consent of CRC Press LLC does not extend to copying for general distribution, for promotion, for creating new works, or for resale. Specific permission must be obtained in writing from CRC Press LLC for such copying. Direct all inquiries to CRC Press LLC, 2000 N.W. Corporate Blvd., Boca Raton, Florida 33431. Trademark Notice: Product or corporate names may be trademarks or registered trademarks, and are used only for identification and explanation, without intent to infringe. Visit the CRC Press Web site at www.crcpress.com © 2004 by CRC Press LLC No claim to original U.S. Government works ISBN 0-203-49219-6 Master e-book ISBN

ISBN 0-203-59054-6 (Adobe e-Reader Format) International Standard Book Number 0-8493-2311-8 (Print Edition) Library of Congress Card Number 2004047812

ACKNOWLEDGMENTS The author gratefully acknowledges the developers of the numerical models described in this monograph—M.L.Gittings of SAIC, R.Weaver, G.Gisler, A.Amsden, C.W.Hirt, F.H.Harlow, B.D. Nichols, L.R.Stein, and D.Butler; and the contributions of J.D. Kershner, R.E.Tangora, A.N.Cox, K.H.Olsen, R.Menikoff, E. M.Kober, B.Gonzales, J.C.Solem, J.G.Hills, R.D.Cowan, B.G. Craig, D.Venable, G.E.Seay, and A.L.Bowman of the Los Alamos National Laboratory. The author wishes to recognize the contributions of Dr. H. Loomis, Dr. W.Adams, Dr. D.Cox, Dr. D.W.Moore, Mr. G. Curtis, Dr. J.Miller, Dr. D.Walker, Dr. S.Furumoto, Dr. L. Spielvogel, Dr. A.Malahoff, Dr. W.Dudley, G.Nabeshima, and S. Lukas of the University of Hawaii. The author wishes to recognize the collaborations with Dr. E. Bernard of the Pacific Marine Environmental Laboratory, Professor George Carrier of Harvard University, Professor Zygmunt Kowalik of the University of Alaska, Dr. Hermann Fritz of the University of Georgia, Professor Valery Kedrinskii of the Institute of Hydrodynamics, Novosibirsk, Russia, Mr. Tom Sokolowski of the West Coast and Alaska Tsunami Warning Center, Mr. Bruce A.Campbell, Mr. Dennis Nottingham of PN&D, Dr. T.S.Murty of Baird & Associates, Dr. David Crawford of Sandia and Dr. Rob T.Sewell of R. T.Sewell Associates. This book was reviewed and determined to be unclassified by the Los Alamos National Laboratory group S-7 and given the designation LA-UR-03–8710.

THE AUTHOR Dr. Charles L.Mader is President of Mader Consulting Co. with offices in Honolulu, Hawaii, Los Alamos, New Mexico and Vail, Colorado. He is an Emeritus Fellow of the Los Alamos National Laboratory. He is the author of Numerical Modeling of Water Waves, Numerical Modeling of Explosives and Propellents, Numerical Modeling of Detonations and four Los Alamos Dynamic Material Properties data volumes. Dr. Mader is a Fellow of the American Institute of Chemists, Editor of the Science of Tsunami Hazards Journal and is listed in Who’s Who in America and in Who’s Who in the World. His web site is http://www.mccohi.com.

INTRODUCTION This second edition of Numerical Modeling of Water Waves describes the technological revolution that has occurred in numerical modeling of water waves during the last decade. A CD-ROM with many of the FORTRAN codes of the numerical methods for solving water wave problems and computer animations of the problems that have been solved using the codes is included with the book. Several PowerPoint presentations describing the modeling results are on the CD-ROM. It will be called the NMWW CDROM in the rest of the book. The objective of this book is to describe the numerical methods for modeling water waves that have been developed primarily at the Los Alamos National Laboratory over the last four decades and to describe some examples of the applications of these methods. Some of the applications of the numerical modeling methods were performed while the author was working at the Joint Institute for Marine and Atmospheric Research at the University of Hawaii, some as Mader Consulting Co. research projects, and the rest as a Fellow of the Los Alamos National Laboratory. Although the two- and three-dimensional numerical methods for modeling water waves had been available in the 1980’s for several decades, they had seldom been used. A major obstacle to their use was the need for access to large and expensive computers. By the 1980’s, inexpensive personal computers were adequate for many applications of these numerical methods. In this book, the basic fluid dynamics associated with water waves are described. The common water wave theories are reviewed in Chapter 1. A computer code called WAVE for personal computers that calculates the wave properties for Airy, third-order Stokes, and Laitone solitary gravity waves is available on the NMWW CD-ROM. The incompressible fluid dynamics model used for shallow water, long waves is described in Chapter 2. A computer code for personal computers using the shallow water model is called SWAN and is available on the NMWW CD-ROM. The SWAN code is used to model the 1946, 1960 and 1964 earthquake generated Hilo, Hawaii tsunamis, the 1964 Crescent City, California tsunami and the 1994 underwater landslide generated Skagway, Alaska tsunami. The two-dimensional incompressible Navier-Stokes model used for solving water wave problems is described in Chapter 3. A computer code for personal computers using the two-dimensional Navier-Stokes model is called ZUNI and is available on the NMWW CD-ROM. The ZUNI code is used to model tsunami wave propagation and flooding. It is also used to model the effect of underwater barriers on tsunami waves. The three-dimensional incompressible Navier-Stokes model for solving water wave problems of any type is described in Chapter 4. The equations used in the computer code called SOLA-3D are described and the FORTRAN code is on the NMWW CD-ROM. The SOLA-3D code is used to model the 1975 Hawaiian tsunami and the 1994 Skagway tsunami by water surface cavities generated by underwater landslides.

It is surprising that most academic and government modelers of water waves have chosen to use shallow water or other incompressible models of limited validity for modeling water waves and not use the incompressible Navier-Stokes model since the first edition of this book was published. The severe limitations of the shallow water model are described in Chapter 5. The generation of water waves by volcanic explosions, conventional or nuclear explosions, projectile and asteroid impacts require the use of the compressible NavierStokes model. The numerical models and codes for solving such problems have recently become available as part of the Accelerated Strategic Computer Initiative program. A computer code for solving one, two and three dimensional compressible problems called SAGE, NOBEL or RAGE is described in Chapter 6 and some of its remarkable capabilities are presented. These include modeling of the KT Chicxulub asteroid impact and the modeling of the largest historical tsunami which occurred July 8, 1958, at Lituya Bay, Alaska. The modeling of the Lituya Bay impact landslide generated tsunami and the flooding to 520 meters altitude is described. A color videotape (VTC-86–4) lecture featuring computer generated films of many of the applications discussed in this book is available from the Los Alamos National Laboratory library. Several web sites have been established that contain additional information and computer movies of the problems described in this book. The major sites are http://www.mccohi.org and http://t14web.lanl.gov/Staff/clm/tsunami.mve/tsunami.html.

CONTENTS

CHAPTER 1— CHAPTER 2— CHAPTER 3— CHAPTER 4— CHAPTER 5—

WATER WAVE THEORY THE SHALLOW WATER MODEL THE TWO-DIMENSIONAL NAVIER-STOKES MODEL THE THREE-DIMENSIONAL NAVIER-STOKES MODEL EVALUATION OF INCOMPRESSIBLE MODELS FOR MODELING WATER WAVES CHAPTER 6— MODELING WAVES USING COMPRESSIBLE MODELS CD-ROM CONTENTS AUTHOR INDEX SUBJECT INDEX

1 24 80 111 134 165 216 223 227

1 WATER WAVE THEORY 1A. The Equations of Fluid Dynamics The equations of fluid dynamics are representations of the laws of conservation of mass, momentum, and energy as applied to a fluid. The general equations are often called the Navier-Stokes equations. In addition, we need an equation of state to describe the properties of the fluids. The equation of state expresses the fact that the pressure is everywhere a function of the density and the energy per unit of mass. For the water wave theories we will describe in this chapter, the fluid will be considered to be incompressible and the initial density of each element of fluid remains constant. The numerical models of fluid dynamics for compressible fluids are described in the Los Alamos monograph Numerical Modeling of Detonations 1 , and in Chapter 6. The equations of fluid dynamics are described by Landau and Lifshitz 2 as a branch of theoretical physics in the classic volume 6 of their Course of Theoretical Physics. The equations of fluid dynamics are treated as a branch of applied mathematics in Batchelor’s 3 monograph entitled An Introduction to Fluid Dynamics. A classic description of the equations of fluid dynamics as a branch of oceanography is Blair Kinsman’s 4 textbook Wind Waves. Any of these exhaustive treatments of the equations of fluid dynamics are recommended if the reader wishes to reaffirm his faith in the efficacy of Newtonian mechanics and the fluid continuum. In this monograph, we will first derive the common water wave description starting with the three-dimensional equations of fluid dynamics for viscous, compressible flow. The Nomenclature gx,gy,gz

acceleration due to gravity in x, y, or z direction

I

internal energy

P

pressure

V

specific volume

Ux

particle velocity in x direction

Uy

particle velocity in y direction

Uz

particle velocity in z direction

ρ

density

q

viscosity

Numerical modeling of water waves qx, qy, qz

viscosity deviators

Sx, Sy, Sz

total deviators

λ, µ

real viscosity coefficients

C

wave velocity

H

maximum wave height

h

wave height

L

wave length

T

wave period

f

frequency

n

wave number

D

depth

2

The three-dimensional partial differential equations for viscous, compressible flow: Eulerian Conservation of Mass

or

1.1

because

and

Water wave theory

3

Viscosity with No Shearing Forces

Sx=qx−q. Sy=qy−q. Sz=qz−q. Eulerian Conservation of Momentum

σ=δijP−Sij.

(1.2)

Numerical modeling of water waves

4

Eulerian Conservation of Energy

where

and

(1.3)

If the flow is incompressible so ρ=ρ0, equation 1.1 becomes

(1.4) equation 1.2 becomes, without viscosity,

(1.5a)

(1.5b)

Water wave theory

5

(1.5c)

and equation 1.3 becomes, without viscosity,

(1.6)

Viscosity terms may be included in equation 1.5 as follows. If shearing forces are included and the fluid is incompressible, then

q=0.0, so equations 1.5 become

(1.7a)

Numerical modeling of water waves

6

(1.7b)

and

(1.7c)

1B. Water Wave Descriptions A useful description of water waves has been given by Bernard Le Mehaute in reference 5. Parts of that report are included in the following description of water waves. The theories for unsteady free surface flow subjected to gravitational forces include motions called water waves. They are also called gravity waves. Other gravity waves include motions in other fluids such as atmospheric motions. A great variety of water waves exists. Water wave motions include storm waves generated by wind in the oceans, flood waves in rivers, seiche or long-period oscillations in harbor basins, tidal bores or moving hydraulic jumps in estuaries, waves generated by a moving ship in a channel, tsunami waves generated by earthquakes, and waves generated by explosions near or under the water. Mathematically, a general solution does not exist for gravity waves and approximations must be made for even simple waves. One of the important problems in water wave theory is to establish the limits of validity of the various solutions that are due to the simplifying assumptions. The mathematical treatments of the water wave motions use all the mathematical physics dealing with linear and nonlinear problems. The main difficulty in the study of water wave motion is that the free surface boundary is unknown. Water wave motions are so varied and complex that any attempt at classification may be misleading. Any definition corresponds to idealized situations that never occur. For example, a purely two-dimensional motion never exists. It is a convenient mathematical concept that is physically best approached in a tank with parallel walls. Boundary layer effects and transverse components still exist although they are small and may be neglected in many applications of the theory. Essentially, two kinds of water waves exist, oscillatory waves and translatory waves. In an oscillatory wave, the transportation of fluid or mass transportation does not occur. The wave motion is then analogous to the transverse oscillation of a rope. A translatory wave involves a transport of fluid in the direction in which the wave travels. For example, a moving hydraulic jump such as a tidal bore is a translatory wave.

Water wave theory

7

An oscillatory wave can be progressive or standing. Consider a disturbance H(x, t) such as a free surface elevation traveling along the OX axis at a velocity C. The characteristics of a progressive wave remain identical for an observer traveling at the same speed and in the same direction as the wave (Figure 1.1). Where h can be expressed as a function of (x−Ct) instead of (x, t), a “steady-state” profile is obtained. h(x−Ct) is the general expression for a progressive wave of steady-state profile traveling in the positive OX direction at a constant velocity C. Where the progressive wave is moving in the opposite direction, its mathematical form is expressed as a function of (x+Ct). The definition of wave velocity C for a non-steady-state profile is unphysical, because each “wave element” travels at its own speed, thereby causing wave deformation.

Fig. 1.1. A Progressive Wave.

The simplest case of a progressive wave is the wave defined by a sine or cosine curve such as

or

Such a wave is called a harmonic wave, where H/2 is the amplitude and H the wave height. The distance between the wave crests is the wave length L, and L=CT, where T is the wave period. The wave number n=2π/L is the number of wave lengths per cycle. The frequency is f=2π/T. Hence, the previous equations can be written

or

Numerical modeling of water waves

8

A standing or stationary wave is characterized by its mathematical description as a product of two independent functions of time and distance, such as

or more generally, F=F1(x)·F2(x). A standing wave can be considered as the superposition of two waves of the same amplitude and period traveling in opposite directions. Where the convective terms are negligible, the standing wave motion is defined by a linear addition of the equations for the two progressive waves.

Fig. 1.2. A Standing Wave.

A standing wave generated by an incident wind wave in relatively shallow water is called a seiche. A seiche is a standing oscillation of long period encountered in lakes and harbor basins. The amplitude at the node is zero, and at the antinode, it is H, as shown in Figure 1.2. Two waves of the same period but different amplitudes traveling in opposite directions can be defined linearly by the summation of A sin(x−Ct)+B sin(x+Ct). It can also be considered as the superposition of a progressive wave with a standing wave and is encountered in front of an obstacle that causes a partial reflection of a wave. The amplitude at the node is N=A+B and at the antinode AN =A−B. The direct measurement of N and AN yields:

and

and the reflection coefficient

Water wave theory

9

. In a translatory wave, water is transported in the direction of the wave travel. Some examples are 1) A tidal bore or moving hydraulic jump 2) Waves generated by the breaking of a dam 3) A surge on a dry bed 4) An undulated, moving hydraulic jump 5) Solitary waves 6) Flood waves in rivers In an Eulerian system of coordinates, a surface wave problem generally involves three unknowns: the free surface elevation (or total water depth), the pressure (generally known at the free surface), and the particle velocity. Since a general method of solution is impossible, a number of simplifying assumptions have been made that apply to a succession of particular cases with varying accuracy. In general, the method of solution used depends on the relative importance of the convective terms. However, instead of dealing with these terms directly, it is more convenient to relate them to more accessible parameters. Three characteristic parameters are used. 1) A typical value of the free surface elevation, such as the wave height H 2) A typical horizontal length such as the wave length L 3) The wave depth D Although the relationships between the convective terms and these three parameters are not simple, their relative values are of considerable help in classifying the water wave theories. As the free surface elevation decreases, the particle velocity also decreases. Thus, when the wave height, H, tends to zero, the convective term, which is related to the square of the particle velocity, is an infinitesimal. Consequently, the convective terms can be neglected and the theory can be linearized. Thus, the three possible parameters to be considered are

The relative importance of the convective terms increases as the values of these three parameters increase. In deep water (small H/D, and small L/D), the most significant parameter is H/L, which is called the wave steepness. In shallow water, the most significant parameter is H/D, which is called the relative height. In intermediate water depth, a significant parameter that also covers the three cases is , and it is called the Ursell parameter. Depending on the problem under consideration and the range of values of the parameters H/L, H/D, and L/D, four mathematical approaches are used. 1) Linearization 2) Power series

Numerical modeling of water waves

10

3) Numerical methods 4) Random functions The simplest cases of water wave theories are the linear wave theories, in which case the convective terms are neglected completely. These theories are valid when H/L, H/D, and L/D are small, i.e., for waves of small amplitude and small wave length in deep water. For the first reason, they are called the “small amplitude wave theory.” It is the infinitesimal wave approximation. The linearization of the basic equation is so suitable to mathematical manipulation that the linear wave theories cover an extreme variety of water wave motions. For example, some phenomena that can be subjected to linear mathematical treatment include the phenomena of wave diffraction and of the waves generated by a moving ship. The solution can be found as a power series in terms of a small quantity by comparison with the other dimensions. This small quantity is H/L for small L/D because in deep water, the most significant parameter is H/L. It is H/D for large L/D because in shallow water the most significant parameter is H/D. In the first case (development in terms of H/L), the first term of the power series is given by application of the linear theory. In the second case, the first term of the series is already a solution of nonlinear equations. The calculation of the successive terms of the series is so cumbersome that these methods are used in a very small number of cases. The most typical case is the progressive periodic wave. In this case, the solution is assumed to be a priori that of a steady-state profile, i.e., a function such as F=f(x−Ct), where C is a constant equal to the wave velocity or phase velocity. The simplification introduced by such an assumption is due to the fact that

and

such that

In such a way, the time derivatives can be eliminated easily and replaced by a space derivative. Typical examples of such treatments include 1) Power series of H/L or the Stokes waves, valid in deep water. The first term of the series is obtained from linear equations and corresponds to the infinitesimal wave

Water wave theory

11

approximation. 2) Power series of H/D: the cnoidal wave or the solitary wave, valid for shallow water. The first terms of the series are obtained as a steady-state solution of already nonlinear equations, but correspond to the shallow water approximation. However, a steady-state profile may not exist as a solution, in which case the method to be used is often a numerical method of calculation where the differentials are replaced by finite differences. This occurs for large values of H/D and L/D, which correspond to the fact that the nonlinear terms such as

are relatively large by comparison with the

linear terms such as . This is the case for long waves in very shallow water. Of course, a numerical method of calculation can be used for solving a linearized system of equations. For example, the relaxation method is used for studying small wave agitation in a basin. Also, an analytical solution of a nonlinear system of equations can be found in some particular cases. Hence, these three methods and the range of application that has been given indicate more of a trend than a general rule. The three previous methods aim at a fully deterministic solution of the water wave problem. Other descriptions of a sea state generally involve the use of random functions. The mathematical operations (such as harmonic analysis) generally imply that the water waves obey linear laws, which are the necessary requirements for assuming the principle of superposition. Such a method loses its validity for describing the sea state in very shallow water (large values of H/D and L/D). In hydrodynamics, the water wave theories are generally classified into two families. They are the “small amplitude wave theories” and the “long wave theories.” The small amplitude wave theories are the linearized theories of the first categories of power series, i.e., the power series in terms of H/L. The long wave theories use the numerical method of solution for the nonlinear long wave equations. These two families include a number of variations and some intermediate cases with some of the characteristics of both families. For example, the cnoidal wave and the solitary wave are considered as particular cases (steady state) of the long wave theories, because they are nonlinear shallow water waves.

1C. Airy Waves The linear theory of Airy in Eulerian coordinates gives the essential characteristics of the wave pattern assuming the wave height is infinitesimal. The free surface is sinusoidal as shown in Figure 1.3, particle paths are elliptic and follow a closed orbit (zero mass transport), and lines of equipressure are also sinusoidal. The terms in (H/L)2 are neglected. The long wave or shallow water theory is the same as the theory of Airy, where it is assumed that D/L is small. As a consequence, the formulas are simplified considerably. The pressure is hydrostatic and the horizontal velocity distribution is uniform. The Airy wave equations used in the WAVE code available on the NMWW CD-ROM are given below.

Numerical modeling of water waves

Fig. 1.3. The Airy Wave.

P=−ρgy+∆P.

where h is wave height from the bottom.

12

Water wave theory

13

Special cases of the Airy wave theory are the deep water (D/L> 0.5) and shallow water (D/L3 even sooner than in the case of the third-order theory, i.e., for larger values of D/L. Consequently, the fifth-order wave theory is less valid than the third-order wave theory for small values of D/L and cannot be used when D/L0.25. In no case should ALPHA exceed 1.0.

4B. Application to Tsunami Wave Formation The first application of the incompressible Navier-Stokes model to tsunami wave formation was by Garcia4 using a two-dimensional arbitrary boundary marker and cell technique to study tsunamis in the vicinity of their source. The method was applied to the Mendocino Escarpment for hypothetical ocean floor displacements of 10 meters in a few seconds to a minute, which might result from a major earthquake on the San Andreas fault. The wave was followed numerically for the first 200 sec. A single hump was first generated that split into two crests moving in opposite directions at a speed slightly less than shallow water wave speed. The waves were slightly dispersive and had a crest elevation above mean water level of about 1 meter and a period of about 1.5 minute. The period of major tsunamis measured in the Pacific Ocean is generally from 10 to 30 minutes. The periods are estimated from tide gauge records assuming that the period remains nearly constant during the shoaling of the wave from the deep ocean. The wave height of major tsunamis is generally given at 1±0.5 meter and is estimated from the tide gauge records assuming various approximate models for the wave height growth and extrapolating back to the deep ocean. The first application of the three-dimensional incompressible Navier-Stokes model to tsunami formation was by Mader, Tangora, and Nichols,5 and by Mader.6 The Hawaiian tsunami of November 29, 1975, has been investigated by Loomis. He described the observed runup heights in reference 7 and a numerical study of the tsunami source in references 8 and 9. The tsunami was generated by an earthquake with a magnitude of 7.2 on the Richter scale, near the Hawaii Volcanoes National Park. Near the source, the first wave was smaller than the second. Coincident with the earthquake was considerable subsidence (up to 3 meters) of the shoreline. Loomis, in reference 8, examined a model of the southeastern coast of Hawaii. The bottom slopes seaward at a ratio of 1:15 until it reaches a constant depth of 6000 meters. The sources examined by Loomis included both initial uplifts and depressions, and he reported that such source motions would not generate the essential features of the tsunami; that is, a second wave larger than the first. In addition to the sources studied by Loomis, a landslide source was modeled in

The three-dimensional navier-stokes model

123

reference 6. The landslide model has been evaluated by Cox in reference 10. He concluded that a landslide could not be distinguished from strictly tectonic displacement by the comparison of arrival times and travel times. The SWAN code described in Chapter 2 was used to solve the shallow water, long wave equations and to examine the tsunami generation problem. The SWAN results confirmed Loomis’s calculated results reported in reference 8. The SOLA-3D code, which solves the three-dimensional incompressible Navier-Stokes equations, was also used to model the tsunami. The Calculated Shallow Water Wave Results The model was identical to that described by Loomis in reference 7. A 40 by 69 rectangular region of 207 kilometers along the coast and 120 kilometers seaward is described using a mesh of 3 kilometers by 3 kilometers. The bottom slopes at a ratio of 1:15 until it reaches a depth of 6 kilometers. The source is a bottom slope subsidence programmed to vary with time. The source is 30 kilometers wide, of which half is included in the calculation and is separated from the other half by a reflective boundary, as shown in Figure 4.1. The calculated wave profile is shown at various locations along the shoreline as a function of time in Figure 4.2 for a source displacement of 3 meters, and in Figure 4.3 for a source displacement of 1 meter, followed by an additional 2 meters displacement 10 minutes later. Surface profiles are shown in Figure 4.4.

Fig. 4.1. Sketch of model used to numerically simulate the tsunami generation.

Numerical modeling of water waves

124

Fig. 4.2. Shoreline wave heights for shallow water wave model resulting from an initial source displacement of 3 meters.

The observed larger second wave can be reproduced by a source that has its displacement change with time. Such a possibility was suggested by Ando, who suggested that the earthquake was a rather

Fig. 4.3. Shoreline wave heights for shallow water wave model resulting from an initial source displacement of 1 meter, followed by an additional 2 meters displacement 10 minutes later.

The three-dimensional navier-stokes model

125

slow rupture lasting 100 sec; however, the displacement change required by the model is of longer duration and includes two fast ruptures. Another source investigated was an undersea landslide. The landslide ocean bottom profile assumed the bottom dropped 3 meters at the shoreline and slid to form the profile shown in Figure 4.5. Landslides are observed to pile up at the bottom one-third of their run. This gives the surface wave profile shown in Figure 4.6. The calculations were performed on the University of Hawaii Harris computer using the Hawaii version of the SWAN code. The shoreline wave heights at various times for the shallow water wave model with the initial water surface displacement of Figure 4.6 are shown in Figure 4.7. Although the wave heights are consistent with the observed behavior of the tsunami, we must check the results with the SOLA code since it was demonstrated in Chapter 3 that the shallow water model is inadequate to describe the waves generated from surface deformations of the water surface.

Numerical modeling of water waves

126

Fig. 4.4. Surface profiles for linear shallow water wave model as a function of time, resulting from an initial shore displacement of 1 meter, followed by an additional 2 meters displacement 10 minutes later.

The three-dimensional navier-stokes model

127

Fig. 4.5. Sketch of the final ocean bottom profile after a landslide for the source region.

Fig. 4.6. Sketch of height of water surface after a landslide on the ocean bottom.

The Calculated Navier-Stokes Results The geometry of the model used to calculate the tsunami is shown in Figure 4.1. The mesh used in the calculation had 20 cells in the x direction, 25 cells in the y direction, and 18 cells in the z direction. The 20 cells in the x direction were 6 kilometers wide. The 18 cells in the z direction starting at the ocean floor were 100 meters high for the first two cells, and 400 meters thereafter. The water depth was 6000 meters and the surface was located at the center of cell 17. The 25 cells in the y direction starting at the source were 3.0 kilometers for the first 5 cells that described the source (15 kilometers wide). The remaining cell widths were 5.75, 6.16, 6.56, 6.97, 7.37, 7.78, 8.18, 8.59, 9.0, 9.4, 9.8, 10.2, 10.6, 11.0, 11.4, 11.8, 12.2, 12.6, 13.0, and 13.5 kilometers, for a total of 206.8 kilometers.

Numerical modeling of water waves

128

Fig. 4.7. Shoreline wave heights for a shallow water wave model resulting from the initial water surface displacement shown in Figure 4.6.

The viscosity coefficient was 2.0 g-sec−1 meter−1 (0.02 poise). The gravity constant gz was −9.8 meters sec−2, and gx and gy were 0.0. The time step for the calculation was 5 sec. The tsunami source was modeled by a 3 meters deep depression or elevation of 90 by 15 kilometers of the water surface, as shown in Figure 4.1. The calculated wave profiles are shown at various locations along the shoreline as a function of time in Figures 4.3, for a source of 3 meters depression of the water surface, The calculation for a 3 meters uplift gave mirror images of the profiles for the 3 meters depression of the water surface. Surface profiles are shown in Figure 4.9. The calculated wave profiles are shown in Figure 4.10 at various locations along the shoreline as a function of time for a landslide source. The observed tsunami wave profile of the 1975 Hawaiian tsunami, near the source of the second wave larger than the first, is not reproduced by a landslide source in an incompressible three-dimensional NavierStokes calculation. This contrasts with results obtained using the shallow water model.

The three-dimensional navier-stokes model

129

Fig. 4.8. Shoreline wave heights for the three-dimensional incompressible Navier-Stokes model as a function of time resulting from an initial source of a 3 meters depression of the water surface.

Numerical modeling of water waves

130

Fig. 4.9. Surface profiles for three-dimensional incompressible NavierStokes model as a function of time resulting from an initial source of a 3 meters depression of the water surface.

The three-dimensional navier-stokes model

131

Fig. 4.10. Shoreline wave heights for a three-dimensional incompressible Navier-Stokes model calculation with the landslide source shown in Figure 4.6.

A source with a 3 meters uplift of the water surface was consistent with the observed tsunami wave profile. These calculations do not support a landslide source for the 1975 Hawaii tsunami. Conclusions The differences between the shallow water and full Navier-Stokes calculations are that the water waves formed in the full Navier-Stokes calculations are deep water waves that move slower than the shallow water waves formed in the shallow water calculations. The nature of the surface collapse is also different, with the collapse occurring throughout the source region in the Navier-Stokes calculations and mostly at the sides in the shallow water calculations. The observed tsunami wave profile of the 1975 Hawaii tsunami, of the second wave being larger than the first wave near the source, is not reproduced by a landslide source, but it is reproduced by a simple uplift or drop of the water surface over the source area. The shallow water approximation is not appropriate for studying waves generated from surface deformations that are small, relative to the water depth. In the next chapter we will investigate the limitations of the shallow water approximation for modeling tsunami waves.

Numerical modeling of water waves

132

4C. The 1994 Skagway Tsunami The Skagway tsunami of November 3, 1994 is described in Section 2G of Chapter 2. A tsunami wave with a period of about 3 minutes and a maximum height of 30 feet occurred in the Taiya Inlet at Skagway, Alaska. The tsunami was a result of a massive underwater landslide of sediment deposited by the Skagway River. The tsunami was modeled in Chapter 2 using the shallow water SWAN code. The Skagway tsunami was also modeled using the SOLA-3D code. A single slide model for the Skagway tsunami, that closely approximated the more realistic 3-slide Landslide model, was used for the three-dimensional SOLA-3D calculation. A SWAN calculation for the same single slide model was also performed for comparision. The SOLA-3D FORTRAN code and input file used to model the Skagway tsunami and the computer movies are on the NMWW CD-ROM in the directory /TSUNAML.MVE/SOLA/SKAGWAY. Since the wave length and period of the resulting tsunami is large compared to the shallow depth of the inlet, the tsunami wave is adequately approximated by a shallow water wave. The tsunami wave profiles generated by the SOLA-3D code and the SWAN code were similar in amplitude and period.

References 1. C.W.Hirt, B.D.Nichols, and N.C.Romero, “SOLA—A Numerical Solution Algorithm for Transient Fluid Flows,” Los Alamos National Laboratory report LA-5852 (1975). 2. C.W.Hirt, B.D.Nichols, and N.C.Romero, “SOLA—A Numerical Solution Algorithm for Transient Fluid Flows—Addendum,” Los Alamos National Laboratory report LA5852, Add. (1976). 3. C.W.Hirt, B.D.Nichols, and L.R.Stein, “SOLA-3D—A Numerical Solution Algorithm for Transient Three-Dimensional Flows,” Los Alamos National Laboratory, Group T3, unpublished internal report (1985). 4. W.J.Garcia, “A Study of Water Waves Generated by Tectonic Displacements,” College of Engineering, University of California at Berkeley report HEL 16–9 (1972). 5. Charles L.Mader, Robert E.Tangora, and B.D.Nichols, “A Model of the 1975 Hawaii Tsunami,” Natural Science of Hazards—The International Journal of the Tsunami Society, pp.C1-C8 (1982). 6. Charles L. Mader, “A Landslide Model for the 1975 Hawaii Tsunami,” Science of Tsunami Hazards, Vol.2, No.2, pp.71–77 (1984). 7. Harold G.Loomis, “The Tsunami of November 29, 1975, in Hawaii,” Hawaii Institute of Geophysics report HIG-75–21 (1975). 8. Harold G.Loomis, “On Defining the Source of the 1975 Tsunami in Hawaii,” JIMAR report to Nuclear Regulatory Commission (1978). 9. M.A.Sklarz, L.Q.Spielvogel, and H.G.Loomis, “Numerical Simulation of the November 29, 1975, Island of Hawaii Tsunami by the Finite Element Method,” Journal of Physical Oceanography, Vol.9, No.5, pp.1022–1031 (1979).

The three-dimensional navier-stokes model

133

10. Doak C.Cox, “Source of the Tsunami Associated with the Kalapana (Hawaii) Earthquake of November 1975,” Hawaii Institute of Geophysics report HIG-80–8 (1980).

5 EVALUATION OF INCOMPRESSIBLE MODELS FOR MODELING WATER WAVES The shallow water code SWAN and the incompressible Navier-Stokes code ZUNI have been used to model the development of a tsunami wave from an initial sea surface displacement, the propagation of a tsunami wave, and the resulting shoaling and flooding. The generation and propagation of a tsunami wave has also been modeled using a linear gravity model. The discrepancies between the shallow water and the more realistic Navier-Stokes and linear gravity models are quantified. These studies were performed by the author at the University of Hawaii Joint Institute for Marine and Atmospheric Research in the late 1980s and were published in references 1, 2 and 3.

5A. Tsunami Wave Generation Initial Surface Displacement Study—SWAN Model A 1 meter high Airy half-wave surface displacement with a width of 45 kilometers in 4550 meters deep water was studied. This approximates within cell resolution a Gaussian wave with a Gaussian break width of 10 kilometers. Calculations were performed using a 250 meters wide mesh of 400 by 4 cells and at 0.5 sec intervals. The wave height in meters as a function of distance is shown in Figure 5.1 at various times. The initial surface displacement separates into two shallow water waves with a height of 0.5 meter and length of 45 kilometers. At maximum height the vertical velocity at the center of the wave is 0.0232 meter/sec. An Airy wave with a 0.5 meter half width, 90 kilometers wave length in 4450 meters of water has a vertical velocity of 0.0235 to 0.0226 meter/sec at maximum height. The wave speed is 208.36 and the group velocity is 202.89 meters/sec. The shallow water approximation used in the SWAN code of constant vertical velocity introduces errors of about 5 percent.

Evaluation of incompressible models for modeling water waves

135

Fig. 5.1. The wave height in meters as a function of distance at various times for a 1 meter high surface displacement with a width of 45 kilometers in 4550 meters deep water for the nonlinear shallow water model using the SWAN code.

The nonlinear feature of the shallow water model included in the SWAN code results in a small trailing wave with an amplitude of less than 0.01 meter.

Numerical modeling of water waves

136

Initial Surface Displacement Study—ZUNI Model An approximately 1 meter high (within cell resolution) Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) in 4550 meter deep water was studied. The calculations for this geometry were performed with 15 cells in the Y or depth direction and 68 cells in the X or distance direction. The cells were rectangles 450 meters high in the Y direction and 4000 meters long in the X direction. The time increment was 3 sec. The water level was placed at 4550 meters or 50 meters up into the eleventh cell. The viscosity coefficient used was 0.02 poise, a value representative of the actual viscosity for water. The viscosity did not significantly affect the results. The wave height in meters as a function of distance is shown in Figure 5.2 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and length of 45 kilometers. At 100 sec two waves have formed. The vertical velocity at the top center of the wave is 0.0231 meter/sec and 0.2158 meter/sec near the bottom of the wave. An Airy wave with a 0.5 meter half width, 90 kilometers wave length in 4450 meters of water has a center vertical velocity of 0.0235 to 0.0226 meter/sec. The Airy wave speed is 208.36 and the group velocity is 202.89 meters/sec. As the wave propagates, the wave height decreases, the slope of the front of the wave becomes less, and small waves form behind the main wave. After the wave has propagated for 500 sec, the wave height has decreased to 0.445 meter, and the vertical velocity at the center of the wave has decreased to 0.0217 meter/sec near the peak of the wave to 0.0203 at the bottom of the wave. A train of waves has developed behind the main waves with maximum negative amplitude of 0.1 meter and positive amplitude of 0.05 meter.

Evaluation of incompressible models for modeling water waves

137

Fig. 5.2. The wave height in meters as a function of distance at various times for a 1 meter high surface displacement with a width of 45 kilometers in 4550 meters deep water for the Navier-Stokes water wave model using the ZUNI code.

Numerical modeling of water waves

138

Initial Surface Displacement Study—LGW Model The analytical methods for solving the linear gravity model were described by Professor George Carrier4,5. The two-dimensional linear gravity wave with a Gaussian, a square wave, or a time dependent Gaussian displacement is solved using Fourier transforms by the LGW code for any time of interest. The wave description is obtained for any uniform depth, density, gravity and Gaussian break width or square wave half width. The wave height, vertical and horizontial velocities, and pressure are calculated at any depths desired. The distance scale is chosen to be a tenth of a Gaussian break width so that a wave half-width is described by about 50 space increments. Since the model is symmetrical about the center of the initial displacement, only half of the wave profile is calculated. The LGW code is included on the NMWW CD-ROM. A 1 meter high Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) in 4550 meters deep water was studied. Since the model is symmetrical about the center of the initial displacement, only half of the wave is modeled analytically. The analytical linear gravity wave calculations were performed using 256 points with a thickness of 1000 meters. The wave profile was calculated for each time selected. The wave height in meters as a function of distance is shown in Figure 5.3 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and length of 45 kilometers. At 100 sec the two waves have formed and the vertical velocity at the center of the wave is 0.0238 meter/sec at the top of the wave to 0.217 meter/sec near the bottom of the wave. An Airy wave with a 0.5 meter half width, 90 kilometers wave length in 4450 meters of water has a vertical velocity of 0.0235 to 0.0226 meter/sec at the center of the wave. The Airy wave speed is 208.36 and the group velocity is 202.89 meters/sec.

Evaluation of incompressible models for modeling water waves

139

Fig. 5.3. The wave height in meters as a function of distance at various times for a 1 meter high Gaussian displacement with a Gaussian break length of 10 kilometers (equivalent to an Airy half wave length of 45 kilometers) in 4550 meters deep water for the linear gravity water wave model using the LGW code.

As the wave propagates, the wave height decreases, the slope of the front of the wave

Numerical modeling of water waves

140

becomes less, and small waves form behind the main wave. After the wave has propagated for 500 sec, the wave height has decreased to 0.4567 meter, and the vertical velocity at the center of the wave has decreased to 0.0216 meter/sec near the peak of the wave to 0.02015 at the bottom of the wave. A train of waves has developed behind the main waves with maximum negative amplitude of 0.1 meter and positive amplitude of 0.03 meter. Comparisions In Table 5.1 are listed for comparision the peak wave height and the associated vertical velocity at the top and bottom of the wave for various times for the nonlinear shallow water, the full Navier-Stokes and the linear gravity wave models. The initial surface displacement was chosen to be as similar as possible for each of the models. The initial surface displacement was chosen to be characteristic of a typical earthquake generated tsunami wave at its source. The surface displacement was one meter high spread over 45 kilometers in 4550 meters of water for a width-to-depth ratio of about 10. The initial surface displacement resulted in two waves- traveling in opposite directions with half the initial height. The nonlinear aspect of the shallow water wave appeared as very small trailing waves. The dispersion aspects of the linear gravity wave and the Navier-Stokes waves resulted in the front of the wave spreading, along with a corresponding decrease in the slope of the front of the wave and the generation of a following train of waves with amplitudes of a tenth or less of the peak wave height after the wave had traveled 105 kilometers in 500 sec (210 meters/sec). The peak wave height of the Navier-Stokes wave was lower than for the linear gravity wave. However, this difference was an order of magnitude less than the difference between either wave and the shallow water wave.

Evaluation of incompressible models for modeling water waves

141

TABLE 5.1 A 1 m High, 48 Kilometers Wide Source

Time sec

Height meters

Bottom Velocity meters/sec

Top Velocity meters/sec

Navier-Stokes Wave -ZUNI Code 0.

1.000

50.

0.0000

0.0000

0.0195

0.0202

100.

0.493

0.0216

0.0231

200.

0.475

0.0216

0.0233

300.

0.464

0.0211

0.0227

400.

0.459

0.0208

0.0224

500.

0.445

0.0203

0.0217

Linear Gravity Wave -LGW Code 0.

1.000

50.

0.0000

0.0000

0.0192

0.0220

100.

0.496

0.0217

0.0238

200.

0.489

0.0214

0.0234

300.

0.478

0.0211

0.0228

400.

0.467

0.0206

0.0222

500.

0.457

0.0201

0.0216

Shallow Water Wave -SWAN Code 0.

1.000

50.

0.0000

0.0000

0.0231

0.0231

100.

0.500

0.0232

0.0232

200.

0.499

0.0232

0.0232

300.

0.500

0.0232

0.0232

400.

0.499

0.0232

0.0232

500.

0.500

0.0232

0.0232

The differences described will become less with increasing wave length-to-depth ratios and become greater for decreasing wave length-to-depth ratios. To quantify this effect we used the linear gravity wave model to examine the propagation of waves of various widths-to-depths ratios.

Numerical modeling of water waves

142

Initial Surface Displacement Width Study-LGW Model The linear gravity wave model was used to investigate the waves formed from initial surface displacements of width-to-depth ratios of 40 to 0.5 (Gaussian break widths of 40 to 0.5 kilometers) in 4550 meters deep water. The wave height after the wave had travel ten times its initial width for depth ratios of 40, 20, 10, 5, 1, 0.5 are shown in Figure 5.4. The Airy half-wave length and half-wave period equivalents for a break width of 5 is 23 kilometers and about 1.8 minutes, a break width of 10 is 46 kilometers and 3.5 minutes, a break width of 20 is 90 kilometers and 7 minutes, and for a break width of 40 is 180 kilometers and 14 minutes. The use of nonlinear shallow water models to describe tsunami waves generated from earthquake generated surface displacements is adequate for tsunamis generated by surface displacements that are at least ten times wider than the depth. The nonlinear shallow water wave becomes less realistic the further it travels from its source and the smaller the width-to-depth ratio. Non-linear shallow water models are adequate for large wave length and long period tsunamis such as generated by the 1964 Alaskan or the 1960 Chile earthquakes where the periods were about 30 minutes and wave length-to-depth ratio in the deep ocean was greater than 80. Tsunami waves generated by earthquakes with small areas and periods of a few minutes will not be realistically described using shallow water models. Either the linear gravity wave model or better the Navier-Stokes model should be used for accurate modeling of tsunamis with small (less than 10) widthto-depth ratios.

Evaluation of incompressible models for modeling water waves

143

Fig. 5.4. The wave height in meters as a function of distance after the wave has traveled ten times its initial width for a 1 meter high Gaussian displacement with a Gaussian break length of 40, 20, 10, 5, 1, and 0.5 kilometers in 4450 meters deep water for the linear gravity wave model using the LGW code.

Generation Conclusions The SWAN code solves the nonlinear, shallow water, long wave equations including the effects of friction and flooding. Both the SWAN and ZUNI code which solves the incompressible Navier-Stokes equations were used to study the development of a tsunami wave from an initial sea surface displacement. The development of a tsunami wave was also modeled using two-dimensional linear gravity wave theory. If the initial displacement is approximately Gaussian and the wave length is very long compared to the depth, similar tsunami waves formed for all three methods. Two tsunami

Numerical modeling of water waves

144

waves traveling in opposite directions formed that were the sum of the original surface displacement. However, dispersion effects resulted in Navier-Stokes and linear gravity waves with decreasing front slopes and amplitudes, and followed by a train of small waves. The shallow water wave has a constant vertical velocity while the Navier-Stokes and linear gravity waves have the more realistic variable vertical velocity. For long wave length tsunamis the constant vertical velocity closely reproduces the velocity characteristics of Navier-Stokes and linear gravity waves, which slowly decrease with depth. With decreasing periods and wave lengths, the discrepancy between the shallow water and the Navier-Stokes and linear gravity waves formed from initial sea surface displacement increases. These studies permit us to determine the parameteric region where the shallow water model can be useful. A summary of that region is presented at the end of this chapter. Sources that involve compressible flow can not be solved with the incompressible models evaluated in this chapter; however they can be modeled using the compressible models described in Chapter 6. Computer movies of the source calculations are in the NMWW CD-ROM TSUNAMLMVE/SOURCE.MVE directory.

5B. Tsunami Wave Propagation The shallow water SWAN code and the incompressible Navier-Stokes ZUNI codes were used to study the propagation of a tsunami wave from an initial sea surface displacement such as the propagation of a tsunami wave from Alaska or the U.S. West Coast to Hawaii. The propagation of a tsunami wave was also modeled using two-dimensional linear wave theory. Tsunami Propagation Study—SWAN Model A 1 meter high Airy half-wave surface displacement with a width of 45 kilometers in 4550 meters deep water was studied. This approximates within cell resolution a Gaussian wave with a Gaussian break width of 10 kilometers. Calculations were performed using a 250 meters wide mesh of 32,000 by 4 cells and at 0.5 sec intervals. The wave height in meters as a function of distance is shown in Figure 5.5 at half hour intervals.

Evaluation of incompressible models for modeling water waves

145

Fig. 5.5. The wave height in meters as a function of distance in meters at half hour intervals for a 1 meter high surface displacement with a half width of 45 kilometers in 4550 meters deep water for the nonlinear shallow water model using the SWAN code.

The initial surface displacement separates into two shallow water waves with a height of 0.5 meter and width of 45 kilometers. The nonlinear feature of the shallow water model included in the SWAN code results in a small trailing wave with an amplitude of less than 0.01 meter. The peak wave amplitude remains nearly constant as the wave propagates. Tsunami Propagation Study—ZUNI Model A 1 meter high approximately (within cell resolution) Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) in 4550 meters deep water was studied as it traveled for 3 hours. The calculations for this geometry were performed with 15 cells in the Y or depth direction and up to 13,600 cells in the X or distance direction. The calculation required 231, 200 cells and over 2.5 million mesh quantities. The cells were rectangles 450 meters high in the Y direction and 4000 meters long in the X direction. The time increment was 3 sec. The water level was placed at 4550 meters or 50 meters up into the eleventh cell. The viscosity coefficient used was 0.02 poise, a value representative of the actual viscosity for water. The viscosity did not significantly affect the results. The wave height in meters as a function of distance is shown in Figure 5.6 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and width of 45 kilometers. At 100 sec the two waves have a

Numerical modeling of water waves

146

vertical velocity at the top center of the wave of 0.0231 meter/sec. Near the bottom of the wave the velocity is 0.02016 meter/sec. An Airy wave with a 0.5 meter half-width, 90 kilometers wave length in 4550 meters of water has a vertical velocity of 0.0235 to 0.0226 meter/sec at the center of the wave. The Airy wave speed is 208.36 meters/sec and the group velocity is 202.89 meters/sec.

Fig. 5.6. The wave height in meters as a function of distance in kilometers at various times for a 1 meter high surface displacement with a width of 45 kilometers in 4550 meters deep water for the NavierStokes water wave model using the ZUNI code.

As the wave propagates, the wave height decreases, the slope of the front of the wave becomes smaller, and a train of small waves forms behind the main wave. After the wave has propagated for 3 hours and over 2.3 megameters, the wave height has decreased to 0.20 meters. A train of waves has developed behind the main wave with maximum negative amplitude of 0.14 meter and positive amplitude of 0.11 meter. Tsunami Propagation Study—LGW Model A 1 meter high Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) propagating for 3 hours in 4550 meters deep water was studied. The wave height in meters as a function of distance is shown in Figure 5.7 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and length of 45 kilometers.

Evaluation of incompressible models for modeling water waves

147

Fig. 5.7. The wave height in meters as a function of distance in kilometers at various times with a width of 45 kilometers in 4550 meters deep water for the linear gravity wave model using the LGW code.

The analytical linear gravity wave calculations were performed using up to 16,384 points. The wave profile was calculated for each time selected. As the wave propagates, the wave height decreases, the slope of the front of the wave becomes smaller, and a train of small waves forms behind the main wave. After the wave has propagated for 3 hours and over 2.3 megameters, the wave height has decreased to 0.22 meter. A train of waves has developed behind the main wave with a maximum negative amplitude of 0.16 meter and positive amplitude of 0.12 meter. The peak height of the linear gravity wave decreases slower and the following wave train heights are higher than for the Navier-Stokes wave. The peak height of the Navier-Stokes wave is about 10 percent less than the linear gravity wave peak height after three hours of travel. This difference is an order of magnitude less than the difference between either wave and the shallow water wave height, which is more than twice as high. The period of the waves in the following wave train are similar; however the amplitudes are slightly higher for the linear gravity wave train. These results support the use of the linear gravity wave model for studying the effect of various wave length-to-depth ratios and of tsunami propagation for longer distances and times. The linear gravity wave model was used to study wave propagation from initial surface displacements of width-to-depth ratios of 40 to 5.0 (Gaussian break widths of 40 to 5

Numerical modeling of water waves

148

kilometers) in 4550 meters deep water. The wave height after the wave had traveled for 0.5, 2.0, 5.0 and 10.0 hours for depth ratios of 40, 20, 10, and 5 are shown in Figure 5.8.

Fig. 5.8. The wave height in meters as a function of distance in meters after the wave has traveled up to 10 hours and 7.6 megameters for a 1 meter high Gaussian displacement with a Gaussian break length of 40, 20, 10, and 5 kilometers in 4550 meters deep water using the LGW code.

The Airy half-wave length and half-wave period equivalents for a break width of 5 is 23 kilometers and about 1.8 minutes, a break width of 10 is 46 kilometers and 3.5 minutes, a break width of 20 is 90 kilometers and 7 minutes, and for a break width of 40 is 180 kilometers and 14 minutes. After traveling for 10 hours and 7.6 megameters, the peak wave amplitude for the 40 break width wave decreased from 0.5 to 0.42 meter, the 20 break width wave amplitude decreased to 0.28 meter, the 10 break width wave amplitude decreased to 0.16 meter, and the 5 break width wave amplitude decreased to 0.08 meter. The addition of nonlinear effects will lower these heights so these are upper bounds on the peak wave amplitude. Propagation Conclusions The nonlinear shallow water wave becomes less realistic the further it travels from its source and the smaller the width-to-depth ratio. Non-linear shallow water models are adequate for large wave length and long period tsunamis such as generated by the 1964 Alaskan or the 1960 Chile earthquakes where the periods were about 30 minutes and

Evaluation of incompressible models for modeling water waves

149

wave length to depth ratio in the deep ocean was greater than 80. Tsunami wave propagation from and generation by earthquakes with small areas and periods of a few minutes are not realistically described using shallow water models. Either the linear gravity wave model or better the Navier-Stokes model should be used for accurate modeling of long distance propagation of tsunamis with small (less than 40) width to depth ratios. Most tsunami waves that have been observed after traveling across the ocean have periods longer than 10 minutes. These studies support the postulate that this is because the shorter wave length tsunamis are so dispersive that as they propagate long distances, their amplitude decreases by an order of magnitude. Since three-dimensional Navier-Stokes solutions to tsunami propagation problems currently have limited practical application, most tsunami flooding studies need to be performed using the shallow water model. These studies permit us to determine the parameteric region where the shallow water model can be useful. A summary of that region is presented at the end of this chapter. The development and application of three-dimensional Navier-Stokes models including a realistic bottom friction treatment will be required for significant improvement in our ability to realistically model tsunami generation, propagation and flooding. The current effort to develop such a capability is described in Chapter 6. Computer movies of the propagation calculations are in the NMWW CD-ROM TSUNAMI.MVE/PRORMVE directory.

5C. Tsunami Wave Flooding The shallow water SWAN code and the Navier-Stokes ZUNI codes were used to study the effect of tsunami wave period, amplitude, bottom slope angle and friction on tsunami shoaling and flooding. The shallow water waves shoal higher, steeper and faster than the Navier-Stokes waves. The differences increase as the periods become shorter and the slopes less steep with large differences for periods less than 500 sec and slopes less than 2 percent.

Numerical modeling of water waves

150

Tsunami Flooding Study—SWAN Model A 3 meters half-height tsunami wave of various periods was propagated in 12 meters deep water and then interacted with slopes of various steepness. The wave height in meters as a function of distance is shown in Figure 5.9 at various times for a tsunami Airy wave with a 900 sec period and a 3 meters half-height traveling 3000 meters in 12 meters deep water before it interacted with a frictionless 1 percent slope. The space resolution in the numerical model grid was 10 meters with 300 cells to the bottom edge of the slope, 120 cells on the slope below and 80 cells above the water surface. The calculations were performed at 0.5 sec intervals. The peak wave height was 6.7 meters, the runup wave height was 6.0 meters and the inundation limit was 600 meters. The steepness of the slope was changed by changing the space resolution and the time step. The wave period was changed assuming that the tsunami was a shallow water Airy wave. The results are shown in Table 5.2 and Figures 5.10–5.12.

Evaluation of incompressible models for modeling water waves

151

Fig. 5.9. The interaction of a 900 second period, 3 meters half-height tsunami in 12 meters of water with a 1 percent slope at various times. The last figure shows the first wave maximum extent of flooding.

Numerical modeling of water waves

152

TABLE 5.2 SWAN Shallow Water Flooding Study

Slope percent

Period seconds

DeChezy coef

Peak Ht. meters

Runup meters

Ht .Flood meters

Initially 3 m Wave at 12 m Depth 4.0

900.

0.0

6.8

6.3

160

2.0

900.

0.0

7.0

6.4

330

1.0

900.

0.0

6.7

6.0

600

0.5

900.

0.0

5.0

4.7

1000

0.25

900.

0.0

3.7

3.2

1300

1.0

1570.

0.0

7.2

6.3

600

1.0

450.

0.0

5.2

4.8

500

1.0

225.

0.0

3.6

3.2

350

2.0

1800.

0.0

6.0

5.6

280

2.0

450.

0.0

6.8

6.0

300

2.0

225.

0.0

5.0

4.8

250

4.0

450.

0.0

7.0

6.4

160

4.0

225.

0.0

6.6

6.0

160

2.0

900.

10.0

4.8

0.5

30

2.0

900.

20.0

5.8

1.8

100

2.0

900.

30.0

6.4

3.0

160

2.0

900.

50.0

6.6

5.2

280

1.0

900.

10.0

4.2

0.4

40

1.0

900.

15.0

5.0

0.8

80

1.0

900.

25.0

6.0

1.6

160

1.0

900.

40.0

6.5

2.8

300

1.0

900.

50.0

6.5

3.8

390

1.0

900.

60.0

6.5

4.4

440

100.

3600.

0.0

4.8

4.8

100.

900.

0.0

6.1

6.1

100.

225.

0.0

6.5

6.5

Evaluation of incompressible models for modeling water waves

153

Initially 1 m Wave at 4550 m Depth 5.0

2666.

0.0

3.5

5.0

1333.

0.0

4.5

5.0

999.

0.0

5.5

5.0

666.

0.0

6.5

5.0

330.

0.0

7.7

The numerical model grid described is satisfactory for the cases studied giving runup and inundation values that are independent of space, time and geometry resolution. The grid model becomes less accurate for steeper slopes and longer wave length waves. The longer wave length waves need longer distances and times for the shoaling to occur; otherwise the calculated runup and inundation will be grid size dependent, decreasing with decreasing distance to the bottom edge of the slope. The effect of the steepness of a frictionless slope on the runup height and inundation limit is shown in Figure 5.10 for a tsunami wave with a 3 meters half-height and a 900 sec period in 12 meters of water. As shown in Figure 5.9, the highest shoaling wave height does not necessarily occur at the front of the wave. The difference between the height at the front of the wave and the peak wave height increases as the slope becomes steeper, with the maximum difference occuring for slopes larger than 1.5 percent. As shown in Figure 5.10, the runup height decreases and the inundation limit increases with decreasing slope angle. Also shown in Figure 5.10 is the predicted runup height for the Bretschneider model without friction. Bretschneider’s model results in smaller runup than SWAN for slopes greater than 1 percent. A popular engineering method for estimating maximum probable tsunami inundation zones is to use historical inputs to provide tsunami wave heights as a function of frequency of occurrence and then to use the Bretschneider6 model to calculate the runup on the shore at selected sites. The Bretschneider model as it is normally used includes surface roughness. The roughness parameters have been calibrated to reproduce observed runups. The Bretschneider model is independent of wave period and the only slope effect is from the friction.

Numerical modeling of water waves

154

Fig. 5.10. The runup wave height and the inundation limit in meters as a function of slope for a 3 meters half-height 900 sec period tsunami wave in 12 meters of water for the nonlinear shallow water model using the SWAN code.

The effect of wave period and frictionless slope steepness on tsunami flooding is shown in Figure 5.11. The maximum runup height is obtained for short period waves interacting with steep slopes, while the maximum inundation occurs for 1000 sec waves on gentler slopes.

Evaluation of incompressible models for modeling water waves

155

Fig. 5.11. The runup wave height and the inundation limit in meters as a function of wave period in sec for a 3 meters half-height 900 sec period tsunami wave in 12 meters of water for the nonlinear shallow water model using the SWAN code.

Numerical modeling of water waves

156

The effect of friction on tsunami flooding of 1 and 2 percent slopes is shown in Figure 5.12. The effect of increasing friction (decreasing DeChezy friction constant) is small on the peak height and becomes larger on the runup height as the slope decreases. In the Bretschneider model, the surface roughness is described using a Manning “n” where 0.04 corresponds to a roughness characteristic of grass or small rocks and 0.03 to a roughness of many trees, large boulders or high grass. The DeChezy friction coefficient is related to Manning “n“by the depth to the 1/6 power. While not directly comparable, for depths in this study a DeChezy friction coefficient of 50 results in about the same friction effect as a Manning “n” in the 0.03–0.04 range. Tsunami Flooding Study—ZUNI Model A 1 meter high tsunami wave of various periods was propagated in 4550 meters deep water and then it interacted with frictionless slopes of various steepness. These results are an extension of those described in Chapter 3. The geometry of the flooding calculation is shown in Figure 5.13 for a 1/15 or 6.66 percent slope. Calculations for this geometry were performed with 15 cells in the Y or depth direction and 68 cells in the X or distance direction. The cells were rectangles 450 meters high in the Y direction and 6750 meters long in the X direction. The time increment was 3 sec. The water level was placed at 4550 meters or 50 meters up into the

Fig. 5.12. The wave height and the inundation limit in meters for a 1 and 2 percent slope as a function of the DeChezy friction coefficient for a 3 meters half-height 900 sec period tsunami wave in 12 meters of water for the nonlinear shallow water model using the SWAN code.

Evaluation of incompressible models for modeling water waves

157

eleventh cell. The viscosity coefficient used was 0.02 poise, a value representative of the actual viscosity for water. The viscosity did not significantly affect the results. The slope in a ZUNI calculation is determined by the diagonal through the cell as discussed in the numerical methods section. A 6.66 percent slope results in a cell aspect ratio of 1 to 15 and a 4.0 percent slope an aspect ratio of 1 to 25. Numerical errors increase with increasing aspect ratio to unacceptable levels for 1 to 50 and greater ratios.

Fig. 5.13. The tsunami flooding geometry for the Navier-Stokes model using the ZUNI code.

Numerical modeling of water waves

158

TABLE 5.3 ZUNI Navier-Stokes Flooding Study

Slope percent

Period seconds

First Peak meters

First Trough meters

2nd Peak meters

2nd Trough meters

Period Study for 1 m Wave at 4550 m Depth 6.66

1333.

+3.2

−3.2

+4.1

−4.1

6.66

900.

+3.0

−4.0

+4.1

6.66

500.

+1.9

−2.9

+2.9

−2.9

6.66

250.

+0.9

−0.5

+0.25

−0.1

Period Study for 1 m Wave at 1011 m Depth 1.48

2500.

+2.8

−3.3

+3.8

1.48

1700.

+2.5

−3.7

+3.6

−3.2

1.48

1333.

+2.1

−3.6

+3.2

−3.3

1.48

900.

+1.5

−2.0

+2.0

1.48

500.

+0.8

−0.3

+0.3

10.0

1333.

+3.0

−3.2

+3.2

6.6

1333.

+3.2

−3.9

+4.1

−4.1

4.0

1333.

+2.9

−4.0

+3.9

−3.8

Slope Study, Slope Width Varied, 4550 m Depth

Slope Study, Depth Varied, Slope Width Constant 13.48

1333.

+3.2

−3.3

+3.5

−3.0

6.66

1333.

+3.2

−3.9

+4.1

−4.1

2.96

1333.

+2.9

−3.9

+3.85

2.22

1333.

+2.8

−3.8

+3.8

−3.7

1.48

1333.

+2.1

−3.6

+3.2

−3.3

0.74

1333.

+1.1

−1.9

+2.8

Evaluation of incompressible models for modeling water waves

159

TABLE 5.4 A 1333 Sec Tsunami Wave with Fixed Slope Width

Slope percent

Period seconds

Depth meters

Vel m/sec

Length km

Slope Width km

13.48

1333.

9100.

299.

398.

68.25

6.66

1333.

4550.

210.

280.

68.25

2.96

1333.

2022.

141.

188.

68.25

2.22

1333.

1517.

122.

163.

68.25

1.48

1333.

1011.

99.

132.

68.25

0.74

1333.

505.

70.

93.

68.25

Fig. 5.14. The runup height as a function of slope for a 1 meter high, 1333 sec period tsunami for the Navier-Stokes model using the ZUNI code.

Another method of changing the slope is to keep the width and aspect ratio constant and change the depth of the calculation. This permits slopes from 13.48 to 0.75 percent to be studied; however, if the period is kept constant, the other wave parameters change as shown in Table 5.4.

Numerical modeling of water waves

160

The slope study was performed using both methods to change the slope with the results shown in Table 5.3 and Figure 5.14. The general features of the calculated runup as a function of slope are the same for the ZUNI and SWAN models as shown in Figures 5.10 and 5.14. The resolution of the ZUNI calculation is inadequate to resolve the peak from the runup height. The ZUNI calculation permits us to realistically calculate wave interactions after the first wave interacts with the slope and to obtain rundown and runup for later waves. This is important because the largest runup both observed and in these calculations is the second or third wave. The interaction of the reflected waves with the later waves often results in the second or third runup being much different than the first runup. The magnitude and direction of the effect depends upon both the slope and the wave period as shown in Table 5.3. To perform these calculations, the source of the wave must be far removed from the slope in order to correctly follow the interaction of the incoming and reflected waves for several waves. A longer wave period needs a larger distance to run and requires more time for the multiple wave interactions. This results in large computing meshes and long running times for the long period tsunamis. Since they also result in the greatest flooding, they are also of the most interest. The effect of wave period was investigated for a 1 meter half-height wave at 4550 meters interacting with a 6.66 percent slope and at 1011 meters depth with a 1.48 percent slope. The results of the study of wave period and slope are shown in Table 5.3 and Figure 5.15. Also shown in Figure 5.15 are the SWAN runup heights as a function of period for the first wave. To obtain adequate resolution in the SWAN calculation, a much smaller mesh is required than in the ZUNI calculation.

Evaluation of incompressible models for modeling water waves

161

Fig. 5.15. The runup height as a function of period and different slopes and the first and second waves for a 1 meter high, 1333 sec period tsunami for the Navier-Stokes model using the ZUNI code. Also shown is the first wave using the SWAN shallow water code.

This is a result of the particle surface following treatment in the ZUNI calculation which resolves the surface in each cell. The SWAN calculations were performed using a

Numerical modeling of water waves

162

780 cell long mesh, a mesh width of 500 meters and a time step of 0.3 sec. The 4500 meters deep bottom was 250 kilometers long (500 cells) to the edge of the 5 percent shelf formed by 90 kilometers (180 cells) below and 50 kilometers (100 cells) above the water surface. Calculations were also performed using a mesh width of 250 meters and a 1560 cell long mesh. As shown in Table 5.3 and Figure 5.15, the shallow water waves flood higher for the first wave than the Navier-Stokes waves. The differences increase with shorter periods and less steep slopes, with large differences for periods less than 500 sec and slopes less than 2 percent. For higher period waves, the second wave is as much as a third larger than the first wave. The second wave is also larger than the shallow water first wave for wave periods greater than 500 sec. Conclusions The effect of tsunami wave period, amplitude, bottom slope angle and friction on tsunami shoaling and flooding has been investigated using a shallow water and a Navier-Stokes model. The study shows higher wave shoaling and flooding for waves interacting with steeper slopes, for waves with longer periods and for waves from deeper water. The shallow water waves shoal higher, steeper and faster than the Navier-Stokes waves. The differences increase as the periods become shorter and slopes less steep with large differences for periods less than 500 sec and slopes less than 2 percent. The interaction of the reflected first wave with the later waves often results in the second or third runup being much different than the first runup. The magnitude and direction of the effect depends upon both the slope and the wave period. For higher period waves, the second wave is as much as a third larger than the first wave. Since the shallow water model distorts the shoaled and reflected wave, the wave interaction and resulting second and third runups are inaccurately modeled. The effect of friction was investigated only for the shallow water model. It is the largest unknown factor in evaluating tsunami flooding. The Navier-Stokes model used does not have a treatment for friction. The development and evaluation of a realistic friction model is an important remaining tsunami flooding problem. A possible solution to the problem is described in Chapter 6. Since three-dimensional Navier-Stokes solutions to tsunami flooding problems currently have limited practical application, many tsunami flooding studies will need to be performed using the shallow water model. The following is a summary of the wave parameters where the shallow water model is valid. SOURCE—wave length must be 10 times longer than the depth. PROPAGATION—wave length must be 40 times longer than deep ocean depth and period needs to be longer than 15 minutes. RUNUP—period needs to be longer than 15 minutes. FLOODING—period needs to be longer than 10 minutes and the slope needs to be greater than 2 percent. Many important water wave problems can not be adequately modeled using the shallow

Evaluation of incompressible models for modeling water waves

163

water model. The incompressible Navier-Stokes model is adequate for many water wave problems that can not be accurately modeled using shallow water models. Currently only an upper limit flooding determination can be obtained without a realistic friction model. As shown in Chapter 6, it is now possible to model with high resolution the interface between the ocean and the sea floor, the mixing and movement of sediment and to include detailed sea floor topography. Future tsunami flooding studies can include friction effects realistically using the modeling methodology described in Chapter 6. One of the more important water wave problems is the determination of Civil Defense evacuation zones for tsunamis. As shown in this chapter the shallow water model is inadequate for the task. The incompressible Navier-Stokes model is unable to describe the initial generation of waves where the fluid-dynamics is compressible such as large earthquakes, impact landslides, volcanic, nuclear or conventional explosions and asteroids. Until recently there was no practical numerical method to solve such problems. Recent advances in the numerical modeling of water waves using the compressible Navier-Stokes model permit us to model almost any water wave generation, propagation and flooding problem. The next chapter describes the compressible model.

References 1. Charles L.Mader, “Numerical Tsunami Flooding Study—I,” Science of Tsunami Hazards, Vol.8, pp.79–96 (1990). 2. Charles L.Mader, Dennis W.Moore, and George F.Carrier, “Numerical Tsunami Source Study—II,” Science of Tsunami Hazards, Vol.11, pp.81–92 (1993). 3. Charles L.Mader, Dennis W.Moore, and George F.Carrier, “Numerical Tsunami Propagation Study—III,” Science of Tsunami Hazards, Vol.11, pp.93–106 (1993). 4. G.F.Carrier and H.P.Greenspan, “Water Waves of Finite Amplitude on a Sloping Beach,” Journal of Fluid Mechanics, Vol.4, pp.97–109 (1958). 5. G.F.Carrier, “The Dynamics of Tsunamis,” Mathematical Problems in the Geophysical Sciences, American Mathematical Society, Vol.1, pp.157–187 (1971). 6. C.L.Bretschneider and P.G.Wybro, “Tsunami Inundation Prediction,” Proceedings of 15th Coastal Engineering Conference (1976).

6 MODELING WAVES USING COMPRESSIBLE MODELS Compressible fluid dynamics is modeled numerically using the Lagrangian and the Eulerian equations of motion. One, two and three-dimensional models and codes have been developed over the last fifty years using various numerical methods for solving the equations of motion and to describe the compressible equation of state. In reference 1 is a complete description of Lagrangian and Eulerian approaches to modeling compressible reactive flows and the computer codes that were developed through the mid 1990’s. The Department of Energy’s program Accelerated Strategic Computing Initiative (ASCI) of the last five years has resulted in major advances in computing technology and in methods for improving the numerical resolution of compressible hydrodynamic calculations. Only recently has it become possible to calculate water wave problems for minutes or even hours using compressible hydrodynamic models that require millions of time steps for each second of flow. Only recently has it become possible to finely resolve a water-air interface and follow the water wave with millimeter resolution in a problem with 40 kilometers of air, 5 kilometers of water and tens of kilometers of ocean crust. The limitations of incompressible flow discussed in the previous chapters no longer must be accepted for modeling water waves. In this chapter the most advanced ASCI computer codes for modeling compressible fluid dynamics will be described. The codes will be used to model water waves generated from impact landslides, explosions, projectile impacts and asteroids. Numerical modeling of water waves has advanced so far so rapidly that it is clearly a technological revolution.

6A. The Three-Dimensional Compressible Model The three-dimensional partial differential equations for nonviscous, nonconducting, compressible fluid flow are The Nomenclature I

internal energy

P

pressure

Ux

velocity in x direction

Uy

velocity in y direction

Numerical modeling of water waves Uz

velocity in z direction

ρ

density

t

time

166

These equations are solved by a high-resolution Godunov differencing scheme using an adaptive grid technique described in reference 2. The solution technique uses Continuous Adaptive Mesh Refinement (CAMR). The decision to refine the grid is made cell-by-cell continuously throughout the calculation. The computing is concentrated on the regions of the problem which require high resolution. Refinement occurs when gradients in physical properties (density, pressure, temperature, material constitution) exceed defined limits, down to a specified minimum cell size for each material. The mesh refinement is shown in Figure 6.1 for a complicated geometry embedded in a material.

Modeling waves using compressible models

167

Fig. 6.1. Adaptive mesh refinement applied to a complicated geometry.

Much larger computational volumes, times and differences in scale can be simulated than possible using previous Eulerian techniques such as those described in reference 1. The original code was called SAGE. A later version with radiation is called RAGE. A recent version with the techniques for modeling reactive flow described in reference 1 is called NOBEL and was used for modeling the Lituya Bay impact landslide generated tsunami and water cavity generation. The codes can describe one-dimensional slab or spherical geometry, two-dimensional slab or cylindrical geometry, and three-dimensional Cartesian geometry. Because modern supercomputing is currently done on clusters of machines containing many identical processors, the parallel implementation of the code is very important. For portability and scalability, the codes use the Message Passing Interface (MPI). Load leveling is accomplished through the use of an adaptive cell pointer list, in which newly created daughter cells are placed immediately after the mother cells. Cells are redistributed among processors at every time step, while keeping mothers and daughters together. If there are a total of M cell and N processors, this technique gives nearly M/N cells per processor. As neighbor cell variables are needed, the MPI gather/scatter routines copy those neighbor variables into local scratch memory. The codes incorporate multiple material equations of state (analytical or SESAME tabular). Every cell can in principle contain a mixture of all the materials in a problem assuming that they are in pressure and temperature equilibrium. As described in reference

Numerical modeling of water waves

168

1, pressure and temperature equilibrium is appropriate only for materials mixed molecularly. The assumption of temperature equilibrium is inappropriate for mixed cells with interfaces between different materials. The errors increase with increasing density differences. While the mixture equations of state described in reference 1 would be more realistic, the problem is minimized by using fine numerical resolution at interfaces. The amount of mass in mixed cells is kept small resulting in small errors being introduced by the temperature equilibrium assumption, The strength is treated using an elastic-plastic model identical to that described in reference 1. Very important for water cavity generation and collapse and the resulting water wave history is the capability to initialize gravity properly, which is included in the code, This results in the initial density and initial pressure changing going from the atmosphere at 40 kilometers altitude down to the ocean surface. Likewise the water density and pressure changes correctly with ocean depth. Some of the remarkable advances in fluid physics using the SAGE code have been the modeling of Richtmyer-Meshkov and shock induced instabilities described in references 3 and 4.

6B. The Lituya Bay Mega-Tsunami Lituya Bay, Alaska is on the northeast shore of the Gulf of Alaska. It is an ice-scoured tidal inlet with a maximum depth of 220 meters and a narrow entrance with a depth of only 10 meters. It is a T-shaped bay, 7 miles long and up to 2 miles wide. The two arms at the head of the bay, Gilbert and Crillon Inlets, are part of a trench along the Fairweather fault. On July 87 1958, a 7.5 Magnitude earthquake occurred along the Fairweather fault with an epicenter near Lituya Bay. A mega-tsunami wave was generated that washed out trees to a maximum altitude of 520 meters at the entrance of Gilbert Inlet. Much of the rest of the shoreline of the bay was denuded by the tsunami from 30 to 200 meters altitude as shown in Figures 6.2– 6.4.

Modeling waves using compressible models

169

Fig. 6.2. The view of Lituya Bay in August 1958 after the earthquake and tsunami. The forest was destroyed to a maximun elevation of 524 meters and a maximum distance of 1000 meters from shoreline. The bay is 7 miles long and up to 2 miles wide. The two arms at the head of the bay are part of the Fairweather fault.

Numerical modeling of water waves

170

Fig. 6.3. Lituya Bay map showing topographic and bathymetric contours, trace of Fairweather fault, 1958 rockslide and trimline of wave runup. Forests destroyed to maximum elevations of 524 meters and 208 meters on north and south shores.

During the last 150 years five giant waves have occurred in Lituya Bay. The previous event occurred on October 27, 1936 which washed out trees to a maximum altitude of 150 meters and was not associated with an earthquake. Don Miller recorded all that was known in 1960 about the giant waves in Lituya Bay in reference 5. The July 9, 1958 earthquake occured at about 10:15 p.m., which is still daylight at Lituya Bay. The weather was clear and the tide was ebbing at about plus 5 feet. Bill and Vivian Swanson were on their boat anchored in Anchorage Cove near the western side of the entrance of Lituya Bay. Their astounding observations are recorded in reference 6 and were as follows: “With the first jolt, I tumbled out of the bunk and looked toward the head of the

Modeling waves using compressible models

171

bay where all the noise was coming from. The mountains were shaking something awful, with slides of rock and snow, but what I noticed mostly was the glacier, the north glacier, the one they call Lituya Glacier. I know you can’t ordinarily see that glacier from where I was anchored. People shake their heads when I tell them I saw it that night. I can’t help it if they don’t believe me. I know the glacier is hidden by the point when you’re in Anchorage Cove, but I know what I saw that night, too. The glacier had risen in the air and moved forward so it was in sight. It must have risen several hundred feet. I don’t mean it was just hanging in the air. It seems to be solid, but it was jumping and shaking like crazy. Big chunks of ice were falling off the face of it and down into the water. That was six miles away and they still looked like big chunks. They came off the glacier like a big load of rocks spilling out of a dump truck. That went on for a little while—it’s hard to tell just how long—and then suddenly the glacier dropped back out of sight and there was a big wall of water going over the point. The wave started for us right after that and I was too busy to tell what else was happening up there.” A 15 meters high wave rushed out of the head of the bay toward the Swansons’ anchored boat. The boat shot upward on the crest of the wave and over the tops of standing spruce trees on the entrance spit of Lituya Bay. Bill Swanson looked down on the trees growing on the spit and said he was more than 25 meters above their tops. The wave crest broke just outside the spit and the boat hit bottom and foundered some distance from the shore. Swanson saw water pouring over the spit, carrying logs and other debris. The Swansons escaped in their skiff to be picked up by another fishing boat 2 hours later. The front of Lituya Glacier on July 10 was a nearly straight, vertical wall almost normal to the trend of the valley. Comparisions with photographs of the glacier taken July 7 indicate that 400 meters of ice had been sheared off of the glacier front. The photograph taken after the event is shown in Figure 6.4.

Numerical modeling of water waves

172

Fig. 6.4. The view of Lituya Bay in August 1958 where the forest was destroyed to a maximun elevation of 524 meters.

After the earthquake there was a fresh scar on the northeast wall of Gilbert Inlet, marking the recent position of a large mass of rock that had plunged down the steep slope into the water. The next day after the earthquake and tsunami, loose rock debris on the fresh scar was still moving at some places, and small masses of rock still were falling from the rock cliffs near the head of the scar. The wave destroyed the forests and even removed the barnacles from the rocks. When the author visited Lituya Bay forty years after the event, almost no vegetation grew on either the slide region or on the cleared area shown in the above photograph. The dimensions of the slide on the slope are accurate but the thickness of the slide mass normal to the slope can only be estimated. The main mass of the slide was a prism of rock that was 730 meters and 900 meters along the slope with a maximum thickness of 90 meters and average thickness of 45 meters normal to the slope. The center of gravity was at about 600 meters altitude. As described in reference 5 this results in an approximate volume of 30 million cubic meters (40 million cubic yards). Miller5 concluded that “the rockslide was the major, if not the sole cause of the 1958 giant wave.” The Swanson observations have not been believed as they indicate that a lot more than a simple landslide occurred. Shallow Water Modeling In reference 7 shallow water modeling was performed using the SWAN nonlinear shallow water code. The generation and propagation of the tsunami wave of July 8, 1958 in

Modeling waves using compressible models

173

Lituya Bay was modeled using a 92.75 by 92.75 meters grid of the topography. The 3 by 6 second land topography was generated from the Rocky Mountain Communication’s CD-ROM compilation of the Defense Mapping Agency (DMA) 1 by 1 degree blocks of 3 arc second elevation data. The sea floor topography was taken from sea floor topographic maps published in reference 5. The grid was 150 by 150 cells and the time step was 0.15 sec. It was concluded that displaced water by a simple landslide or an earthquake along the Fairweather fault at the head of the bay could not result in the observed 550 meters high runup. The water in the inlet with the width of the landslide and between the landslide and the 520 meter high runup was sufficient to cover the runup region to a 100 meters height. In reference 7 it was shown that this high water layer was sufficient to form a wave that will reproduce the observed flooding of the bay beyond the inlet. It was concluded that a landslide impact model was required similar to that for asteroid generated waves. In 1999 it appeared that the numerical technology required would not become available for many decades. Physical Modeling During the Summer of 2000, Hermann Fritz8 conducted experiments that reproduced the 1958 Lituya Bay event. The 1:675 scale laboratory model of Lituya Bay shown in Figure 6.5 was built at VAW at the Swiss Federal Institute of Technology at Zurich, Switzerland. The laboratory experiments indicated that the 1958 Lituya Bay 524 meters runup on the spur ridge of Gilbert Inlet could have been caused by a landslide impact. The study was reported in reference 8. A novel pneumatic landslide generator was used to generate a high-speed granular slide with a controlled impact velocity and shape. A granular slide with the density and volume given by Miller5 was impacted with a mean velocity of 110 meters/sec. It generated a large air cavity and an extremely nonlinear wave with a maximum scaled height of about 160 meters which ran up to a scaled elevation of 530 meters above mean sea level.

Numerical modeling of water waves

174

Fig. 6.5. Experimental setup with pneumatic installation and measurement systems that included laser distance sensors (LDS), capacitance wave gauges (CWG) and particle image velocimetry (PIV).

The initial geometry used by Fritz8 in his experimental modeling in reference 8 is shown in Figure 6.6.

Fig. 6.6. Cross section of Gibert Inlet along slide axis used in the physical model.

The rockslide was simulated using 4 mm diameter Barium Sulfate grains with a density

Modeling waves using compressible models

175

of 1.6 g/cc and a void fraction of 39 percent. To what extent the rockslide broke up before impacting Gilbert Inlet is unknown. Scanned slide impact profiles were used to calculate the mean slide impact velocity of 110 meters/sec. The slide profile showed a gentle increase in slide thickness with time to scaled maximum of 134 meters and a fast decay back to zero. The scaled maximum slide thickness of 134 meters is 1.4 times the pre-motion slide thickness of 92 meters estimated roughly by Miller in reference 5. The three phases—granular material, water and air—are clearly separated along distinct border lines before flow reattachment occurs. Flow reattachment traps a large volume of air in the back of the rockslide, which leads to large bubble formation, bubble breakup and massive mixing. The slide is deformed due to impact and deflection at the channel bottom. A sequence of twelve instantaneous velocity vector fields computed the PIV data are shown in Figure 6.7. The sequence starts at 0.76 sec after rockslide impact and continues with a time step of 1.73 sec.

Numerical modeling of water waves

176

Fig. 6.7. Photosequence of the Fritz granular impact slide experiment. The time interval is 1.73 sec with the first image at 0.76 sec after impact.

Compressible Navier-Stokes Modeling The Lituya Bay impact landslide generated tsunami was modeled in reference 9 with the Navier-Stokes AMR (Adaptive Mesh Refinement) Eulerian compressible hydrodynamic

Modeling waves using compressible models

177

code NOBEL including the effects of gravity. The geometry of the calculation is shown in Figure 6.8 superimposed on the actual topography.

Fig. 6.8. Gilbert Inlet illustration showing rockslide dimensions, impact site and wave runup to 524 meters on the spur ridge directly opposite to the rockslide impact. Direction of view is North. The front of the Lituya Glacier is set to the 1958 post slide position. The preslide position extended 400 meters closer to the entrance of Gilbert Inlet.

The calculated density profiles are shown in Figures 6.9 and 6.10 for a rockslide moving with a resultant velocity of 110 meters/sec (X and Y component velocities of 77.8 meters/sec). The rockslide had an area of 21,000 square meters and was basalt with a density of 2.868 g/cc. The initial slide profile was triangular with its apex near the back of the slide to reproduce the Fritz slide profile shape described previously. The initial water depth was 120 meters and the length was 1.4 kilometers. The calculated maximum wave height in the bay was about 250 meters above sea level which ran up to 580 meters which is to be compared to the observed 524 meters. Such a wave could lift the glacier ice and result in the spectacular behavior observed by Swanson.

Numerical modeling of water waves

178

Fig. 6.9. Calculated density profiles of the Lituya impact landslide generated tsunami. The times are 0, 6, 10, 16, 22, 24, 26, 30, 36, 40, 48 and 58 sec.

A PowerPoint presentation with Fritz’s experimental movies and computer animations is available on the NMWW CD-ROM in the NOBEL/LITUYA directory.

Modeling waves using compressible models

179

Numerical modeling of water waves

180

Fig. 6.10. Calculated density profiles of the Lituya impact landslide generated tsunami. The slide material is shown as a dark region running across the ocean floor. At later times it comes to rest mostly in the middle of the bay. This illustrates how friction effects can be modeled realistically.

Conclusions The mega-tsunami that occurred on July 8, 1958 in Lituya Bay washed out trees to a maximum altitude of 520 meters at the entrance of Gilbert Inlet. Much of the rest of the shoreline of the bay was denuded by the tsunami from 30 to 200 meters altitude. The Lituya Bay impact generated tsunami was modeled in reference 9 with the AMR Eulerian compressible hydrodynamic code called NOBEL/SAGE. The capability now exists to evaluate the potential impact landslide tsunami hazards for vunerable regions of the world.

Modeling waves using compressible models

181

6C. Water Cavity Generation Introduction In the mid 1960’s, B.G.Craig10 at the Los Alamos National Laboratory performed experiments designed to characterize the formation of water waves from explosives detonated near the water surface. He reported observing the formation of ejecta water jets above and jets or “roots” below the expanding gas cavity. This was unexpected and a scientific mystery which has remained unsolved until it was finally modeled using the NOBEL code in December of 2002. In the early 1980’s, the hypervelocity impact (1.25 to 6 kilometers/sec) of projectiles into water was studied at the University of Arizona by Gault and Sonett11. They observed quite different behavior of the water cavity as it expanded when the atmospheric pressure was reduced from one to a tenth atmosphere. Above about a third of an atmosphere, a jet of water formed above the expanding cavity and a jet or “root” emerged below the bottom of the cavity. In the mid 1980’s, similar results were observed by Kedrinskii12 at the Institute of Hydrodynamics in Novosibirsk, Russia, who created cavities in water by sending large electrical currents through small lengths of small diameter Gold wires (bridge wires) causing the Gold to vaporize. The “exploding bridge wire” is a common method used to initiate propagating detonation in explosives. He observed water ejecta jets and roots forming for normal atmospheric pressure and not for reduced pressures. Thus the earlier Craig observations were not caused by some unique feature of generation of the cavity by an explosion. The process of cavity generation by projectiles or explosives in the ocean surface and the resulting complicated fluid flows has been an important unsolved problem for over 50 years. The prediction of water waves generated by large-yield explosions and asteroid impacts has been based on extrapolation of empirical correlations of small-yield experimental data or numerical modeling assuming incompressible flow, which does not reproduce the above experimental observations. The NOBEL code has been used to model the experimental geometries of Sonett and of Craig. The experimental observations were reproduced as the atmospheric pressure was varied as described in reference 13. Projectile and Exploding Bridge Wire Generated Cavities In the early 1980’s, experiments were being performed at the University of Arizona to simulate asteroid impacts in the ocean. The hypervelocity impact (1.25 to 6 kilometers/sec) of various solid spherical projectiles (Pyrex or Aluminum) with water was performed by Gault and Sonett11. Their observations were similar to those previously observed by Craig10. While the water cavity was expanding, an ejecta jet was formed at the axis above the water plume and a jet or “root” emerged along the axis below the cavity. The water cavity appeared to close and descend into deeper water. To improve the photographic resolution and reduce the light from the air shock, Sonett repeated his impact experiments under reduced atmospheric pressure. Much to

Numerical modeling of water waves

182

everyone’s surprise, when the pressure was reduced from one to a tenth atmosphere, the ejecta jet and the root did not occur and the water cavity expanded and collapsed upward toward the surface. This was what had been expected to occur in both the earlier Craig experiments and the projectile impacts. It became evident that the atmospheric pressure and the pressure differences inside and outside the water plume above the water surface was the cause of the formation of the jet, the root, the cavity closure and descent into deeper water. Figure 6.11 shows the Gault and Sonett11 results for a 0.25 cm diameter aluminum projectile moving at 1.8 kilometers/sec impacting water at one atmosphere (760 mm), and at 16 mm air pressure. A water stem and jet occurs at one atmosphere and not at low pressure. Figure 6.12 shows the Gault and Sonett results for a 0.635 cm diameter aluminum projectile moving at 2.5 kilometers/sec impacting water at one atmosphere (760 mm) and a 0.3175 cm diameter Pyrex projectile moving at 2.32 kilometers/sec impacting water at 16 mm air pressure. A water stem and jet occurs at one atmosphere and not at low pressure. Professor Kedrinskii12 at the Russian Institute for Hydrodynamics was also studying the generation of water cavities from exploding bridge wires. He was observing the formation of ejecta jets and roots as the water cavity expanded similar to those observed by Craig using explosives and by Gault and Sonett using projectiles. After we showed him the effect of reduced atmospheric pressure, he proceeded to repeat his exploding bridge wire experiments under reduced pressure. He observed that the jets and roots did not form when the atmospheric pressure was reduced to 0.2 atmosphere. Figure 6.13 shows the Kedrinskii results for an exploding bridge wire in water at one atmosphere and at 0.2 atmosphere air pressure. Compressible Navier-Stokes Modeling The projectile impact and explosive generated water cavity was modeled with the recently developed full Navier-Stokes AMR (Adaptive Mesh Refinement) Eulerian compressible hydrodynamic code NOBEL described earlier. The continuous adaptive mesh refinement permits the following of shocks and contact discontinuities with a very fine grid while using a coarse grid in smooth flow regions. It can resolve the water plume and the pressure gradients across the water plume and follow the generation of the water ejecta jet and root. Figure 6.14 shows the calculated density profiles for a 0.25 cm diameter aluminum projectile moving at 2.0 kilometers/sec impacting water at five atmosphere air pressure. The water plume collapses at the axis creating a jet moving upward and downward. The jet passes down through the cavity, penetrating the bottom of the cavity at the axis forming the stem. The flow results in the cavity descending down into the water. Figure 6.15 shows the calculated density profiles for a 0.25 cm diameter aluminum projectile moving at 2.0 kilometers/sec impacting water at one atmosphere air pressure. Figure 6.16 shows the calculated density profiles for a 0.25 cm diameter aluminum projectile moving at 2.0 kilometers/sec impacting water at 0.1 atmosphere air pressure. The tip of the water plume continues to expand in contrast to what is observed at

Modeling waves using compressible models

183

atmospheric pressures higher than 200 mm.

Fig. 6.11. The Gault and Sonett experimental results for a 0.635 cm diameter aluminum projectile moving at 1.8 kilometers/sec impacting water.

Numerical modeling of water waves

184

Fig. 6.12. The Gault and Sonett experimental results for a 0.25 cm diameter aluminum projectile moving at 2.5 kilometers/sec at 760 mm air pressure. A 0.3175 diameter Pyrex projectile moving at 2.32 kilometers/sec in 16 mm of air is shown in the bottom frame.

Modeling waves using compressible models

185

Fig. 6.13. The Kedrinskii results for an exploding bridge wire in water at 760 mm and 150 mm air pressure.

Numerical modeling of water waves

186

Fig. 6.14. The density profiles for a 0.25 cm diameter Aluminum projectile moving at 2.0 kilometers/sec impacting water at 5 atmosphere air pressure. The times are 0, 5, 10, 15, 30, and 70 milliseconds. The graphs are 100 cm wide by 120 cm tall, with 50 cm of water.

Modeling waves using compressible models

187

Fig. 6.15. The density profiles for a 0.25 cm diameter Aluminum projectile moving at 2.0 kilometers/sec impacting water at 1 atmosphere air pressure. The times are 0, 25, 75, and 125 milliseconds. The graphs are 100 cm wide and 120 cm tall, with 50 cm of water and 70 cm air.

Numerical modeling of water waves

188

Fig. 6.16. The density profiles for a 0.25 cm diameter Aluminum projectile moving at 2.0 kilometers/sec impacting water at 76 mm (0.1 atmosphere) air pressure. The times are 0, 5, 20, and 50 milliseconds. The graphs are 80 cm wide by 100 cm tall with 40 cm of water and 60 cm air.

Explosive Generated Cavities In reference 14, detailed one-dimensional compressible hydrodynamic modeling is described for explosives detonated in deep water. Agreement was obtained with the experimentally observed explosive cavity maximum radius and the period of the

Modeling waves using compressible models

189

oscillation. It was concluded that the detonation product equation of state over the required range of 1 megabar to 0.01 atmosphere was adequate for accurately determining the energy partition between detonation products and the water. It was also concluded that the equations of state for water and detonation products were sufficiently accurate that they could be realiably used in multidimensional studies of water cavity formation and the resulting water wave generation in the region of the “upper critical depth”. The “upper critical depth” is the experimentally observed location of a explosive charge relative to the initial water surface that results in the maximum water wave height. It occurs when the explosive charge is approximately two-thirds submerged. The observed wave height at the upper critical depth is twice that observed for completely submerged explosive charges at the “lower critical depth.” If the waves formed are shallow water waves capable of forming tsunamis, then the upper critical depth phenomenon would be important in evaluating the magnitude of a tsunami event from other than tectonic events.

Fig. 6.17. The scaled wave height as a function of scaled explosive charge depth. The “upper critical depth” is the explosive charge depth when the maximum wave height occurs which is approximately two-thirds submerged. The second smaller increase in wave height is the “lower critical depth.”

The water wave amplitude as a function of the depth the explosive is immersed in water is shown in Figure 6.17. The scaled amplitude is AR/W and the scaled depth is D/W where A is maximum wave amplitude at a distance R from the explosive charge of weight W. The “upper critical depth” is at the first wave height maximum which occurs when an explosive charge is located at the water surface. The second smaller increase in wave

Numerical modeling of water waves

190

height is at the “lower critical depth” which is about half the upper critical height but results in longer wave length water waves. Data are included for explosives with weights of 0.017 to 175 kilograms. The Craig10 experimental results are shown with a large x. During the study of the upper critical depth phenomenon in the 1960’s evidence of complicated and unexpected fluid flows during water cavity formation was generated by B.G.Craig and described in references 1, 10 and 14. A sphere of explosive consisting of a 0.635 cm radius XTX 8003 (80/20 PETN/Silicon Binder at 1.5 g/cc) explosive and a 0.635 cm radius PBX-9404 (94/6 HMX/binder at 1.84 g/cc) explosive was detonated at its center. The sphere was submerged at various depths in water. PHERMEX15 radiographs and photographs were taken with framing and movie cameras. The radiographs are shown in Figure 6.18. The cavity, water ejecta and water surface profiles shown in the PHERMEX radiography in Figure 6.18 were closely approximated by the compressible hydrodynamic modeling described in reference 14 using the 2DE code and in Figure 6.19 using the NOBEL code.

Modeling waves using compressible models

191

Fig. 6.18. Dynamic radiographs of a 2.54 cm diameter PBX-9404 explosive sphere detonated at its center and half submerged in water at 1 atmosphere air pressure. The times are 15.8, 26.3 and 61.3 microseconds. The sketch shows the prominent features of the radiographs with the water shock dashed.

Numerical modeling of water waves

192

Fig. 6.19. NOBEL density profiles for a 2.54 cm diameter PBX-9404 explosive sphere detonated at its center and half submerged in water at 1 atmospheric air pressure. The times are 15.0, 26.0 and 60.0 microseconds for comparision with the radiographs in Figure 6.18. The width is 16 cm and the height is 24 cm, of which 16 cm is water.

At later times, while the water cavity was expanding the upper ejecta plume collapsed and converged on the axis generating an upward water ejecta jet on the axis above the

Modeling waves using compressible models

193

water plume and a downward water jet which generated a root on the axis below the bottom of the cavity. These results were not anticipated and neither was the observation that the water cavity proceeded to close at its top and descend down into deeper water. At first it was assumed that there was something unique about the explosive source that was resulting in these remarkable observations. The reactive compressible hydrodynamic numerical models available were unable to reproduce the experimental observations or suggest any possible physical mechanisms unique to explosives. As described previously, the different behavior of the water cavity as it expanded when the atmospheric pressure was reduced from one atmosphere to less than a third of an atmosphere is independent of the method used to generate the cavity, such as an exploding bridge wire or a hypervelocity projectile impact. These remarkable experimental observations resisted all modeling attempts for over 25 years. The numerical simulations could not describe the thin water ejecta plumes formed above the cavity or the interaction with the atmosphere on the outside of the ejecta plume and the pressure inside the expanding cavity and plume. Figure 6.20 shows the calculated water profiles for a 0.25 cm diameter PBX-9404 explosive sphere detonated at its center half submerged in water at one atmosphere air pressure. All the projectile, exploding bridge wire and explosive experimental observations were reproduced as the atmospheric pressure was varied. At sufficiently high atmospheric pressure, the difference between the pressure outside the ejecta plume and the decreasing pressure inside the water plume and cavity as it expanded, resulted in the ejecta plume converging and colliding at the axis. A jet of water formed and proceeded above and back into the bubble cavity along the axis. The jet proceedes back through the bubble cavity penetrating the bottom of the cavity and forms the root observed experimentally. The complicated cavity collapse and resulting descent into deeper water was also numerically modeled.

Numerical modeling of water waves

194

Fig. 6.20. The density profiles for a 2.54 cm diameter PBX-9404 explosive sphere detonated at its center and half submerged in water at 1 atmosphere air pressure. The times are 0, 25, 75, and 275 milliseconds. The graphs are 100 cm wide and 120 cm tall, with 50 cm of water and 70 cm of air.

Explosive Generated Water Wave Craig10 measured the wave amplitude as a function of time for the first few seconds at a

Modeling waves using compressible models

195

distance of 4 rneters from a 2.54 cm diameter PBX-9404 explosive sphere initiated at its center in 3 meters of water. He included mass markers in the water. Mass markers located 1 meter below the water surface and markers located 0.5 meter below the surface and 1 meter from the explosive showed no appreciable movement compared with those located nearer the surface or explosive charge. These results showed that the wave formed was not a shallow water wave. The experimental and calculated wave parameters are summarized in Table 6.1. The parameters are given after 4 meters of travel from the center of the explosive charge. The wave parameters for the Airy wave were calculated using the WAVE code described in Chapter 1 for a depth of 3 meters and the experimentally observed wave length. Since the group velocity is almost exactly half the wave velocity, the Airy wave is a deep water wave. The shallow water results are from reference 14. A small wave from the initial cavity formation is followed by a larger negative and then a positive wave resulting from the cavity collapse. Only the second wave parameters are given in the table.

TABLE 6.1 Calculated and Experimental Wave Parameters

Experimental Wave Velocity (m/sec) Amplitude (cm) Wave Length (m) Period (sec)

Airy Wave

Shallow Water

NOBEL

2.50±0. 2

2.41

5.42

2.50±0.10

0.8

1.0

10.1

0.6

3.75

3.75(input)

1.0

3.75

1.5

1.55

0.18

1.5±0.1

Group Velocity (m/sec)

1.21

The wave gauge was close to the edge of the water tank, which resulted in reflected waves which perturbed the subsequent wave measurements. Since the wave gauge was located close (0.69 meter) to the side of the tank, the reflections from the first small wave perturbed the second wave, which probably explains the larger than calculated amplitude, as the calculated wave was unperturbed by any boundary. Conclusions In the late 1960’s and early 1970’s, B.G.Craig at the Los Alamos National Laboratory reported observing the formation of ejecta jets and roots from cavities generated by small spherical explosives detonated near the water surface while the gas cavity was expanding. The hypervelocity impact (1.25 to 6 kilometers/sec) of projectiles into water was studied at the University of Arizona in the early 1980’s by Gault and Sonett. They observed quite different behavior of the water cavity as it expanded when the atmospheric pressure was reduced from one to a tenth atmosphere. Above about a third of

Numerical modeling of water waves

196

an atmosphere, a jet of water formed above the expanding cavity and a root developed below the bottom of the bubble cavity. They did not occur for atmospheric pressures below a third of an atmosphere. Similar results were observed in the middle 1980’s by Kedrinskii at the Institute of Hydrodynamics in Novosibirsk, Russia when the water cavity was generated by exploding bridge wires, with jets and roots forming for normal atmospheric pressure and not for reduced pressures. During the last decade, a compressible Eulerian hydrodynamic code called SAGE has been under development by the Los Alamos National Laboratory and Science Applications International (SAIC) which has continuous adaptive mesh refinement (AMR) for following shocks and contact discontinuities with a very fine grid. A version of the SAGE code that models explosives, called NOBEL, has been used to model the experimental geometries of Sonett and of Craig. The experimental observations were reproduced as the atmospheric pressure was varied. When the atmospheric pressure was increased, the difference between the pressure outside the ejecta plume above the water cavity and the decreasing pressure inside the water plume and cavity as it expanded resulted in the ejecta plume converging and colliding at the axis forming a jet of water proceeding above and back into the bubble cavity along the axis. The jet proceeding back through the bubble cavity penetrated the bottom of the cavity and formed the root observed experimentally. The complicated cavity collapse was numerically modeled. Now that a code is available that can describe the experimentally observed features of projectile interaction with the ocean, we have a tool that can be used to evaluate impact landslide, projectile or asteroid interactions with the ocean, and the resulting generation of tsunami waves. A PowerPoint presentation with Craig’s and Sonett’s experimental movies and computer animations is available on the NMWW CD-ROM in the NOBEL/CAVITY directory.

6D. Asteroid Generated Tsunamis Introduction Two- and three-dimensional compressible hydrodynamic modeling using the SAGE code has been performed to study the generation of tsunamis by asteroid impacts with the ocean16,17. The goal was to determine the characteristics of impact generated tsunami waves as a function of the size and energy of the asteroid. The evaluation in the popular press of the potential threats from modest sized earth crossing asteroids often overestimates the tsunami threat from small (100 to 250 meters in diameter) asteroids. The studies reported are based on faulty shallow water modeling, from which they conclude that the resulting short (less than 5 minutes) period tsunamis are a major threat throughout an ocean basin. The news media are often delighted to promote such imaginary ocean basin wide “mega-tsunami” threats from small asteroid impacts, underwater landslides and volcanic

Modeling waves using compressible models

197

or other explosions. Part of the problem is that until recently there has not been a method for realistically modeling such threats. That changed with the development of the NOBEL/SAGE/RAGE compressible AMR code described in this chapter. While tsunamis up to a kilometer in initial height are generated by impactors of a kilometer diameter, they do not propagate as long period tsunamis such as generated by large earthquakes. Asteroid impact generated waves are highly complex in form and interact strongly with shocks propagating through the water and ocean crust. They decay more rapidly than the inverse of the distance from the impact point. A panel of tsunami scientists of The Tsunami Society have evaluated the proposed tsunami hazards from events that generate initially high amplitude short period tsunamis. Their evaluation is available on The Tsunami Society web site at http://www.sthjournal.org/media.htm. The Tsunami Society study concludes the following: No mega-tsunami has occurred in either the Atlantic or Pacific oceans in recorded history. The colossal collapse of Krakatau or Santorin generated catastropic waves in the immediate areas but hazardous waves did not propagate to distance shores. Carefully performed numerical and experimental model experiments on such events verify that the relatively short period waves from these small, though intense, occurances do not travel as do tsunami waves from a major earthquake. The 1958 Lituya Bay mega-tsunami described earlier in this chapter was the largest wave in recorded history, but it barely registered at nearby tide gauges. It had a short period and wavelength, which resulted in the wave rapidly decaying as It traveled out into the deep ocean. Two-Dimensional Vertical Impact Asteroid Generated Tsunamis A series of calculations in two-dimensions (cylindrical symmetry) for asteroids impacting the ocean vertically at 20 kilometers/sec was performed using the ASCI BlueMountain machine at the Los Alamos National Laboratory. These simulations were designed to follow the passage of an asteroid through the atmosphere, its impact with the ocean, the cavity generation, and subsequent re-collapse and generation of the tsunami. The results of these studies are described in reference 16, The parameter study included different asteroid masses. Stony and iron bodies of diameters 250 meters, 500 meters and 1000 meters were used. The kinetic energies of the impacts ranged from 1.3 Gigatons to 195 Gigatons of equivalent TNT yield. The projectile moved through an exponential atmosphere into a 5 kilometers deep ocean. The impactors were composed of mockup mantle material (dunite with a density of 3.32 g/cc) or iron (density of 7.8 g/cc). as a mockup for nickel-iron asteroids, For these projectiles, analytical Mie-Gruneisen equations of state were used. An elastic-plastic strength model was used for the crust and asteroid materials. Water was described using the SAIC/Pactech tables. An example montage from the two-dimensional parameter study is shown in Figure 6.21.

Numerical modeling of water waves

198

Fig. 6.21. A 1 kilometer iron vertical impactor craters the basalt crust, excavates a cavity in the ocean with a diameter of 25 kilometers.

Modeling waves using compressible models

199

The density plot shortly after the collapse of the transient water crater is shown in Figure 6.22.

Fig. 6.22. The wave train density plots for a 1 kilometer iron asteroid impacting 5 kilometers of water and a basalt crust. The complex wave train is a result of the interaction of multiple shocks propagating through the water and the basalt crust.

The complex wave train is a result of the reflections and interaction of multiple shocks propagating through the water and basalt crust. A complicated time-dependent flow occurs which is best studied using the PowerPoint presentation and computer movies on

Numerical modeling of water waves

200

the NMWW CD-ROM in the directory /NOBEL/ASTWAVE/. It becomes obvious that our previous attempts to model the tsunami wave generation by asteroids using incompressible models were inadequate for the task. A tabular summary of the parameter study is presented in Table 6.2. Listed are the impact characteristics of the asteroid (composition, diameter, density, mass, velocity and kinetic energy) and the measured characteristics of the impact (maximum depth and diameter of the cavity, quantity of water displaced, time of maximum cavity, maximum jet and jet rebound, tsunami wave length and velocity).

TABLE 6.2 Asteroid Generated Tsunamis

Asteroid Material

Dunite

Iron

Dunite

Iron

Dunite

Iron

Asteroid Diameter

250 m

250 m

500 m

500 m

1000 m

1000 m

Asteroid Density

3.32 g/cc

7.81 g/cc

3.32 g/cc

7.81 g/cc

3.32 g/cc

7.81 g/cc

Asteroid Mass

2.7e13 g

6.4e13 g

2.2e14 g

5.1e14 g

1.7e15 g

4.1e15 g

Asteroid Velocity

20 km/s

20 km/s

20 km/s

20 km/s

20 km/s

20 km/s

Kinetic Energy

1.3 GT

3.1 GT

10.4 GT

24.4 GT

83 GT

195 GT

Max Cavity Dia

4.4 km

5.2 km

10.0 km

12.6 km

18.6 km

25.2 km

Max Cavity Depth

2.9km

4.3 km

4.5 km

5.7 km

6.6 km

9.7 km

Water Displacement

4.4e16 g

9.1e16 g

3.5e17 g

7.1e17 g

1.8e18 g

4.8e18 g

Time of Max Cavity

13.5 s

16.0 s

22.5 s

28.0 s

28.5 s

33.0 s

Time of Max Jet

54.5 s

65.0 s

96.5 s

111 s

128 s

142 S

Time of Rebound

100.5 S

118.5 S

137.5 S

162 s

187 s

218 s

9 km

12 km

17 km

20 km

23 km

27 km

120 m/s

140 m/s

150 m/s

160 m/s

170 m/s

175 m/s

Tsunami Wavelength Tsunami Velocity

The shallow water tsunami velocity in 5 kilometers deep water is 221 meters/sec. The amount of water displaced during the formation of the cavity was found to scale linearly with the kinetic energy of the asteroid, as illustrated in Figure 6.23. A fraction of this displaced mass is actually vaporized during the explosive phase of the encounter, while the rest is pushed aside by the pressure of the vapor to form the crown and rim of the transient cavity.

Modeling waves using compressible models

201

Fig. 6.23. The mass of water displaced in the initial cavity formation scales with the asteroid kinetic energy. The squares are the results from the parameter study tabulated in Table 6.2. The solid line shows direct proportionality. A fraction (≈5–10 percent) of this mass is vaporized at impact.

The tsunami amplitude evolves in a complex manner, eventually decaying faster than the expected recipocal of the distance of propagation from the impact point as shown in Figure 6.24. The wave trains are highly complex because of the multiple shock reflections and interactions involving the seafloor. The complexity of the wave breaking associated with the reflection and interaction of shocks results in a very dynamic maximum wave height. Realistic seafloor topography would undoubtedly influence the development of the water wave.

Numerical modeling of water waves

202

Fig. 6.24. The tsunami amplitude has a very complex early time profile. It eventually decreases with distance faster than 1/(Radius). The legend identifies the points associated with individual cases. The notation signifies the asteroid composition (“Dn” for dunite and “Fe” for iron) and the diameter in meters. A “tr” indicates calculated tracer heights and a “Isq” indicates a least square fitted line to the tracer heights.

For an ocean of 5 kilometers depth, the shallow water velocity is 221 meters/sec. In Figure 6.25 are shown the wave crest positions as a function of time for the simulations in the parameter study, along with constant velocity lines of 150 and 221 meters/sec. The wave velocities are substantially lower than the shallow water limit. The tsunami wave length is found to scale roughly with the 1/4 power of the asteroid kinetic energy. The reason for this is that the wavelength is determined by the cavity-jetrebound cycle. The time scale for this varies as where jh is the mean jet height and g is gravity. The mean jet height varies as the square root of the asteroid kinetic energy.

Modeling waves using compressible models

203

Fig. 6.25. The tsunami wave crest positions as a function of time for the parameter study cases. The notation is the same as for Figure 6.24. The solid lines of constant velocity illustrate that these waves are substantially slower than the shallow water prediction.

Three-Dimensional Oblique Impact Asteroid Generated Tsunamis Three-dimensional simulations of a 1 kilometer diameter iron asteroid impacting the ocean at a 45 degree angle at 20 kilometers/sec were performed using the ASCI White machine. The calculations used up to 1200 processors and several weeks of computing time. Up to 200,000,000 computational cells were used, with a total computation time of 1,300,000 cpu-hours. The computational volume was a rectangular box 200 kilometers long in the direction of the asteroid trajectory, 100 kilometers wide and 60 kilometers tall. The height was divided into 42 kilometers of atmosphere, 5 kilometers of ocean water, 7 kilometers basalt crust and 6 kilometers of mantle material. The asteroid was started at a point 30 kilometers above the surface of the water as shown in Figure 6.26. The atmosphere was a standard exponential atmosphere, so the air initially surrounding the asteroid was only 1.5 percent of the sea level density.

Numerical modeling of water waves

204

Fig. 6.26. A 1 kilometer diameter iron impactor at an angle of 45 degrees craters a 5 kilometers deep ocean and the basalt crust. Twodimension slices in the vertical plane containing the asteroid trajectory are shown. The initial asymmetry disappears with time.

During the 2.1 sec of the asteroid’s atmospheric passage at about Mach 60, a strong bow air shock develops, heating the air to temperatures of 1 electron volt (11,600 K). Less than 1 percent of the asteroid’s kinetic energy is dissipated in the atmospheric passage. The water is much more effective at slowing the asteroid, and essentially all of its kinetic energy is absorbed by the ocean and seafloor within 0.7 sec. The water immediately surrounding the trajectory is vaporized. The rapid expansion of the vapor cloud evacuates a cavity in the water that eventually expands to a diameter of 15 kilometers. This initial cavity is asymmetric because of the inclined trajectory of the asteroid. The splash, or crown, is markedly higher on the side opposite the incoming trajectory, as shown in Figure 6.27.

Modeling waves using compressible models

205

Fig. 6.27. A 1 kilometer diameter iron asteroid impacted the ocean from the right at 45 degrees. A perspective cutaway view from the same run as in Figure 6.26 at the time of maximum cavity (27.5 sec after impact) is shown. The density profiles are shown. The cut passes through the ocean to the basalt ocean floor. The diameter of the cavity is about 25 kilometers.

The maximum height of the crown on the downstream side is nearly 30 kilometers at 70 sec after impact. The collapse of the bulk of the crown makes a “rim wave” or precursor tsunami that propagates outward, somewhat higher on the downstream side. The higher portion of the crown breaks up into droplets that fall back into the water giving this precursor tsunami a very uneven and asymmetric profile. The rapid dissipation of the asteroid’s kinetic energy is very much like an explosion. Shocks propagate outward from the cavity in the water, in the basalt crust and in the mantle beneath. Multiple reflections of shocks and acoustic waves between the material interfaces complicate the dynamics. The hot vapor from the initial cavity expands into the atmosphere, mainly downstream because of the momentum of the asteroid. When the pressure of the vapor in the cavity has diminished sufficiently, at about 35 sec after the impact, water begins to fill the cavity from the bottom. This filling has a high degree of symmetry because of the uniform gravity responsible for the water pressure. An asymmetric fill could result from non uniform seafloor topography, but that case was not modeled. The filling water converges on the center of the cavity and the implosion produces another series of shock waves, and a jet that rises vertically in the atmosphere to a height in excess of 20 kilometers at a time of 115 sec after impact, as shown in Figure 6.28. It is the collapse of this central vertical jet that produces the principal tsunami wave shown in Figure 6.29.

Numerical modeling of water waves

206

Fig. 6.28. Similar to Figure 6.27 but 115 sec after impact, which is the time of the formation of the central vertical jet. The collapse of the crown has produced a circular rim that is propagating outward, but the principal tsunami wave will be produced by the collapse of the central vertical jet. The crown splash has collapsed, pockmarking the surface of the first wave.

The evolution of this wave was modeled out to a time of 400 sec after impact. The inclined impact eventually produces a tsunami that is nearly circularly symmetric at late times, as shown in Figure 6.29. The tsunami has an initial height in excess of 1 kilometer, and decreases to 100 meters at a distance of 40 kilometers from the initial impact. Its propagation speed is 175 meters/sec, which is much less than the shallow water speed of 221 meters/sec.

Modeling waves using compressible models

207

Fig. 6.29. Similar to Figure 6.27 and 6.28 but at 270 sec. The central jet has collapsed and formed the tsunami wave. The pock-marked precursor wave and the smoother principal wave are evident. The principal wave is about 1.5 kilometers in initial amplitude and moves with a speed of 175 meters/sec.

The evolution of this wave was modeled out to a time of 400 sec after impact. The inclined impact eventually produces a tsunami that is nearly circularly symmetric at late times as shown in Figure 6.29. The tsunami has an initial height in excess of 1 kilometer, and decreases to 100 meters at a distance of 40 kilometers from the initial impact. Its propagation speed is 175 meters/sec, which is much less than the shallow water speed of 221 meters/sec. The 45 degree angle chosen for the three-dimensional simulation is the most probable angle for impacts. Impacts at 30 degrees and 60 degrees were also modeled. Significant differences with impact angle were found for the KT impactor described in reference 17 and in the next section of this chapter. Conclusions All asteroid impacts, including the oblique ones, produce a large underwater cavity with nearly vertical walls followed by a collapse starting from the bottom and subsequent vertical jetting. Substantial amounts of water are vaporized and lofted into the atmosphere. In the larger impacts, significant amounts of the crustal and even mantle material are lofted as well. Tsunamis up to a kilometer in initial height are generated by the collapse of the vertical jet. These waves are initially very complex in form and interact strongly with shocks propagating through the water and the crust. The tsunami waves were followed out to 100 kilometers from the point of impact. Their periods and wavelengths show them to be intermediate type waves and not shallow water waves. At

Numerical modeling of water waves

208

great distances, the waves decay as the inverse of the distance from the impact point, ignoring sea floor topography. For all impactors smaller than about 2 kilometers diameter, the impacting body is highly fragmented and it remains lofted into the stratosphere with the water vapor and crustal material, hence very little trace of the impacting body will be found for most oceanic impacts. In the oblique impacts, the initial asymmetry of the transient crater and crown does not persist beyond a tsunami propagation length of 50 kilometers.

6E. KT Chicxulub Event Introduction The impact that created the Chicxulub crater in Mexico’s Yucatan Peninsula is widely believed to be responsible for the mass extinctions at the end of the Cretceous period. This event has been extensively studied and modeled during the few years since its discovery. The stratigraphy at Chicxulub involves rocks like calcite and anhydrite that are highly volatile at pressures reached during impact. The volatility of the Chicxulub strata made the event very dangerous to the megafauna of the late Cretaceous. On a geological time scale, impacts of asteroids and comets with the earth are a relatively frequent occurrence, causing significant disturbances to biological communities and strongly perturbing the course of evolution. Most famous among known catastrophic impacts is the one that marked the end of the Cretceous period and the dominance of the dinosaurs. It is now widely accepted that the world sequence of mass extinctions at the Cretaceous-Tertiary (KT) boundary 65 million years ago was directly caused by the collision of an asteroid or comet with the earth. Evidence for this includes the large (200 kilometers diameter) buried impact structure at Chicxulub, Yucatan, Mexico, the world wide distributed Iridium layer at the KT boundary, and tsunami deposits well inland in North America, all dated to the same epoch as the extinction event. As described in reference 17, the KT impactor was an asteroid of diameter roughly 10 kilometers. Its impact was oblique, either from the SE at 30 degrees to the horizontal or from the SW at 60 degrees. The asteroid encountered layers of water, anhydrite, gypsum, and calcium carbonate, which resulted in the lofting of many hundreds of cubic kilometers of these materials into the stratosphere, where they resided for many years and produced a global climate cooling that was fatal to many large animal species on earth. Three-Dimensional KT Impact Model The Chicxulub impact structure was discovered by scientists from Petroleos Mexicanos (Pemex), the Mexican national oil company18, 19. It is believed to be the location of an impact that was responsible for the mass extinction at the end of the Cretaceous period, as proposed by Alvarez20 et al. on the basis of the anomaly in abundances of Iridium and other Platinum group elements in the KT boundary bedding plane. Paleogeographic data suggest that the crater site, which presently straddles the Yucatan

Modeling waves using compressible models

209

coastline, was submerged at the end of the Cretaceous period on the continental shelf. The substrate consisted of fossilized coral reef over continental crust. In the Gisler, et al.17 simulation, the multilayered target was 300 meters of water, 3 kilometers of calcite, 30 kilometers of granite and 18 kilometers of mantle material. Above the target was an exponential atmosphere up to 106 kilometers altitude. The asteroid’s plunge was numerically started at the 106 kilometers altitude. Threedimensional simulations were performed with impact angles of 30, 45 and 60 degrees to the horizontal. The computational domain extended 256 kilometers by 128 kilometers and simulated the half space of 256 by 256 kilometers. The calculations were performed on the new ASCI Q computer at Los Alamos, a cluster of ES45 Alpha boxes from HP/Compaq. The problems used 1024 processors and a total of about 1 million cpu hours. The adaptive mesh included up to a third of a billion computation cells. Three prominent features of the simulation are illustrated by the 45 degree impact shown in Figure 6.30.

Fig. 6.30. The density at seven sec after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. The isosurface, at density 0.005 g/cc is chosen to show everything denser than air. The scale is set by the back boundary which is 256 kilometers long. The maximum height of the uplifted material is 50 kilometers.

After impact, billions of tons of very hot material are lofted into the atmosphere. A water plume directed away from the impact angle carries much of the horizontal component of the asteroid’s momentum in the downrange direction. This material, consisting of vaporized fragments of the projectile mixed with the target, is extremely hot, and will ignite vegetation many hundreds of kilometers away from the impact site. The highly turbulent and energetic ejecta plume is directed predominatly upward as

Numerical modeling of water waves

210

shown in Figure 6.31.

Fig. 6.31. The density at 42 sec after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. The isosurface, at density 0.005 g/cc is chosen to show everything denser than air. The scale is set by the back boundary which is 256 kilometers long. The maximum height of the uplifted material is 50 kilometers.

Some material is projected into orbits that terminate far outside the computational volume. The dissipation of the kinetic energy, some 300 teratons TNT equivalent, produces a stupendous explosion that melts, vaporizes and ejects a substantial volume of water and granite. Ballistic trajectories carry some of this material back to earth in a conical debris curtain that gradually moves away from the crater lip and deposits a blanket of ejecta around the forming crater as shown in Figure 6.32. The debris curtain has separated from the rim of the still forming crater as material in the curtain falls to earth. The debris from the curtain is deposited in a blanket of ejecta that is asymmetric around the crater with more in the downrange than in the uprange direction. The distribution of material in the ejecta blanket was used as a diagnostic to determine the direction and angle of impact of the asteroid. Some material is projected into orbits that terminate far outside the computational volume.

Modeling waves using compressible models

211

Fig. 6.32. The density at two minutes after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. The isosurface, at density 0.005 g/cc is chosen to show everything denser than air. The scale is set by the back boundary which is 256 kilometers long. The debris curtain has separated from the crater.

The blanket of ejecta is found to be strongly asymmetrical around the crater, with the uprange portion much thinner than the rest. This is a result of the coupling of the horizontal component of the asteroid’s momentum to the debris, and to the ionized and shocked atmosphere in the asteroid’s wake.

Numerical modeling of water waves

212

Fig. 6.33. The density at 122 sec after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. A cut through the center of the impact cavity shows that no water is present from the center of the crater out to the debris curtain.

As shown in Figure 6.33, the ocean water is mostly vaporized in the crater region and beyond the crater to the debris curtain which has separated from the rim of the crustal crater. The initial shallow 300 meters deep water where the asteroid impact occurred would not support the generation of major long period tsunami waves. However, the ocean water is driven by strong shock waves and the expanding curtain out into surrounding deep ocean. This would result in large, long period tsunamis that would propagate around the world and leave the tsunami deposits well inland in North America, as observed during the KT epoch. Spectacular computer movies of the KT impact and a PowerPoint presentation produced by Dr. Galen Gisler are on the NMWW CD-ROM in the directory /NOBEL/KTIMPACT/. Dr. Gisler performed the KT asteroid impact calculations as part of the Los Alamos Crestone Project validation studies. The NOBEL/SAGE/RAGE hydrodynamic codes and the 3-D graphics were developed by the Crestone Project team. The Chicxulub PowerPoint presentation has been shown to many members of the legislative and executive branches of the U.S. government and to Los Alamos visiting dignitaries including Prince William.

Modeling waves using compressible models

213

If life on earth is to continue we must develop the capability of deflecting asteroids before they collide with the earth. Solem21 has described the options available for deflection of asteroids and comets. He concludes “the most cost effective and the ONLY currently available means of asteroid disruption (deflection or pulverization) is a nuclear explosive.” So the Los Alamos nuclear weapon program that has developed most of the water wave modeling technologies described in this book as part of its nuclear weapon effects studies also has developed the nuclear technology that is both a major threat to life on earth and that which makes it possible to defend the earth from the extinction of life by asteroids and comets.

References 1. Charles L.Mader, Numerical Modeling of Explosives and Propellants, CRC Press, Boca Raton, Florida (1998). 2. M.L.Gittings, “1992 SAIC’s Adaptive Grid Eulerian Code,” Defense Nuclear Agency Numerical Methods Symposium, pp. 28–30 (1992). 3. R.L.Holmes, G.Dimonte, B.Fryxell, M.L.Gittings, J.W. Grove, M.Schneider, D.H.Sharp, A.L.Velikovich, R.P.Weaver and Q.Zhang, “Richtmyer-Meshkov Instability Growth: Experiment, Simulation and Theory,” Journal of Fluid Mechanics, Vol. 9, pp.55–79 (1999). 4. R.M.Baltrusaitis, M.L.Gittings, R.P.Weaver, R.F.Benjamin and J.M.Budzinski, “Simulation of Shock-Generated Instabilities,” Physics of Fluids, Vol.8, pp.2471–2483 (1996). 5. Don J.Miller, “Giant Waves in Lituya Bay, Alaska” Geological Survey Professional Paper 354-C, U.S. Government Printing Office, Washington, D.C. (1960). 6. Frances E.Calwell, Land of The Ocean Mists-The Wild Ocean Coast West of Glacier Bay, Alaska Northwest Publishing Company, Edmonds, Washington (1986). 7. Charles L.Mader, “Modeling the 1958 Lituya Bay Tsunami,” Science of Tsunami Hazards, Vol.17, pp.57–67 (1999). 8. Hermann M.Fritz, Willi H.Hager and Hans-Erwin Minor, “Lituya Bay Case: Rockslide Impact and Wave Runup,” Science of Tsunami Hazards, Vol.19, pp.3–22 (2001). 9. Charles L.Mader, “Modeling the 1958 Lituya Bay Mega-Tsunami, II,” Science of Tsunami Hazards, Vol.20, pp.241–250 (2002). 10. Bobby G.Craig, “Experimental Observations of Underwater Detonations Near the Water Surface,” Los Alamos Scientific Laboratory report LA-5548-MS (1972). 11. Donald E.Gault and Charles P.Sonett, “Laboratory Simulation of Pelagic Asteroid Impact: Atmospheric Injection, Benthic Topography, and the Surface Wave Radiation Field,” Geological Society of America, Special Paper 190, pp.69–92 (1982). 12. Valery Kedrinskii, private communication (1985). 13. Charles L.Mader and Michael L.Gittings, “Dynamics of Water Cavity Generation,” Science of Tsunami Hazards, Vol.21, pp. 91–118 (2003). 14. Charles L.Mader, Numerical Modeling of Detonations, University of California Press, Berkeley, California (1979). 15. Charles L.Mader, Timothy R.Neal and Richard D.Dick, LASL Phermex Data, Volume I, II, and III, University of California Press, Berkeley, California (1980). Available at http://lib-www.lanl.gov/ladcdmp/phl.pdf. 16. Galen Gisler, Robert Weaver, Michael L.Gittings and Charles Mader, “Two- and

Numerical modeling of water waves

214

Three-Dimensional Simulations of Asteroid Ocean Impacts,” Science of Tsunami Hazards, Vol.21, pp.119–134 (2003). 17. Galen Gisler, Robert Weaver, Michael L.Gittings and Charles Mader, “Two- and Three-Dimensional Asteroid Impact Simulations,” Computers in Science and Engineering (2004). 18. A.R.Hildebrand, G.T.Penfield, D.A.Kring, M.Pilkington, Z.A.Camargo, S.B.Jacobsen and W.V.Boynton, “Chicxulub Crater: A Possible Cretaceous/Tertiary Boundary Impact Crater on the Yucatan Peninsula, Mexico.” Geology, Vol.19, pp.867–871 (1991). 19. V.L.Sharpton, G.B.Dalrymple, L.E.Marin, G.Ryder, B. C.Schuraytz and J.UrruitaFucugauchi, “New Links Between the Chicxulub Impact Structure and the Cretaceous/Tertiary Boundary,” Nature, Vol.359, pp.819–821 (1992). 20. L.W.Alvarez, W.Alvarez, F.Asaro and H.V.Michel, “Extraterrestrial Cause for the Cretaceous/Teritary Extinction,” Science, Vol.208, pp.1095–1108 (1980). 21. Johndale C. Solem, “Comet and Asteroid Hazards: Threat and Mitigation,” Science of Tsunami Hazards, Vol.17, pp.141–154 (1999).

CD-ROM CONTENTS Tsunami Animations A collection of tsunami animations performed using the SWAN code described in the Numerical Modeling of Water Waves—Second Edition, by Dr. Charles L.Mader. A few calculations were performed using the full Navier-Stokes ZUNI or SOLA codes. The tsunami animations are archived at http://tl4web.lanl.gov/Staff/clm/tsunami.mve/tsunami.htm. To view animation type or click—MOVIE—To study individual cases type—MM TAPEXX where XX is file number. Follow instructions on screen. 1960.MVE—May 23, 1960 tsunami generation in Chile, propagation across the Pacific Ocean, and indundation of Hilo, Hawaii. Described in “Modeling Hilo, Hawaii Tsunami Inundation,” Science of Tsunami Hazards, Vol.9, pp.85–94 (1991), and Scientific Computing and Automation, June issue, pp.19–23 (1993). 1964.MVE—Tsunami of April 1, 1964 generation in Gulf of Alaska, propagation across the Pacific Ocean, and inundation of Crescent City, California. See “Tsunami Inundation Model Study of Eureka and Crescent City, California,” NOAA Tech. Memo. ERL PMEL-105 (1994). 60THEAT.MVE—The interaction of the tsunami of May 23, 1960 with the Hilo, Hawaii Theater. Described in PACON 1993. 90HILO.MVE—1990 Hilo topography and buildings inundated by a 1960 tsunami wave. See also HOTEL.MVE. 2ATAST.MVE—The inundation of the U.S. East Coast by a 100 meters, 2000 sec tsunami wave that could be generated by an asteroid. 10NYAST.MVE—The inundation of the U.S. East Coast by a wave from the incompressible collapse of a 10 kilometers radius cavity, 3000 meters deep and a 100 kilometers radius cavity in the Atlantic Ocean off New York City. AIMPACT.MVE—An impact cavity collapse and tsunami generation study using shallow water and full Navier-Stokes models. Described in “Modeling Asteroid Impact and Tsunami,” Science of Tsunami Hazards, Vol.16, pp.21–30 (1998). ATLAST.MVE—A tsunami in the Atlantic Ocean generated by the incompressible collapse of a cavity 150 kilometers wide and 3500 meters deep. AUSAST.MVE—Interaction of a tsunami with Australia from a Hawaii landslide generated tsunami and from a cavity collapse generated tsunami. Described in “Modeling of Tsunami Propagation Directed at Wave Erosion on Southeastern Australia Coast 105,000 Years Ago,” Science of Tsunami Hazards, Vol.13, pp.45–52 (1995). BBAY.MVE—A study of the vulnerability of Berau Bay, Indonesia to tsunamis. CASCAD.MVE—Inundation of U.S. West Coast by a tsunami from the Cascadia fault. ECAST.MVE—The inundation of the U.S. East Coast by a tsunami generated by the incompressible collapse of a 150 kilometers wide, 3000 meters deep cavity. See also

CD-Rom contents

217

NYAST.MVE ELTAST.MVE—Described in “Modeling the Eltanin Asteroid Tsunami,” Science of Tsunami Hazards, Vol.16, pp.17–20 (1998). EURAST.MVE—The inundation of Europe by a 100 meters high and 2000 sec period tsunami. EUREKA.MVE—The Eureka, California tsunami of April 25, 1992. See “Tsunami Inundation Model Study of Eureka and Crescent City, California,” NOAA Tech. Memo. ERL PMEL-105 (1994). GUS.MVE—The Furumoto sources for the Hawaiian tsunamis of 1946, 1957, 1964 and 1965. Part of a source modeling project for Dr. A.Furumoto, Hawaii Civil Defense Tsunami Advisor. HIAST.MVE—The inundation of the Hawaiian Islands by a 100 meters high, 2000 sec period tsunami wave. Described in “Asteroid Tsunami Inundation of Hawaii,” Science of Tsunami Hazards, Vol. 14, pp.85–88 (1996). HILAND.MVE—The tsunami generated by a landslide off the Kona coast of the island of Hawaii about 105 Ka years ago. Described in “Modeling the 105 Ka Landslide Lanai Tsunami,” Science of Tsunami Hazards, Vol.12, pp.33–38 (1994). HKAI.MVE—Inundation of Hawaii Kai, Hawaii by a typical off shore 3 meters high, 1500 sec tsunami wave. HOTEL.MVE—The interaction of a May 23, 1960 tsunami wave with current Hilo, Hawaii tourist hotels. See also 90HILO.MVE. HUMBOL.MVE—Tsunami inundation of Humboldt Bay, California by an offshore maximum expectable 10 meters high, 2000 sec tsunami wave. See “Tsunami Inundation Model Study of Eureka and Crescent City, California,” NOAA Tech. Memo. ERL PMEL-105 (1994). ICEAST.MVE—The inundation of Iceland by a 100 meters high and 2000 sec period tsunami. INDIA.MVE—Tsunami in the Indian Ocean generated by the incompressible collapse of a cavity 38 kilometers wide and 4000 meters deep. INDONES.MVE—Indonesia tsunami of December 12, 1992. JAPAST.MVE—The inundation of Tokyo, Japan by a tsunami generated by a incompressible cavity collapse. Described in “Asteroid Tsunami Inundation of Japan,” Science of Tsunami Hazards, Vol.16, pp.11–16 (1998). KAIAKA.MVE—Tsunami inundation of Kaiaka Bay, Oahu, Hawaii by the 1952 tsunami. KBAY.MVE—Tsunami inundation of Kaneohe Bay, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami and by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave. KONA.MVE—Tsunami inundation of Kona, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami wave. KURIL.MVE—The tsunami of October 1994 generated off the Kuril islands of Japan. LAAST.MVE—Inundation of Los Angeles, California by a 100 meters high, 2000 sec period tsunami wave. LAPALMA.MVE—Modeling the proposed La Palma landslide tsunami. Published in “Modeling the La Palma Landslide Tsunami,” Science of Tsunami Hazards, Vol.19,

CD-Rom contents

218

pp.160–180 (2001). LAUP.MVE—The April 1, 1946 tsunami inundation of Laupahoehoe, Hawaii. LITUYA.MVE—The July 8, 1958 mega-tsunami at Lituya Bay, Alaska with inundations up to 520 meters. Described in “Modeling the 1958 Lituya Bay Mega-Tsunami,” Science of Tsunami Hazards, Vol.17, pp.57–67 (1999). The Lituya Bay impact landslide generation of the tsunami is described in Chapter 6 and in Science of Tsunami Hazards, Vol.20, pp.241–250 (2002). LISBON.MVE—Modeling the 1755 Lisbon tsunami generation and propagation across the Atlantic Ocean to the Caribbean. Science of Tsunami Hazards, Vol.19, pp.93–98 (2001). LOIHI.MVE—A study using the ZUNI full Navier-Stokes code of the tsunami wave generation and propagation from the collapse of the Loihi, Hawaii summit in August, 1996. M9CALIF.MVE—An M9 earthquake generated tsunami interacting with the Oregon and California coasts. NIC.MVE—The tsunami generated off the coast of Nicaragua in 1992. Described in “Modeling the 1992 Nicaragua Tsunami,” Science of Tsunami Hazards, Vol.11, pp.107– 110 (1993). NYAST.MVE—The inundation of the U.S. East Coast by the incompressible collapse of a 100 kilometers radius 3000 meters deep cavity. Another tsunami wave had a height of 100 meters and a 2000 sec period. See also 10NYAST.MVE. ORAST.MVE—A 100 meters high, 2000 sec period tsunami interacting with the Oregon coast. OREGM9.MVE—An M9 earthquake generated tsunami interacting with the Oregon coast. PACAST.MVE—Tsunami in the middle of the Pacific formed from the incompressible collapse of a cavity 150 kilometers wide and 4500 meters deep. PROP.MVE—Described in “Numerical Tsunami Propagation Study,” Science of Tsunami Hazards, Vol.11, pp.93–106(1993) and in Chapter 5 of Numerical Modeling of Water Waves—Second Edition. SANDY.MVE—Tsunami inundation of Sandy Beach region of Oahu, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami and by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave. SANFAST.MVE—Inundation of San Francisco, California by a 100 meters high, 2000 sec tsunami wave. SKAGWAY.MVE—The landslide generated tsunami of November 3, 1994 at Skagway, Alaska. The Skagway modeling is described in “Modeling the 1994 Skagway Tsunami,” Science of Tsunami Hazards, Vol.15, pp.41–48 (1997). See also SOLA.MVE. SMSFAST.MVE—Inundation of San Francisco by a tsunami wave generated by the incompressible collapse of a 20 kilometers wide, 3000 meters deep cavity. SOLA.MVE—Three-dimensional, full Navier-Stokes modeling using the MCC SOLA code of the November 3, 1994 Skagway, Alaska tsunami. See also SKAGWAY.MVE. SOURCE.MVE—Described in “Numerical Tsunami Source Study,” Science of Tsunami Hazards, Vol.11, pp.81–92 (1993) and in Chapter 5 of Numerical Modeling of Water Waves—Second Edition.

CD-Rom contents

219

VSLIDE.MVE—A landslide generated tsunami from the Chain of Craters road region of the island of Hawaii. WAIANAE.MVE—The inundation of the leeward side of Oahu, Hawaii by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave. WAIPIO.MVE—The interaction of the May 23, 1960 tsunami with the Waipio, Hawaii region. The 50 foot inundation is the largest recorded in Hawaii. WALKER.MVE—An evaluation of the vulnerability of Hawaii to tsunamis generated south of Honolulu, either along the Kona Coast or in the Tonga trench. Modeling requested by Dr. D.Walker, Oahu Civil Defence Tsunami Advisor. WINDWARD.MVE—Tsunami inundation of the Windward side of Oahu, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami and by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave.

NOBEL PowerPoint Presentations A collection of PowerPoint presentations describing water wave studies performed using the compressible hydrodynamic code NOBEL. The studies are described in Chapter 6 of Numerical Modeling of Water Waves—Second Edition and archived at http://t14web.lanl.gov/Staff/clm/tsunami.mve/tsunami.htm. The PowerPoint presentations may be viewed using PPVIEW in /NOBEL/PPRESENT/. LITUYA—The July 8, 1958 Lituya Bay, Alaska impact landslide tsunami generation. A mega-tsunami was generated that reached an altitude of 520 meters. Laboratory experiments and numerical modeling results are presented. Described in “Modeling the 1958 Lituya Bay Mega-Tsunami, II,” Science of Tsunami Hazards, Vol. 20, pp.241–250 (2002). CAVITY—The generation of cavities in water by projectile impacts and by explosives is described both experimentally and using compressible hydrodynamic models. Described in “Dynamics of Water Cavity Generation,” Science of Tsunami Hazards, Vol.21, pp. 91–118 (2003). ASTWAVE—The generation of tsunamis by the impact of a 0.25 to 1 kilometer diameter asteroid at 20 kilometers/sec with 5 kilometers of ocean and 5 kilometers of basalt is modeled using compressible hydrodynamics in two and three dimensions. Described in “Two- and Three-Dimensional Simulations of Asteroid Ocean Impacts,” Science of Tsunami Hazards, Vol.21, pp.119–134 (2003). KTIMPACT—The KT Chicxulub asteroid impact event is modeled using the threedimensional compressible Navier-Stokes model. Described in “Two and ThreeDimensional Asteroid Impact Simulations,” Computers in Science and Engineering (2004). CD-ROM CODE DIRECTORIES WAVE—The WAVE code described in Chapter 1 solves the equations for Airy, third-

CD-Rom contents

220

order Stokes and Laitone solitary gravity waves. The directory contains the FORTRAN source code, the executable code for DOS or Windows and WAVE.PDF which describes the code. SWAN—The shallow water SWAN code described in Chapter 2 solves the long wave, shallow water, nonlinear equations of fluid flow. The directory contains the FORTRAN source and executable codes which generate a graphics file that may be processed using the programs included. It also includes a description of the input to the code in the file SWAN.PDF. Examples and topographic files are furnished. ZUNI—The incompressible Navier-Stokes ZUNI code described in Chapter 3 solves the incompressible, viscous fluid flows with a free surface using the Navier-Stokes equations. A detailed description of the computer program and its input file is included in the file ZUNI.PDF. The FORTRAN source and the executable codes are included. SOLA—The incompressible three-dimensional Navier-Stokes ZUNI code described in Chapter 4 solves the incompressible viscous fluid flows with a free surface using the Navier-Stokes equations. The FORTRAN source and the executable codes are included. The Skagway 1994 tsunami is used as an example. LGW—The Carrier linear gravity wave LGW code described in Chapter 5 uses analytical methods for solving the linear gravity model. The FORTRAN source and executable codes are included. Examples of Gaussian tsunamis described in Chapter 5 are furnished. TIDE—A classic computer program for calculating tides with the FORTRAN source and executable codes furnished.

Science of Tsunami Hazards Directory All the Science of Tsunami Hazards journals through 2003 in PDF format may be searched using Adobe Acrobat 4.0 or higher. Issues of the journal are archived at http://epubs.lanl.gov/tsunami and current issues are available at http://www.sthjournal.org. Dir=TS211.PDF, TS212.PDF, TS213.PDF, TS214.PDF * Volume 21 (No.1), (No.2), (No.3), (No.4) 2003 Dir=TS201.PDF, TS202.PDF, TS203.PDF, TS204.PDF, TS205.PDF * Volume 20 (No.1), (No.2), (No.3), (No.4), (No. 5), 2002 DIR=TS191.PDF, TS192.PDF, TS193.PDF * Volume 19 (No.1) (No.2) (No.3), 2001 DIR=TS181.PDF, TS182.PDF * Volume 18 (No.1) (No.2), 2000 DIR=TS171.PDF, TS172.PDF, TS173.PDF * Volume 17 (No.1) (No.2) (No.3), 1999 DIR=00394718.PDF * Volume 16 (No.1), 1998 DIR=00394719.PDF, 00394720.PDF

CD-Rom contents

221

* Volume 15 (No.1) (No.2), 1997 DIR=00394721.PDF, 00394722.PDF, 00394723.PDF * Volume 14 (No.3) (No.2) (No.1), 1996 DIR=00394724.PDF * Volume 13 (No.1), 1995 DIR=00394725.PDF, 00394726.PDF * Volume 12 (No.2) (No.1), 1994 DIRECTORY=00394727.PDF, 00394728.PDF * Volume 11 (No.2) (No.1), 1993 DIRECTORY=00394729.PDF * Volume 10 (No.1), 1992 DIRECTORY=00394730.PDF, 00394731.PDF * Volume 9 (No.2) (No.1), 1991 DIRECTORY=00394732.PDF, 00394733.PDF * Volume 8 (No.2) (No.1), 1990 DIRECTORY=00394734.PDF, 00394735.PDF * Volume 7 (No.2) (No.1), 1989 DIRECTORY=00394736.PDF * Volume 6 (No.1), 1988 DIRECTORY=00394737.PDF, 00394738.PDF * Volume 5 (No.2) (No.1), 1987 DIRECTORY=00394739.PDF, 00394740.PDF, 00394741.PDF * Volume 4 (No.3) (No.2) (No.1), 1986 DIRECTORY=00394742.PDF * Volume 3 (No.1), 1985 DIRECTORY=00394743.PDF, 00394744.PDF * Volume 2 (No.2) (No.1), 1984 DIRECTORY=00394745.PDF * Volume 1 (No.1), 1982 CNMWW.PDF is a searchable PDF file of the book Numerical Modeling of Water Waves—Second Edition with many figures in color.

AUTHOR INDEX Numbers in italic type indicate pages where references are listed in full.

Alvarez, L.W., 224, 231 Alvarez, W., 224, 231 Amsden, A.A., 86, 118 Asaro, F., 223, 231 Baltrusaitis, R.M., 180, 229 Batchelor, G.K., viii, 21 Benjamin, R.F. 180, 229 Bernard, E., 56, 82 Bornhold, B.D., 65, 82 Boynton, W.V., 223, 230 Bretschneider, C.L., 166, 177 Budzinski, J.M., 180, 229 Butler, H.L., 39, 82 Calwell, F.E., 183, 229 Camargo, Z.A., 223, 230 Campbell, B., 65, 70, 78, 83 Carrier, G.F., 146, 151, 177 Chan, R.K.C., 86, 92–3, 95, 118 Chubarov, L.B., 26, 30, 80 Cox, D.C., 135, 146 Craig, B.G., 110, 114, 119, 194, 196, 198, 203, 209–10, 230 Curtis, G.D., 26, 31, 48, 56, 81, 82 Dalrymple, G.B., 223, 231 Daly, B.J., 85, 93, 118 Dick, R.D., 203, 230 Dimonte, G., 180, 229 Divoky, D., 39, 44, 46, 82 Dronkers, J.J., 22, 79 Fritz, H.M., 186, 190, 230 Fromm, J.E., 93, 95, 118 Fryxell, B., 180, 229 Fuchs, R.A., 105, 107, 119 Fucugauchi, J.U., 223, 231 Furumoto, A.S., 34, 81

Author index

224

Garcia, A.W., 39, 82 Garcia, W.J., 134, 145 Gault, D.E., 194, 210, 230 Gisler, G., 211, 212, 221, 224, 228, 230 Gittings, M.L., 114, 119, 179, 180, 193, 211, 212, 221, 223, 224, 229, 230 Green, C.K., 34, 81 Greenspan, H.P., 150, 177 Grove, J.W., 180, 229 Hadi, S., 30, 80 Hager, W.H., 186, 230 Hansen, W., 22, 24, 79 Harlow, F.H., 85, 93, 118 Hildebrand, A.R., 223, 230 Hirt, C.W., 86, 92, 94, 118, 119, 145 Holmes, R.L., 180, 229 Houston, J.R., 39, 82 Hwang, L., 39, 44, 46, 82 Jacobsen, S.B., 223, 230 Johnson, J.W., 105, 108, 119 Kedrinskii, V., 194, 195, 199, 210, 230 Kinsman, B., viii, 21 Kowalik, Z., 26, 30, 47, 76, 80, 82 Kring, D.A., 223, 230 Kulikov, E.A., 65, 82 Landau, L.D., viii, 21 Lander, J.F., 64, 82 Le Mehaute, B., 3, 21 Lee, J.J., 64, 66, 76, 77, 82 Leenderstse, J.J., 80 Lifshitz, E.M., viii, 21 Loomis, H.G., 135, 146 Lukas, S., 30, 80 Mader, C.L., 21, 26, 29, 30, 31, 48, 56, 80, 81, 82, 88, 93, 109, 118, 119, 135, 145, 146, 177, 179, 185, 189, 193, 202–4, 208, 211, 212, 221, 223, 224, 229, 230 Madsen, O.S., 97, 100, 118 Marchuk, A.G., 26, 30, 80 Marin, L.E., 223, 231 Mei, C.C., 97, 100, 118 Michel, H.V., 223, 231 Miller, D.J., 183, 186, 229 Minor, H.E., 186, 230

Author index Moore, D.W., 146, 177 Morison, J.R., 105, 107, 119 Murty, T.S., 26, 30, 80 Nabeshima, G., 48, 82 Neal, T.R., 203, 230 Nichols, B.D., 86, 92, 94, 119, 135, 145 Nottingham, D., 65, 69, 70, 78, 83 Pelinovsky, E.N., 26, 30, 80 Penfield, G.T., 223, 230 Peregrine, D.H., 97, 118 Petroff, C., 64, 66, 76, 77, 83 Pilkington, M., 223, 230 Plafker, G., 39, 44, 58, 82 Proudman, J., 22, 79 Rabinovich, A.B., 65, 82 Raichlen, F., 66, 77, 83 Ramming, H.G., 26, 30, 80 Romero, N.C., 119, 145 Ryder, G., 223, 231 Satake, K., 56, 82 Savage, J.C., 44, 82 Schneider, M., 180, 229 Schuraytz, C., 223, 231 Shannon, J.P., 85, 86, 93, 118 Sharp, D.H., 180, 229 Sharpton, V.L., 223, 231 Shokin, I.I., 26, 30, 80 Sklarz, M.A., 134, 146 Sokolowski, T.J., 47, 82 Solem, J.C., 229, 231 Sonett, C.P., 194, 197, 198, 210, 230 Spielvogel, L.Q., 134, 146 Stein, L.R., 119, 145 Street, R.L., 86, 92–3, 96, 118 Tangora, R.E., 135, 145 Thomson, R.E., 65, 82 Uusitalo, S., 24, 80 Van Dorn, W.G., 39, 59, 81 Velikovich, A.L., 180, 229 Vitousek, M., 30, 81

225

Author index Watts, P, 64, 66, 76, 77, 83 Weaver, R.P., 180, 211, 212, 223, 224, 229, 230 Welander, P., 24, 79 Welch, J.E., 85, 93, 118 Whalin, R.W., 39, 82 Whitmore, P.M., 47, 82 Wiegel, R.L., 13, 14, 21 Wybro, P.G., 166, 177 Zhang, Q., 180, 229

226

SUBJECT INDEX 1946, April 1 Tsunami, 31–46, 57–61, 232, 233 1958, July 8, Lituya Bay Tsunami, 181–93 1960 Tsunami Type, 52–3, 232 1960, May 23 Tsunami, 31, 32, 41–8, 52, 81, 155, 162, 231, 235 1964, March 28 Tsunami, 31, 32, 38–46, 56–62, 79–81, 155, 162, 231, 232 1975, November 29, Hawaii Tsunami, 134–42 1994, November 3, Skagway Tsunami, 63–78, 144 Accuracy, 7, 26, 110, 125–6, 134 Airy Wave, viii, ix, 9–12, 20, 29, 94, 146, 148, 150, 158, 159, 163 AMR, 190, 193, 196, 210, 211 ASCI, 177, 212, 218, 224 Asteroid, 186, 194, 210–2, 214–6, 222–37 Barrier, 105–8 Beach, 80, 86, 118, 177, 235 Berau Bay, 232 BlueMountain, 212 Bore, 3, 6, 32, 42, 46, 48 Bottom Motion, 25, 27 Boundary Layer, 4 Bridgewire, 194, 206 Brigham Victory Ship, 39 Bubble Cavity, 206, 209, 210 Bubble Radius, 110 Building, 32, 43, 47–52 Cartesian, 21, 84, 120, 124 Chicxulub, 223, 229, 230, 237 Comet, 223, 229 Compressible Flow, vii, ix, 110, 156 Compressible Model, 176, 178–80 Conservation of Energy, viii, 1 Conservation of Momentum, viii, x, 85, 124 Continental Slope, 30, 93, 95–102 Continuative Boundary, 131, 132 Continuum Boundary, 29 Convergence Criterion, 94, 111 Coriolis, 20, 22, 27, 29, 31, 32 Crestone Project, 228

Subject index

228

DeChezy, 24, 27, 36, 44, 169 Deep Water Wave, 113, 115, 116, 143, 208 Detonation, 110, 115, 194, 203 Dockslide, 68, 69, 75 Donor Cell, 124–6, 134 Earth Rotation, 22 Earthquake, 32, 33, 39, 41, 44, 46, 57–8, 62, 82, 135–7, 142, 153, 155, 162, 176, 181, 182, 185, 211 Estuaries, 3, 30, 80 Eulerian, ix, x, 1, 7, 10, 85, 177, 180, 190, 193, 196, 210, 229 Explosion, 3, 110, 114, 176, 178, 194, 211, 220, 226 Explosive Generated, 195, 202, 208 Extinction of Life, 229 Flood Wave, 3, 6 Flooding, vii, x, 2, 3, 25, 26, 29, 31–62, 156, 162–75 FORTRAN, 237 Free Slip, 85, 130 Gravity, viii, 21, 23, 93, 111, 141, 181, 185, 190, 217, 220 Gravity Waves, 3, 20, 93, 104 Hawaiian Islands, 32, 33, 40, 44, 46, 233 Hilo Bay, 31, 32, 36, 38, 39, 43, 46, 52, 55 Hilo Harbor, 32, 37, 38, 39, 43, 46, 81 Hilo Pier No.1, 38 Hilo Theater, 46–52, 231 Hotel, 32, 52, 55, 232 Hydrostatic, 10, 23, 26, 94, 111 Impact Landslide, 180, 190, 193, 210, 233, 236 Impact Landslide Experiment, 186–9 Incompressible Flow, 85, 118, 130, 177, 194 Indonesia, 30, 232 Intermediate Wave, 7, 9, 222 Jet, 115, 194–5, 207, 210, 213, 217, 221, 222 JIMAR, 81, 146 Kinetic Energy, 114, 215–6, 219, 226 KT Impact, 221, 223, 227, 228, 236 Laitone Solitary Waves, 12, 13, 14, 20, 237 Landslide, 30, 65, 68–82, 135–45, 176, 178, 180, 186–93, 210, 211, 232–6 Linear Gravity Waves, 115, 146, 151–62, 237 Lituya Bay, 180–7, 193, 212, 229, 233, 236

Subject index

229

Loihi, 115–6 Long Wave Paradox, 25, 26 Los Alamos, viii, 224, 228, 229 Los Alamos National Laboratory, 118, 145, 194, 209, 212, 231 MAC, 85, 89, 92–3, 118, 119, 124, 133 Marker and Cell, 85, 119, 134 Mega-tsunami, 181, 193, 211, 212, 233, 236 Momentum, 89, 123, 124, 129, 131–2, 134, 225, 227 Musi Upang, 30, 80 National Academy of Science, 81 Navier-Stokes, vii, 30, 31, 47, 57, 78, 83, 85, 87, 93, 110113, 116, 119, 120, 134, 135, 140–2, 146, 153, 155, 156, 159–63, 170, 174–6, 189, 196, 234–7 News Media, 211 NMWW CD-ROM, 11, 13, 14, 20, 26, 27, 31, 52, 62, 75, 78 114, 115, 144, 151, 162, 191, 210, 215, 227 No-slip, 130, 133 NOAA, 32, 39, 57, 66, 70, 82, 232, 235 NOBEL, 110, 180, 190, 193, 196, 205, 209–11, 215, 228, 236, 237 Nomenclature, viii, 27, 178 NTHMP, 78 Nuclear, 176, 229 Oahu, 234–5 Obstacle, 6, 86, 130, 133 Oscillatory Wave, 3 PARN Dock, 64, 68, 72, 74, 77 PBX-9404 Explosive 110, 203–8 PEMEX, 224 PHERMEX, 204, 230 Piston Boundary, 29, 51, 98, 100, 103 Potential Energy, 7 Potential Function, 4, 86, 90 PowerPoint Presentation, 115, 191, 211, 215, 228, 236 Progressive Waves, 5 Projectile, 38, 194–5, 210, 212 Q Computer, 224 Radiograph, 204–5 RAGE, 180, 211, 229 Reflective Boundary, 29, 136 Reynolds Number, 90 Rigid Wall, 89–91, 130 Root, 115, 194–6, 207, 210, 217

Subject index

230

Runup, 34, 86, 93, 96, 100, 103, 163, 165–75, 230 SAGE, 180, 193, 210, 211, 231 SAIC, 210, 212, 229, 231 Scotch Cap Lighthouse, 33, 34 Sea Floor Displacement, 30 Sea of Japan, 30 Seiche, 3, 6 Shallow Water, 6–11, 15, 21, 22, 26, 30–46, 57, 61, 66, 78, 80, 85, 93–5, 99–114, 135–59, 161–2, 165, 166, 167, 169, 174–5, 180, 185, 193, 202, 208, 211, 216, 217, 220–2, 231 Shoaling, 77, 78, 116, 119, 120, 122, 125, 126, 131, 133, 135, 137, 144, 231, 234, 235, 237 Shock, 181, 195, 205, 216, 219, 220, 227, 229 Skagway, 64–6, 71, 76, 78, 82 SMAC, 88–90, 92, 94, 97 SOLA-3D, 78, 117, 119, 120, 123, 125, 126, 131, 133, 136, 138, 145, 231, 235, 237 Source, 26, 31–46, 58, 62, 72, 82, 115, 134–46, 153–7, 162, 173, 176, 207, 236, 237 Stability, 94, 133, 134 Standing Wave, 5, 6 Stationary Wave, 5 Steady State, 4, 8, 9 Stem, 195 Stokes Wave, 14, 15, 20, 237 Storm Wave, 3 Submerged Barrier, 105 Surface Marker Particles, 85, 92 SWAN, 22, 26, 29, 30, 32, 46, 57, 62, 67, 76, 78, 81, 100, 110–2, 135, 138, 145, 146–8, 154, 156–8, 163, 165–7, 170, 172–4, 186, 231, 237 Three-dimensional, vii, ix, 17, 78, 116, 119, 134, 135, 141–4, 162, 176, 178, 211, 218, 224, 230, 235, 237 Tide Gauge, 34, 38, 64, 65, 66, 68, 74–7, 83, 117 Tide Generating Force, 22 Time Step, 27, 29, 47, 68, 93, 110, 116, 122, 128, 133, 134, 163, 174, 180, 186, 188 Translatory Waves, 3 Transmission Coefficient, 106, 108 Tsunami Museum, 42 Two-dimensional, 20, 25, 30, 83, 93, 116, 134, 151, 157, 180, 212 Underwater Barrier, 107 Unimak Island, 33, 34 Upper Critical Depth, 109, 114, 122, 203 Viscosity, ix, x, 1, 2, 24, 94, 111, 121, 134, 141, 149, 159, 169 Volcanic, 115, 176, 211 Vorticity, 86, 90, 93 Waianae Harbor, 30, 80, 236

Subject index

231

Water Cavity, viii, 15, 112, 119, 180, 193–5, 202, 206, 209, 230, 231, 236 WAVE Code, 10, 13, 14, 17, 20, 80, 208, 237 Wave Energy, 30, 104 Wave Length, viii, 5, 7, 8, 13, 14, 20, 23, 95, 98, 99, 102, 104, 114, 115, 144, 154, 156, 162, 166, 175, 203, 208, 214, 216 Wave Number, viii, 5 Wave Period, viii, 5, 33, 38, 75, 95, 105, 115, 154, 162, 163, 166–7, 172, 174 White Machine, 218 Wind Waves, 24 XTX-8003 Explosive, 203 ZUNI, 85, 93–4, 100, 110, 113, 116–7, 146, 149, 153, 156–9, 163, 169–74, 231, 234, 237

NUMERICAL MODELING of WATER WAVES

SECOND EDITION

NUMERICAL MODELING of WATER WAVES Charles L.Mader

CRC PRESS Boca Raton London New York Washington, D.C.

This edition published in the Taylor & Francis e-Library, 2005. To purchase your own copy of this or any of Taylor & Francis or Routledge’s collection of thousands of eBooks please go to www.eBookstore.tandf.co.uk. Library of Congress Cataloging-in-Publication Data Mader, Charles L. Numerical modeling of water waves/Charles L.Mader.—2nd ed. p. cm. Includes bibliographical references and index. ISBN 0-8493-2311-8 (alk. paper) 1. Water waves—Mathematical model. 2. hydrodynamics. I. Title. QA927.M29 2004 532′.593–dc22 2004047812 This book contains information obtained from authentic and highly regarded sources. Reprinted material is quoted with permission, and sources are indicated. A wide variety of references are listed. Reasonable efforts have been made to publish reliable data and information, but the author and the publisher cannot assume responsibility for the validity of all materials or for the consequences of their use. Neither this book nor any part may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, microfilming, and recording, or by any information storage or retrieval system, without prior permission in writing from the publisher. The consent of CRC Press LLC does not extend to copying for general distribution, for promotion, for creating new works, or for resale. Specific permission must be obtained in writing from CRC Press LLC for such copying. Direct all inquiries to CRC Press LLC, 2000 N.W. Corporate Blvd., Boca Raton, Florida 33431. Trademark Notice: Product or corporate names may be trademarks or registered trademarks, and are used only for identification and explanation, without intent to infringe. Visit the CRC Press Web site at www.crcpress.com © 2004 by CRC Press LLC No claim to original U.S. Government works ISBN 0-203-49219-6 Master e-book ISBN

ISBN 0-203-59054-6 (Adobe e-Reader Format) International Standard Book Number 0-8493-2311-8 (Print Edition) Library of Congress Card Number 2004047812

ACKNOWLEDGMENTS The author gratefully acknowledges the developers of the numerical models described in this monograph—M.L.Gittings of SAIC, R.Weaver, G.Gisler, A.Amsden, C.W.Hirt, F.H.Harlow, B.D. Nichols, L.R.Stein, and D.Butler; and the contributions of J.D. Kershner, R.E.Tangora, A.N.Cox, K.H.Olsen, R.Menikoff, E. M.Kober, B.Gonzales, J.C.Solem, J.G.Hills, R.D.Cowan, B.G. Craig, D.Venable, G.E.Seay, and A.L.Bowman of the Los Alamos National Laboratory. The author wishes to recognize the contributions of Dr. H. Loomis, Dr. W.Adams, Dr. D.Cox, Dr. D.W.Moore, Mr. G. Curtis, Dr. J.Miller, Dr. D.Walker, Dr. S.Furumoto, Dr. L. Spielvogel, Dr. A.Malahoff, Dr. W.Dudley, G.Nabeshima, and S. Lukas of the University of Hawaii. The author wishes to recognize the collaborations with Dr. E. Bernard of the Pacific Marine Environmental Laboratory, Professor George Carrier of Harvard University, Professor Zygmunt Kowalik of the University of Alaska, Dr. Hermann Fritz of the University of Georgia, Professor Valery Kedrinskii of the Institute of Hydrodynamics, Novosibirsk, Russia, Mr. Tom Sokolowski of the West Coast and Alaska Tsunami Warning Center, Mr. Bruce A.Campbell, Mr. Dennis Nottingham of PN&D, Dr. T.S.Murty of Baird & Associates, Dr. David Crawford of Sandia and Dr. Rob T.Sewell of R. T.Sewell Associates. This book was reviewed and determined to be unclassified by the Los Alamos National Laboratory group S-7 and given the designation LA-UR-03–8710.

THE AUTHOR Dr. Charles L.Mader is President of Mader Consulting Co. with offices in Honolulu, Hawaii, Los Alamos, New Mexico and Vail, Colorado. He is an Emeritus Fellow of the Los Alamos National Laboratory. He is the author of Numerical Modeling of Water Waves, Numerical Modeling of Explosives and Propellents, Numerical Modeling of Detonations and four Los Alamos Dynamic Material Properties data volumes. Dr. Mader is a Fellow of the American Institute of Chemists, Editor of the Science of Tsunami Hazards Journal and is listed in Who’s Who in America and in Who’s Who in the World. His web site is http://www.mccohi.com.

INTRODUCTION This second edition of Numerical Modeling of Water Waves describes the technological revolution that has occurred in numerical modeling of water waves during the last decade. A CD-ROM with many of the FORTRAN codes of the numerical methods for solving water wave problems and computer animations of the problems that have been solved using the codes is included with the book. Several PowerPoint presentations describing the modeling results are on the CD-ROM. It will be called the NMWW CDROM in the rest of the book. The objective of this book is to describe the numerical methods for modeling water waves that have been developed primarily at the Los Alamos National Laboratory over the last four decades and to describe some examples of the applications of these methods. Some of the applications of the numerical modeling methods were performed while the author was working at the Joint Institute for Marine and Atmospheric Research at the University of Hawaii, some as Mader Consulting Co. research projects, and the rest as a Fellow of the Los Alamos National Laboratory. Although the two- and three-dimensional numerical methods for modeling water waves had been available in the 1980’s for several decades, they had seldom been used. A major obstacle to their use was the need for access to large and expensive computers. By the 1980’s, inexpensive personal computers were adequate for many applications of these numerical methods. In this book, the basic fluid dynamics associated with water waves are described. The common water wave theories are reviewed in Chapter 1. A computer code called WAVE for personal computers that calculates the wave properties for Airy, third-order Stokes, and Laitone solitary gravity waves is available on the NMWW CD-ROM. The incompressible fluid dynamics model used for shallow water, long waves is described in Chapter 2. A computer code for personal computers using the shallow water model is called SWAN and is available on the NMWW CD-ROM. The SWAN code is used to model the 1946, 1960 and 1964 earthquake generated Hilo, Hawaii tsunamis, the 1964 Crescent City, California tsunami and the 1994 underwater landslide generated Skagway, Alaska tsunami. The two-dimensional incompressible Navier-Stokes model used for solving water wave problems is described in Chapter 3. A computer code for personal computers using the two-dimensional Navier-Stokes model is called ZUNI and is available on the NMWW CD-ROM. The ZUNI code is used to model tsunami wave propagation and flooding. It is also used to model the effect of underwater barriers on tsunami waves. The three-dimensional incompressible Navier-Stokes model for solving water wave problems of any type is described in Chapter 4. The equations used in the computer code called SOLA-3D are described and the FORTRAN code is on the NMWW CD-ROM. The SOLA-3D code is used to model the 1975 Hawaiian tsunami and the 1994 Skagway tsunami by water surface cavities generated by underwater landslides.

It is surprising that most academic and government modelers of water waves have chosen to use shallow water or other incompressible models of limited validity for modeling water waves and not use the incompressible Navier-Stokes model since the first edition of this book was published. The severe limitations of the shallow water model are described in Chapter 5. The generation of water waves by volcanic explosions, conventional or nuclear explosions, projectile and asteroid impacts require the use of the compressible NavierStokes model. The numerical models and codes for solving such problems have recently become available as part of the Accelerated Strategic Computer Initiative program. A computer code for solving one, two and three dimensional compressible problems called SAGE, NOBEL or RAGE is described in Chapter 6 and some of its remarkable capabilities are presented. These include modeling of the KT Chicxulub asteroid impact and the modeling of the largest historical tsunami which occurred July 8, 1958, at Lituya Bay, Alaska. The modeling of the Lituya Bay impact landslide generated tsunami and the flooding to 520 meters altitude is described. A color videotape (VTC-86–4) lecture featuring computer generated films of many of the applications discussed in this book is available from the Los Alamos National Laboratory library. Several web sites have been established that contain additional information and computer movies of the problems described in this book. The major sites are http://www.mccohi.org and http://t14web.lanl.gov/Staff/clm/tsunami.mve/tsunami.html.

CONTENTS

CHAPTER 1— CHAPTER 2— CHAPTER 3— CHAPTER 4— CHAPTER 5—

WATER WAVE THEORY THE SHALLOW WATER MODEL THE TWO-DIMENSIONAL NAVIER-STOKES MODEL THE THREE-DIMENSIONAL NAVIER-STOKES MODEL EVALUATION OF INCOMPRESSIBLE MODELS FOR MODELING WATER WAVES CHAPTER 6— MODELING WAVES USING COMPRESSIBLE MODELS CD-ROM CONTENTS AUTHOR INDEX SUBJECT INDEX

1 24 80 111 134 165 216 223 227

1 WATER WAVE THEORY 1A. The Equations of Fluid Dynamics The equations of fluid dynamics are representations of the laws of conservation of mass, momentum, and energy as applied to a fluid. The general equations are often called the Navier-Stokes equations. In addition, we need an equation of state to describe the properties of the fluids. The equation of state expresses the fact that the pressure is everywhere a function of the density and the energy per unit of mass. For the water wave theories we will describe in this chapter, the fluid will be considered to be incompressible and the initial density of each element of fluid remains constant. The numerical models of fluid dynamics for compressible fluids are described in the Los Alamos monograph Numerical Modeling of Detonations 1 , and in Chapter 6. The equations of fluid dynamics are described by Landau and Lifshitz 2 as a branch of theoretical physics in the classic volume 6 of their Course of Theoretical Physics. The equations of fluid dynamics are treated as a branch of applied mathematics in Batchelor’s 3 monograph entitled An Introduction to Fluid Dynamics. A classic description of the equations of fluid dynamics as a branch of oceanography is Blair Kinsman’s 4 textbook Wind Waves. Any of these exhaustive treatments of the equations of fluid dynamics are recommended if the reader wishes to reaffirm his faith in the efficacy of Newtonian mechanics and the fluid continuum. In this monograph, we will first derive the common water wave description starting with the three-dimensional equations of fluid dynamics for viscous, compressible flow. The Nomenclature gx,gy,gz

acceleration due to gravity in x, y, or z direction

I

internal energy

P

pressure

V

specific volume

Ux

particle velocity in x direction

Uy

particle velocity in y direction

Uz

particle velocity in z direction

ρ

density

q

viscosity

Numerical modeling of water waves qx, qy, qz

viscosity deviators

Sx, Sy, Sz

total deviators

λ, µ

real viscosity coefficients

C

wave velocity

H

maximum wave height

h

wave height

L

wave length

T

wave period

f

frequency

n

wave number

D

depth

2

The three-dimensional partial differential equations for viscous, compressible flow: Eulerian Conservation of Mass

or

1.1

because

and

Water wave theory

3

Viscosity with No Shearing Forces

Sx=qx−q. Sy=qy−q. Sz=qz−q. Eulerian Conservation of Momentum

σ=δijP−Sij.

(1.2)

Numerical modeling of water waves

4

Eulerian Conservation of Energy

where

and

(1.3)

If the flow is incompressible so ρ=ρ0, equation 1.1 becomes

(1.4) equation 1.2 becomes, without viscosity,

(1.5a)

(1.5b)

Water wave theory

5

(1.5c)

and equation 1.3 becomes, without viscosity,

(1.6)

Viscosity terms may be included in equation 1.5 as follows. If shearing forces are included and the fluid is incompressible, then

q=0.0, so equations 1.5 become

(1.7a)

Numerical modeling of water waves

6

(1.7b)

and

(1.7c)

1B. Water Wave Descriptions A useful description of water waves has been given by Bernard Le Mehaute in reference 5. Parts of that report are included in the following description of water waves. The theories for unsteady free surface flow subjected to gravitational forces include motions called water waves. They are also called gravity waves. Other gravity waves include motions in other fluids such as atmospheric motions. A great variety of water waves exists. Water wave motions include storm waves generated by wind in the oceans, flood waves in rivers, seiche or long-period oscillations in harbor basins, tidal bores or moving hydraulic jumps in estuaries, waves generated by a moving ship in a channel, tsunami waves generated by earthquakes, and waves generated by explosions near or under the water. Mathematically, a general solution does not exist for gravity waves and approximations must be made for even simple waves. One of the important problems in water wave theory is to establish the limits of validity of the various solutions that are due to the simplifying assumptions. The mathematical treatments of the water wave motions use all the mathematical physics dealing with linear and nonlinear problems. The main difficulty in the study of water wave motion is that the free surface boundary is unknown. Water wave motions are so varied and complex that any attempt at classification may be misleading. Any definition corresponds to idealized situations that never occur. For example, a purely two-dimensional motion never exists. It is a convenient mathematical concept that is physically best approached in a tank with parallel walls. Boundary layer effects and transverse components still exist although they are small and may be neglected in many applications of the theory. Essentially, two kinds of water waves exist, oscillatory waves and translatory waves. In an oscillatory wave, the transportation of fluid or mass transportation does not occur. The wave motion is then analogous to the transverse oscillation of a rope. A translatory wave involves a transport of fluid in the direction in which the wave travels. For example, a moving hydraulic jump such as a tidal bore is a translatory wave.

Water wave theory

7

An oscillatory wave can be progressive or standing. Consider a disturbance H(x, t) such as a free surface elevation traveling along the OX axis at a velocity C. The characteristics of a progressive wave remain identical for an observer traveling at the same speed and in the same direction as the wave (Figure 1.1). Where h can be expressed as a function of (x−Ct) instead of (x, t), a “steady-state” profile is obtained. h(x−Ct) is the general expression for a progressive wave of steady-state profile traveling in the positive OX direction at a constant velocity C. Where the progressive wave is moving in the opposite direction, its mathematical form is expressed as a function of (x+Ct). The definition of wave velocity C for a non-steady-state profile is unphysical, because each “wave element” travels at its own speed, thereby causing wave deformation.

Fig. 1.1. A Progressive Wave.

The simplest case of a progressive wave is the wave defined by a sine or cosine curve such as

or

Such a wave is called a harmonic wave, where H/2 is the amplitude and H the wave height. The distance between the wave crests is the wave length L, and L=CT, where T is the wave period. The wave number n=2π/L is the number of wave lengths per cycle. The frequency is f=2π/T. Hence, the previous equations can be written

or

Numerical modeling of water waves

8

A standing or stationary wave is characterized by its mathematical description as a product of two independent functions of time and distance, such as

or more generally, F=F1(x)·F2(x). A standing wave can be considered as the superposition of two waves of the same amplitude and period traveling in opposite directions. Where the convective terms are negligible, the standing wave motion is defined by a linear addition of the equations for the two progressive waves.

Fig. 1.2. A Standing Wave.

A standing wave generated by an incident wind wave in relatively shallow water is called a seiche. A seiche is a standing oscillation of long period encountered in lakes and harbor basins. The amplitude at the node is zero, and at the antinode, it is H, as shown in Figure 1.2. Two waves of the same period but different amplitudes traveling in opposite directions can be defined linearly by the summation of A sin(x−Ct)+B sin(x+Ct). It can also be considered as the superposition of a progressive wave with a standing wave and is encountered in front of an obstacle that causes a partial reflection of a wave. The amplitude at the node is N=A+B and at the antinode AN =A−B. The direct measurement of N and AN yields:

and

and the reflection coefficient

Water wave theory

9

. In a translatory wave, water is transported in the direction of the wave travel. Some examples are 1) A tidal bore or moving hydraulic jump 2) Waves generated by the breaking of a dam 3) A surge on a dry bed 4) An undulated, moving hydraulic jump 5) Solitary waves 6) Flood waves in rivers In an Eulerian system of coordinates, a surface wave problem generally involves three unknowns: the free surface elevation (or total water depth), the pressure (generally known at the free surface), and the particle velocity. Since a general method of solution is impossible, a number of simplifying assumptions have been made that apply to a succession of particular cases with varying accuracy. In general, the method of solution used depends on the relative importance of the convective terms. However, instead of dealing with these terms directly, it is more convenient to relate them to more accessible parameters. Three characteristic parameters are used. 1) A typical value of the free surface elevation, such as the wave height H 2) A typical horizontal length such as the wave length L 3) The wave depth D Although the relationships between the convective terms and these three parameters are not simple, their relative values are of considerable help in classifying the water wave theories. As the free surface elevation decreases, the particle velocity also decreases. Thus, when the wave height, H, tends to zero, the convective term, which is related to the square of the particle velocity, is an infinitesimal. Consequently, the convective terms can be neglected and the theory can be linearized. Thus, the three possible parameters to be considered are

The relative importance of the convective terms increases as the values of these three parameters increase. In deep water (small H/D, and small L/D), the most significant parameter is H/L, which is called the wave steepness. In shallow water, the most significant parameter is H/D, which is called the relative height. In intermediate water depth, a significant parameter that also covers the three cases is , and it is called the Ursell parameter. Depending on the problem under consideration and the range of values of the parameters H/L, H/D, and L/D, four mathematical approaches are used. 1) Linearization 2) Power series

Numerical modeling of water waves

10

3) Numerical methods 4) Random functions The simplest cases of water wave theories are the linear wave theories, in which case the convective terms are neglected completely. These theories are valid when H/L, H/D, and L/D are small, i.e., for waves of small amplitude and small wave length in deep water. For the first reason, they are called the “small amplitude wave theory.” It is the infinitesimal wave approximation. The linearization of the basic equation is so suitable to mathematical manipulation that the linear wave theories cover an extreme variety of water wave motions. For example, some phenomena that can be subjected to linear mathematical treatment include the phenomena of wave diffraction and of the waves generated by a moving ship. The solution can be found as a power series in terms of a small quantity by comparison with the other dimensions. This small quantity is H/L for small L/D because in deep water, the most significant parameter is H/L. It is H/D for large L/D because in shallow water the most significant parameter is H/D. In the first case (development in terms of H/L), the first term of the power series is given by application of the linear theory. In the second case, the first term of the series is already a solution of nonlinear equations. The calculation of the successive terms of the series is so cumbersome that these methods are used in a very small number of cases. The most typical case is the progressive periodic wave. In this case, the solution is assumed to be a priori that of a steady-state profile, i.e., a function such as F=f(x−Ct), where C is a constant equal to the wave velocity or phase velocity. The simplification introduced by such an assumption is due to the fact that

and

such that

In such a way, the time derivatives can be eliminated easily and replaced by a space derivative. Typical examples of such treatments include 1) Power series of H/L or the Stokes waves, valid in deep water. The first term of the series is obtained from linear equations and corresponds to the infinitesimal wave

Water wave theory

11

approximation. 2) Power series of H/D: the cnoidal wave or the solitary wave, valid for shallow water. The first terms of the series are obtained as a steady-state solution of already nonlinear equations, but correspond to the shallow water approximation. However, a steady-state profile may not exist as a solution, in which case the method to be used is often a numerical method of calculation where the differentials are replaced by finite differences. This occurs for large values of H/D and L/D, which correspond to the fact that the nonlinear terms such as

are relatively large by comparison with the

linear terms such as . This is the case for long waves in very shallow water. Of course, a numerical method of calculation can be used for solving a linearized system of equations. For example, the relaxation method is used for studying small wave agitation in a basin. Also, an analytical solution of a nonlinear system of equations can be found in some particular cases. Hence, these three methods and the range of application that has been given indicate more of a trend than a general rule. The three previous methods aim at a fully deterministic solution of the water wave problem. Other descriptions of a sea state generally involve the use of random functions. The mathematical operations (such as harmonic analysis) generally imply that the water waves obey linear laws, which are the necessary requirements for assuming the principle of superposition. Such a method loses its validity for describing the sea state in very shallow water (large values of H/D and L/D). In hydrodynamics, the water wave theories are generally classified into two families. They are the “small amplitude wave theories” and the “long wave theories.” The small amplitude wave theories are the linearized theories of the first categories of power series, i.e., the power series in terms of H/L. The long wave theories use the numerical method of solution for the nonlinear long wave equations. These two families include a number of variations and some intermediate cases with some of the characteristics of both families. For example, the cnoidal wave and the solitary wave are considered as particular cases (steady state) of the long wave theories, because they are nonlinear shallow water waves.

1C. Airy Waves The linear theory of Airy in Eulerian coordinates gives the essential characteristics of the wave pattern assuming the wave height is infinitesimal. The free surface is sinusoidal as shown in Figure 1.3, particle paths are elliptic and follow a closed orbit (zero mass transport), and lines of equipressure are also sinusoidal. The terms in (H/L)2 are neglected. The long wave or shallow water theory is the same as the theory of Airy, where it is assumed that D/L is small. As a consequence, the formulas are simplified considerably. The pressure is hydrostatic and the horizontal velocity distribution is uniform. The Airy wave equations used in the WAVE code available on the NMWW CD-ROM are given below.

Numerical modeling of water waves

Fig. 1.3. The Airy Wave.

P=−ρgy+∆P.

where h is wave height from the bottom.

12

Water wave theory

13

Special cases of the Airy wave theory are the deep water (D/L> 0.5) and shallow water (D/L3 even sooner than in the case of the third-order theory, i.e., for larger values of D/L. Consequently, the fifth-order wave theory is less valid than the third-order wave theory for small values of D/L and cannot be used when D/L0.25. In no case should ALPHA exceed 1.0.

4B. Application to Tsunami Wave Formation The first application of the incompressible Navier-Stokes model to tsunami wave formation was by Garcia4 using a two-dimensional arbitrary boundary marker and cell technique to study tsunamis in the vicinity of their source. The method was applied to the Mendocino Escarpment for hypothetical ocean floor displacements of 10 meters in a few seconds to a minute, which might result from a major earthquake on the San Andreas fault. The wave was followed numerically for the first 200 sec. A single hump was first generated that split into two crests moving in opposite directions at a speed slightly less than shallow water wave speed. The waves were slightly dispersive and had a crest elevation above mean water level of about 1 meter and a period of about 1.5 minute. The period of major tsunamis measured in the Pacific Ocean is generally from 10 to 30 minutes. The periods are estimated from tide gauge records assuming that the period remains nearly constant during the shoaling of the wave from the deep ocean. The wave height of major tsunamis is generally given at 1±0.5 meter and is estimated from the tide gauge records assuming various approximate models for the wave height growth and extrapolating back to the deep ocean. The first application of the three-dimensional incompressible Navier-Stokes model to tsunami formation was by Mader, Tangora, and Nichols,5 and by Mader.6 The Hawaiian tsunami of November 29, 1975, has been investigated by Loomis. He described the observed runup heights in reference 7 and a numerical study of the tsunami source in references 8 and 9. The tsunami was generated by an earthquake with a magnitude of 7.2 on the Richter scale, near the Hawaii Volcanoes National Park. Near the source, the first wave was smaller than the second. Coincident with the earthquake was considerable subsidence (up to 3 meters) of the shoreline. Loomis, in reference 8, examined a model of the southeastern coast of Hawaii. The bottom slopes seaward at a ratio of 1:15 until it reaches a constant depth of 6000 meters. The sources examined by Loomis included both initial uplifts and depressions, and he reported that such source motions would not generate the essential features of the tsunami; that is, a second wave larger than the first. In addition to the sources studied by Loomis, a landslide source was modeled in

The three-dimensional navier-stokes model

123

reference 6. The landslide model has been evaluated by Cox in reference 10. He concluded that a landslide could not be distinguished from strictly tectonic displacement by the comparison of arrival times and travel times. The SWAN code described in Chapter 2 was used to solve the shallow water, long wave equations and to examine the tsunami generation problem. The SWAN results confirmed Loomis’s calculated results reported in reference 8. The SOLA-3D code, which solves the three-dimensional incompressible Navier-Stokes equations, was also used to model the tsunami. The Calculated Shallow Water Wave Results The model was identical to that described by Loomis in reference 7. A 40 by 69 rectangular region of 207 kilometers along the coast and 120 kilometers seaward is described using a mesh of 3 kilometers by 3 kilometers. The bottom slopes at a ratio of 1:15 until it reaches a depth of 6 kilometers. The source is a bottom slope subsidence programmed to vary with time. The source is 30 kilometers wide, of which half is included in the calculation and is separated from the other half by a reflective boundary, as shown in Figure 4.1. The calculated wave profile is shown at various locations along the shoreline as a function of time in Figure 4.2 for a source displacement of 3 meters, and in Figure 4.3 for a source displacement of 1 meter, followed by an additional 2 meters displacement 10 minutes later. Surface profiles are shown in Figure 4.4.

Fig. 4.1. Sketch of model used to numerically simulate the tsunami generation.

Numerical modeling of water waves

124

Fig. 4.2. Shoreline wave heights for shallow water wave model resulting from an initial source displacement of 3 meters.

The observed larger second wave can be reproduced by a source that has its displacement change with time. Such a possibility was suggested by Ando, who suggested that the earthquake was a rather

Fig. 4.3. Shoreline wave heights for shallow water wave model resulting from an initial source displacement of 1 meter, followed by an additional 2 meters displacement 10 minutes later.

The three-dimensional navier-stokes model

125

slow rupture lasting 100 sec; however, the displacement change required by the model is of longer duration and includes two fast ruptures. Another source investigated was an undersea landslide. The landslide ocean bottom profile assumed the bottom dropped 3 meters at the shoreline and slid to form the profile shown in Figure 4.5. Landslides are observed to pile up at the bottom one-third of their run. This gives the surface wave profile shown in Figure 4.6. The calculations were performed on the University of Hawaii Harris computer using the Hawaii version of the SWAN code. The shoreline wave heights at various times for the shallow water wave model with the initial water surface displacement of Figure 4.6 are shown in Figure 4.7. Although the wave heights are consistent with the observed behavior of the tsunami, we must check the results with the SOLA code since it was demonstrated in Chapter 3 that the shallow water model is inadequate to describe the waves generated from surface deformations of the water surface.

Numerical modeling of water waves

126

Fig. 4.4. Surface profiles for linear shallow water wave model as a function of time, resulting from an initial shore displacement of 1 meter, followed by an additional 2 meters displacement 10 minutes later.

The three-dimensional navier-stokes model

127

Fig. 4.5. Sketch of the final ocean bottom profile after a landslide for the source region.

Fig. 4.6. Sketch of height of water surface after a landslide on the ocean bottom.

The Calculated Navier-Stokes Results The geometry of the model used to calculate the tsunami is shown in Figure 4.1. The mesh used in the calculation had 20 cells in the x direction, 25 cells in the y direction, and 18 cells in the z direction. The 20 cells in the x direction were 6 kilometers wide. The 18 cells in the z direction starting at the ocean floor were 100 meters high for the first two cells, and 400 meters thereafter. The water depth was 6000 meters and the surface was located at the center of cell 17. The 25 cells in the y direction starting at the source were 3.0 kilometers for the first 5 cells that described the source (15 kilometers wide). The remaining cell widths were 5.75, 6.16, 6.56, 6.97, 7.37, 7.78, 8.18, 8.59, 9.0, 9.4, 9.8, 10.2, 10.6, 11.0, 11.4, 11.8, 12.2, 12.6, 13.0, and 13.5 kilometers, for a total of 206.8 kilometers.

Numerical modeling of water waves

128

Fig. 4.7. Shoreline wave heights for a shallow water wave model resulting from the initial water surface displacement shown in Figure 4.6.

The viscosity coefficient was 2.0 g-sec−1 meter−1 (0.02 poise). The gravity constant gz was −9.8 meters sec−2, and gx and gy were 0.0. The time step for the calculation was 5 sec. The tsunami source was modeled by a 3 meters deep depression or elevation of 90 by 15 kilometers of the water surface, as shown in Figure 4.1. The calculated wave profiles are shown at various locations along the shoreline as a function of time in Figures 4.3, for a source of 3 meters depression of the water surface, The calculation for a 3 meters uplift gave mirror images of the profiles for the 3 meters depression of the water surface. Surface profiles are shown in Figure 4.9. The calculated wave profiles are shown in Figure 4.10 at various locations along the shoreline as a function of time for a landslide source. The observed tsunami wave profile of the 1975 Hawaiian tsunami, near the source of the second wave larger than the first, is not reproduced by a landslide source in an incompressible three-dimensional NavierStokes calculation. This contrasts with results obtained using the shallow water model.

The three-dimensional navier-stokes model

129

Fig. 4.8. Shoreline wave heights for the three-dimensional incompressible Navier-Stokes model as a function of time resulting from an initial source of a 3 meters depression of the water surface.

Numerical modeling of water waves

130

Fig. 4.9. Surface profiles for three-dimensional incompressible NavierStokes model as a function of time resulting from an initial source of a 3 meters depression of the water surface.

The three-dimensional navier-stokes model

131

Fig. 4.10. Shoreline wave heights for a three-dimensional incompressible Navier-Stokes model calculation with the landslide source shown in Figure 4.6.

A source with a 3 meters uplift of the water surface was consistent with the observed tsunami wave profile. These calculations do not support a landslide source for the 1975 Hawaii tsunami. Conclusions The differences between the shallow water and full Navier-Stokes calculations are that the water waves formed in the full Navier-Stokes calculations are deep water waves that move slower than the shallow water waves formed in the shallow water calculations. The nature of the surface collapse is also different, with the collapse occurring throughout the source region in the Navier-Stokes calculations and mostly at the sides in the shallow water calculations. The observed tsunami wave profile of the 1975 Hawaii tsunami, of the second wave being larger than the first wave near the source, is not reproduced by a landslide source, but it is reproduced by a simple uplift or drop of the water surface over the source area. The shallow water approximation is not appropriate for studying waves generated from surface deformations that are small, relative to the water depth. In the next chapter we will investigate the limitations of the shallow water approximation for modeling tsunami waves.

Numerical modeling of water waves

132

4C. The 1994 Skagway Tsunami The Skagway tsunami of November 3, 1994 is described in Section 2G of Chapter 2. A tsunami wave with a period of about 3 minutes and a maximum height of 30 feet occurred in the Taiya Inlet at Skagway, Alaska. The tsunami was a result of a massive underwater landslide of sediment deposited by the Skagway River. The tsunami was modeled in Chapter 2 using the shallow water SWAN code. The Skagway tsunami was also modeled using the SOLA-3D code. A single slide model for the Skagway tsunami, that closely approximated the more realistic 3-slide Landslide model, was used for the three-dimensional SOLA-3D calculation. A SWAN calculation for the same single slide model was also performed for comparision. The SOLA-3D FORTRAN code and input file used to model the Skagway tsunami and the computer movies are on the NMWW CD-ROM in the directory /TSUNAML.MVE/SOLA/SKAGWAY. Since the wave length and period of the resulting tsunami is large compared to the shallow depth of the inlet, the tsunami wave is adequately approximated by a shallow water wave. The tsunami wave profiles generated by the SOLA-3D code and the SWAN code were similar in amplitude and period.

References 1. C.W.Hirt, B.D.Nichols, and N.C.Romero, “SOLA—A Numerical Solution Algorithm for Transient Fluid Flows,” Los Alamos National Laboratory report LA-5852 (1975). 2. C.W.Hirt, B.D.Nichols, and N.C.Romero, “SOLA—A Numerical Solution Algorithm for Transient Fluid Flows—Addendum,” Los Alamos National Laboratory report LA5852, Add. (1976). 3. C.W.Hirt, B.D.Nichols, and L.R.Stein, “SOLA-3D—A Numerical Solution Algorithm for Transient Three-Dimensional Flows,” Los Alamos National Laboratory, Group T3, unpublished internal report (1985). 4. W.J.Garcia, “A Study of Water Waves Generated by Tectonic Displacements,” College of Engineering, University of California at Berkeley report HEL 16–9 (1972). 5. Charles L.Mader, Robert E.Tangora, and B.D.Nichols, “A Model of the 1975 Hawaii Tsunami,” Natural Science of Hazards—The International Journal of the Tsunami Society, pp.C1-C8 (1982). 6. Charles L. Mader, “A Landslide Model for the 1975 Hawaii Tsunami,” Science of Tsunami Hazards, Vol.2, No.2, pp.71–77 (1984). 7. Harold G.Loomis, “The Tsunami of November 29, 1975, in Hawaii,” Hawaii Institute of Geophysics report HIG-75–21 (1975). 8. Harold G.Loomis, “On Defining the Source of the 1975 Tsunami in Hawaii,” JIMAR report to Nuclear Regulatory Commission (1978). 9. M.A.Sklarz, L.Q.Spielvogel, and H.G.Loomis, “Numerical Simulation of the November 29, 1975, Island of Hawaii Tsunami by the Finite Element Method,” Journal of Physical Oceanography, Vol.9, No.5, pp.1022–1031 (1979).

The three-dimensional navier-stokes model

133

10. Doak C.Cox, “Source of the Tsunami Associated with the Kalapana (Hawaii) Earthquake of November 1975,” Hawaii Institute of Geophysics report HIG-80–8 (1980).

5 EVALUATION OF INCOMPRESSIBLE MODELS FOR MODELING WATER WAVES The shallow water code SWAN and the incompressible Navier-Stokes code ZUNI have been used to model the development of a tsunami wave from an initial sea surface displacement, the propagation of a tsunami wave, and the resulting shoaling and flooding. The generation and propagation of a tsunami wave has also been modeled using a linear gravity model. The discrepancies between the shallow water and the more realistic Navier-Stokes and linear gravity models are quantified. These studies were performed by the author at the University of Hawaii Joint Institute for Marine and Atmospheric Research in the late 1980s and were published in references 1, 2 and 3.

5A. Tsunami Wave Generation Initial Surface Displacement Study—SWAN Model A 1 meter high Airy half-wave surface displacement with a width of 45 kilometers in 4550 meters deep water was studied. This approximates within cell resolution a Gaussian wave with a Gaussian break width of 10 kilometers. Calculations were performed using a 250 meters wide mesh of 400 by 4 cells and at 0.5 sec intervals. The wave height in meters as a function of distance is shown in Figure 5.1 at various times. The initial surface displacement separates into two shallow water waves with a height of 0.5 meter and length of 45 kilometers. At maximum height the vertical velocity at the center of the wave is 0.0232 meter/sec. An Airy wave with a 0.5 meter half width, 90 kilometers wave length in 4450 meters of water has a vertical velocity of 0.0235 to 0.0226 meter/sec at maximum height. The wave speed is 208.36 and the group velocity is 202.89 meters/sec. The shallow water approximation used in the SWAN code of constant vertical velocity introduces errors of about 5 percent.

Evaluation of incompressible models for modeling water waves

135

Fig. 5.1. The wave height in meters as a function of distance at various times for a 1 meter high surface displacement with a width of 45 kilometers in 4550 meters deep water for the nonlinear shallow water model using the SWAN code.

The nonlinear feature of the shallow water model included in the SWAN code results in a small trailing wave with an amplitude of less than 0.01 meter.

Numerical modeling of water waves

136

Initial Surface Displacement Study—ZUNI Model An approximately 1 meter high (within cell resolution) Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) in 4550 meter deep water was studied. The calculations for this geometry were performed with 15 cells in the Y or depth direction and 68 cells in the X or distance direction. The cells were rectangles 450 meters high in the Y direction and 4000 meters long in the X direction. The time increment was 3 sec. The water level was placed at 4550 meters or 50 meters up into the eleventh cell. The viscosity coefficient used was 0.02 poise, a value representative of the actual viscosity for water. The viscosity did not significantly affect the results. The wave height in meters as a function of distance is shown in Figure 5.2 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and length of 45 kilometers. At 100 sec two waves have formed. The vertical velocity at the top center of the wave is 0.0231 meter/sec and 0.2158 meter/sec near the bottom of the wave. An Airy wave with a 0.5 meter half width, 90 kilometers wave length in 4450 meters of water has a center vertical velocity of 0.0235 to 0.0226 meter/sec. The Airy wave speed is 208.36 and the group velocity is 202.89 meters/sec. As the wave propagates, the wave height decreases, the slope of the front of the wave becomes less, and small waves form behind the main wave. After the wave has propagated for 500 sec, the wave height has decreased to 0.445 meter, and the vertical velocity at the center of the wave has decreased to 0.0217 meter/sec near the peak of the wave to 0.0203 at the bottom of the wave. A train of waves has developed behind the main waves with maximum negative amplitude of 0.1 meter and positive amplitude of 0.05 meter.

Evaluation of incompressible models for modeling water waves

137

Fig. 5.2. The wave height in meters as a function of distance at various times for a 1 meter high surface displacement with a width of 45 kilometers in 4550 meters deep water for the Navier-Stokes water wave model using the ZUNI code.

Numerical modeling of water waves

138

Initial Surface Displacement Study—LGW Model The analytical methods for solving the linear gravity model were described by Professor George Carrier4,5. The two-dimensional linear gravity wave with a Gaussian, a square wave, or a time dependent Gaussian displacement is solved using Fourier transforms by the LGW code for any time of interest. The wave description is obtained for any uniform depth, density, gravity and Gaussian break width or square wave half width. The wave height, vertical and horizontial velocities, and pressure are calculated at any depths desired. The distance scale is chosen to be a tenth of a Gaussian break width so that a wave half-width is described by about 50 space increments. Since the model is symmetrical about the center of the initial displacement, only half of the wave profile is calculated. The LGW code is included on the NMWW CD-ROM. A 1 meter high Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) in 4550 meters deep water was studied. Since the model is symmetrical about the center of the initial displacement, only half of the wave is modeled analytically. The analytical linear gravity wave calculations were performed using 256 points with a thickness of 1000 meters. The wave profile was calculated for each time selected. The wave height in meters as a function of distance is shown in Figure 5.3 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and length of 45 kilometers. At 100 sec the two waves have formed and the vertical velocity at the center of the wave is 0.0238 meter/sec at the top of the wave to 0.217 meter/sec near the bottom of the wave. An Airy wave with a 0.5 meter half width, 90 kilometers wave length in 4450 meters of water has a vertical velocity of 0.0235 to 0.0226 meter/sec at the center of the wave. The Airy wave speed is 208.36 and the group velocity is 202.89 meters/sec.

Evaluation of incompressible models for modeling water waves

139

Fig. 5.3. The wave height in meters as a function of distance at various times for a 1 meter high Gaussian displacement with a Gaussian break length of 10 kilometers (equivalent to an Airy half wave length of 45 kilometers) in 4550 meters deep water for the linear gravity water wave model using the LGW code.

As the wave propagates, the wave height decreases, the slope of the front of the wave

Numerical modeling of water waves

140

becomes less, and small waves form behind the main wave. After the wave has propagated for 500 sec, the wave height has decreased to 0.4567 meter, and the vertical velocity at the center of the wave has decreased to 0.0216 meter/sec near the peak of the wave to 0.02015 at the bottom of the wave. A train of waves has developed behind the main waves with maximum negative amplitude of 0.1 meter and positive amplitude of 0.03 meter. Comparisions In Table 5.1 are listed for comparision the peak wave height and the associated vertical velocity at the top and bottom of the wave for various times for the nonlinear shallow water, the full Navier-Stokes and the linear gravity wave models. The initial surface displacement was chosen to be as similar as possible for each of the models. The initial surface displacement was chosen to be characteristic of a typical earthquake generated tsunami wave at its source. The surface displacement was one meter high spread over 45 kilometers in 4550 meters of water for a width-to-depth ratio of about 10. The initial surface displacement resulted in two waves- traveling in opposite directions with half the initial height. The nonlinear aspect of the shallow water wave appeared as very small trailing waves. The dispersion aspects of the linear gravity wave and the Navier-Stokes waves resulted in the front of the wave spreading, along with a corresponding decrease in the slope of the front of the wave and the generation of a following train of waves with amplitudes of a tenth or less of the peak wave height after the wave had traveled 105 kilometers in 500 sec (210 meters/sec). The peak wave height of the Navier-Stokes wave was lower than for the linear gravity wave. However, this difference was an order of magnitude less than the difference between either wave and the shallow water wave.

Evaluation of incompressible models for modeling water waves

141

TABLE 5.1 A 1 m High, 48 Kilometers Wide Source

Time sec

Height meters

Bottom Velocity meters/sec

Top Velocity meters/sec

Navier-Stokes Wave -ZUNI Code 0.

1.000

50.

0.0000

0.0000

0.0195

0.0202

100.

0.493

0.0216

0.0231

200.

0.475

0.0216

0.0233

300.

0.464

0.0211

0.0227

400.

0.459

0.0208

0.0224

500.

0.445

0.0203

0.0217

Linear Gravity Wave -LGW Code 0.

1.000

50.

0.0000

0.0000

0.0192

0.0220

100.

0.496

0.0217

0.0238

200.

0.489

0.0214

0.0234

300.

0.478

0.0211

0.0228

400.

0.467

0.0206

0.0222

500.

0.457

0.0201

0.0216

Shallow Water Wave -SWAN Code 0.

1.000

50.

0.0000

0.0000

0.0231

0.0231

100.

0.500

0.0232

0.0232

200.

0.499

0.0232

0.0232

300.

0.500

0.0232

0.0232

400.

0.499

0.0232

0.0232

500.

0.500

0.0232

0.0232

The differences described will become less with increasing wave length-to-depth ratios and become greater for decreasing wave length-to-depth ratios. To quantify this effect we used the linear gravity wave model to examine the propagation of waves of various widths-to-depths ratios.

Numerical modeling of water waves

142

Initial Surface Displacement Width Study-LGW Model The linear gravity wave model was used to investigate the waves formed from initial surface displacements of width-to-depth ratios of 40 to 0.5 (Gaussian break widths of 40 to 0.5 kilometers) in 4550 meters deep water. The wave height after the wave had travel ten times its initial width for depth ratios of 40, 20, 10, 5, 1, 0.5 are shown in Figure 5.4. The Airy half-wave length and half-wave period equivalents for a break width of 5 is 23 kilometers and about 1.8 minutes, a break width of 10 is 46 kilometers and 3.5 minutes, a break width of 20 is 90 kilometers and 7 minutes, and for a break width of 40 is 180 kilometers and 14 minutes. The use of nonlinear shallow water models to describe tsunami waves generated from earthquake generated surface displacements is adequate for tsunamis generated by surface displacements that are at least ten times wider than the depth. The nonlinear shallow water wave becomes less realistic the further it travels from its source and the smaller the width-to-depth ratio. Non-linear shallow water models are adequate for large wave length and long period tsunamis such as generated by the 1964 Alaskan or the 1960 Chile earthquakes where the periods were about 30 minutes and wave length-to-depth ratio in the deep ocean was greater than 80. Tsunami waves generated by earthquakes with small areas and periods of a few minutes will not be realistically described using shallow water models. Either the linear gravity wave model or better the Navier-Stokes model should be used for accurate modeling of tsunamis with small (less than 10) widthto-depth ratios.

Evaluation of incompressible models for modeling water waves

143

Fig. 5.4. The wave height in meters as a function of distance after the wave has traveled ten times its initial width for a 1 meter high Gaussian displacement with a Gaussian break length of 40, 20, 10, 5, 1, and 0.5 kilometers in 4450 meters deep water for the linear gravity wave model using the LGW code.

Generation Conclusions The SWAN code solves the nonlinear, shallow water, long wave equations including the effects of friction and flooding. Both the SWAN and ZUNI code which solves the incompressible Navier-Stokes equations were used to study the development of a tsunami wave from an initial sea surface displacement. The development of a tsunami wave was also modeled using two-dimensional linear gravity wave theory. If the initial displacement is approximately Gaussian and the wave length is very long compared to the depth, similar tsunami waves formed for all three methods. Two tsunami

Numerical modeling of water waves

144

waves traveling in opposite directions formed that were the sum of the original surface displacement. However, dispersion effects resulted in Navier-Stokes and linear gravity waves with decreasing front slopes and amplitudes, and followed by a train of small waves. The shallow water wave has a constant vertical velocity while the Navier-Stokes and linear gravity waves have the more realistic variable vertical velocity. For long wave length tsunamis the constant vertical velocity closely reproduces the velocity characteristics of Navier-Stokes and linear gravity waves, which slowly decrease with depth. With decreasing periods and wave lengths, the discrepancy between the shallow water and the Navier-Stokes and linear gravity waves formed from initial sea surface displacement increases. These studies permit us to determine the parameteric region where the shallow water model can be useful. A summary of that region is presented at the end of this chapter. Sources that involve compressible flow can not be solved with the incompressible models evaluated in this chapter; however they can be modeled using the compressible models described in Chapter 6. Computer movies of the source calculations are in the NMWW CD-ROM TSUNAMLMVE/SOURCE.MVE directory.

5B. Tsunami Wave Propagation The shallow water SWAN code and the incompressible Navier-Stokes ZUNI codes were used to study the propagation of a tsunami wave from an initial sea surface displacement such as the propagation of a tsunami wave from Alaska or the U.S. West Coast to Hawaii. The propagation of a tsunami wave was also modeled using two-dimensional linear wave theory. Tsunami Propagation Study—SWAN Model A 1 meter high Airy half-wave surface displacement with a width of 45 kilometers in 4550 meters deep water was studied. This approximates within cell resolution a Gaussian wave with a Gaussian break width of 10 kilometers. Calculations were performed using a 250 meters wide mesh of 32,000 by 4 cells and at 0.5 sec intervals. The wave height in meters as a function of distance is shown in Figure 5.5 at half hour intervals.

Evaluation of incompressible models for modeling water waves

145

Fig. 5.5. The wave height in meters as a function of distance in meters at half hour intervals for a 1 meter high surface displacement with a half width of 45 kilometers in 4550 meters deep water for the nonlinear shallow water model using the SWAN code.

The initial surface displacement separates into two shallow water waves with a height of 0.5 meter and width of 45 kilometers. The nonlinear feature of the shallow water model included in the SWAN code results in a small trailing wave with an amplitude of less than 0.01 meter. The peak wave amplitude remains nearly constant as the wave propagates. Tsunami Propagation Study—ZUNI Model A 1 meter high approximately (within cell resolution) Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) in 4550 meters deep water was studied as it traveled for 3 hours. The calculations for this geometry were performed with 15 cells in the Y or depth direction and up to 13,600 cells in the X or distance direction. The calculation required 231, 200 cells and over 2.5 million mesh quantities. The cells were rectangles 450 meters high in the Y direction and 4000 meters long in the X direction. The time increment was 3 sec. The water level was placed at 4550 meters or 50 meters up into the eleventh cell. The viscosity coefficient used was 0.02 poise, a value representative of the actual viscosity for water. The viscosity did not significantly affect the results. The wave height in meters as a function of distance is shown in Figure 5.6 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and width of 45 kilometers. At 100 sec the two waves have a

Numerical modeling of water waves

146

vertical velocity at the top center of the wave of 0.0231 meter/sec. Near the bottom of the wave the velocity is 0.02016 meter/sec. An Airy wave with a 0.5 meter half-width, 90 kilometers wave length in 4550 meters of water has a vertical velocity of 0.0235 to 0.0226 meter/sec at the center of the wave. The Airy wave speed is 208.36 meters/sec and the group velocity is 202.89 meters/sec.

Fig. 5.6. The wave height in meters as a function of distance in kilometers at various times for a 1 meter high surface displacement with a width of 45 kilometers in 4550 meters deep water for the NavierStokes water wave model using the ZUNI code.

As the wave propagates, the wave height decreases, the slope of the front of the wave becomes smaller, and a train of small waves forms behind the main wave. After the wave has propagated for 3 hours and over 2.3 megameters, the wave height has decreased to 0.20 meters. A train of waves has developed behind the main wave with maximum negative amplitude of 0.14 meter and positive amplitude of 0.11 meter. Tsunami Propagation Study—LGW Model A 1 meter high Gaussian surface water displacement with a Gaussian break width of 10 kilometers (which is equivalent to an Airy wave with a half-wave length of 45 kilometers) propagating for 3 hours in 4550 meters deep water was studied. The wave height in meters as a function of distance is shown in Figure 5.7 at various times. The initial surface displacement separates into two nearly shallow water waves with a height of 0.5 meter and length of 45 kilometers.

Evaluation of incompressible models for modeling water waves

147

Fig. 5.7. The wave height in meters as a function of distance in kilometers at various times with a width of 45 kilometers in 4550 meters deep water for the linear gravity wave model using the LGW code.

The analytical linear gravity wave calculations were performed using up to 16,384 points. The wave profile was calculated for each time selected. As the wave propagates, the wave height decreases, the slope of the front of the wave becomes smaller, and a train of small waves forms behind the main wave. After the wave has propagated for 3 hours and over 2.3 megameters, the wave height has decreased to 0.22 meter. A train of waves has developed behind the main wave with a maximum negative amplitude of 0.16 meter and positive amplitude of 0.12 meter. The peak height of the linear gravity wave decreases slower and the following wave train heights are higher than for the Navier-Stokes wave. The peak height of the Navier-Stokes wave is about 10 percent less than the linear gravity wave peak height after three hours of travel. This difference is an order of magnitude less than the difference between either wave and the shallow water wave height, which is more than twice as high. The period of the waves in the following wave train are similar; however the amplitudes are slightly higher for the linear gravity wave train. These results support the use of the linear gravity wave model for studying the effect of various wave length-to-depth ratios and of tsunami propagation for longer distances and times. The linear gravity wave model was used to study wave propagation from initial surface displacements of width-to-depth ratios of 40 to 5.0 (Gaussian break widths of 40 to 5

Numerical modeling of water waves

148

kilometers) in 4550 meters deep water. The wave height after the wave had traveled for 0.5, 2.0, 5.0 and 10.0 hours for depth ratios of 40, 20, 10, and 5 are shown in Figure 5.8.

Fig. 5.8. The wave height in meters as a function of distance in meters after the wave has traveled up to 10 hours and 7.6 megameters for a 1 meter high Gaussian displacement with a Gaussian break length of 40, 20, 10, and 5 kilometers in 4550 meters deep water using the LGW code.

The Airy half-wave length and half-wave period equivalents for a break width of 5 is 23 kilometers and about 1.8 minutes, a break width of 10 is 46 kilometers and 3.5 minutes, a break width of 20 is 90 kilometers and 7 minutes, and for a break width of 40 is 180 kilometers and 14 minutes. After traveling for 10 hours and 7.6 megameters, the peak wave amplitude for the 40 break width wave decreased from 0.5 to 0.42 meter, the 20 break width wave amplitude decreased to 0.28 meter, the 10 break width wave amplitude decreased to 0.16 meter, and the 5 break width wave amplitude decreased to 0.08 meter. The addition of nonlinear effects will lower these heights so these are upper bounds on the peak wave amplitude. Propagation Conclusions The nonlinear shallow water wave becomes less realistic the further it travels from its source and the smaller the width-to-depth ratio. Non-linear shallow water models are adequate for large wave length and long period tsunamis such as generated by the 1964 Alaskan or the 1960 Chile earthquakes where the periods were about 30 minutes and

Evaluation of incompressible models for modeling water waves

149

wave length to depth ratio in the deep ocean was greater than 80. Tsunami wave propagation from and generation by earthquakes with small areas and periods of a few minutes are not realistically described using shallow water models. Either the linear gravity wave model or better the Navier-Stokes model should be used for accurate modeling of long distance propagation of tsunamis with small (less than 40) width to depth ratios. Most tsunami waves that have been observed after traveling across the ocean have periods longer than 10 minutes. These studies support the postulate that this is because the shorter wave length tsunamis are so dispersive that as they propagate long distances, their amplitude decreases by an order of magnitude. Since three-dimensional Navier-Stokes solutions to tsunami propagation problems currently have limited practical application, most tsunami flooding studies need to be performed using the shallow water model. These studies permit us to determine the parameteric region where the shallow water model can be useful. A summary of that region is presented at the end of this chapter. The development and application of three-dimensional Navier-Stokes models including a realistic bottom friction treatment will be required for significant improvement in our ability to realistically model tsunami generation, propagation and flooding. The current effort to develop such a capability is described in Chapter 6. Computer movies of the propagation calculations are in the NMWW CD-ROM TSUNAMI.MVE/PRORMVE directory.

5C. Tsunami Wave Flooding The shallow water SWAN code and the Navier-Stokes ZUNI codes were used to study the effect of tsunami wave period, amplitude, bottom slope angle and friction on tsunami shoaling and flooding. The shallow water waves shoal higher, steeper and faster than the Navier-Stokes waves. The differences increase as the periods become shorter and the slopes less steep with large differences for periods less than 500 sec and slopes less than 2 percent.

Numerical modeling of water waves

150

Tsunami Flooding Study—SWAN Model A 3 meters half-height tsunami wave of various periods was propagated in 12 meters deep water and then interacted with slopes of various steepness. The wave height in meters as a function of distance is shown in Figure 5.9 at various times for a tsunami Airy wave with a 900 sec period and a 3 meters half-height traveling 3000 meters in 12 meters deep water before it interacted with a frictionless 1 percent slope. The space resolution in the numerical model grid was 10 meters with 300 cells to the bottom edge of the slope, 120 cells on the slope below and 80 cells above the water surface. The calculations were performed at 0.5 sec intervals. The peak wave height was 6.7 meters, the runup wave height was 6.0 meters and the inundation limit was 600 meters. The steepness of the slope was changed by changing the space resolution and the time step. The wave period was changed assuming that the tsunami was a shallow water Airy wave. The results are shown in Table 5.2 and Figures 5.10–5.12.

Evaluation of incompressible models for modeling water waves

151

Fig. 5.9. The interaction of a 900 second period, 3 meters half-height tsunami in 12 meters of water with a 1 percent slope at various times. The last figure shows the first wave maximum extent of flooding.

Numerical modeling of water waves

152

TABLE 5.2 SWAN Shallow Water Flooding Study

Slope percent

Period seconds

DeChezy coef

Peak Ht. meters

Runup meters

Ht .Flood meters

Initially 3 m Wave at 12 m Depth 4.0

900.

0.0

6.8

6.3

160

2.0

900.

0.0

7.0

6.4

330

1.0

900.

0.0

6.7

6.0

600

0.5

900.

0.0

5.0

4.7

1000

0.25

900.

0.0

3.7

3.2

1300

1.0

1570.

0.0

7.2

6.3

600

1.0

450.

0.0

5.2

4.8

500

1.0

225.

0.0

3.6

3.2

350

2.0

1800.

0.0

6.0

5.6

280

2.0

450.

0.0

6.8

6.0

300

2.0

225.

0.0

5.0

4.8

250

4.0

450.

0.0

7.0

6.4

160

4.0

225.

0.0

6.6

6.0

160

2.0

900.

10.0

4.8

0.5

30

2.0

900.

20.0

5.8

1.8

100

2.0

900.

30.0

6.4

3.0

160

2.0

900.

50.0

6.6

5.2

280

1.0

900.

10.0

4.2

0.4

40

1.0

900.

15.0

5.0

0.8

80

1.0

900.

25.0

6.0

1.6

160

1.0

900.

40.0

6.5

2.8

300

1.0

900.

50.0

6.5

3.8

390

1.0

900.

60.0

6.5

4.4

440

100.

3600.

0.0

4.8

4.8

100.

900.

0.0

6.1

6.1

100.

225.

0.0

6.5

6.5

Evaluation of incompressible models for modeling water waves

153

Initially 1 m Wave at 4550 m Depth 5.0

2666.

0.0

3.5

5.0

1333.

0.0

4.5

5.0

999.

0.0

5.5

5.0

666.

0.0

6.5

5.0

330.

0.0

7.7

The numerical model grid described is satisfactory for the cases studied giving runup and inundation values that are independent of space, time and geometry resolution. The grid model becomes less accurate for steeper slopes and longer wave length waves. The longer wave length waves need longer distances and times for the shoaling to occur; otherwise the calculated runup and inundation will be grid size dependent, decreasing with decreasing distance to the bottom edge of the slope. The effect of the steepness of a frictionless slope on the runup height and inundation limit is shown in Figure 5.10 for a tsunami wave with a 3 meters half-height and a 900 sec period in 12 meters of water. As shown in Figure 5.9, the highest shoaling wave height does not necessarily occur at the front of the wave. The difference between the height at the front of the wave and the peak wave height increases as the slope becomes steeper, with the maximum difference occuring for slopes larger than 1.5 percent. As shown in Figure 5.10, the runup height decreases and the inundation limit increases with decreasing slope angle. Also shown in Figure 5.10 is the predicted runup height for the Bretschneider model without friction. Bretschneider’s model results in smaller runup than SWAN for slopes greater than 1 percent. A popular engineering method for estimating maximum probable tsunami inundation zones is to use historical inputs to provide tsunami wave heights as a function of frequency of occurrence and then to use the Bretschneider6 model to calculate the runup on the shore at selected sites. The Bretschneider model as it is normally used includes surface roughness. The roughness parameters have been calibrated to reproduce observed runups. The Bretschneider model is independent of wave period and the only slope effect is from the friction.

Numerical modeling of water waves

154

Fig. 5.10. The runup wave height and the inundation limit in meters as a function of slope for a 3 meters half-height 900 sec period tsunami wave in 12 meters of water for the nonlinear shallow water model using the SWAN code.

The effect of wave period and frictionless slope steepness on tsunami flooding is shown in Figure 5.11. The maximum runup height is obtained for short period waves interacting with steep slopes, while the maximum inundation occurs for 1000 sec waves on gentler slopes.

Evaluation of incompressible models for modeling water waves

155

Fig. 5.11. The runup wave height and the inundation limit in meters as a function of wave period in sec for a 3 meters half-height 900 sec period tsunami wave in 12 meters of water for the nonlinear shallow water model using the SWAN code.

Numerical modeling of water waves

156

The effect of friction on tsunami flooding of 1 and 2 percent slopes is shown in Figure 5.12. The effect of increasing friction (decreasing DeChezy friction constant) is small on the peak height and becomes larger on the runup height as the slope decreases. In the Bretschneider model, the surface roughness is described using a Manning “n” where 0.04 corresponds to a roughness characteristic of grass or small rocks and 0.03 to a roughness of many trees, large boulders or high grass. The DeChezy friction coefficient is related to Manning “n“by the depth to the 1/6 power. While not directly comparable, for depths in this study a DeChezy friction coefficient of 50 results in about the same friction effect as a Manning “n” in the 0.03–0.04 range. Tsunami Flooding Study—ZUNI Model A 1 meter high tsunami wave of various periods was propagated in 4550 meters deep water and then it interacted with frictionless slopes of various steepness. These results are an extension of those described in Chapter 3. The geometry of the flooding calculation is shown in Figure 5.13 for a 1/15 or 6.66 percent slope. Calculations for this geometry were performed with 15 cells in the Y or depth direction and 68 cells in the X or distance direction. The cells were rectangles 450 meters high in the Y direction and 6750 meters long in the X direction. The time increment was 3 sec. The water level was placed at 4550 meters or 50 meters up into the

Fig. 5.12. The wave height and the inundation limit in meters for a 1 and 2 percent slope as a function of the DeChezy friction coefficient for a 3 meters half-height 900 sec period tsunami wave in 12 meters of water for the nonlinear shallow water model using the SWAN code.

Evaluation of incompressible models for modeling water waves

157

eleventh cell. The viscosity coefficient used was 0.02 poise, a value representative of the actual viscosity for water. The viscosity did not significantly affect the results. The slope in a ZUNI calculation is determined by the diagonal through the cell as discussed in the numerical methods section. A 6.66 percent slope results in a cell aspect ratio of 1 to 15 and a 4.0 percent slope an aspect ratio of 1 to 25. Numerical errors increase with increasing aspect ratio to unacceptable levels for 1 to 50 and greater ratios.

Fig. 5.13. The tsunami flooding geometry for the Navier-Stokes model using the ZUNI code.

Numerical modeling of water waves

158

TABLE 5.3 ZUNI Navier-Stokes Flooding Study

Slope percent

Period seconds

First Peak meters

First Trough meters

2nd Peak meters

2nd Trough meters

Period Study for 1 m Wave at 4550 m Depth 6.66

1333.

+3.2

−3.2

+4.1

−4.1

6.66

900.

+3.0

−4.0

+4.1

6.66

500.

+1.9

−2.9

+2.9

−2.9

6.66

250.

+0.9

−0.5

+0.25

−0.1

Period Study for 1 m Wave at 1011 m Depth 1.48

2500.

+2.8

−3.3

+3.8

1.48

1700.

+2.5

−3.7

+3.6

−3.2

1.48

1333.

+2.1

−3.6

+3.2

−3.3

1.48

900.

+1.5

−2.0

+2.0

1.48

500.

+0.8

−0.3

+0.3

10.0

1333.

+3.0

−3.2

+3.2

6.6

1333.

+3.2

−3.9

+4.1

−4.1

4.0

1333.

+2.9

−4.0

+3.9

−3.8

Slope Study, Slope Width Varied, 4550 m Depth

Slope Study, Depth Varied, Slope Width Constant 13.48

1333.

+3.2

−3.3

+3.5

−3.0

6.66

1333.

+3.2

−3.9

+4.1

−4.1

2.96

1333.

+2.9

−3.9

+3.85

2.22

1333.

+2.8

−3.8

+3.8

−3.7

1.48

1333.

+2.1

−3.6

+3.2

−3.3

0.74

1333.

+1.1

−1.9

+2.8

Evaluation of incompressible models for modeling water waves

159

TABLE 5.4 A 1333 Sec Tsunami Wave with Fixed Slope Width

Slope percent

Period seconds

Depth meters

Vel m/sec

Length km

Slope Width km

13.48

1333.

9100.

299.

398.

68.25

6.66

1333.

4550.

210.

280.

68.25

2.96

1333.

2022.

141.

188.

68.25

2.22

1333.

1517.

122.

163.

68.25

1.48

1333.

1011.

99.

132.

68.25

0.74

1333.

505.

70.

93.

68.25

Fig. 5.14. The runup height as a function of slope for a 1 meter high, 1333 sec period tsunami for the Navier-Stokes model using the ZUNI code.

Another method of changing the slope is to keep the width and aspect ratio constant and change the depth of the calculation. This permits slopes from 13.48 to 0.75 percent to be studied; however, if the period is kept constant, the other wave parameters change as shown in Table 5.4.

Numerical modeling of water waves

160

The slope study was performed using both methods to change the slope with the results shown in Table 5.3 and Figure 5.14. The general features of the calculated runup as a function of slope are the same for the ZUNI and SWAN models as shown in Figures 5.10 and 5.14. The resolution of the ZUNI calculation is inadequate to resolve the peak from the runup height. The ZUNI calculation permits us to realistically calculate wave interactions after the first wave interacts with the slope and to obtain rundown and runup for later waves. This is important because the largest runup both observed and in these calculations is the second or third wave. The interaction of the reflected waves with the later waves often results in the second or third runup being much different than the first runup. The magnitude and direction of the effect depends upon both the slope and the wave period as shown in Table 5.3. To perform these calculations, the source of the wave must be far removed from the slope in order to correctly follow the interaction of the incoming and reflected waves for several waves. A longer wave period needs a larger distance to run and requires more time for the multiple wave interactions. This results in large computing meshes and long running times for the long period tsunamis. Since they also result in the greatest flooding, they are also of the most interest. The effect of wave period was investigated for a 1 meter half-height wave at 4550 meters interacting with a 6.66 percent slope and at 1011 meters depth with a 1.48 percent slope. The results of the study of wave period and slope are shown in Table 5.3 and Figure 5.15. Also shown in Figure 5.15 are the SWAN runup heights as a function of period for the first wave. To obtain adequate resolution in the SWAN calculation, a much smaller mesh is required than in the ZUNI calculation.

Evaluation of incompressible models for modeling water waves

161

Fig. 5.15. The runup height as a function of period and different slopes and the first and second waves for a 1 meter high, 1333 sec period tsunami for the Navier-Stokes model using the ZUNI code. Also shown is the first wave using the SWAN shallow water code.

This is a result of the particle surface following treatment in the ZUNI calculation which resolves the surface in each cell. The SWAN calculations were performed using a

Numerical modeling of water waves

162

780 cell long mesh, a mesh width of 500 meters and a time step of 0.3 sec. The 4500 meters deep bottom was 250 kilometers long (500 cells) to the edge of the 5 percent shelf formed by 90 kilometers (180 cells) below and 50 kilometers (100 cells) above the water surface. Calculations were also performed using a mesh width of 250 meters and a 1560 cell long mesh. As shown in Table 5.3 and Figure 5.15, the shallow water waves flood higher for the first wave than the Navier-Stokes waves. The differences increase with shorter periods and less steep slopes, with large differences for periods less than 500 sec and slopes less than 2 percent. For higher period waves, the second wave is as much as a third larger than the first wave. The second wave is also larger than the shallow water first wave for wave periods greater than 500 sec. Conclusions The effect of tsunami wave period, amplitude, bottom slope angle and friction on tsunami shoaling and flooding has been investigated using a shallow water and a Navier-Stokes model. The study shows higher wave shoaling and flooding for waves interacting with steeper slopes, for waves with longer periods and for waves from deeper water. The shallow water waves shoal higher, steeper and faster than the Navier-Stokes waves. The differences increase as the periods become shorter and slopes less steep with large differences for periods less than 500 sec and slopes less than 2 percent. The interaction of the reflected first wave with the later waves often results in the second or third runup being much different than the first runup. The magnitude and direction of the effect depends upon both the slope and the wave period. For higher period waves, the second wave is as much as a third larger than the first wave. Since the shallow water model distorts the shoaled and reflected wave, the wave interaction and resulting second and third runups are inaccurately modeled. The effect of friction was investigated only for the shallow water model. It is the largest unknown factor in evaluating tsunami flooding. The Navier-Stokes model used does not have a treatment for friction. The development and evaluation of a realistic friction model is an important remaining tsunami flooding problem. A possible solution to the problem is described in Chapter 6. Since three-dimensional Navier-Stokes solutions to tsunami flooding problems currently have limited practical application, many tsunami flooding studies will need to be performed using the shallow water model. The following is a summary of the wave parameters where the shallow water model is valid. SOURCE—wave length must be 10 times longer than the depth. PROPAGATION—wave length must be 40 times longer than deep ocean depth and period needs to be longer than 15 minutes. RUNUP—period needs to be longer than 15 minutes. FLOODING—period needs to be longer than 10 minutes and the slope needs to be greater than 2 percent. Many important water wave problems can not be adequately modeled using the shallow

Evaluation of incompressible models for modeling water waves

163

water model. The incompressible Navier-Stokes model is adequate for many water wave problems that can not be accurately modeled using shallow water models. Currently only an upper limit flooding determination can be obtained without a realistic friction model. As shown in Chapter 6, it is now possible to model with high resolution the interface between the ocean and the sea floor, the mixing and movement of sediment and to include detailed sea floor topography. Future tsunami flooding studies can include friction effects realistically using the modeling methodology described in Chapter 6. One of the more important water wave problems is the determination of Civil Defense evacuation zones for tsunamis. As shown in this chapter the shallow water model is inadequate for the task. The incompressible Navier-Stokes model is unable to describe the initial generation of waves where the fluid-dynamics is compressible such as large earthquakes, impact landslides, volcanic, nuclear or conventional explosions and asteroids. Until recently there was no practical numerical method to solve such problems. Recent advances in the numerical modeling of water waves using the compressible Navier-Stokes model permit us to model almost any water wave generation, propagation and flooding problem. The next chapter describes the compressible model.

References 1. Charles L.Mader, “Numerical Tsunami Flooding Study—I,” Science of Tsunami Hazards, Vol.8, pp.79–96 (1990). 2. Charles L.Mader, Dennis W.Moore, and George F.Carrier, “Numerical Tsunami Source Study—II,” Science of Tsunami Hazards, Vol.11, pp.81–92 (1993). 3. Charles L.Mader, Dennis W.Moore, and George F.Carrier, “Numerical Tsunami Propagation Study—III,” Science of Tsunami Hazards, Vol.11, pp.93–106 (1993). 4. G.F.Carrier and H.P.Greenspan, “Water Waves of Finite Amplitude on a Sloping Beach,” Journal of Fluid Mechanics, Vol.4, pp.97–109 (1958). 5. G.F.Carrier, “The Dynamics of Tsunamis,” Mathematical Problems in the Geophysical Sciences, American Mathematical Society, Vol.1, pp.157–187 (1971). 6. C.L.Bretschneider and P.G.Wybro, “Tsunami Inundation Prediction,” Proceedings of 15th Coastal Engineering Conference (1976).

6 MODELING WAVES USING COMPRESSIBLE MODELS Compressible fluid dynamics is modeled numerically using the Lagrangian and the Eulerian equations of motion. One, two and three-dimensional models and codes have been developed over the last fifty years using various numerical methods for solving the equations of motion and to describe the compressible equation of state. In reference 1 is a complete description of Lagrangian and Eulerian approaches to modeling compressible reactive flows and the computer codes that were developed through the mid 1990’s. The Department of Energy’s program Accelerated Strategic Computing Initiative (ASCI) of the last five years has resulted in major advances in computing technology and in methods for improving the numerical resolution of compressible hydrodynamic calculations. Only recently has it become possible to calculate water wave problems for minutes or even hours using compressible hydrodynamic models that require millions of time steps for each second of flow. Only recently has it become possible to finely resolve a water-air interface and follow the water wave with millimeter resolution in a problem with 40 kilometers of air, 5 kilometers of water and tens of kilometers of ocean crust. The limitations of incompressible flow discussed in the previous chapters no longer must be accepted for modeling water waves. In this chapter the most advanced ASCI computer codes for modeling compressible fluid dynamics will be described. The codes will be used to model water waves generated from impact landslides, explosions, projectile impacts and asteroids. Numerical modeling of water waves has advanced so far so rapidly that it is clearly a technological revolution.

6A. The Three-Dimensional Compressible Model The three-dimensional partial differential equations for nonviscous, nonconducting, compressible fluid flow are The Nomenclature I

internal energy

P

pressure

Ux

velocity in x direction

Uy

velocity in y direction

Numerical modeling of water waves Uz

velocity in z direction

ρ

density

t

time

166

These equations are solved by a high-resolution Godunov differencing scheme using an adaptive grid technique described in reference 2. The solution technique uses Continuous Adaptive Mesh Refinement (CAMR). The decision to refine the grid is made cell-by-cell continuously throughout the calculation. The computing is concentrated on the regions of the problem which require high resolution. Refinement occurs when gradients in physical properties (density, pressure, temperature, material constitution) exceed defined limits, down to a specified minimum cell size for each material. The mesh refinement is shown in Figure 6.1 for a complicated geometry embedded in a material.

Modeling waves using compressible models

167

Fig. 6.1. Adaptive mesh refinement applied to a complicated geometry.

Much larger computational volumes, times and differences in scale can be simulated than possible using previous Eulerian techniques such as those described in reference 1. The original code was called SAGE. A later version with radiation is called RAGE. A recent version with the techniques for modeling reactive flow described in reference 1 is called NOBEL and was used for modeling the Lituya Bay impact landslide generated tsunami and water cavity generation. The codes can describe one-dimensional slab or spherical geometry, two-dimensional slab or cylindrical geometry, and three-dimensional Cartesian geometry. Because modern supercomputing is currently done on clusters of machines containing many identical processors, the parallel implementation of the code is very important. For portability and scalability, the codes use the Message Passing Interface (MPI). Load leveling is accomplished through the use of an adaptive cell pointer list, in which newly created daughter cells are placed immediately after the mother cells. Cells are redistributed among processors at every time step, while keeping mothers and daughters together. If there are a total of M cell and N processors, this technique gives nearly M/N cells per processor. As neighbor cell variables are needed, the MPI gather/scatter routines copy those neighbor variables into local scratch memory. The codes incorporate multiple material equations of state (analytical or SESAME tabular). Every cell can in principle contain a mixture of all the materials in a problem assuming that they are in pressure and temperature equilibrium. As described in reference

Numerical modeling of water waves

168

1, pressure and temperature equilibrium is appropriate only for materials mixed molecularly. The assumption of temperature equilibrium is inappropriate for mixed cells with interfaces between different materials. The errors increase with increasing density differences. While the mixture equations of state described in reference 1 would be more realistic, the problem is minimized by using fine numerical resolution at interfaces. The amount of mass in mixed cells is kept small resulting in small errors being introduced by the temperature equilibrium assumption, The strength is treated using an elastic-plastic model identical to that described in reference 1. Very important for water cavity generation and collapse and the resulting water wave history is the capability to initialize gravity properly, which is included in the code, This results in the initial density and initial pressure changing going from the atmosphere at 40 kilometers altitude down to the ocean surface. Likewise the water density and pressure changes correctly with ocean depth. Some of the remarkable advances in fluid physics using the SAGE code have been the modeling of Richtmyer-Meshkov and shock induced instabilities described in references 3 and 4.

6B. The Lituya Bay Mega-Tsunami Lituya Bay, Alaska is on the northeast shore of the Gulf of Alaska. It is an ice-scoured tidal inlet with a maximum depth of 220 meters and a narrow entrance with a depth of only 10 meters. It is a T-shaped bay, 7 miles long and up to 2 miles wide. The two arms at the head of the bay, Gilbert and Crillon Inlets, are part of a trench along the Fairweather fault. On July 87 1958, a 7.5 Magnitude earthquake occurred along the Fairweather fault with an epicenter near Lituya Bay. A mega-tsunami wave was generated that washed out trees to a maximum altitude of 520 meters at the entrance of Gilbert Inlet. Much of the rest of the shoreline of the bay was denuded by the tsunami from 30 to 200 meters altitude as shown in Figures 6.2– 6.4.

Modeling waves using compressible models

169

Fig. 6.2. The view of Lituya Bay in August 1958 after the earthquake and tsunami. The forest was destroyed to a maximun elevation of 524 meters and a maximum distance of 1000 meters from shoreline. The bay is 7 miles long and up to 2 miles wide. The two arms at the head of the bay are part of the Fairweather fault.

Numerical modeling of water waves

170

Fig. 6.3. Lituya Bay map showing topographic and bathymetric contours, trace of Fairweather fault, 1958 rockslide and trimline of wave runup. Forests destroyed to maximum elevations of 524 meters and 208 meters on north and south shores.

During the last 150 years five giant waves have occurred in Lituya Bay. The previous event occurred on October 27, 1936 which washed out trees to a maximum altitude of 150 meters and was not associated with an earthquake. Don Miller recorded all that was known in 1960 about the giant waves in Lituya Bay in reference 5. The July 9, 1958 earthquake occured at about 10:15 p.m., which is still daylight at Lituya Bay. The weather was clear and the tide was ebbing at about plus 5 feet. Bill and Vivian Swanson were on their boat anchored in Anchorage Cove near the western side of the entrance of Lituya Bay. Their astounding observations are recorded in reference 6 and were as follows: “With the first jolt, I tumbled out of the bunk and looked toward the head of the

Modeling waves using compressible models

171

bay where all the noise was coming from. The mountains were shaking something awful, with slides of rock and snow, but what I noticed mostly was the glacier, the north glacier, the one they call Lituya Glacier. I know you can’t ordinarily see that glacier from where I was anchored. People shake their heads when I tell them I saw it that night. I can’t help it if they don’t believe me. I know the glacier is hidden by the point when you’re in Anchorage Cove, but I know what I saw that night, too. The glacier had risen in the air and moved forward so it was in sight. It must have risen several hundred feet. I don’t mean it was just hanging in the air. It seems to be solid, but it was jumping and shaking like crazy. Big chunks of ice were falling off the face of it and down into the water. That was six miles away and they still looked like big chunks. They came off the glacier like a big load of rocks spilling out of a dump truck. That went on for a little while—it’s hard to tell just how long—and then suddenly the glacier dropped back out of sight and there was a big wall of water going over the point. The wave started for us right after that and I was too busy to tell what else was happening up there.” A 15 meters high wave rushed out of the head of the bay toward the Swansons’ anchored boat. The boat shot upward on the crest of the wave and over the tops of standing spruce trees on the entrance spit of Lituya Bay. Bill Swanson looked down on the trees growing on the spit and said he was more than 25 meters above their tops. The wave crest broke just outside the spit and the boat hit bottom and foundered some distance from the shore. Swanson saw water pouring over the spit, carrying logs and other debris. The Swansons escaped in their skiff to be picked up by another fishing boat 2 hours later. The front of Lituya Glacier on July 10 was a nearly straight, vertical wall almost normal to the trend of the valley. Comparisions with photographs of the glacier taken July 7 indicate that 400 meters of ice had been sheared off of the glacier front. The photograph taken after the event is shown in Figure 6.4.

Numerical modeling of water waves

172

Fig. 6.4. The view of Lituya Bay in August 1958 where the forest was destroyed to a maximun elevation of 524 meters.

After the earthquake there was a fresh scar on the northeast wall of Gilbert Inlet, marking the recent position of a large mass of rock that had plunged down the steep slope into the water. The next day after the earthquake and tsunami, loose rock debris on the fresh scar was still moving at some places, and small masses of rock still were falling from the rock cliffs near the head of the scar. The wave destroyed the forests and even removed the barnacles from the rocks. When the author visited Lituya Bay forty years after the event, almost no vegetation grew on either the slide region or on the cleared area shown in the above photograph. The dimensions of the slide on the slope are accurate but the thickness of the slide mass normal to the slope can only be estimated. The main mass of the slide was a prism of rock that was 730 meters and 900 meters along the slope with a maximum thickness of 90 meters and average thickness of 45 meters normal to the slope. The center of gravity was at about 600 meters altitude. As described in reference 5 this results in an approximate volume of 30 million cubic meters (40 million cubic yards). Miller5 concluded that “the rockslide was the major, if not the sole cause of the 1958 giant wave.” The Swanson observations have not been believed as they indicate that a lot more than a simple landslide occurred. Shallow Water Modeling In reference 7 shallow water modeling was performed using the SWAN nonlinear shallow water code. The generation and propagation of the tsunami wave of July 8, 1958 in

Modeling waves using compressible models

173

Lituya Bay was modeled using a 92.75 by 92.75 meters grid of the topography. The 3 by 6 second land topography was generated from the Rocky Mountain Communication’s CD-ROM compilation of the Defense Mapping Agency (DMA) 1 by 1 degree blocks of 3 arc second elevation data. The sea floor topography was taken from sea floor topographic maps published in reference 5. The grid was 150 by 150 cells and the time step was 0.15 sec. It was concluded that displaced water by a simple landslide or an earthquake along the Fairweather fault at the head of the bay could not result in the observed 550 meters high runup. The water in the inlet with the width of the landslide and between the landslide and the 520 meter high runup was sufficient to cover the runup region to a 100 meters height. In reference 7 it was shown that this high water layer was sufficient to form a wave that will reproduce the observed flooding of the bay beyond the inlet. It was concluded that a landslide impact model was required similar to that for asteroid generated waves. In 1999 it appeared that the numerical technology required would not become available for many decades. Physical Modeling During the Summer of 2000, Hermann Fritz8 conducted experiments that reproduced the 1958 Lituya Bay event. The 1:675 scale laboratory model of Lituya Bay shown in Figure 6.5 was built at VAW at the Swiss Federal Institute of Technology at Zurich, Switzerland. The laboratory experiments indicated that the 1958 Lituya Bay 524 meters runup on the spur ridge of Gilbert Inlet could have been caused by a landslide impact. The study was reported in reference 8. A novel pneumatic landslide generator was used to generate a high-speed granular slide with a controlled impact velocity and shape. A granular slide with the density and volume given by Miller5 was impacted with a mean velocity of 110 meters/sec. It generated a large air cavity and an extremely nonlinear wave with a maximum scaled height of about 160 meters which ran up to a scaled elevation of 530 meters above mean sea level.

Numerical modeling of water waves

174

Fig. 6.5. Experimental setup with pneumatic installation and measurement systems that included laser distance sensors (LDS), capacitance wave gauges (CWG) and particle image velocimetry (PIV).

The initial geometry used by Fritz8 in his experimental modeling in reference 8 is shown in Figure 6.6.

Fig. 6.6. Cross section of Gibert Inlet along slide axis used in the physical model.

The rockslide was simulated using 4 mm diameter Barium Sulfate grains with a density

Modeling waves using compressible models

175

of 1.6 g/cc and a void fraction of 39 percent. To what extent the rockslide broke up before impacting Gilbert Inlet is unknown. Scanned slide impact profiles were used to calculate the mean slide impact velocity of 110 meters/sec. The slide profile showed a gentle increase in slide thickness with time to scaled maximum of 134 meters and a fast decay back to zero. The scaled maximum slide thickness of 134 meters is 1.4 times the pre-motion slide thickness of 92 meters estimated roughly by Miller in reference 5. The three phases—granular material, water and air—are clearly separated along distinct border lines before flow reattachment occurs. Flow reattachment traps a large volume of air in the back of the rockslide, which leads to large bubble formation, bubble breakup and massive mixing. The slide is deformed due to impact and deflection at the channel bottom. A sequence of twelve instantaneous velocity vector fields computed the PIV data are shown in Figure 6.7. The sequence starts at 0.76 sec after rockslide impact and continues with a time step of 1.73 sec.

Numerical modeling of water waves

176

Fig. 6.7. Photosequence of the Fritz granular impact slide experiment. The time interval is 1.73 sec with the first image at 0.76 sec after impact.

Compressible Navier-Stokes Modeling The Lituya Bay impact landslide generated tsunami was modeled in reference 9 with the Navier-Stokes AMR (Adaptive Mesh Refinement) Eulerian compressible hydrodynamic

Modeling waves using compressible models

177

code NOBEL including the effects of gravity. The geometry of the calculation is shown in Figure 6.8 superimposed on the actual topography.

Fig. 6.8. Gilbert Inlet illustration showing rockslide dimensions, impact site and wave runup to 524 meters on the spur ridge directly opposite to the rockslide impact. Direction of view is North. The front of the Lituya Glacier is set to the 1958 post slide position. The preslide position extended 400 meters closer to the entrance of Gilbert Inlet.

The calculated density profiles are shown in Figures 6.9 and 6.10 for a rockslide moving with a resultant velocity of 110 meters/sec (X and Y component velocities of 77.8 meters/sec). The rockslide had an area of 21,000 square meters and was basalt with a density of 2.868 g/cc. The initial slide profile was triangular with its apex near the back of the slide to reproduce the Fritz slide profile shape described previously. The initial water depth was 120 meters and the length was 1.4 kilometers. The calculated maximum wave height in the bay was about 250 meters above sea level which ran up to 580 meters which is to be compared to the observed 524 meters. Such a wave could lift the glacier ice and result in the spectacular behavior observed by Swanson.

Numerical modeling of water waves

178

Fig. 6.9. Calculated density profiles of the Lituya impact landslide generated tsunami. The times are 0, 6, 10, 16, 22, 24, 26, 30, 36, 40, 48 and 58 sec.

A PowerPoint presentation with Fritz’s experimental movies and computer animations is available on the NMWW CD-ROM in the NOBEL/LITUYA directory.

Modeling waves using compressible models

179

Numerical modeling of water waves

180

Fig. 6.10. Calculated density profiles of the Lituya impact landslide generated tsunami. The slide material is shown as a dark region running across the ocean floor. At later times it comes to rest mostly in the middle of the bay. This illustrates how friction effects can be modeled realistically.

Conclusions The mega-tsunami that occurred on July 8, 1958 in Lituya Bay washed out trees to a maximum altitude of 520 meters at the entrance of Gilbert Inlet. Much of the rest of the shoreline of the bay was denuded by the tsunami from 30 to 200 meters altitude. The Lituya Bay impact generated tsunami was modeled in reference 9 with the AMR Eulerian compressible hydrodynamic code called NOBEL/SAGE. The capability now exists to evaluate the potential impact landslide tsunami hazards for vunerable regions of the world.

Modeling waves using compressible models

181

6C. Water Cavity Generation Introduction In the mid 1960’s, B.G.Craig10 at the Los Alamos National Laboratory performed experiments designed to characterize the formation of water waves from explosives detonated near the water surface. He reported observing the formation of ejecta water jets above and jets or “roots” below the expanding gas cavity. This was unexpected and a scientific mystery which has remained unsolved until it was finally modeled using the NOBEL code in December of 2002. In the early 1980’s, the hypervelocity impact (1.25 to 6 kilometers/sec) of projectiles into water was studied at the University of Arizona by Gault and Sonett11. They observed quite different behavior of the water cavity as it expanded when the atmospheric pressure was reduced from one to a tenth atmosphere. Above about a third of an atmosphere, a jet of water formed above the expanding cavity and a jet or “root” emerged below the bottom of the cavity. In the mid 1980’s, similar results were observed by Kedrinskii12 at the Institute of Hydrodynamics in Novosibirsk, Russia, who created cavities in water by sending large electrical currents through small lengths of small diameter Gold wires (bridge wires) causing the Gold to vaporize. The “exploding bridge wire” is a common method used to initiate propagating detonation in explosives. He observed water ejecta jets and roots forming for normal atmospheric pressure and not for reduced pressures. Thus the earlier Craig observations were not caused by some unique feature of generation of the cavity by an explosion. The process of cavity generation by projectiles or explosives in the ocean surface and the resulting complicated fluid flows has been an important unsolved problem for over 50 years. The prediction of water waves generated by large-yield explosions and asteroid impacts has been based on extrapolation of empirical correlations of small-yield experimental data or numerical modeling assuming incompressible flow, which does not reproduce the above experimental observations. The NOBEL code has been used to model the experimental geometries of Sonett and of Craig. The experimental observations were reproduced as the atmospheric pressure was varied as described in reference 13. Projectile and Exploding Bridge Wire Generated Cavities In the early 1980’s, experiments were being performed at the University of Arizona to simulate asteroid impacts in the ocean. The hypervelocity impact (1.25 to 6 kilometers/sec) of various solid spherical projectiles (Pyrex or Aluminum) with water was performed by Gault and Sonett11. Their observations were similar to those previously observed by Craig10. While the water cavity was expanding, an ejecta jet was formed at the axis above the water plume and a jet or “root” emerged along the axis below the cavity. The water cavity appeared to close and descend into deeper water. To improve the photographic resolution and reduce the light from the air shock, Sonett repeated his impact experiments under reduced atmospheric pressure. Much to

Numerical modeling of water waves

182

everyone’s surprise, when the pressure was reduced from one to a tenth atmosphere, the ejecta jet and the root did not occur and the water cavity expanded and collapsed upward toward the surface. This was what had been expected to occur in both the earlier Craig experiments and the projectile impacts. It became evident that the atmospheric pressure and the pressure differences inside and outside the water plume above the water surface was the cause of the formation of the jet, the root, the cavity closure and descent into deeper water. Figure 6.11 shows the Gault and Sonett11 results for a 0.25 cm diameter aluminum projectile moving at 1.8 kilometers/sec impacting water at one atmosphere (760 mm), and at 16 mm air pressure. A water stem and jet occurs at one atmosphere and not at low pressure. Figure 6.12 shows the Gault and Sonett results for a 0.635 cm diameter aluminum projectile moving at 2.5 kilometers/sec impacting water at one atmosphere (760 mm) and a 0.3175 cm diameter Pyrex projectile moving at 2.32 kilometers/sec impacting water at 16 mm air pressure. A water stem and jet occurs at one atmosphere and not at low pressure. Professor Kedrinskii12 at the Russian Institute for Hydrodynamics was also studying the generation of water cavities from exploding bridge wires. He was observing the formation of ejecta jets and roots as the water cavity expanded similar to those observed by Craig using explosives and by Gault and Sonett using projectiles. After we showed him the effect of reduced atmospheric pressure, he proceeded to repeat his exploding bridge wire experiments under reduced pressure. He observed that the jets and roots did not form when the atmospheric pressure was reduced to 0.2 atmosphere. Figure 6.13 shows the Kedrinskii results for an exploding bridge wire in water at one atmosphere and at 0.2 atmosphere air pressure. Compressible Navier-Stokes Modeling The projectile impact and explosive generated water cavity was modeled with the recently developed full Navier-Stokes AMR (Adaptive Mesh Refinement) Eulerian compressible hydrodynamic code NOBEL described earlier. The continuous adaptive mesh refinement permits the following of shocks and contact discontinuities with a very fine grid while using a coarse grid in smooth flow regions. It can resolve the water plume and the pressure gradients across the water plume and follow the generation of the water ejecta jet and root. Figure 6.14 shows the calculated density profiles for a 0.25 cm diameter aluminum projectile moving at 2.0 kilometers/sec impacting water at five atmosphere air pressure. The water plume collapses at the axis creating a jet moving upward and downward. The jet passes down through the cavity, penetrating the bottom of the cavity at the axis forming the stem. The flow results in the cavity descending down into the water. Figure 6.15 shows the calculated density profiles for a 0.25 cm diameter aluminum projectile moving at 2.0 kilometers/sec impacting water at one atmosphere air pressure. Figure 6.16 shows the calculated density profiles for a 0.25 cm diameter aluminum projectile moving at 2.0 kilometers/sec impacting water at 0.1 atmosphere air pressure. The tip of the water plume continues to expand in contrast to what is observed at

Modeling waves using compressible models

183

atmospheric pressures higher than 200 mm.

Fig. 6.11. The Gault and Sonett experimental results for a 0.635 cm diameter aluminum projectile moving at 1.8 kilometers/sec impacting water.

Numerical modeling of water waves

184

Fig. 6.12. The Gault and Sonett experimental results for a 0.25 cm diameter aluminum projectile moving at 2.5 kilometers/sec at 760 mm air pressure. A 0.3175 diameter Pyrex projectile moving at 2.32 kilometers/sec in 16 mm of air is shown in the bottom frame.

Modeling waves using compressible models

185

Fig. 6.13. The Kedrinskii results for an exploding bridge wire in water at 760 mm and 150 mm air pressure.

Numerical modeling of water waves

186

Fig. 6.14. The density profiles for a 0.25 cm diameter Aluminum projectile moving at 2.0 kilometers/sec impacting water at 5 atmosphere air pressure. The times are 0, 5, 10, 15, 30, and 70 milliseconds. The graphs are 100 cm wide by 120 cm tall, with 50 cm of water.

Modeling waves using compressible models

187

Fig. 6.15. The density profiles for a 0.25 cm diameter Aluminum projectile moving at 2.0 kilometers/sec impacting water at 1 atmosphere air pressure. The times are 0, 25, 75, and 125 milliseconds. The graphs are 100 cm wide and 120 cm tall, with 50 cm of water and 70 cm air.

Numerical modeling of water waves

188

Fig. 6.16. The density profiles for a 0.25 cm diameter Aluminum projectile moving at 2.0 kilometers/sec impacting water at 76 mm (0.1 atmosphere) air pressure. The times are 0, 5, 20, and 50 milliseconds. The graphs are 80 cm wide by 100 cm tall with 40 cm of water and 60 cm air.

Explosive Generated Cavities In reference 14, detailed one-dimensional compressible hydrodynamic modeling is described for explosives detonated in deep water. Agreement was obtained with the experimentally observed explosive cavity maximum radius and the period of the

Modeling waves using compressible models

189

oscillation. It was concluded that the detonation product equation of state over the required range of 1 megabar to 0.01 atmosphere was adequate for accurately determining the energy partition between detonation products and the water. It was also concluded that the equations of state for water and detonation products were sufficiently accurate that they could be realiably used in multidimensional studies of water cavity formation and the resulting water wave generation in the region of the “upper critical depth”. The “upper critical depth” is the experimentally observed location of a explosive charge relative to the initial water surface that results in the maximum water wave height. It occurs when the explosive charge is approximately two-thirds submerged. The observed wave height at the upper critical depth is twice that observed for completely submerged explosive charges at the “lower critical depth.” If the waves formed are shallow water waves capable of forming tsunamis, then the upper critical depth phenomenon would be important in evaluating the magnitude of a tsunami event from other than tectonic events.

Fig. 6.17. The scaled wave height as a function of scaled explosive charge depth. The “upper critical depth” is the explosive charge depth when the maximum wave height occurs which is approximately two-thirds submerged. The second smaller increase in wave height is the “lower critical depth.”

The water wave amplitude as a function of the depth the explosive is immersed in water is shown in Figure 6.17. The scaled amplitude is AR/W and the scaled depth is D/W where A is maximum wave amplitude at a distance R from the explosive charge of weight W. The “upper critical depth” is at the first wave height maximum which occurs when an explosive charge is located at the water surface. The second smaller increase in wave

Numerical modeling of water waves

190

height is at the “lower critical depth” which is about half the upper critical height but results in longer wave length water waves. Data are included for explosives with weights of 0.017 to 175 kilograms. The Craig10 experimental results are shown with a large x. During the study of the upper critical depth phenomenon in the 1960’s evidence of complicated and unexpected fluid flows during water cavity formation was generated by B.G.Craig and described in references 1, 10 and 14. A sphere of explosive consisting of a 0.635 cm radius XTX 8003 (80/20 PETN/Silicon Binder at 1.5 g/cc) explosive and a 0.635 cm radius PBX-9404 (94/6 HMX/binder at 1.84 g/cc) explosive was detonated at its center. The sphere was submerged at various depths in water. PHERMEX15 radiographs and photographs were taken with framing and movie cameras. The radiographs are shown in Figure 6.18. The cavity, water ejecta and water surface profiles shown in the PHERMEX radiography in Figure 6.18 were closely approximated by the compressible hydrodynamic modeling described in reference 14 using the 2DE code and in Figure 6.19 using the NOBEL code.

Modeling waves using compressible models

191

Fig. 6.18. Dynamic radiographs of a 2.54 cm diameter PBX-9404 explosive sphere detonated at its center and half submerged in water at 1 atmosphere air pressure. The times are 15.8, 26.3 and 61.3 microseconds. The sketch shows the prominent features of the radiographs with the water shock dashed.

Numerical modeling of water waves

192

Fig. 6.19. NOBEL density profiles for a 2.54 cm diameter PBX-9404 explosive sphere detonated at its center and half submerged in water at 1 atmospheric air pressure. The times are 15.0, 26.0 and 60.0 microseconds for comparision with the radiographs in Figure 6.18. The width is 16 cm and the height is 24 cm, of which 16 cm is water.

At later times, while the water cavity was expanding the upper ejecta plume collapsed and converged on the axis generating an upward water ejecta jet on the axis above the

Modeling waves using compressible models

193

water plume and a downward water jet which generated a root on the axis below the bottom of the cavity. These results were not anticipated and neither was the observation that the water cavity proceeded to close at its top and descend down into deeper water. At first it was assumed that there was something unique about the explosive source that was resulting in these remarkable observations. The reactive compressible hydrodynamic numerical models available were unable to reproduce the experimental observations or suggest any possible physical mechanisms unique to explosives. As described previously, the different behavior of the water cavity as it expanded when the atmospheric pressure was reduced from one atmosphere to less than a third of an atmosphere is independent of the method used to generate the cavity, such as an exploding bridge wire or a hypervelocity projectile impact. These remarkable experimental observations resisted all modeling attempts for over 25 years. The numerical simulations could not describe the thin water ejecta plumes formed above the cavity or the interaction with the atmosphere on the outside of the ejecta plume and the pressure inside the expanding cavity and plume. Figure 6.20 shows the calculated water profiles for a 0.25 cm diameter PBX-9404 explosive sphere detonated at its center half submerged in water at one atmosphere air pressure. All the projectile, exploding bridge wire and explosive experimental observations were reproduced as the atmospheric pressure was varied. At sufficiently high atmospheric pressure, the difference between the pressure outside the ejecta plume and the decreasing pressure inside the water plume and cavity as it expanded, resulted in the ejecta plume converging and colliding at the axis. A jet of water formed and proceeded above and back into the bubble cavity along the axis. The jet proceedes back through the bubble cavity penetrating the bottom of the cavity and forms the root observed experimentally. The complicated cavity collapse and resulting descent into deeper water was also numerically modeled.

Numerical modeling of water waves

194

Fig. 6.20. The density profiles for a 2.54 cm diameter PBX-9404 explosive sphere detonated at its center and half submerged in water at 1 atmosphere air pressure. The times are 0, 25, 75, and 275 milliseconds. The graphs are 100 cm wide and 120 cm tall, with 50 cm of water and 70 cm of air.

Explosive Generated Water Wave Craig10 measured the wave amplitude as a function of time for the first few seconds at a

Modeling waves using compressible models

195

distance of 4 rneters from a 2.54 cm diameter PBX-9404 explosive sphere initiated at its center in 3 meters of water. He included mass markers in the water. Mass markers located 1 meter below the water surface and markers located 0.5 meter below the surface and 1 meter from the explosive showed no appreciable movement compared with those located nearer the surface or explosive charge. These results showed that the wave formed was not a shallow water wave. The experimental and calculated wave parameters are summarized in Table 6.1. The parameters are given after 4 meters of travel from the center of the explosive charge. The wave parameters for the Airy wave were calculated using the WAVE code described in Chapter 1 for a depth of 3 meters and the experimentally observed wave length. Since the group velocity is almost exactly half the wave velocity, the Airy wave is a deep water wave. The shallow water results are from reference 14. A small wave from the initial cavity formation is followed by a larger negative and then a positive wave resulting from the cavity collapse. Only the second wave parameters are given in the table.

TABLE 6.1 Calculated and Experimental Wave Parameters

Experimental Wave Velocity (m/sec) Amplitude (cm) Wave Length (m) Period (sec)

Airy Wave

Shallow Water

NOBEL

2.50±0. 2

2.41

5.42

2.50±0.10

0.8

1.0

10.1

0.6

3.75

3.75(input)

1.0

3.75

1.5

1.55

0.18

1.5±0.1

Group Velocity (m/sec)

1.21

The wave gauge was close to the edge of the water tank, which resulted in reflected waves which perturbed the subsequent wave measurements. Since the wave gauge was located close (0.69 meter) to the side of the tank, the reflections from the first small wave perturbed the second wave, which probably explains the larger than calculated amplitude, as the calculated wave was unperturbed by any boundary. Conclusions In the late 1960’s and early 1970’s, B.G.Craig at the Los Alamos National Laboratory reported observing the formation of ejecta jets and roots from cavities generated by small spherical explosives detonated near the water surface while the gas cavity was expanding. The hypervelocity impact (1.25 to 6 kilometers/sec) of projectiles into water was studied at the University of Arizona in the early 1980’s by Gault and Sonett. They observed quite different behavior of the water cavity as it expanded when the atmospheric pressure was reduced from one to a tenth atmosphere. Above about a third of

Numerical modeling of water waves

196

an atmosphere, a jet of water formed above the expanding cavity and a root developed below the bottom of the bubble cavity. They did not occur for atmospheric pressures below a third of an atmosphere. Similar results were observed in the middle 1980’s by Kedrinskii at the Institute of Hydrodynamics in Novosibirsk, Russia when the water cavity was generated by exploding bridge wires, with jets and roots forming for normal atmospheric pressure and not for reduced pressures. During the last decade, a compressible Eulerian hydrodynamic code called SAGE has been under development by the Los Alamos National Laboratory and Science Applications International (SAIC) which has continuous adaptive mesh refinement (AMR) for following shocks and contact discontinuities with a very fine grid. A version of the SAGE code that models explosives, called NOBEL, has been used to model the experimental geometries of Sonett and of Craig. The experimental observations were reproduced as the atmospheric pressure was varied. When the atmospheric pressure was increased, the difference between the pressure outside the ejecta plume above the water cavity and the decreasing pressure inside the water plume and cavity as it expanded resulted in the ejecta plume converging and colliding at the axis forming a jet of water proceeding above and back into the bubble cavity along the axis. The jet proceeding back through the bubble cavity penetrated the bottom of the cavity and formed the root observed experimentally. The complicated cavity collapse was numerically modeled. Now that a code is available that can describe the experimentally observed features of projectile interaction with the ocean, we have a tool that can be used to evaluate impact landslide, projectile or asteroid interactions with the ocean, and the resulting generation of tsunami waves. A PowerPoint presentation with Craig’s and Sonett’s experimental movies and computer animations is available on the NMWW CD-ROM in the NOBEL/CAVITY directory.

6D. Asteroid Generated Tsunamis Introduction Two- and three-dimensional compressible hydrodynamic modeling using the SAGE code has been performed to study the generation of tsunamis by asteroid impacts with the ocean16,17. The goal was to determine the characteristics of impact generated tsunami waves as a function of the size and energy of the asteroid. The evaluation in the popular press of the potential threats from modest sized earth crossing asteroids often overestimates the tsunami threat from small (100 to 250 meters in diameter) asteroids. The studies reported are based on faulty shallow water modeling, from which they conclude that the resulting short (less than 5 minutes) period tsunamis are a major threat throughout an ocean basin. The news media are often delighted to promote such imaginary ocean basin wide “mega-tsunami” threats from small asteroid impacts, underwater landslides and volcanic

Modeling waves using compressible models

197

or other explosions. Part of the problem is that until recently there has not been a method for realistically modeling such threats. That changed with the development of the NOBEL/SAGE/RAGE compressible AMR code described in this chapter. While tsunamis up to a kilometer in initial height are generated by impactors of a kilometer diameter, they do not propagate as long period tsunamis such as generated by large earthquakes. Asteroid impact generated waves are highly complex in form and interact strongly with shocks propagating through the water and ocean crust. They decay more rapidly than the inverse of the distance from the impact point. A panel of tsunami scientists of The Tsunami Society have evaluated the proposed tsunami hazards from events that generate initially high amplitude short period tsunamis. Their evaluation is available on The Tsunami Society web site at http://www.sthjournal.org/media.htm. The Tsunami Society study concludes the following: No mega-tsunami has occurred in either the Atlantic or Pacific oceans in recorded history. The colossal collapse of Krakatau or Santorin generated catastropic waves in the immediate areas but hazardous waves did not propagate to distance shores. Carefully performed numerical and experimental model experiments on such events verify that the relatively short period waves from these small, though intense, occurances do not travel as do tsunami waves from a major earthquake. The 1958 Lituya Bay mega-tsunami described earlier in this chapter was the largest wave in recorded history, but it barely registered at nearby tide gauges. It had a short period and wavelength, which resulted in the wave rapidly decaying as It traveled out into the deep ocean. Two-Dimensional Vertical Impact Asteroid Generated Tsunamis A series of calculations in two-dimensions (cylindrical symmetry) for asteroids impacting the ocean vertically at 20 kilometers/sec was performed using the ASCI BlueMountain machine at the Los Alamos National Laboratory. These simulations were designed to follow the passage of an asteroid through the atmosphere, its impact with the ocean, the cavity generation, and subsequent re-collapse and generation of the tsunami. The results of these studies are described in reference 16, The parameter study included different asteroid masses. Stony and iron bodies of diameters 250 meters, 500 meters and 1000 meters were used. The kinetic energies of the impacts ranged from 1.3 Gigatons to 195 Gigatons of equivalent TNT yield. The projectile moved through an exponential atmosphere into a 5 kilometers deep ocean. The impactors were composed of mockup mantle material (dunite with a density of 3.32 g/cc) or iron (density of 7.8 g/cc). as a mockup for nickel-iron asteroids, For these projectiles, analytical Mie-Gruneisen equations of state were used. An elastic-plastic strength model was used for the crust and asteroid materials. Water was described using the SAIC/Pactech tables. An example montage from the two-dimensional parameter study is shown in Figure 6.21.

Numerical modeling of water waves

198

Fig. 6.21. A 1 kilometer iron vertical impactor craters the basalt crust, excavates a cavity in the ocean with a diameter of 25 kilometers.

Modeling waves using compressible models

199

The density plot shortly after the collapse of the transient water crater is shown in Figure 6.22.

Fig. 6.22. The wave train density plots for a 1 kilometer iron asteroid impacting 5 kilometers of water and a basalt crust. The complex wave train is a result of the interaction of multiple shocks propagating through the water and the basalt crust.

The complex wave train is a result of the reflections and interaction of multiple shocks propagating through the water and basalt crust. A complicated time-dependent flow occurs which is best studied using the PowerPoint presentation and computer movies on

Numerical modeling of water waves

200

the NMWW CD-ROM in the directory /NOBEL/ASTWAVE/. It becomes obvious that our previous attempts to model the tsunami wave generation by asteroids using incompressible models were inadequate for the task. A tabular summary of the parameter study is presented in Table 6.2. Listed are the impact characteristics of the asteroid (composition, diameter, density, mass, velocity and kinetic energy) and the measured characteristics of the impact (maximum depth and diameter of the cavity, quantity of water displaced, time of maximum cavity, maximum jet and jet rebound, tsunami wave length and velocity).

TABLE 6.2 Asteroid Generated Tsunamis

Asteroid Material

Dunite

Iron

Dunite

Iron

Dunite

Iron

Asteroid Diameter

250 m

250 m

500 m

500 m

1000 m

1000 m

Asteroid Density

3.32 g/cc

7.81 g/cc

3.32 g/cc

7.81 g/cc

3.32 g/cc

7.81 g/cc

Asteroid Mass

2.7e13 g

6.4e13 g

2.2e14 g

5.1e14 g

1.7e15 g

4.1e15 g

Asteroid Velocity

20 km/s

20 km/s

20 km/s

20 km/s

20 km/s

20 km/s

Kinetic Energy

1.3 GT

3.1 GT

10.4 GT

24.4 GT

83 GT

195 GT

Max Cavity Dia

4.4 km

5.2 km

10.0 km

12.6 km

18.6 km

25.2 km

Max Cavity Depth

2.9km

4.3 km

4.5 km

5.7 km

6.6 km

9.7 km

Water Displacement

4.4e16 g

9.1e16 g

3.5e17 g

7.1e17 g

1.8e18 g

4.8e18 g

Time of Max Cavity

13.5 s

16.0 s

22.5 s

28.0 s

28.5 s

33.0 s

Time of Max Jet

54.5 s

65.0 s

96.5 s

111 s

128 s

142 S

Time of Rebound

100.5 S

118.5 S

137.5 S

162 s

187 s

218 s

9 km

12 km

17 km

20 km

23 km

27 km

120 m/s

140 m/s

150 m/s

160 m/s

170 m/s

175 m/s

Tsunami Wavelength Tsunami Velocity

The shallow water tsunami velocity in 5 kilometers deep water is 221 meters/sec. The amount of water displaced during the formation of the cavity was found to scale linearly with the kinetic energy of the asteroid, as illustrated in Figure 6.23. A fraction of this displaced mass is actually vaporized during the explosive phase of the encounter, while the rest is pushed aside by the pressure of the vapor to form the crown and rim of the transient cavity.

Modeling waves using compressible models

201

Fig. 6.23. The mass of water displaced in the initial cavity formation scales with the asteroid kinetic energy. The squares are the results from the parameter study tabulated in Table 6.2. The solid line shows direct proportionality. A fraction (≈5–10 percent) of this mass is vaporized at impact.

The tsunami amplitude evolves in a complex manner, eventually decaying faster than the expected recipocal of the distance of propagation from the impact point as shown in Figure 6.24. The wave trains are highly complex because of the multiple shock reflections and interactions involving the seafloor. The complexity of the wave breaking associated with the reflection and interaction of shocks results in a very dynamic maximum wave height. Realistic seafloor topography would undoubtedly influence the development of the water wave.

Numerical modeling of water waves

202

Fig. 6.24. The tsunami amplitude has a very complex early time profile. It eventually decreases with distance faster than 1/(Radius). The legend identifies the points associated with individual cases. The notation signifies the asteroid composition (“Dn” for dunite and “Fe” for iron) and the diameter in meters. A “tr” indicates calculated tracer heights and a “Isq” indicates a least square fitted line to the tracer heights.

For an ocean of 5 kilometers depth, the shallow water velocity is 221 meters/sec. In Figure 6.25 are shown the wave crest positions as a function of time for the simulations in the parameter study, along with constant velocity lines of 150 and 221 meters/sec. The wave velocities are substantially lower than the shallow water limit. The tsunami wave length is found to scale roughly with the 1/4 power of the asteroid kinetic energy. The reason for this is that the wavelength is determined by the cavity-jetrebound cycle. The time scale for this varies as where jh is the mean jet height and g is gravity. The mean jet height varies as the square root of the asteroid kinetic energy.

Modeling waves using compressible models

203

Fig. 6.25. The tsunami wave crest positions as a function of time for the parameter study cases. The notation is the same as for Figure 6.24. The solid lines of constant velocity illustrate that these waves are substantially slower than the shallow water prediction.

Three-Dimensional Oblique Impact Asteroid Generated Tsunamis Three-dimensional simulations of a 1 kilometer diameter iron asteroid impacting the ocean at a 45 degree angle at 20 kilometers/sec were performed using the ASCI White machine. The calculations used up to 1200 processors and several weeks of computing time. Up to 200,000,000 computational cells were used, with a total computation time of 1,300,000 cpu-hours. The computational volume was a rectangular box 200 kilometers long in the direction of the asteroid trajectory, 100 kilometers wide and 60 kilometers tall. The height was divided into 42 kilometers of atmosphere, 5 kilometers of ocean water, 7 kilometers basalt crust and 6 kilometers of mantle material. The asteroid was started at a point 30 kilometers above the surface of the water as shown in Figure 6.26. The atmosphere was a standard exponential atmosphere, so the air initially surrounding the asteroid was only 1.5 percent of the sea level density.

Numerical modeling of water waves

204

Fig. 6.26. A 1 kilometer diameter iron impactor at an angle of 45 degrees craters a 5 kilometers deep ocean and the basalt crust. Twodimension slices in the vertical plane containing the asteroid trajectory are shown. The initial asymmetry disappears with time.

During the 2.1 sec of the asteroid’s atmospheric passage at about Mach 60, a strong bow air shock develops, heating the air to temperatures of 1 electron volt (11,600 K). Less than 1 percent of the asteroid’s kinetic energy is dissipated in the atmospheric passage. The water is much more effective at slowing the asteroid, and essentially all of its kinetic energy is absorbed by the ocean and seafloor within 0.7 sec. The water immediately surrounding the trajectory is vaporized. The rapid expansion of the vapor cloud evacuates a cavity in the water that eventually expands to a diameter of 15 kilometers. This initial cavity is asymmetric because of the inclined trajectory of the asteroid. The splash, or crown, is markedly higher on the side opposite the incoming trajectory, as shown in Figure 6.27.

Modeling waves using compressible models

205

Fig. 6.27. A 1 kilometer diameter iron asteroid impacted the ocean from the right at 45 degrees. A perspective cutaway view from the same run as in Figure 6.26 at the time of maximum cavity (27.5 sec after impact) is shown. The density profiles are shown. The cut passes through the ocean to the basalt ocean floor. The diameter of the cavity is about 25 kilometers.

The maximum height of the crown on the downstream side is nearly 30 kilometers at 70 sec after impact. The collapse of the bulk of the crown makes a “rim wave” or precursor tsunami that propagates outward, somewhat higher on the downstream side. The higher portion of the crown breaks up into droplets that fall back into the water giving this precursor tsunami a very uneven and asymmetric profile. The rapid dissipation of the asteroid’s kinetic energy is very much like an explosion. Shocks propagate outward from the cavity in the water, in the basalt crust and in the mantle beneath. Multiple reflections of shocks and acoustic waves between the material interfaces complicate the dynamics. The hot vapor from the initial cavity expands into the atmosphere, mainly downstream because of the momentum of the asteroid. When the pressure of the vapor in the cavity has diminished sufficiently, at about 35 sec after the impact, water begins to fill the cavity from the bottom. This filling has a high degree of symmetry because of the uniform gravity responsible for the water pressure. An asymmetric fill could result from non uniform seafloor topography, but that case was not modeled. The filling water converges on the center of the cavity and the implosion produces another series of shock waves, and a jet that rises vertically in the atmosphere to a height in excess of 20 kilometers at a time of 115 sec after impact, as shown in Figure 6.28. It is the collapse of this central vertical jet that produces the principal tsunami wave shown in Figure 6.29.

Numerical modeling of water waves

206

Fig. 6.28. Similar to Figure 6.27 but 115 sec after impact, which is the time of the formation of the central vertical jet. The collapse of the crown has produced a circular rim that is propagating outward, but the principal tsunami wave will be produced by the collapse of the central vertical jet. The crown splash has collapsed, pockmarking the surface of the first wave.

The evolution of this wave was modeled out to a time of 400 sec after impact. The inclined impact eventually produces a tsunami that is nearly circularly symmetric at late times, as shown in Figure 6.29. The tsunami has an initial height in excess of 1 kilometer, and decreases to 100 meters at a distance of 40 kilometers from the initial impact. Its propagation speed is 175 meters/sec, which is much less than the shallow water speed of 221 meters/sec.

Modeling waves using compressible models

207

Fig. 6.29. Similar to Figure 6.27 and 6.28 but at 270 sec. The central jet has collapsed and formed the tsunami wave. The pock-marked precursor wave and the smoother principal wave are evident. The principal wave is about 1.5 kilometers in initial amplitude and moves with a speed of 175 meters/sec.

The evolution of this wave was modeled out to a time of 400 sec after impact. The inclined impact eventually produces a tsunami that is nearly circularly symmetric at late times as shown in Figure 6.29. The tsunami has an initial height in excess of 1 kilometer, and decreases to 100 meters at a distance of 40 kilometers from the initial impact. Its propagation speed is 175 meters/sec, which is much less than the shallow water speed of 221 meters/sec. The 45 degree angle chosen for the three-dimensional simulation is the most probable angle for impacts. Impacts at 30 degrees and 60 degrees were also modeled. Significant differences with impact angle were found for the KT impactor described in reference 17 and in the next section of this chapter. Conclusions All asteroid impacts, including the oblique ones, produce a large underwater cavity with nearly vertical walls followed by a collapse starting from the bottom and subsequent vertical jetting. Substantial amounts of water are vaporized and lofted into the atmosphere. In the larger impacts, significant amounts of the crustal and even mantle material are lofted as well. Tsunamis up to a kilometer in initial height are generated by the collapse of the vertical jet. These waves are initially very complex in form and interact strongly with shocks propagating through the water and the crust. The tsunami waves were followed out to 100 kilometers from the point of impact. Their periods and wavelengths show them to be intermediate type waves and not shallow water waves. At

Numerical modeling of water waves

208

great distances, the waves decay as the inverse of the distance from the impact point, ignoring sea floor topography. For all impactors smaller than about 2 kilometers diameter, the impacting body is highly fragmented and it remains lofted into the stratosphere with the water vapor and crustal material, hence very little trace of the impacting body will be found for most oceanic impacts. In the oblique impacts, the initial asymmetry of the transient crater and crown does not persist beyond a tsunami propagation length of 50 kilometers.

6E. KT Chicxulub Event Introduction The impact that created the Chicxulub crater in Mexico’s Yucatan Peninsula is widely believed to be responsible for the mass extinctions at the end of the Cretceous period. This event has been extensively studied and modeled during the few years since its discovery. The stratigraphy at Chicxulub involves rocks like calcite and anhydrite that are highly volatile at pressures reached during impact. The volatility of the Chicxulub strata made the event very dangerous to the megafauna of the late Cretaceous. On a geological time scale, impacts of asteroids and comets with the earth are a relatively frequent occurrence, causing significant disturbances to biological communities and strongly perturbing the course of evolution. Most famous among known catastrophic impacts is the one that marked the end of the Cretceous period and the dominance of the dinosaurs. It is now widely accepted that the world sequence of mass extinctions at the Cretaceous-Tertiary (KT) boundary 65 million years ago was directly caused by the collision of an asteroid or comet with the earth. Evidence for this includes the large (200 kilometers diameter) buried impact structure at Chicxulub, Yucatan, Mexico, the world wide distributed Iridium layer at the KT boundary, and tsunami deposits well inland in North America, all dated to the same epoch as the extinction event. As described in reference 17, the KT impactor was an asteroid of diameter roughly 10 kilometers. Its impact was oblique, either from the SE at 30 degrees to the horizontal or from the SW at 60 degrees. The asteroid encountered layers of water, anhydrite, gypsum, and calcium carbonate, which resulted in the lofting of many hundreds of cubic kilometers of these materials into the stratosphere, where they resided for many years and produced a global climate cooling that was fatal to many large animal species on earth. Three-Dimensional KT Impact Model The Chicxulub impact structure was discovered by scientists from Petroleos Mexicanos (Pemex), the Mexican national oil company18, 19. It is believed to be the location of an impact that was responsible for the mass extinction at the end of the Cretaceous period, as proposed by Alvarez20 et al. on the basis of the anomaly in abundances of Iridium and other Platinum group elements in the KT boundary bedding plane. Paleogeographic data suggest that the crater site, which presently straddles the Yucatan

Modeling waves using compressible models

209

coastline, was submerged at the end of the Cretaceous period on the continental shelf. The substrate consisted of fossilized coral reef over continental crust. In the Gisler, et al.17 simulation, the multilayered target was 300 meters of water, 3 kilometers of calcite, 30 kilometers of granite and 18 kilometers of mantle material. Above the target was an exponential atmosphere up to 106 kilometers altitude. The asteroid’s plunge was numerically started at the 106 kilometers altitude. Threedimensional simulations were performed with impact angles of 30, 45 and 60 degrees to the horizontal. The computational domain extended 256 kilometers by 128 kilometers and simulated the half space of 256 by 256 kilometers. The calculations were performed on the new ASCI Q computer at Los Alamos, a cluster of ES45 Alpha boxes from HP/Compaq. The problems used 1024 processors and a total of about 1 million cpu hours. The adaptive mesh included up to a third of a billion computation cells. Three prominent features of the simulation are illustrated by the 45 degree impact shown in Figure 6.30.

Fig. 6.30. The density at seven sec after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. The isosurface, at density 0.005 g/cc is chosen to show everything denser than air. The scale is set by the back boundary which is 256 kilometers long. The maximum height of the uplifted material is 50 kilometers.

After impact, billions of tons of very hot material are lofted into the atmosphere. A water plume directed away from the impact angle carries much of the horizontal component of the asteroid’s momentum in the downrange direction. This material, consisting of vaporized fragments of the projectile mixed with the target, is extremely hot, and will ignite vegetation many hundreds of kilometers away from the impact site. The highly turbulent and energetic ejecta plume is directed predominatly upward as

Numerical modeling of water waves

210

shown in Figure 6.31.

Fig. 6.31. The density at 42 sec after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. The isosurface, at density 0.005 g/cc is chosen to show everything denser than air. The scale is set by the back boundary which is 256 kilometers long. The maximum height of the uplifted material is 50 kilometers.

Some material is projected into orbits that terminate far outside the computational volume. The dissipation of the kinetic energy, some 300 teratons TNT equivalent, produces a stupendous explosion that melts, vaporizes and ejects a substantial volume of water and granite. Ballistic trajectories carry some of this material back to earth in a conical debris curtain that gradually moves away from the crater lip and deposits a blanket of ejecta around the forming crater as shown in Figure 6.32. The debris curtain has separated from the rim of the still forming crater as material in the curtain falls to earth. The debris from the curtain is deposited in a blanket of ejecta that is asymmetric around the crater with more in the downrange than in the uprange direction. The distribution of material in the ejecta blanket was used as a diagnostic to determine the direction and angle of impact of the asteroid. Some material is projected into orbits that terminate far outside the computational volume.

Modeling waves using compressible models

211

Fig. 6.32. The density at two minutes after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. The isosurface, at density 0.005 g/cc is chosen to show everything denser than air. The scale is set by the back boundary which is 256 kilometers long. The debris curtain has separated from the crater.

The blanket of ejecta is found to be strongly asymmetrical around the crater, with the uprange portion much thinner than the rest. This is a result of the coupling of the horizontal component of the asteroid’s momentum to the debris, and to the ionized and shocked atmosphere in the asteroid’s wake.

Numerical modeling of water waves

212

Fig. 6.33. The density at 122 sec after a 10 kilometers diameter granite asteroid impacts the earth at 45 degrees. A cut through the center of the impact cavity shows that no water is present from the center of the crater out to the debris curtain.

As shown in Figure 6.33, the ocean water is mostly vaporized in the crater region and beyond the crater to the debris curtain which has separated from the rim of the crustal crater. The initial shallow 300 meters deep water where the asteroid impact occurred would not support the generation of major long period tsunami waves. However, the ocean water is driven by strong shock waves and the expanding curtain out into surrounding deep ocean. This would result in large, long period tsunamis that would propagate around the world and leave the tsunami deposits well inland in North America, as observed during the KT epoch. Spectacular computer movies of the KT impact and a PowerPoint presentation produced by Dr. Galen Gisler are on the NMWW CD-ROM in the directory /NOBEL/KTIMPACT/. Dr. Gisler performed the KT asteroid impact calculations as part of the Los Alamos Crestone Project validation studies. The NOBEL/SAGE/RAGE hydrodynamic codes and the 3-D graphics were developed by the Crestone Project team. The Chicxulub PowerPoint presentation has been shown to many members of the legislative and executive branches of the U.S. government and to Los Alamos visiting dignitaries including Prince William.

Modeling waves using compressible models

213

If life on earth is to continue we must develop the capability of deflecting asteroids before they collide with the earth. Solem21 has described the options available for deflection of asteroids and comets. He concludes “the most cost effective and the ONLY currently available means of asteroid disruption (deflection or pulverization) is a nuclear explosive.” So the Los Alamos nuclear weapon program that has developed most of the water wave modeling technologies described in this book as part of its nuclear weapon effects studies also has developed the nuclear technology that is both a major threat to life on earth and that which makes it possible to defend the earth from the extinction of life by asteroids and comets.

References 1. Charles L.Mader, Numerical Modeling of Explosives and Propellants, CRC Press, Boca Raton, Florida (1998). 2. M.L.Gittings, “1992 SAIC’s Adaptive Grid Eulerian Code,” Defense Nuclear Agency Numerical Methods Symposium, pp. 28–30 (1992). 3. R.L.Holmes, G.Dimonte, B.Fryxell, M.L.Gittings, J.W. Grove, M.Schneider, D.H.Sharp, A.L.Velikovich, R.P.Weaver and Q.Zhang, “Richtmyer-Meshkov Instability Growth: Experiment, Simulation and Theory,” Journal of Fluid Mechanics, Vol. 9, pp.55–79 (1999). 4. R.M.Baltrusaitis, M.L.Gittings, R.P.Weaver, R.F.Benjamin and J.M.Budzinski, “Simulation of Shock-Generated Instabilities,” Physics of Fluids, Vol.8, pp.2471–2483 (1996). 5. Don J.Miller, “Giant Waves in Lituya Bay, Alaska” Geological Survey Professional Paper 354-C, U.S. Government Printing Office, Washington, D.C. (1960). 6. Frances E.Calwell, Land of The Ocean Mists-The Wild Ocean Coast West of Glacier Bay, Alaska Northwest Publishing Company, Edmonds, Washington (1986). 7. Charles L.Mader, “Modeling the 1958 Lituya Bay Tsunami,” Science of Tsunami Hazards, Vol.17, pp.57–67 (1999). 8. Hermann M.Fritz, Willi H.Hager and Hans-Erwin Minor, “Lituya Bay Case: Rockslide Impact and Wave Runup,” Science of Tsunami Hazards, Vol.19, pp.3–22 (2001). 9. Charles L.Mader, “Modeling the 1958 Lituya Bay Mega-Tsunami, II,” Science of Tsunami Hazards, Vol.20, pp.241–250 (2002). 10. Bobby G.Craig, “Experimental Observations of Underwater Detonations Near the Water Surface,” Los Alamos Scientific Laboratory report LA-5548-MS (1972). 11. Donald E.Gault and Charles P.Sonett, “Laboratory Simulation of Pelagic Asteroid Impact: Atmospheric Injection, Benthic Topography, and the Surface Wave Radiation Field,” Geological Society of America, Special Paper 190, pp.69–92 (1982). 12. Valery Kedrinskii, private communication (1985). 13. Charles L.Mader and Michael L.Gittings, “Dynamics of Water Cavity Generation,” Science of Tsunami Hazards, Vol.21, pp. 91–118 (2003). 14. Charles L.Mader, Numerical Modeling of Detonations, University of California Press, Berkeley, California (1979). 15. Charles L.Mader, Timothy R.Neal and Richard D.Dick, LASL Phermex Data, Volume I, II, and III, University of California Press, Berkeley, California (1980). Available at http://lib-www.lanl.gov/ladcdmp/phl.pdf. 16. Galen Gisler, Robert Weaver, Michael L.Gittings and Charles Mader, “Two- and

Numerical modeling of water waves

214

Three-Dimensional Simulations of Asteroid Ocean Impacts,” Science of Tsunami Hazards, Vol.21, pp.119–134 (2003). 17. Galen Gisler, Robert Weaver, Michael L.Gittings and Charles Mader, “Two- and Three-Dimensional Asteroid Impact Simulations,” Computers in Science and Engineering (2004). 18. A.R.Hildebrand, G.T.Penfield, D.A.Kring, M.Pilkington, Z.A.Camargo, S.B.Jacobsen and W.V.Boynton, “Chicxulub Crater: A Possible Cretaceous/Tertiary Boundary Impact Crater on the Yucatan Peninsula, Mexico.” Geology, Vol.19, pp.867–871 (1991). 19. V.L.Sharpton, G.B.Dalrymple, L.E.Marin, G.Ryder, B. C.Schuraytz and J.UrruitaFucugauchi, “New Links Between the Chicxulub Impact Structure and the Cretaceous/Tertiary Boundary,” Nature, Vol.359, pp.819–821 (1992). 20. L.W.Alvarez, W.Alvarez, F.Asaro and H.V.Michel, “Extraterrestrial Cause for the Cretaceous/Teritary Extinction,” Science, Vol.208, pp.1095–1108 (1980). 21. Johndale C. Solem, “Comet and Asteroid Hazards: Threat and Mitigation,” Science of Tsunami Hazards, Vol.17, pp.141–154 (1999).

CD-ROM CONTENTS Tsunami Animations A collection of tsunami animations performed using the SWAN code described in the Numerical Modeling of Water Waves—Second Edition, by Dr. Charles L.Mader. A few calculations were performed using the full Navier-Stokes ZUNI or SOLA codes. The tsunami animations are archived at http://tl4web.lanl.gov/Staff/clm/tsunami.mve/tsunami.htm. To view animation type or click—MOVIE—To study individual cases type—MM TAPEXX where XX is file number. Follow instructions on screen. 1960.MVE—May 23, 1960 tsunami generation in Chile, propagation across the Pacific Ocean, and indundation of Hilo, Hawaii. Described in “Modeling Hilo, Hawaii Tsunami Inundation,” Science of Tsunami Hazards, Vol.9, pp.85–94 (1991), and Scientific Computing and Automation, June issue, pp.19–23 (1993). 1964.MVE—Tsunami of April 1, 1964 generation in Gulf of Alaska, propagation across the Pacific Ocean, and inundation of Crescent City, California. See “Tsunami Inundation Model Study of Eureka and Crescent City, California,” NOAA Tech. Memo. ERL PMEL-105 (1994). 60THEAT.MVE—The interaction of the tsunami of May 23, 1960 with the Hilo, Hawaii Theater. Described in PACON 1993. 90HILO.MVE—1990 Hilo topography and buildings inundated by a 1960 tsunami wave. See also HOTEL.MVE. 2ATAST.MVE—The inundation of the U.S. East Coast by a 100 meters, 2000 sec tsunami wave that could be generated by an asteroid. 10NYAST.MVE—The inundation of the U.S. East Coast by a wave from the incompressible collapse of a 10 kilometers radius cavity, 3000 meters deep and a 100 kilometers radius cavity in the Atlantic Ocean off New York City. AIMPACT.MVE—An impact cavity collapse and tsunami generation study using shallow water and full Navier-Stokes models. Described in “Modeling Asteroid Impact and Tsunami,” Science of Tsunami Hazards, Vol.16, pp.21–30 (1998). ATLAST.MVE—A tsunami in the Atlantic Ocean generated by the incompressible collapse of a cavity 150 kilometers wide and 3500 meters deep. AUSAST.MVE—Interaction of a tsunami with Australia from a Hawaii landslide generated tsunami and from a cavity collapse generated tsunami. Described in “Modeling of Tsunami Propagation Directed at Wave Erosion on Southeastern Australia Coast 105,000 Years Ago,” Science of Tsunami Hazards, Vol.13, pp.45–52 (1995). BBAY.MVE—A study of the vulnerability of Berau Bay, Indonesia to tsunamis. CASCAD.MVE—Inundation of U.S. West Coast by a tsunami from the Cascadia fault. ECAST.MVE—The inundation of the U.S. East Coast by a tsunami generated by the incompressible collapse of a 150 kilometers wide, 3000 meters deep cavity. See also

CD-Rom contents

217

NYAST.MVE ELTAST.MVE—Described in “Modeling the Eltanin Asteroid Tsunami,” Science of Tsunami Hazards, Vol.16, pp.17–20 (1998). EURAST.MVE—The inundation of Europe by a 100 meters high and 2000 sec period tsunami. EUREKA.MVE—The Eureka, California tsunami of April 25, 1992. See “Tsunami Inundation Model Study of Eureka and Crescent City, California,” NOAA Tech. Memo. ERL PMEL-105 (1994). GUS.MVE—The Furumoto sources for the Hawaiian tsunamis of 1946, 1957, 1964 and 1965. Part of a source modeling project for Dr. A.Furumoto, Hawaii Civil Defense Tsunami Advisor. HIAST.MVE—The inundation of the Hawaiian Islands by a 100 meters high, 2000 sec period tsunami wave. Described in “Asteroid Tsunami Inundation of Hawaii,” Science of Tsunami Hazards, Vol. 14, pp.85–88 (1996). HILAND.MVE—The tsunami generated by a landslide off the Kona coast of the island of Hawaii about 105 Ka years ago. Described in “Modeling the 105 Ka Landslide Lanai Tsunami,” Science of Tsunami Hazards, Vol.12, pp.33–38 (1994). HKAI.MVE—Inundation of Hawaii Kai, Hawaii by a typical off shore 3 meters high, 1500 sec tsunami wave. HOTEL.MVE—The interaction of a May 23, 1960 tsunami wave with current Hilo, Hawaii tourist hotels. See also 90HILO.MVE. HUMBOL.MVE—Tsunami inundation of Humboldt Bay, California by an offshore maximum expectable 10 meters high, 2000 sec tsunami wave. See “Tsunami Inundation Model Study of Eureka and Crescent City, California,” NOAA Tech. Memo. ERL PMEL-105 (1994). ICEAST.MVE—The inundation of Iceland by a 100 meters high and 2000 sec period tsunami. INDIA.MVE—Tsunami in the Indian Ocean generated by the incompressible collapse of a cavity 38 kilometers wide and 4000 meters deep. INDONES.MVE—Indonesia tsunami of December 12, 1992. JAPAST.MVE—The inundation of Tokyo, Japan by a tsunami generated by a incompressible cavity collapse. Described in “Asteroid Tsunami Inundation of Japan,” Science of Tsunami Hazards, Vol.16, pp.11–16 (1998). KAIAKA.MVE—Tsunami inundation of Kaiaka Bay, Oahu, Hawaii by the 1952 tsunami. KBAY.MVE—Tsunami inundation of Kaneohe Bay, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami and by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave. KONA.MVE—Tsunami inundation of Kona, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami wave. KURIL.MVE—The tsunami of October 1994 generated off the Kuril islands of Japan. LAAST.MVE—Inundation of Los Angeles, California by a 100 meters high, 2000 sec period tsunami wave. LAPALMA.MVE—Modeling the proposed La Palma landslide tsunami. Published in “Modeling the La Palma Landslide Tsunami,” Science of Tsunami Hazards, Vol.19,

CD-Rom contents

218

pp.160–180 (2001). LAUP.MVE—The April 1, 1946 tsunami inundation of Laupahoehoe, Hawaii. LITUYA.MVE—The July 8, 1958 mega-tsunami at Lituya Bay, Alaska with inundations up to 520 meters. Described in “Modeling the 1958 Lituya Bay Mega-Tsunami,” Science of Tsunami Hazards, Vol.17, pp.57–67 (1999). The Lituya Bay impact landslide generation of the tsunami is described in Chapter 6 and in Science of Tsunami Hazards, Vol.20, pp.241–250 (2002). LISBON.MVE—Modeling the 1755 Lisbon tsunami generation and propagation across the Atlantic Ocean to the Caribbean. Science of Tsunami Hazards, Vol.19, pp.93–98 (2001). LOIHI.MVE—A study using the ZUNI full Navier-Stokes code of the tsunami wave generation and propagation from the collapse of the Loihi, Hawaii summit in August, 1996. M9CALIF.MVE—An M9 earthquake generated tsunami interacting with the Oregon and California coasts. NIC.MVE—The tsunami generated off the coast of Nicaragua in 1992. Described in “Modeling the 1992 Nicaragua Tsunami,” Science of Tsunami Hazards, Vol.11, pp.107– 110 (1993). NYAST.MVE—The inundation of the U.S. East Coast by the incompressible collapse of a 100 kilometers radius 3000 meters deep cavity. Another tsunami wave had a height of 100 meters and a 2000 sec period. See also 10NYAST.MVE. ORAST.MVE—A 100 meters high, 2000 sec period tsunami interacting with the Oregon coast. OREGM9.MVE—An M9 earthquake generated tsunami interacting with the Oregon coast. PACAST.MVE—Tsunami in the middle of the Pacific formed from the incompressible collapse of a cavity 150 kilometers wide and 4500 meters deep. PROP.MVE—Described in “Numerical Tsunami Propagation Study,” Science of Tsunami Hazards, Vol.11, pp.93–106(1993) and in Chapter 5 of Numerical Modeling of Water Waves—Second Edition. SANDY.MVE—Tsunami inundation of Sandy Beach region of Oahu, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami and by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave. SANFAST.MVE—Inundation of San Francisco, California by a 100 meters high, 2000 sec tsunami wave. SKAGWAY.MVE—The landslide generated tsunami of November 3, 1994 at Skagway, Alaska. The Skagway modeling is described in “Modeling the 1994 Skagway Tsunami,” Science of Tsunami Hazards, Vol.15, pp.41–48 (1997). See also SOLA.MVE. SMSFAST.MVE—Inundation of San Francisco by a tsunami wave generated by the incompressible collapse of a 20 kilometers wide, 3000 meters deep cavity. SOLA.MVE—Three-dimensional, full Navier-Stokes modeling using the MCC SOLA code of the November 3, 1994 Skagway, Alaska tsunami. See also SKAGWAY.MVE. SOURCE.MVE—Described in “Numerical Tsunami Source Study,” Science of Tsunami Hazards, Vol.11, pp.81–92 (1993) and in Chapter 5 of Numerical Modeling of Water Waves—Second Edition.

CD-Rom contents

219

VSLIDE.MVE—A landslide generated tsunami from the Chain of Craters road region of the island of Hawaii. WAIANAE.MVE—The inundation of the leeward side of Oahu, Hawaii by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave. WAIPIO.MVE—The interaction of the May 23, 1960 tsunami with the Waipio, Hawaii region. The 50 foot inundation is the largest recorded in Hawaii. WALKER.MVE—An evaluation of the vulnerability of Hawaii to tsunamis generated south of Honolulu, either along the Kona Coast or in the Tonga trench. Modeling requested by Dr. D.Walker, Oahu Civil Defence Tsunami Advisor. WINDWARD.MVE—Tsunami inundation of the Windward side of Oahu, Hawaii by a typical offshore 3 meters high, 2000 sec tsunami and by a maximum expectable offshore 10 meters high, 2000 sec tsunami wave.

NOBEL PowerPoint Presentations A collection of PowerPoint presentations describing water wave studies performed using the compressible hydrodynamic code NOBEL. The studies are described in Chapter 6 of Numerical Modeling of Water Waves—Second Edition and archived at http://t14web.lanl.gov/Staff/clm/tsunami.mve/tsunami.htm. The PowerPoint presentations may be viewed using PPVIEW in /NOBEL/PPRESENT/. LITUYA—The July 8, 1958 Lituya Bay, Alaska impact landslide tsunami generation. A mega-tsunami was generated that reached an altitude of 520 meters. Laboratory experiments and numerical modeling results are presented. Described in “Modeling the 1958 Lituya Bay Mega-Tsunami, II,” Science of Tsunami Hazards, Vol. 20, pp.241–250 (2002). CAVITY—The generation of cavities in water by projectile impacts and by explosives is described both experimentally and using compressible hydrodynamic models. Described in “Dynamics of Water Cavity Generation,” Science of Tsunami Hazards, Vol.21, pp. 91–118 (2003). ASTWAVE—The generation of tsunamis by the impact of a 0.25 to 1 kilometer diameter asteroid at 20 kilometers/sec with 5 kilometers of ocean and 5 kilometers of basalt is modeled using compressible hydrodynamics in two and three dimensions. Described in “Two- and Three-Dimensional Simulations of Asteroid Ocean Impacts,” Science of Tsunami Hazards, Vol.21, pp.119–134 (2003). KTIMPACT—The KT Chicxulub asteroid impact event is modeled using the threedimensional compressible Navier-Stokes model. Described in “Two and ThreeDimensional Asteroid Impact Simulations,” Computers in Science and Engineering (2004). CD-ROM CODE DIRECTORIES WAVE—The WAVE code described in Chapter 1 solves the equations for Airy, third-

CD-Rom contents

220

order Stokes and Laitone solitary gravity waves. The directory contains the FORTRAN source code, the executable code for DOS or Windows and WAVE.PDF which describes the code. SWAN—The shallow water SWAN code described in Chapter 2 solves the long wave, shallow water, nonlinear equations of fluid flow. The directory contains the FORTRAN source and executable codes which generate a graphics file that may be processed using the programs included. It also includes a description of the input to the code in the file SWAN.PDF. Examples and topographic files are furnished. ZUNI—The incompressible Navier-Stokes ZUNI code described in Chapter 3 solves the incompressible, viscous fluid flows with a free surface using the Navier-Stokes equations. A detailed description of the computer program and its input file is included in the file ZUNI.PDF. The FORTRAN source and the executable codes are included. SOLA—The incompressible three-dimensional Navier-Stokes ZUNI code described in Chapter 4 solves the incompressible viscous fluid flows with a free surface using the Navier-Stokes equations. The FORTRAN source and the executable codes are included. The Skagway 1994 tsunami is used as an example. LGW—The Carrier linear gravity wave LGW code described in Chapter 5 uses analytical methods for solving the linear gravity model. The FORTRAN source and executable codes are included. Examples of Gaussian tsunamis described in Chapter 5 are furnished. TIDE—A classic computer program for calculating tides with the FORTRAN source and executable codes furnished.

Science of Tsunami Hazards Directory All the Science of Tsunami Hazards journals through 2003 in PDF format may be searched using Adobe Acrobat 4.0 or higher. Issues of the journal are archived at http://epubs.lanl.gov/tsunami and current issues are available at http://www.sthjournal.org. Dir=TS211.PDF, TS212.PDF, TS213.PDF, TS214.PDF * Volume 21 (No.1), (No.2), (No.3), (No.4) 2003 Dir=TS201.PDF, TS202.PDF, TS203.PDF, TS204.PDF, TS205.PDF * Volume 20 (No.1), (No.2), (No.3), (No.4), (No. 5), 2002 DIR=TS191.PDF, TS192.PDF, TS193.PDF * Volume 19 (No.1) (No.2) (No.3), 2001 DIR=TS181.PDF, TS182.PDF * Volume 18 (No.1) (No.2), 2000 DIR=TS171.PDF, TS172.PDF, TS173.PDF * Volume 17 (No.1) (No.2) (No.3), 1999 DIR=00394718.PDF * Volume 16 (No.1), 1998 DIR=00394719.PDF, 00394720.PDF

CD-Rom contents

221

* Volume 15 (No.1) (No.2), 1997 DIR=00394721.PDF, 00394722.PDF, 00394723.PDF * Volume 14 (No.3) (No.2) (No.1), 1996 DIR=00394724.PDF * Volume 13 (No.1), 1995 DIR=00394725.PDF, 00394726.PDF * Volume 12 (No.2) (No.1), 1994 DIRECTORY=00394727.PDF, 00394728.PDF * Volume 11 (No.2) (No.1), 1993 DIRECTORY=00394729.PDF * Volume 10 (No.1), 1992 DIRECTORY=00394730.PDF, 00394731.PDF * Volume 9 (No.2) (No.1), 1991 DIRECTORY=00394732.PDF, 00394733.PDF * Volume 8 (No.2) (No.1), 1990 DIRECTORY=00394734.PDF, 00394735.PDF * Volume 7 (No.2) (No.1), 1989 DIRECTORY=00394736.PDF * Volume 6 (No.1), 1988 DIRECTORY=00394737.PDF, 00394738.PDF * Volume 5 (No.2) (No.1), 1987 DIRECTORY=00394739.PDF, 00394740.PDF, 00394741.PDF * Volume 4 (No.3) (No.2) (No.1), 1986 DIRECTORY=00394742.PDF * Volume 3 (No.1), 1985 DIRECTORY=00394743.PDF, 00394744.PDF * Volume 2 (No.2) (No.1), 1984 DIRECTORY=00394745.PDF * Volume 1 (No.1), 1982 CNMWW.PDF is a searchable PDF file of the book Numerical Modeling of Water Waves—Second Edition with many figures in color.

AUTHOR INDEX Numbers in italic type indicate pages where references are listed in full.

Alvarez, L.W., 224, 231 Alvarez, W., 224, 231 Amsden, A.A., 86, 118 Asaro, F., 223, 231 Baltrusaitis, R.M., 180, 229 Batchelor, G.K., viii, 21 Benjamin, R.F. 180, 229 Bernard, E., 56, 82 Bornhold, B.D., 65, 82 Boynton, W.V., 223, 230 Bretschneider, C.L., 166, 177 Budzinski, J.M., 180, 229 Butler, H.L., 39, 82 Calwell, F.E., 183, 229 Camargo, Z.A., 223, 230 Campbell, B., 65, 70, 78, 83 Carrier, G.F., 146, 151, 177 Chan, R.K.C., 86, 92–3, 95, 118 Chubarov, L.B., 26, 30, 80 Cox, D.C., 135, 146 Craig, B.G., 110, 114, 119, 194, 196, 198, 203, 209–10, 230 Curtis, G.D., 26, 31, 48, 56, 81, 82 Dalrymple, G.B., 223, 231 Daly, B.J., 85, 93, 118 Dick, R.D., 203, 230 Dimonte, G., 180, 229 Divoky, D., 39, 44, 46, 82 Dronkers, J.J., 22, 79 Fritz, H.M., 186, 190, 230 Fromm, J.E., 93, 95, 118 Fryxell, B., 180, 229 Fuchs, R.A., 105, 107, 119 Fucugauchi, J.U., 223, 231 Furumoto, A.S., 34, 81

Author index

224

Garcia, A.W., 39, 82 Garcia, W.J., 134, 145 Gault, D.E., 194, 210, 230 Gisler, G., 211, 212, 221, 224, 228, 230 Gittings, M.L., 114, 119, 179, 180, 193, 211, 212, 221, 223, 224, 229, 230 Green, C.K., 34, 81 Greenspan, H.P., 150, 177 Grove, J.W., 180, 229 Hadi, S., 30, 80 Hager, W.H., 186, 230 Hansen, W., 22, 24, 79 Harlow, F.H., 85, 93, 118 Hildebrand, A.R., 223, 230 Hirt, C.W., 86, 92, 94, 118, 119, 145 Holmes, R.L., 180, 229 Houston, J.R., 39, 82 Hwang, L., 39, 44, 46, 82 Jacobsen, S.B., 223, 230 Johnson, J.W., 105, 108, 119 Kedrinskii, V., 194, 195, 199, 210, 230 Kinsman, B., viii, 21 Kowalik, Z., 26, 30, 47, 76, 80, 82 Kring, D.A., 223, 230 Kulikov, E.A., 65, 82 Landau, L.D., viii, 21 Lander, J.F., 64, 82 Le Mehaute, B., 3, 21 Lee, J.J., 64, 66, 76, 77, 82 Leenderstse, J.J., 80 Lifshitz, E.M., viii, 21 Loomis, H.G., 135, 146 Lukas, S., 30, 80 Mader, C.L., 21, 26, 29, 30, 31, 48, 56, 80, 81, 82, 88, 93, 109, 118, 119, 135, 145, 146, 177, 179, 185, 189, 193, 202–4, 208, 211, 212, 221, 223, 224, 229, 230 Madsen, O.S., 97, 100, 118 Marchuk, A.G., 26, 30, 80 Marin, L.E., 223, 231 Mei, C.C., 97, 100, 118 Michel, H.V., 223, 231 Miller, D.J., 183, 186, 229 Minor, H.E., 186, 230

Author index Moore, D.W., 146, 177 Morison, J.R., 105, 107, 119 Murty, T.S., 26, 30, 80 Nabeshima, G., 48, 82 Neal, T.R., 203, 230 Nichols, B.D., 86, 92, 94, 119, 135, 145 Nottingham, D., 65, 69, 70, 78, 83 Pelinovsky, E.N., 26, 30, 80 Penfield, G.T., 223, 230 Peregrine, D.H., 97, 118 Petroff, C., 64, 66, 76, 77, 83 Pilkington, M., 223, 230 Plafker, G., 39, 44, 58, 82 Proudman, J., 22, 79 Rabinovich, A.B., 65, 82 Raichlen, F., 66, 77, 83 Ramming, H.G., 26, 30, 80 Romero, N.C., 119, 145 Ryder, G., 223, 231 Satake, K., 56, 82 Savage, J.C., 44, 82 Schneider, M., 180, 229 Schuraytz, C., 223, 231 Shannon, J.P., 85, 86, 93, 118 Sharp, D.H., 180, 229 Sharpton, V.L., 223, 231 Shokin, I.I., 26, 30, 80 Sklarz, M.A., 134, 146 Sokolowski, T.J., 47, 82 Solem, J.C., 229, 231 Sonett, C.P., 194, 197, 198, 210, 230 Spielvogel, L.Q., 134, 146 Stein, L.R., 119, 145 Street, R.L., 86, 92–3, 96, 118 Tangora, R.E., 135, 145 Thomson, R.E., 65, 82 Uusitalo, S., 24, 80 Van Dorn, W.G., 39, 59, 81 Velikovich, A.L., 180, 229 Vitousek, M., 30, 81

225

Author index Watts, P, 64, 66, 76, 77, 83 Weaver, R.P., 180, 211, 212, 223, 224, 229, 230 Welander, P., 24, 79 Welch, J.E., 85, 93, 118 Whalin, R.W., 39, 82 Whitmore, P.M., 47, 82 Wiegel, R.L., 13, 14, 21 Wybro, P.G., 166, 177 Zhang, Q., 180, 229

226

SUBJECT INDEX 1946, April 1 Tsunami, 31–46, 57–61, 232, 233 1958, July 8, Lituya Bay Tsunami, 181–93 1960 Tsunami Type, 52–3, 232 1960, May 23 Tsunami, 31, 32, 41–8, 52, 81, 155, 162, 231, 235 1964, March 28 Tsunami, 31, 32, 38–46, 56–62, 79–81, 155, 162, 231, 232 1975, November 29, Hawaii Tsunami, 134–42 1994, November 3, Skagway Tsunami, 63–78, 144 Accuracy, 7, 26, 110, 125–6, 134 Airy Wave, viii, ix, 9–12, 20, 29, 94, 146, 148, 150, 158, 159, 163 AMR, 190, 193, 196, 210, 211 ASCI, 177, 212, 218, 224 Asteroid, 186, 194, 210–2, 214–6, 222–37 Barrier, 105–8 Beach, 80, 86, 118, 177, 235 Berau Bay, 232 BlueMountain, 212 Bore, 3, 6, 32, 42, 46, 48 Bottom Motion, 25, 27 Boundary Layer, 4 Bridgewire, 194, 206 Brigham Victory Ship, 39 Bubble Cavity, 206, 209, 210 Bubble Radius, 110 Building, 32, 43, 47–52 Cartesian, 21, 84, 120, 124 Chicxulub, 223, 229, 230, 237 Comet, 223, 229 Compressible Flow, vii, ix, 110, 156 Compressible Model, 176, 178–80 Conservation of Energy, viii, 1 Conservation of Momentum, viii, x, 85, 124 Continental Slope, 30, 93, 95–102 Continuative Boundary, 131, 132 Continuum Boundary, 29 Convergence Criterion, 94, 111 Coriolis, 20, 22, 27, 29, 31, 32 Crestone Project, 228

Subject index

228

DeChezy, 24, 27, 36, 44, 169 Deep Water Wave, 113, 115, 116, 143, 208 Detonation, 110, 115, 194, 203 Dockslide, 68, 69, 75 Donor Cell, 124–6, 134 Earth Rotation, 22 Earthquake, 32, 33, 39, 41, 44, 46, 57–8, 62, 82, 135–7, 142, 153, 155, 162, 176, 181, 182, 185, 211 Estuaries, 3, 30, 80 Eulerian, ix, x, 1, 7, 10, 85, 177, 180, 190, 193, 196, 210, 229 Explosion, 3, 110, 114, 176, 178, 194, 211, 220, 226 Explosive Generated, 195, 202, 208 Extinction of Life, 229 Flood Wave, 3, 6 Flooding, vii, x, 2, 3, 25, 26, 29, 31–62, 156, 162–75 FORTRAN, 237 Free Slip, 85, 130 Gravity, viii, 21, 23, 93, 111, 141, 181, 185, 190, 217, 220 Gravity Waves, 3, 20, 93, 104 Hawaiian Islands, 32, 33, 40, 44, 46, 233 Hilo Bay, 31, 32, 36, 38, 39, 43, 46, 52, 55 Hilo Harbor, 32, 37, 38, 39, 43, 46, 81 Hilo Pier No.1, 38 Hilo Theater, 46–52, 231 Hotel, 32, 52, 55, 232 Hydrostatic, 10, 23, 26, 94, 111 Impact Landslide, 180, 190, 193, 210, 233, 236 Impact Landslide Experiment, 186–9 Incompressible Flow, 85, 118, 130, 177, 194 Indonesia, 30, 232 Intermediate Wave, 7, 9, 222 Jet, 115, 194–5, 207, 210, 213, 217, 221, 222 JIMAR, 81, 146 Kinetic Energy, 114, 215–6, 219, 226 KT Impact, 221, 223, 227, 228, 236 Laitone Solitary Waves, 12, 13, 14, 20, 237 Landslide, 30, 65, 68–82, 135–45, 176, 178, 180, 186–93, 210, 211, 232–6 Linear Gravity Waves, 115, 146, 151–62, 237 Lituya Bay, 180–7, 193, 212, 229, 233, 236

Subject index

229

Loihi, 115–6 Long Wave Paradox, 25, 26 Los Alamos, viii, 224, 228, 229 Los Alamos National Laboratory, 118, 145, 194, 209, 212, 231 MAC, 85, 89, 92–3, 118, 119, 124, 133 Marker and Cell, 85, 119, 134 Mega-tsunami, 181, 193, 211, 212, 233, 236 Momentum, 89, 123, 124, 129, 131–2, 134, 225, 227 Musi Upang, 30, 80 National Academy of Science, 81 Navier-Stokes, vii, 30, 31, 47, 57, 78, 83, 85, 87, 93, 110113, 116, 119, 120, 134, 135, 140–2, 146, 153, 155, 156, 159–63, 170, 174–6, 189, 196, 234–7 News Media, 211 NMWW CD-ROM, 11, 13, 14, 20, 26, 27, 31, 52, 62, 75, 78 114, 115, 144, 151, 162, 191, 210, 215, 227 No-slip, 130, 133 NOAA, 32, 39, 57, 66, 70, 82, 232, 235 NOBEL, 110, 180, 190, 193, 196, 205, 209–11, 215, 228, 236, 237 Nomenclature, viii, 27, 178 NTHMP, 78 Nuclear, 176, 229 Oahu, 234–5 Obstacle, 6, 86, 130, 133 Oscillatory Wave, 3 PARN Dock, 64, 68, 72, 74, 77 PBX-9404 Explosive 110, 203–8 PEMEX, 224 PHERMEX, 204, 230 Piston Boundary, 29, 51, 98, 100, 103 Potential Energy, 7 Potential Function, 4, 86, 90 PowerPoint Presentation, 115, 191, 211, 215, 228, 236 Progressive Waves, 5 Projectile, 38, 194–5, 210, 212 Q Computer, 224 Radiograph, 204–5 RAGE, 180, 211, 229 Reflective Boundary, 29, 136 Reynolds Number, 90 Rigid Wall, 89–91, 130 Root, 115, 194–6, 207, 210, 217

Subject index

230

Runup, 34, 86, 93, 96, 100, 103, 163, 165–75, 230 SAGE, 180, 193, 210, 211, 231 SAIC, 210, 212, 229, 231 Scotch Cap Lighthouse, 33, 34 Sea Floor Displacement, 30 Sea of Japan, 30 Seiche, 3, 6 Shallow Water, 6–11, 15, 21, 22, 26, 30–46, 57, 61, 66, 78, 80, 85, 93–5, 99–114, 135–59, 161–2, 165, 166, 167, 169, 174–5, 180, 185, 193, 202, 208, 211, 216, 217, 220–2, 231 Shoaling, 77, 78, 116, 119, 120, 122, 125, 126, 131, 133, 135, 137, 144, 231, 234, 235, 237 Shock, 181, 195, 205, 216, 219, 220, 227, 229 Skagway, 64–6, 71, 76, 78, 82 SMAC, 88–90, 92, 94, 97 SOLA-3D, 78, 117, 119, 120, 123, 125, 126, 131, 133, 136, 138, 145, 231, 235, 237 Source, 26, 31–46, 58, 62, 72, 82, 115, 134–46, 153–7, 162, 173, 176, 207, 236, 237 Stability, 94, 133, 134 Standing Wave, 5, 6 Stationary Wave, 5 Steady State, 4, 8, 9 Stem, 195 Stokes Wave, 14, 15, 20, 237 Storm Wave, 3 Submerged Barrier, 105 Surface Marker Particles, 85, 92 SWAN, 22, 26, 29, 30, 32, 46, 57, 62, 67, 76, 78, 81, 100, 110–2, 135, 138, 145, 146–8, 154, 156–8, 163, 165–7, 170, 172–4, 186, 231, 237 Three-dimensional, vii, ix, 17, 78, 116, 119, 134, 135, 141–4, 162, 176, 178, 211, 218, 224, 230, 235, 237 Tide Gauge, 34, 38, 64, 65, 66, 68, 74–7, 83, 117 Tide Generating Force, 22 Time Step, 27, 29, 47, 68, 93, 110, 116, 122, 128, 133, 134, 163, 174, 180, 186, 188 Translatory Waves, 3 Transmission Coefficient, 106, 108 Tsunami Museum, 42 Two-dimensional, 20, 25, 30, 83, 93, 116, 134, 151, 157, 180, 212 Underwater Barrier, 107 Unimak Island, 33, 34 Upper Critical Depth, 109, 114, 122, 203 Viscosity, ix, x, 1, 2, 24, 94, 111, 121, 134, 141, 149, 159, 169 Volcanic, 115, 176, 211 Vorticity, 86, 90, 93 Waianae Harbor, 30, 80, 236

Subject index

231

Water Cavity, viii, 15, 112, 119, 180, 193–5, 202, 206, 209, 230, 231, 236 WAVE Code, 10, 13, 14, 17, 20, 80, 208, 237 Wave Energy, 30, 104 Wave Length, viii, 5, 7, 8, 13, 14, 20, 23, 95, 98, 99, 102, 104, 114, 115, 144, 154, 156, 162, 166, 175, 203, 208, 214, 216 Wave Number, viii, 5 Wave Period, viii, 5, 33, 38, 75, 95, 105, 115, 154, 162, 163, 166–7, 172, 174 White Machine, 218 Wind Waves, 24 XTX-8003 Explosive, 203 ZUNI, 85, 93–4, 100, 110, 113, 116–7, 146, 149, 153, 156–9, 163, 169–74, 231, 234, 237

Our partners will collect data and use cookies for ad personalization and measurement. Learn how we and our ad partner Google, collect and use data. Agree & close