Title: Numerical calculation of the response of coastal waters to storm systems-with application to hurricane Camille of August 17-22, 1969
CITATION PDF VIEWER THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00097632/00001
 Material Information
Title: Numerical calculation of the response of coastal waters to storm systems-with application to hurricane Camille of August 17-22, 1969
Physical Description: xvi, 153 leaves. : illus. ; 28 cm.
Language: English
Creator: Pearce, Bryan Rowell, 1944-
Publication Date: 1972
Copyright Date: 1972
 Subjects
Subject: Storm surges   ( lcsh )
Hurricane protection   ( lcsh )
Hurricanes -- Mathematical models   ( lcsh )
Civil and Coastal Engineering thesis Ph. D   ( lcsh )
Dissertations, Academic -- Civil and Coastal Engineering -- UF   ( lcsh )
Genre: bibliography   ( marcgt )
non-fiction   ( marcgt )
 Notes
Thesis: Thesis--University of Florida.
Bibliography: Bibliography: leaves 150-152.
Additional Physical Form: Also available on World Wide Web
General Note: Typescript.
General Note: Vita.
 Record Information
Bibliographic ID: UF00097632
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: alephbibnum - 000869248
notis - AEG6273
oclc - 014258470

Downloads

This item has the following downloads:

PDF ( 5 MBs ) ( PDF )


Full Text













NU14ERICAL CALCULATION OF THE RESPONSE OF COASTAL ',ATERS
TO STOPI SYSTEMS WITH APPLICATION TO HURRICANE CA,1IlLLE
OF AUGUST 17-22, 1969





By



Bryan Rowell Pearce









A DISSERTATION PRESENTED TO THE GRADUATE
COUNCIL OF THE UNIVERSITY OF FLORIDA lIN PARTIAL
FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PH ILOSOPHY


UNIVERSITY OF FLORIDA
1972














DEDICATION


This dissertation is dedicated to my son,

Nicholas James, born August 14, 1972.












ACKNOWLEDGEMENTS


The author wishes to thank Dr. R. G. Dean, who initiated

and supervised this investigation, and Dr. B. A. Christensen.

The hours of discussion and field trips with Dr. Dean have provided

a unique educational experience and contributed to a better under-

standing of the subject. The helpful criticism and conferences

with Dr. Christensen were invaluable and provided the extrospective

view that is necessary to any extended work.

Also to be thanked are Susan Phillips, Elizabeth Maltz and

Marilyn Morrison who typed the manuscript and Denise Frank who

drafted the figures.

The study was financed by the National Science Foundation.

All computations were carried out using the University of Florida

IBM 360/65 computer.












TABLE OF CONTENTS


Page

AC IOWLEDGEIMENTS iii

LIST OF FiGURES vi

LIST OF SYMBOLS x

ABSTRACT xv

CHAPTER

1. INTRODUCTION 1

2. STOPRM SURGE CALCULATIONS THEORY 7

a. Volume Transport Equations 7

b. Assumptions and Boundary Conditions 13

c. The Hurricane Model 20

d. Analytical Treatment Steady State,
Without Bottom Friction 29

i. Problem Definition 29

ii. Uniform Wind Stress Acting Over a
Shelf of Uniform Depth 30

iii. Uniform Wind Stress Acting Over a
Shelf of Uniform Slope 34

e. Analytical Treatment Moving Wind Stress
and Pressure Field 38

i. Problem Fonnulation 38

ii. Solution 40







Page


3. THE NUMERICAL MODEL 56

a. One-Dimensional Numerical Model 56

b. Development in Two Dimensions 64

c. Effects of Discretization and Checks of
the Numerical Model 73

4. RESULTS OF NUMERICAL CALCULATIONS 87

a. Introduction 87

b. Numerical Calculations and Comparison
with Available Data 94

5. FIELD DATA AND INSTRUMENTATION 118

6. SUMMARY AND CONCLUSIONS 132

APPENDIX

A. DERIVATION OF GOVERNING EQUATIONS FOR SURGE
HYDRODYNAMICS 134

B. TRANSFOPSATION OF GOVERNI;IG DIFFERENTIAL
EQUATIONS INTO FINITE DIFFERENCE FORM 140

1. Transformation Without Convective Terms 140

2. Effect of Convective terms 146

REFERENCES 146

BIOGRAPHICAL SKETCH 153













LIST OF FIGURES


Figure Page

2-1 PERSPECTIVE OF NEARSHORE ZONE SHOWING COORDINATE
SYSTEM AND JNOMENICLATURE 8

2-2 COMPARISON OF VAN DOPN 'S k WITH OCEANIC WIND STRESS
MEASUREMENTS 15

2-3 TRANSIENT WATER VELOCITY DISTRIBUTION IN A CLOSED
CHANNEL DUE TO IMPULSIVE WIND 18

2-4 DEFINITION SKETCH FOR STORM TRAVELING AT CONSTANT
VELOCITY 22

2-5 STORM TRACK FOR HURRICANE CAMILLE, AUGUST 16-18,
1969 23

2-6a PRESSURE CHARACTERISTICS OF HUPRICANE CAMILLE AS
A FUNCTION OF TIME 24

2-6b RADIUS TO MAXIMUM WINDS AND FORWARD SPEED OF
HURRICANE CAMILLE AS FUNCTIONS OF TIME 24

2-6c MAXI.MUIM WIND SPEED IN HURRICANE CAMILLE AS A
FUNCTION OF TIME 25

2-7a COMPARISON OF 30 FOOT SURFACE ISOTACHS 8/17/1200 26

2-7b COMPARISON OF 30 FOOT SURFACE ISOTACHS 8/17/1800 27

2-7c COMPARISON OF 30 FOOT SURFACE ISOTACHS 8/18/0000 28

2-8 NOTATION FOR ONE-DIMENSIONAL PROBLEM 30

2-9 ONE-DIMENSIONAL SHELF WITH UNIFORM DEPTH 31

2-10 SURGE HEIGHT VERSUS FETCH/DEPTH FOR UNIFORM DEPTH 33

2-11 ONE-DIMIENSIONAL SHELF WITH UNIFORMLY SLOPING BOTTOM 34

2-12 SOLUTIONS TO UNIFORM SLOPE CASE AT THE COASTLINE 36

2-13 IDEALIZED SHELF BATHYM.ETRY 40

2-14 IMPULSIVE FORCING FUNCTION MOVING AT CONSTANT SPEED 43






Figure Page

2-15a ARBITRARY FORCING FUNCTION 44

2-15b TRIANGULAR FORCING FUNCTION 44

2-16 DEFINITION SKETCH FOR SHEAR STRESS DISTRIBUTION 47

2-17 RESPONSE TO TRIANGULAR SHEAR STRESS DISTRIBUTION 51

2-18 DIMENSIONLESS SURGE HEIGHT MAXIMA AS A FUNCTION
OF DIMEISIONLESS STORM WIDTH AND VELOCITY 52

2-19a DEFINITION SKETCH FOR NUMERICAL CALCULATIONS 54

2-19b FINITE DIFFERENCE MODEL OF TRIANGULAR SHEAR
STRESS DISTRIBUTION 54

2-20 COMPARISON OF NUMERICAL AND ANALYTICAL SOLUTIONS 55

3-1 DEFINITION SKETCH FOR ONE-DIMENSIONAL NUMERICAL
CALCULATIONS, SECTION NORMAL TO THE COAST 58

3-2 MOVEMENT OF STORM 61

3-3 GRID ELEMENT DESCRIPTION FOR TWO-DIMENSIONAL
SCHEME 66

3-4a COMPONENTS OF q (i,j,n) 68

3-4b COMPONENTS OF qx(i,j,n) 68

3-5 SCHEMA OF GRID SYSTEM AND BOUNDARY CONDITIONS
EMPLOYED IN TWO-DIMENSIONAL CALCULATIONS 69

3-6 FILTERS USED TO CALCULATE WATER VELOCITIES 71

3-7a COMPARISON OF MAXIMUM SURGE HEIGHTS CALCULATED
BY THE ONE- AND TWO-DIMENSIONAL MODELS 74

3-7b PROFILES OF STORM CHARACTERISTICS USED IN
COMPARE I SON 74

3-7c Y-COMPONENT OF WIND VELOCITY 75

3-8 COMPARISON OF THREE ONE-DIMENSIONAL MODELS 77

3-9 IDEALIZED SHELF RESPONSE PROBLEM WITH UNIFORM
DEPTH 79

3-10a ISOLINES OF DIM'ENSIONLESS STORM TIDE, ANALYTICAL
SOLUTION 81





Figure Page

3-10b ISOLINES OF DIMEN SIONLESS TRANSPORT IN
x-DIRECTIO;, ANALYTICAL SOLUTION 82

3-10c ISOLINES OF D1;;ENSIOILESS TRANSPORT IN
y-DIRECTIOil, ANALYTICAL SOLUTION 83

3-11 DEFINITION SKETCH FO, ;NU:MIERICAL MODEL 84

3-12 EFFECT OF LATERAL BOUNDARY PROXIMITY AND
RELATIVE GRID SIZE ON MODEL RESPONSE 86

4-1 COASTAL REGION MODELED 88

4-2a LARGE GRID MODEL 89

4-2b DEPTHS USED IN LARGE GRID MODEL 90

4-2c SMALL GRID MODEL 92

4-2d DEPTHS USED IN SMALL GRID MODEL 93

4-3a CALCULATED STORM SURGE ELEVATIONS 8/17/0920 95

4-3b CALCULATED STORM SURGE ELEVATIONS 8/17/1640 96

4-3c CALCULATED STORM SURGE ELEVATIONS 8/17/2000 97

4-3d CALCULATED STORM SURGE ELEVATIONS 8/17/2200 98

4-3e CALCULATED STORM SURGE ELEVATIONS 8/17/2340 99

4-4 COMPARISON OF OBSERVED AND CALCULATED SURGE
HEIGHT MAXIMA 100

4-5 SENSITIVITY OF SURGE CALCULATIONS TO VARIATIONS
IN f AND Xs 101

4-6 SCHEME USED TO EXTRAPOLATE SURGE HEIGHTS TO COAST 102

4-7 COMPARISON OF EXTRAPOLATED AND NON-EXTRAPOLATED
SURGE HEIGHTS AT BILOXI BAY, MISSISSIPPI 104

4-8a COMPARISON OF OBSERVED AND PREDICTED SURGE AT
DAUPHIN ISLAND, ALABAMA LARGE GRID MODEL 105

4-8b COMPARISON OF OBSERVED AND PREDICTED SURGE AT
DAUPHIN ISLAND, ALABA'A SMALL GRID MODEL 107

4-9a WATER VELOCITY AND WIND VELOCITY RECORDS FOR
HURRICANE CAMILLE AT EGLI;I AFB, FLORIDA 108


vii i







Figure Page

4-9b COASTAL TOPOGRAPHY AT SITE OF WATER VELOCITY
MEASUREMENTS 109

4-10 COMPARISON OF CALCULATED AND MEASURED WATER
VELOCITIES 110

4-11 LONGSHORE CURRENT AND BEACH PROFILE MEASURED
AT MUNICIPAL PIER, PANAMA CITY, FLORIDA 112

4-12 COMPARISON OF MEASURED AND CALCULATED WIND
SPEEDS AT SITE OF WATER VELOCITY MEASUREMENTS 113

4-13 COMPARISON OF PREDICTED AND ACTUAL FLOODING
ON MISSISSIPPI COAST 115

4-14 COMPARISON OF PREDICTED AND MEASURED FLOODED
AREAS PER UNIT LENGTH ALONG COAST 117

5-1 SUMMARY OF CRITICAL WARNINGS 120

5-2 HINDCAST WAVE FIELDS FOR HURRICANE CAMILLE 122

5-3 NUMBER OF TIMES DESTRUCTION WAS CAUSED BY
TROPICAL STORMS, 1901-1955 124

5-4 PROPOSED FIXED CONTINUOUSLY MAINTAINED
INSTRUMENTATION SITES 125

5-5 SCHEMATIC OF RECORDING UNDERWATER TIDE GAUGE 126

5-6 RECORDING UNDERWATER TIDE GAUGE 127

5-7 TIDE MEASUREMENT-UNDERWATER GAUGE 128

5-8a CURRENT METER DEPLOYED 130

5-8b PHOTOGRAPHIC RECORDING MECHANISM OF CURRENT
METER 130











LIST OF SYMBOLS


an Fourier coefficients

Bn Fourier coefficients for forcing function

c wave celerity, Jgh

c1,c2 constant coefficients used in moving wind stress development

C storm or forcing function velocity

CDT central daylight time

c* C/c

D total depth

f Darcy-leisbach friction factor

fL linear friction factor
F,F' forcing function

Fx,F resistance coefficients for x and y equations of motion

g gravitational acceleration

G;.;T Greenwich ;.lean Time

h water depth, referenced to Mean Sea Level (';.SL)

h water depth at shelf break

hi water depth at coastline

H11,12 definitions used to simplify moving wind stress development

HF hindcast wave height

i indicial functional parameter

k wind stress coefficient

z shelf width






m bottom slope

MSL mean sea level

n indicial functional parameter

p pressure

pO surface pressure at center of hurricane

p surface pressure far from storm

pn surface pressure

q water transport

q water transport component normal to boundary

qx,qy water transport components in x and y directions
r distance from center of hurricane

R radius to maximum winds

t time variable

t* dimensionless time variable, t/T

t' dummy variable

t initial time

tl,t2,t3 time that zeros of idealized forcing function cross shelf break

T function of time; 4z/C

TF hindcast wave periods

u component of water velocity in x direction

u time average of u

u' instantaneous fluctuation of u

U wind velocity at 30 feet above the surface

Uy x and y components of U

Ucr U critical, in Van Dorn's expression for -.,ind stress coefficient







Uc cyclostrophic wind

U geostrophic wind

UG gradient wind

v component of water velocity in y direction

v time average of v

v' instantaneous fluctuation of v

V forward velocity of storm

V' V sine

w component of water velocity in z direction; wind. field breadth

w time average of w

w' instantaneous fluctuation of w

w* lateral dimension of shelf incorporated into numerical model

x,xl horizontal coordinate, normal to coast

x*,x' dimensionless horizontal coordinates, x/Q

xc coordinate of coastline

X function of x; dimensionless depth, x/xB

Xs starting point of storm as distance from landfall

y,X2 horizontal coordinate, parallel to shore

y' dimensionless horizontal coordinate, y/1

z,X3 vertical coordinate

aB angle between breaking waves and shoreline

an coefficient in moving wind stress analytical development
a* dimensionless parameter, w/w*

a** dimensionless parameter, c/w*

B constant, T /og
X


xii






"'V Dg

6 impulse function

Ap central pressure index

Lt time increment

AX,Ay grid element dimensions

AT shear stress increment

c dimensionless parameter, mh/l

o, dimensionless parameter, mho /

c dimensionless parameter, mhl/B

n surge height, referenced to MSL

n' dimensionless surge height

e angle between point of observation and line of storm movement

shear stress parameter

\ width of idealized storm

\* dimensionless storm width, .\/

p molecular viscosity

( mn/E; dummy variable

o density of seawater

oa density of air

0,on cnn/2z

o',an' Cnn/2Z

T shear stress

T surface shear stress

T components of shear stress
xTy,
1x and y components of surface shear stress
xq,y


xiii







Tb x and y conoonents of bottom shear stress
x,y
T' Reynolds stress

( latitude

w the earth's angular velocity

V gradient operator

v2 Laplacian operator










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



NUMERICAL CALCULATION OF THE RESPONSE OF COASTAL WATERS
TO STORM SYSTE;IS WITH APPLICATION TO HURRICANE CA'ILLE
OF AUGUST 17-22, 1969


By

Bryan Row.ell Pearce

August, 1972

Chairman: Robert G. Dean
Cochairman: Bent A. Christensen
Hajor Department: Civil and Coastal Engineering


A numerical model is employed to describe the water motion in a

coastal region associated with the passage of a hurricane or severe

storm. The model is two-dimensional, employs the vertically integrated

or "tidal" equations of motion, and is used to describe the specific

case of Hurricane Camille of August 1'69. Two models are employed, a

large grid with 16 nautical mile grid elements and a small grid model

with 6 nautical mile grid elements. The results of the two models were

not significantly different. The maximum storm surge or maximum

water level above mean sea level was calculated and found to agree to

within twenty per cent with maximum surge heights determined from

high water marks. Calculated surge height was found to be insensitive

to bottom friction coefficients varying from .005 to .02. inclusion






of th? nonlinear convective terms affected the calculated surge

height maximum for Hurricane Camille by less than two per cent and

affected the surge in any grid element adjacent to the coast by less

than five per cent. Two-dimensional storm surge plots at different

times are presented.

Hurricane generated currents were calculated and compared to data

taken off the Florida coast. It is concluded that more actual current

data are necessary before hurricane generated currents can be cal-

culated with confidence. The current calculations were found to be

sensitive to bottom friction and subject to "wind-up".

In addition to numerical calculations, basic analytical cases

are covered, including the response of a shelf of uniform depth to

a triangular wind stress distribution moving with constant velocity.


xvi











CHAPTER 1

INTRODUCTION


The problem of predicting the response of coastal waters to

extreme wind systems is becoming increasingly important as the coastal

regions of this country are developed for residential, recreational

and industrial uses. Disasters such as the ones occurring in recent

years at Galveston as a result of Hurricane Carla and along the

Mississippi coast as a result of Hurricane Camille demonstrate that

action should be taken to minimize the tragic consequences of such

events. Prediction of the vulnerability of particular areas to flood-

ing in similar storms would be of considerable interest to local

residents and governments as well as insurance underwriters. Moreover,

real time forecasting of storm surge elevation, with confidence, for

an approaching storm would be of great value in planning and imple-

menting evacuation and other measures to minimize damage. This is

clearly illustrated by Hurricane Agnes (1972). Agnes was a hurricane

of large areal extent but relatively low intensity. Little attention

was paid to surge forecasting, and as a result the residents of some

areas on the west coast of Florida were surprised by a storm tide 3 to

8 feet above normal. Damage in the Tampa-St. Petersburg area was

estimated at about $20 million. Undoubtedly some of this could have

been prevented if adequate real time surge forecasts had been avail-

able.




2


The design of permanent facilities at or near the shoreline

must include consideration of the maximum and minimum tides that may

occur as well as the wave action. In the operation of nuclear power

plants, for example, one consideration is the availability of

sufficient cooling water to ensure an orderly shut down prior to a low

tide occurence below the level necessary for plant operation. Thus,

a proper study of the situation would establish a maximum elevation at

which the inlet facilities could be placed.

Various approaches are used in solving this type of problem.

The purely analytical approaches are limited to idealized wind

stresses, pressure distributions, and topography. Lamb [1] considers

the case of "free" long waves in an infinite canal of uniform depth.

In addition to gravity, small disturbing forces are considered to be

acting on the fluid, which vary only by a small fraction of their value

within distances comparable to the depth. This leads to an in-

homogeneous wave equation similar to that of Section 3-e of the text.

Lamb then finds a solution for a pressure distribution that travels

with unchanged form at constant velocity. The presence of boundary or

terminal conditions greatly complicates the solution. The case of a

finite shelf is considered in the text.

Reid [2] presents the approximate response of a sloping shelf

due to a triangular shear stress distribution traveling directly on-

shore. The results, obtained at the shoreline, are based upon the

linear one-dimensional wave equation. The solutions are obtained

using a graphical technique and the method of characteristics. He

considers many different cases of fetch lengths and storm speed and

summarizes his results in graphical and tabular form.







Further analytical work in the area of storm surges has been

and will be done, yet the complexities of the actual forcing functions

and coastal topographies lead one to realize that another approach is

necessary, if realistic solutions to actual problems are to be

obtained. Hansen [3,4] is one of the pioneers of storm surge cal-

culations. He calculated surges arising from storms occurring in the

North Sea, using the linearized equations of motion. That is, he

ignored surge height with respect to depth and ignored the convective

terms. He laid the foundation, however, for the present work in

numerical surge calculations. Many contemporary investigators still

use, in some form, the staggered integration scheme he devised.

Platzman [5] analyzes the case of a pressure distribution or

squall line moving over Lake Michigan causing a surge of more than

six feet on the Chicago lake front. The volume transport or ver-

tically integrated equations of motion were used in central difference

form in conjunction with a two-dimensional bathymetry and a moving

pressure disturbance. Platzman's calculated surge is about one-

half of the observed surge and he concludes that this is probably due

to the omission of wind stress.

Miyazaki [6] treats the Gulf of Mexico as a closed basin,

assuming the surge height to be zero at the mouth between Yucatan,

Cuba and Florida. This procedure eliminates the difficult problem

of specifying lateral boundary conditions in open water. The grid

he uses is variable and finest at points of interest in the coastal

regions. He includes bottom stress in a form given by Reid [7] which

is a function of both flux, q, and of surface shear stress, Tn

Marinos and Woodward [8] calculate storm surges using the








assumption of bathystrophic flow (currents parallel to the bottom

contours). They arrive at a "quasi-one-dimensional" scheme which

considers only one-dimensional depth profiles but includes the

Coriolis effect. The bathystrophic procedure was originally developed

by Freeman, Baer and Jung [9].

Reid and Bodine [10] present a detailed numerical model for

the prediction of surges in Galveston Bay, Texas. The bay is treated

as a nearly closed system which, again, eliminates the problem of

open water, lateral boundary conditions. Measured tidal data at

the entrances to the bay are used as an input to the model and flow

through these areas is allowed. The model includes bottom friction

and allows for over-topping and flooding of adjacent low-lying areas.

Harris and Jelesnianski [11] and Jelesnianski [12,13,14,15]

have presented a series of papers concerning storm surge calculations.

The numerical technique follows that of Platzman. In the paper

"Numerical Computations of Storm Surges with Bottom Stress," Jel-

esnianski [14] includes bottom stress in his calculation, using an

eddy viscosity and a slip coefficient. Jelesnianski notes that for

a storm not moving too slowly, the bottom friction term is not im-

portant and may be omitted, but for slow storms or storms moving

parallel to the coast bottom friction must be included. In the paper

"Bottom Stress Time-History in Linearized Equations of Motion for

Storm Surges," Jelesnianski [15] abandons the slip coefficient and

maintains the eddy viscosity with the assumption that it is constant

in both the vertical and horizontal. Good agreement is obtained

between surge height data and the calculations. He notes, however, that

"currents are only of casual interest compared to the slope or







height variations." In general this is true, yet, it is becoming

apparent that the currents are very important in terms of beach

erosion and forces on various ocean structures. The development of a

more generally applicable dissipation mechanism would decrease the

sensitivity of the current calculations thus decreasing the necessity

for "calibrating" a model which may mask errors and deficiencies in

the model by appropriately selecting bottom friction coefficients.

The volume transport equations are formulated in terms of

flux components and do not describe the velocity distribution near the

bottom, which is related to the bottom friction. Jelesnianski [15]

deals with the problem of bottom friction by providing bottom shear

stress as a convolution integral of the surface shear stress function.

In this technique he linearizes the equation of motion by making the

assumption that the surge height is very small compared to the depth,

n < < h. In doing so he precludes flooding in low-lying areas and

calculations for very shallow waters. In this report it is antici-

pated that the modeling of flooding in low-lying areas and surge

calculations in very shallow water will both be desirable, thus a

quadratic bottom friction function which allows the use of the non-

linear equations is used. This is equivalent to making the assumption

that the vertical velocity profile is maintained in a near steady

state condition. Future development of the model will rely on measure-

ments to provide insight into the processes involved.

Surprisingly little work has been done in measuring hurricane

phenomena. The storm surge data available are primarily from high

water marks. Tidal records are predominantly from enclosed bays and

rivers making comparison to the numerical model difficult. The only








known current data for Hurricane Camille were taken by Murray [16,17]

on the Florida Panhandle. These currents, measured in the surf zone,

were due to wave-induced mass transport as well as wind shear stress and

water slope. Clearly more data of this type are needed.

The difficulty of reaching offshore areas in a hurricane makes

manned observation of these phenomena nearly impossible, so alternate

measures must be found. The appropriate instrumentation must be

either in place continuously waiting for a storm or it must be placed

in the path of the storm prior to the development of adverse sea

conditions in that area. In the first case, the instrumentation sites

must be periodically serviced with the knowledge that the wait for a

hurricane may be a long one. In both cases, there is the problem of

instrument recovery. An inexpensive underwater gage and integrating

current meter have been developed and are described in the report.











CHAPTER 2

STOPRI SURGE CALCULATIONS THEORY

a. Volume Transport Equations


The equations of motion and continuity as used in the develop-

ment of the storm surge model are developed qualitatively in this

section with the details of the integration over depth presented in

Appendix A. The coordinate system and nomenclature are presented in

Figure 2-1. The origin of the coordinate system is arbitrarily chosen

at the shelfbreak and is approximately centered in the region of

interest. The x axis is directed toward shore, y is oriented

approximately parallel to the coast and z is positive upwards. The

Navier-Stokes or momentum equations for fluid motion are the governing

dynamic equations:


-+ u + -u w-au 2w(sin)v -1 a + 72u (2-1a)



v+ uv + vv + + 2w(sin)u 2v (2-1b)
at ax ay az P ay p


aw aw dw )w
+ U-%- + v~ + = - + V2w -g (2-1c)
at ax ay az p ,z P

where u, v and w are the components of the water velocity in the x, y

and z directions respectively, t is time, p is the water density, p is

pressure, u is molecular viscosity, u is the angular velocity of the















Elemental Water Column
Showing Shear Stresses
and Surface Pressure, P7


Mean Sea Level
(MSL)


FIGURE 2-I


PERSPECTIVE OF NEARSHORE
SHOWING COORDINATE SYSTEM

NOMENCLATURE


A z
Y A


Origin


Bottom


ZONE

AND







Earth, } is the latitude, v2 is the Laplacian operator v2 -=
dX-
a2 a2
+ + and the fluid has been assumed incompressible. The
continuity equation is:

au + av aw (2-d)
ax ay az-1d)


Equations (2-la) through (2-1c) portray explicitly only the laminar

processes. In this development, the turbulent effects are dominant,

especially the turbulent or Reynolds stresses. Although it is

impossible to describe the turbulent motion in detail, it is useful

to separate the effects of the mean and fluctuating velocity components.

By describing the variables in Equations (2-1) as the sum of a time

mean value and a fluctuation and then taking a time mean, the turbulent

equations of motion result. Thus setting u = u + u', v = v + v'

w = w + w', and p = p + p' leads to the x-equation of motion


au + au 2u au 1 ap- + ,!2 a'r+ au'v' + au'w') (2-2)
-+ u -+ v + w ..-.u-(.+ + (2-2)
at ax ay az p ax p ax ay az

where
T
UI I u'iu dt


and uj, u2 and u. indicate u', v' and w' respectively. The terms

in parenthesis are the turbulent stresses and T is a period large

compared to the period of a turbulent fluctuation but small compared

to any change in u or v. It is assumed here that terms of the form
7 v2 u are negligible when compared to the turbulent stresses. Thus

the turbulent equations to be used are:








au au + vu au- 2w(sin)v B- 2 1 xx + xz)
at ax + y az P ax ax 2y a-3


av + v V v +p _l a=v aT'
-+ + v- + w-+ 2w(sin,)u 1-P a + y v + z) (2-3b)
at ax 3y 3z p ay P ax ay aZ

a 'I a -' a Iz
aw aw W a w_ 1 1 xz -z z
+ u-- + v- + az z g +4 x + + ) (2-3c)
at ax ay TZ p aZ p ax ay 3Z


where it is now understood that u, v, w, p, T' indicate time mean
xixj
values. Tix the Reynolds stress, is represented by -puu.1 u where
;ixj x1J
pui-u =PU'.. This is explained in detail in reference [18]. The
I3j 3 1
continuity equation in turbulent form is unchanged


au + av aw (2-3d)
ax ay az

Frequently a Boussinesq approximation is used for the expression

in parenthesis in Equation (2-2) which maintains the original form of

Equation (2-1). The equations remain intractable either in the form

(2-1)or(2-3). Since the velocity gradients in the vertical are antici-

pated to be, in most cases, an order of magnitude larger than the

velocity gradients in the horizontal, the lateral stresses are

neglected. Thus it is assumed that T' = T' = 0. The normal stress
aT' xy yx
terms of the form xixi are assumed small compared to the surface
axi 1 apn
pressure term - and the hydrostatic pressure. To

further simplify Equations (2-3) the assumption is made that the

vertical accelerations are small, or as a consequence that the pressure

distribution in the vertical is hydrostatic. Thus,


p(x,y,z,t) = Pr(x,y) + y(n(x,y,t)-z)


(2-4)







where n(x,y,t) is the free surface displacement from mean sea level

(MSL), p (x,y) is the pressure at the free surface, and y = pg is the

specific weight of water.

The Equations (2-3) are now integrated over the vertical from
-h to n where -h is the z coordinate of the bottom. This procedure

together with Equation (2-4) and the free surface and bottom boundary

conditions leads to the vertically integrated or tidal equations of

motion. The details carried out in Appendix A are straightforward

and with approximations to be defined later lead to the following,

using Leibnitz's rule,


qx qx qx + !L sqx D n gD-+1 ) (2-Sa)
at D ax D ay y p a aX p nx (2-5a)
x x


qy + x + + 2w(sing)q gD + 1 T ) (2-5b)
at D ax D aByx p y pny b


aqx anq
x a + 0 (2-5c)
ax ay at

where the subscripts, n, and b indicate the subscripted quantities are

to be evaluated at the surface and bottom respectively. D, qx and q
are defined as follows:

D(x,y) = h(x,y) + n(x,y,t)

1n(x,y,t)
qx(x,y,t) = u(x,y,t) dz
-h(x,y)
rln(x,y,t)
q (x,yt) = v(x,y,t) dz
-h(x,y)







It should be noted also that the convective terms of the form
q aq
x as written are only approximate; the actual terms resulting

from the integration are of the form

n
ay J uv dz,
a -h


but are written in the differential form in analogy to Equations (2-1).

The convective terms will be neglected in the further development of

the model until considered again in Appendix B. In summary, the

equations of motion in the final form to be used are:


aq xD apri Dr, 1
aq -2(sin.)qy - gD a +- (T -T ) (2-6a)
at y pax ax p x



3qy + 2~(sin.t.)q D -a gD 1 + 1 ( -Tb ) (2-6b)
at x 0 ay ;y p n y


The form of the continuity equation is:


aq aq
a ay n 0 (2-6c)
ax ay at













b. Assumptions and Boundary Conditions

The equations of motion were developed in the preceding

section. To specify a particular problem, various conditions

must be imposed on the boundaries of the model. On the water

surface a forcing function consisting of a normal pressure and a

shear stress is applied. At the bottom, a shear stress or friction

dependent on the water velocity is imposed. In addition conditions

on the surge height and water transport are necessary at the areal

boundaries of the model.

The shear stress at the free surface is calculated by the

method proposed by Van Dorn [19]. The shear stress is based on the

wind velocity and is given as follows:

T p klUlUx
f x
(2-7)
T = kpklUU
ny


where U = U + U and U and U are the x and y components of the
V y x y
wind velocity, respectively, at a reference elevation of approximately

10 meters. The constant k has been determined by Van Dorn as the

following function of velocity:

k U < Ucr
k = (2-8)

S+ 2(l Ucr/U)2 ; U Ucr







where

kI = 1.1 x 10-6


k2 = 2.5 x 10-6


Ucr = 14 knots = 23.6 feet/second


There are uncertainties in applying Van Dorn's results to hurricane

winds. Van Dorn's data are generally for wind speeds low in comparison

to the wind speeds of a large hurricane. Although there is no proof

that the results are directly applicable, they are used here because

the associated measurements showed little scatter and the results

are reasonably consistent with other data sources. Wu [20] has pub-

lished more recent information concerning the wind stress coefficient,

k. Figure 2-2 summarizes available data. Van Dorn's k, calculated

from Equation (2-3), is also plotted in Figure 2-2 to provide a com-

parison. Since the hurricane model to be discussed requires a shear

stress field as input, a better shear stress model, if it becomes

available, can be easily incorporated into the surge model without

changes in concept or approach. The wind velocity input will be

discussed in Section 2-c.

The treatment of the bottom friction, rb Tb is a much
x y
more difficult problem. In making storm surge calculations the

simplest assumption is to consider b = Tb = o. This is employed
x y
in many developments of storm surge calculations, especially the

earlier ones. Under appropriate conditions this is a justifiable

assumption. For the case of a rapidly moving storm traveling per-

pendicular or nearly perpendicular to the shore it is a good assump-















.003

0
e o.oe 0

/
/
.002 0


0 oLegend
/ Vo--- n Dorn's k

S.001 Oceanic Measurements
S0
) 0 Note:

.c- Oceanic Measurements

From Reference 20]
S ------Ii


5
Wind


10
Velocity,


15 20 25
U, (Meters/Second)


FIGURE 2-2


COMPARISON OF VAN

WITH OCEANIC WIND

MEASUREMENTS


DORN'S k

STRESS








tion since there is little time for "wind-up" and the friction

effects are small due to low water velocities [13].

A better assumption in some cases is that of a linear or

linearizedd" bottom stress of the form fLqx, f qy where fL is a

linearized friction factor. This is the type of friction used, for

example, by Dean and Pearce [21]. A linearizedd" friction term is

generally necessary if friction is required and any type of analytical

solution is desired. The bottom shear stress employed in this model

is the Darcy-Weisbach expression and is the most accurate of the

three,
pflulux pflulu
b 8 Tb 8
x y

or in terms of volume transport

Pflqlqx Pfjqlqy
Tb = 8D ; Tb y 8DP
x y

where u = JUx2 + Uy2 q = /q2 + q 2 and f is the quadratic

bottom friction factor. In incorporating this steady state rela-

tionship into the model it is assumed that the law is valid for the

propagation of long waves of the period of interest. This is

equivalent to assuming that the velocity variation over depth is

nearly the same for slowly varying flows as for steady state.

Although this technique has limitations, especially in deep water,

it is the best method available that can be used in conjunction with

the nonlinear equations to allow for flooding and to permit cal-

culations in very shallow water. Another important consideration is

that in deep water the friction is relatively less important than in








shallow water since the water velocities are in general much smaller.

As an illustrative example, the response of a water column to an

impulsive wind is depicted in Figure 2-3.

Returning to the problem specification, the boundary conditions

are illustrated in Figure 3-5. At the seaward boundary, the wind

setup is assumed zero as a result of considering an infinite depth

at that point. The barometric tide still exists at the seaward

boundary and is given by n(x,y) = 1.15(p -p (x,y)) where pc designates

the barometric pressure at a distance large compared to the storm

dimension, where r, is in feet, and the pressure p is in inches of

mercury.

At any land-sea boundary the flux normal to that boundary is

zero, q = 0, where the subscript n indicates the component normal

to that boundary.

Most difficult from a physical standpoint are the lateral

boundary conditions. They are, however, necessary to carry out the

calculations. Since no real physical limitations exist, conditions

must be chosen that least disturb the model and therefore allow a

realistic solution to be obtained. A further condition is then

imposed that the boundaries are at a distance large compared to a

characteristic dimension of the storm. This means, in effect, that

the boundary conditions on the lateral boundaries become a second

order effect if they are far away. The boundary conditions used at

the lateral boundaries of this model are: the flow in the x direction

is equal to zero, qx = 0, and the surge height n is set equal to the

barometric tide at the lateral boundaries.




18





Direction of Wind Velocity -



7.--
_!1-- ----- _--- -_-_-J----
-- .

it i"7 2"f<''


0.9-



0.8-



0.7-



0.6-


0.5-



0.4-



0.3-


0.2-



0.1-


I =0

I5

10

20

30

40

50

Notes:

Wind

Depth


0

Woter Velocity,


FIGURE 2-3


sec.

sec.

sec.

sec.

sec.

sec.

sec.


eociy 0; t<0
velocity 9.8 ft/sec, t 0 -

= 0.48 Feet


0.1

u, Feet/Second


0.2


TRANSIENT WATER VELOCITY

DISTRIBUTION IN A CLOSED CHANNEL

DUE TO IMPULSIVE WIND (FROM

REFERENCE [22])


I .'I r
(
I
I
I
1
I
.1





A







:I
I
i
I
I


I


ii: '
Y


01
-0.


I







The initial conditions are simply a quiescent system. The surge

heights can be initially zero or in barometric equilibrium: this makes

little difference since the storm is initially far from the region

of interest. Initial barometric equilibrium is used in this study.

More importantly the storm is assumed to start at a distance offshore

large enough that the effects due to model start-up are negligible.

The effects of various starting points can be seen in Figures 4-5 and

4-10, where the surge height is relatively insensitive but the currents

exhibit wind-up.













c. The Hurricane Model


The actual storm surge model uses T T and p as input.
x y
Since the measured wind shear stress is in general not available,

Van Dorn's technique is used to convert U and U to T and T as
x y
discussed in Section 2-b. Actual wind fields may be applied or an

idealized model used to represent the storm. In this section an

analytical model is used to represent wind velocity and pressure fields.

Since the actual surge model is the primary objective of this study,

the hurricane model will be presented and used without further comment,

except for a short note on the numerical technique in Chapter 3.

The hurricane representation used is not unique and a different

one could be substituted into the numerical model. The hurricane

model used is that presented by Wilson [23]. The important parameters

of an idealized hurricane are: the central pressure index CPI or

equivalently the central pressure anomaly Ap = pm p0 (where pm is

the normal pressure or the pressure at a large distance from the

hurricane center and po is the central pressure of the hurricane), the

distance from the hurricane center to the point of maximum winds R,

and the forward velocity V of the hurricane. These parameters Ap,

R, and V are reported by the National Weather Service when character-

izing a storm, thereby providing a method to approximate a hurricane.

The isobars are assumed circular and the pressure distribution at
the surface is taken as pR/r
the surface is taken as p = p + Apep The streamlines are
n o






assumed circular for a stationary storm. The equations presented

include the effect of a constant translational velocity, V,

superimposed on the storm as shown in Figure 2-4. Using Wilson's

notation:

UG = Uc( J2 + 1 -

where
V U
Y 1/2(- + U)
c g

and

Uc = Ap R e-R/r
P r

Pa
Ap R -R/r

U =
g 2w(sinT)


V' = V(sine)

The subscript G refers to gradient wind, c to cyclostrophic

wind, and g to geostrophic wind. pa is the density of air, V the

forward velocity of the storm, e is defined in Figure 2-4 and 2w(sing)

is the Coriolis parameter as used in Equations (2-1). The gradient

wind UG is the desired parameter. These winds are calculated at a

point high enough above the surface to allow the assumption of zero

friction. To reduce these to the 10 meter elevation a reduction factor

is needed. The factor depends on latitude, divergence of wind

direction from the isobars and other factors. In this study a

constant factor of 0.83 was used. See references [24,25].















UG = Gradient Wind Velocity








SV = Forward Velocity
Storm Center o Stor
of Storm








FIGURE 2-4 DEFINITION SKETCH FOR STORM TRAVELING AT CONSTANT VELOCITY


The winds of Hurricane Camille are evaluated using this method.

The hurricane position as a function of time is given in Figure 2-5.

The central pressure p and the peripheral pressure pm, are given as

functions of time in Figure 2-6a. Figure 2-6b depicts radius to

maximum winds R, and forward speed of the storm V, as functions of

time. Figure 2-6c gives maximum sustained wind speeds in knots as

a function of time. These data are from the National Weather Service

[26,27]. Figures 2-7a through 2-7c present comparisons, at different

times, of the surface wind fields for Hurricane Camille as developed

by the National Weather Service [27] with the wind fields calculated

by the above model.




























Notes: I
Times are GMT
Dotes ore in A


FIGURE 2-5 STORM TRACK FOR HURRICANE
CAMILLE, AUGUST 16-18, 1969









Legend:


940




o
S920



90

0-900

o
c
(-


--.---- --.------- 0------------ -- ----0- -


0600


1200 1800
Time (GMT)


FIGURE 2-60


PRESSURE CHARACTERISTICS OF

HURRICANE CAMILLE AS A

FUNCTION OF TIME

(DATA FROM REFERENCE [26])


-I 20- |---
2 odius to Maximum Winds
E _,- 0odiu 0 Mx-- -u m Wiuds
Er


o| C I0. -
"^ i Forward Velocity 15

*0 / 5 /
-Y
o / \ -


10 10
0000 0600 1200 1800 0000 0600
Aug 17 Time (GMT) Aug 18

FIGURE 2-6b RADIUS TO MAXIMUM WINDS AND

FORWARD SPEED OF HURRICANE

CAMILLE AS FUNCTIONS OF TIME

(DATA FROM REFERENCE [26])


---- Centrol Pressure

0-------0 Peripheral Pressure
Nole:
Storr Landfall at 0430 GMT, Aug 18


-_ ,


oou
0000
Aug 17


)


0000
Aug 18


1010
0600





L.-
o




030-





Q)
,-
L-
020a
I-

0)


1l

















o)




S125
U,

'o
c








Aug 171
E 120 2
D
E
0


115
0000
Aug 17








FIGURE 2-6c


0600 1200 1800 OC
Time (GMT) A


MAXIMUM WIND SPEED IN

CAMILLE AS A FUNCTION

(DATA FROM REFERENCE


)00 0600
ug 18








HURRICANE

OF TIME

[26])
















cn
'.0






4--














0
cn

C>
C)

















cJ1
4-


U
L
0


0



0



C.J
0











rL.










O.J


- O














C)
I i4
J












.I U

-
- ,a
0












1 u





U-e
>
> rn





27








cn


I>







"-/ "
o -. co
z CD
01
m 0O






C) CO



if


C-,
OL
I--I /2







t/ C




.O / '
, t,,, l- ,,

0 (D
\ 0"---O"\ o o0 /

C O




co o L-)
n,,
00
9 I / \\ (
r o\ co
s.
In \0 \ 1 \*



0 0 0

0 I I.
-0 0 0 CL..









\ 'o^ \fe\ --- -- ^


mL L
0 0

0 0


o II
LI L1



I
e '



EJ I
00













---C--!












d. Analytical Treatment-Steady State, Without Bottom Friction


i. Problem Definition

Having defined the problem and established equations of motion

and boundary conditions, it remains to find a solution for a particular

storm and affected coastal region. Before presenting the numerical

approach necessary to attain a solution of complicated prototype con-

ditions and geometry, several simplified analytical approaches will be

considered. These will serve three basic purposes. First, they

allow highly simplified calculations for basic engineering purposes.

Second, they allow intuitive investigations into the physical processes

involved. Third, they produce solutions with which to check the

numerical approaches.

Figure 2-8 is a definition sketch for this section. The present

treatment includes only motions in the x-direction; Equations (2-5a)

and (2-5c) therefore reduce to Equations (2-9a) and (2-9b)

respectively.


aq q aq D ap
ax +x x D gD a n + 1 ( ) (2-9a)
5at D -x p xx ;x p n b
x x

x +an -0 (2-9b)
ax at






















Bottom


FIGURE 2-8 NOTATION FOR ONE-DIMENSIONAL PROBLEM


The solution is now independent of y. Two cases will be considered,

the uniform depth case h(x) = h = constant and the case of h(x) = h

mx or uniform slope. In both cases it will be assumed that the wind

stress on the surface is uniform. For the steady state case, the

case of no net transport in the x-direction, q = 0 and therefore the

shear stress ,b = 0, regardless of f, the bottom friction coefficient.
v.
(Note: it is assumed here that there is zero velocity over the entire
q aq
depth.) For this problem the D- (-x) term is zero and the atmospheric
D dX
pressure gradient --- is assumed small in relation to the other terms.
ax

ii. Uniform wind stress Acting Over a Shelf of Uniform Depth

In the first example h(x) is uniform (=h) as shown in Figure


2-9.








z




9^ ----- 1--------L------- -:--,--7- P. X
Origin

h x= x






CD


FIGURE 2-9 ONE-DIMENSIONAL SHELF WITH UNIFORM DEPTH


The further assumption is made that the shear stress is uniform and

constant at the water surface z = n. The problem can be stated as

follows:
T
a2n x
-gD --- + 0
ax p

aq
qx at 0


Since D = h + n


pg(h+n) an T (2-10)
ax n


Rewriting Equation (2-10) leads to

T
aX 7 -x (2-11)
ax pg(h+[r)


At x 0, the boundary condition for r is determined by noting that

seaward of that point h m and therefore









an 0 as h -,
ax

Thus n = constant as h -, <, and n = 0 is chosen as the boundary con-

dition at the shelf break, x = 0. Solving for n yields

T
r(h + n/2) = x
pg

or

n = -h + h2 + 2x
pg


Evaluating r using Van Dorn's constants, completes the solution.
nX
An even simpler procedure is frequently the basis for storm surge

calculations. Assuming that the storm surge, n, is small compared to

the undisturbed depth (i.e. h + n -, h) leads to the solution

T
n-x
pg-

Designating xc as the fetch or shelf width and h as depth leads to

the storm surge at the coastline x = xc

T
nx Fetch
x=Xc pg Depth


Figure 2-10 provides a simplified procedure for making these

calculations and is based on the latter solution. As an example,

assume a depth of 50 feet for a 50 mile shelf width (or fetch) and

a wind speed of 50 knots, then the storm surge calculated from

Figure 2-10 is 3.23 feet.







































s-

r. 0
0


LT


(\J
O
O


(1-aj) L '10SP!H a6jnS




34




iii. Uniform Wind Stress Acting Over a Shelf of Uniform Slope

As a slightly more realistic situation, the same assumptions

are made as in the previous case except now a sloping bottom

(or h(x) = ho mx) is considered as shown in Figure 2-11.


T"q
77x


X/ I
xKxK


FIGURE 2-11


ONE-DIMENSIONAL SHELF WITH UNIFORMLY SLOPING BOTTOM


From Equation (2-11)


(h + rj) ar- x constant =
ax og

and noting that h(x) = h mx
0

aD (h+n) an + (h+) (h+n) m (h+n)
ax ax ax ax


DdD
thus dx = Timd

at x = 0 to obtain


Integrate and apply the boundary condition n = 0


1 4
= ln (1 - )


(2-12)








where
mn


mh


mh


mhi
=o





To simplify the calculations obtain n at the coastline x = xc

Equation (2-12) then reduces to

1 -
S= In ( ) (2-13)


Solutions to Equation (2-13) at x = xc are provided in Figure 2-12.

As an example consider a shelf with xc = 50 nautical miles, ho = 80

feet, and hI = 20 feet, this shelf has the same width and the same

average depth as in case ii. For the same 50 knot wind, n at the coast

is 3.53 feet. Heuristically, one would expect a larger n for the

case of the sloping shelf.

Consider the location on the sloping shelf which has the same

depth as the uniform depth shelf. For a point x towards shore

there is a change in depth ah and conversely for Ax away from shore.

Thus

+ constant constant 'h h 2 2
Ar, [ -- +- (- 0 (Ah)3]


constant constant [ +h Ah2h)3]
n h.h. (F) + 0 (Ah) ]
hF -fah FI
















w
I--
S0
N ,0 o



0\ 0







\ \n
0
C\%J 0




C)





-"-J


0
LLi






)\ (/)
z
S I0
5 X \ \
-JD



0 \





4------ CL
L-




37


and the total over 2Ax is

n = An + An- = constant [2 + 2(Ah + 0 (Ah)4]


It can be seen, then, that a shelf of uniform depth will, in all

cases, produce a smaller surge than a shelf of constant slope with

the same average depth.













e. Analytical Treatment Moving Wind Stress and Pressure Field


i. Problem Formulation

The purpose of this section is to contribute to the development of

a simple procedure for the analysis of the dynamic effects involved

in storm surge calculations. This analysis uses the one-dimensional

tidal equations and the development will be for an arbitrary forcing

function such as a hurricane. Also the bottom profile must be in an

analytical form simple enough to allow a solution to the equations of

motion. Conceptually, the approach is as follows: find the storm

surge n(x,t) resulting from an impulsive forcing function representing

the wind shear stress, pressure gradient or both. If the impulse

occurred at time t' and position x', the resulting surge at (x,t)

could be written n(x,t) = G(x,t; x',t'); this is sometimes referred

to as the Green's function. Having found the general form of the

solution G(x,t; x',t') for a simple offshore bathymetry one can then

construct the solution for an arbitrary forcing function. This may

be obtained by formal integration or the results may be obtained

numerically.

The development in this section will be for the special case of

a shelf of finite width and uniform depth as shown in Figure 2-13

and a triangular shear stress distribution as shown in Figure 2-15

For the equations of motion and continuity, consider Equations

(2-5) where qx is abbreviated as q, nx as in, and bx as Tb,
\ x




39




9 + gD an D aP- b) (2-14a)
at D ax ax pax p (2-14a)


S+ = 0 (2-14b)


The convective acceleration terms will be neglected and the assumption
n < < h is made, and since D(x,t) = h + n(x,t) = h, Equations (2-14a)

and (2-14b) reduce to


aq(x,t) + gh an(x,t) h apn 1
at + gh ,ax p (x,t) + (n(x,t) Tb(x,t))

aq(x,t) + an(x,t) = 0
ax at

Ignoring the bottom friction term Tb and letting (1-T (x,t) -
h ap (x,t) n
- ax -) = F'(x,t) we have:
p ax

a+ gh F= F'(x,t) (2-15a)

aqt an
ax a 0 (2-15b)


The variables are described in Figure 2-13, where q is the volume
transport per unit width. Differentiating Equation (2-15a) with

respect to x and Equation (2-15b) with respect to t and subtracting,

produces the inhomogeneous wave equation


2n gh 2r- aF'(,t) F(x,t) (2-16)
Tt= -TX 7 x










T (x, t)
7-



Origin ]

h
q(x,t)



S-II== Il l_ i _ _


CO


FIGURE 2-13 IDEALIZED SHELF BATHYMETRY

ii. Solution

First a solution to the homogeneous problem F(x,t) = 0 is

obtained. Using separation of variables, let n = X(x)T(t)

For the spacial part of the problem, we obtain

gh = a2X (2-17)

where an2 is a constant as yet undetermined. The solution to

Equation (2-17) is

X = cl sin -- x + c cos x (2-18)


Referring to Figure 2-13 note that the depth is infinite at x = 0,

and the surge must be zero at that point. At x = z, there can be no

flow through the wall and this leads to the condition = 0 at

x = 9..







The condition n = 0 at x = 0 determines

C2 0 and I lan (cos n ) T(t) = 0
ax x=Z fgh 'gh

leading to -= (2n + 1) The general case for the solution

becomes a sum of the solutions in Equation (2-18). Thus the solution

is taken to be in the form


n E a (t) sin nx (2-19)
n=1,3,5"*

Referring to the homogeneous form of Equation (2-16), obtain with

Equation (2-19) the following:

32a (t) nx C2 n2 a innj- =
S {at sin + C2 -) a (t) sin 2 = 0
n=1,3,5 -


From orthogonality, obtain

a "(t) + n2 a (t) = 0

where on = 1 and c = g

The two homogeneous solutions are of the form a (t) = sin ant,
Cnj
cos ont. Writing F(x,t) = Z B (t) sin "- and B (t) =
n n=1,3,5 n
SF(x,t) sin 2x dx for odd n, and using the method of variation of
parameters, arrive at the following for the inhomogeneous equation;
with the lower limit determined by the fact that a (0) = 0.

sin oat o t c)s o t C(t
an (t) B=(t') cosant'dt' n -n B(t )sinot'dt'
0 0







t B (t )
a (t) = B sinan(t-t')dt'




r a(t) sin n = B (t')sina (t-t')dt'sin 2
n= ,3,5 n=1,3,5 n 0

and finally

n(x,t) = Z F(x',t')sino (t-t')sin n dt'dx'sin (2-20)
n=1,3,5 n. O 0



Letting F(x',t') = ,6(x'-x ) 6(t'-t ) one then has an impulse function

located at x = x and t = to where n represents the Green's function
or response to this impulse function, and : is a constant representing
the strength of the instantaneous surface shear stress function


n(x,t) = G(x,t; Xo ,t) =


C O. 2 z tx n x
n=l n= 1,3,5 0


G(x,t;x ,to) = ? sine, (t-to) ssin sin x (2-21)
S n=1,3,5 n

It is now possible to construct a solution by the superposition of an

appropriate series of impulses in space and time.
For some hurricane motions, a more convenient form of the Green's
function is obtained by using an impulse function which travels with
constant velocity and remains unchanged during its traverse over the




43


shelf. Assuming the forcing function to always start at x = 0, to

travel at constant velocity C, and that the time that it crosses

x = 0 is arbitrary at to, the impulse function is:


F(x,t) = K<(x C(t-t )) 0 < x < z


Without loss of generality let t = t + to which yields


F(x,t) = K6(x Ct) 0 < x < R


Thus for the case of a suddenly increased shear stress as shown in
T
Figure 2-14, K where T is the amount of the increase in
p 0
shear stress.





P
-C-0


P
C

x



-I ar -r 8 (x-Ct) = F(x, t)
P ax P


0C

Sx


FIGURE 2-14 IMPULSIVE FORCING FUNCTION MOVING AT CONSTANT SPEED







For the case of an arbitrary shear stress distribution as shown in
Figure 2-15a, locally LT = x0
ax


FIGURE 2-15a


ARBITRARY FORCING FUNCTION


For a triangular shear stress distribution as shown in Figure 2-15b


I a.
P ad

+O


r I


FIGURE 2-15b TRIANGULAR FORCING FUNCTION


I~


>-C


------;-------------


.--r"'/p


-KO --------







1 r
r: becomes - Ux, and the following is obtained.
p dX


F(x',t') = 1 - %(x'-Ct') dx = AT 6(x'-Ct')
p dX p

The resulting Green's Function can then be written in integral form as

follows.

J 2 ar j -Ct't' nx' n7Tx
G(x,t;O,O)= (x'-Ct)sin(t-t)sin-dt'dx'sinx
n=1,3,5 n p 0n


Note that this is really the solution for an incremental element having

a shear stress jump of AT. Since 6(x'-Ct') 1 0 only when x' = Ct'


G(x,t;0) =


t
2AT nnCt' nux
E 2A sinn(t-t')sin 2 dt' sin 2
n=1,3,5 Pn 0


Letting an = = nV n = o


and integrating yields


2AT t n rx
G(x,t;O) = po(Z [o' sinot a sin'ot] sin 2x
n=1,3,5

Note that this is the solution for a solitary continuous impulse

function traveling at speed C; at time t > */C the forcing function

crosses the shelf, and the solution becomes periodic in time.

Letting c* = C/c, c = /gh and o' = o c* one arrives at


G(x,t;) = +2T [c*sinat-sin(ac*t)] sin nx
n=1,3,5 pz(c*zl

with the constraints t < r/C and c* / 1. When t > z/C









G(x,t;0) =


2 fAT k/C

n=1,3,5 p
o


sina(t-t') sin n "Ct' d sin n7x
22 21


where a' i/C = nT
2


n7T
cos = 0

n-l
sin = (-1) 2
2


G(x,t; 0) thus becomes


G(x,t;0) =


S 2Airsin(nnx/2)[+cs t
n=I*3,5 )P(c*-) +c*sn
n=1,3,5


n-1
- coso(t-i/C) (-1) 2 ]


a'-
with the constraints t > z/C and c*2 / I. If AT is written as r dx'
dX
and it is noted that x = Ct then AT = C -1-dt. Letting t = t t

yields G(x,t; t ) where to represents the time at which *(x-C(t-t ))

crosses the point x = 0. The desired result is then the sum or integral

of the above solutions


n(x,t) = I


G(x,t; t ) dt


The shear stress profile of a hurricane can be approximated by a

triangular function [2]. The derivative of this is two rectangular

pulses as indicated in Figures 2-15b. The response can be found

using Equation (2-22) and remembering that G(x,t; t ) is represented

by two functions depending on whether t is less than or greater than

E/C. The solution is straightforward but unfortunately tedious.


n = 1,3,5


(2-22)








The solution is simplified by putting n into a dimensionless form.

Letting


T = 4k/c

= aT = 2nTr

X* = X/Z

t* = t/T

allows n' to be represented as Equation (2-23).


l'(x*,t*) = -n(xt) (2n)3
(2K 'T3/z)


(2-23)


The appropriate variables and constants for the triangular stress

distribution are illustrated in Figure 2-16, where tI has been

arbitrarily set equal to zero for convenience,


C aT
P ax

Ko






-K0


tI
.--- X


FIGURE 2-16


DEFIINITIO;, SKETCH FOR SHEAR STRESS DISTRIBUTION


and are defined as follows:

C 3a
K ax
p ax




48




ti = 0 tl* = 0


t2 -2C


St2 8c*


t3 C 3* 4c*


X+


K = -r ;


I =
o


0 < t < -


T < t Co 2C
2C 'o < C


<' = ; < t
C 0

The solution must be presented as six different forms, each for a

different time span. They are as follows:

(1) In the first case the storm has completely passed by the

position, x = -, thus

,1+e **+l
t > C t* >-


and the solution for n is


n'(x*,t*) =
n-l to=t3
sin(n7x/29) [sin((t-t )-./C) (-) 2 +c*coso (t-t)]
n=1,3,5 ) n



By defining the contents of the brackets as H (to), the response becomes


n'(x*,t*) = r sin(n x/2) [H 3) H ) 2H )
n=l1,3,5 nF (c- l) [HI(t3) + Hl(tl) -2HI(t2)1







where tl, t2, t3 represent values of to and are defined in Figure 2-16.

(2) In case 2 the storm is completely on the shelf 0 < x < , thus


A/C < t < R/C T I

and n can be written as follows

n'(x*,t*) =
t=t3

sin(ncx/2) [+c*cos( n(t-to))- L cos(c*.n(t-to))]
n=1,3,5
t=tl

Defining H2(to) as the contents of the brackets, as above, leads to


n'(x*,t*) = sin(nnx/2) [H
n=,3,5 nJ(c*z-l) [H2(t3)+H2(tl) 2H2(t2)]
n=1,3,5 nc 1

The time periods from 0 < t < A/C and ,/C < t < have not yet been

covered since they introduce the additional complication of having

only part of the forcing function in the region of interest. The

solutions for these remaining time periods are:

(3)
z + A/2 < + A 1 + \*/2 1 + I *
C t C 4c* < < 4c*

n' (x*,t*) =

sin(nnx/22)
z n(nc ) [H2(t3)-H2(t-/C)+H1(t-/C)-2H (t2)+H (tl)
n=1,3,5

(4)
2 + A/2 1 1 + .0/2
r/C < t < ? + 2 1< <1 + /2
t- C c









sin(nux/2z) [H (t)
n in x/2)l [H2(t3)-2H2(t2) + H2(t-t/C)-Hl(t-z/C)+H(tl)]
n=1,3,5,...


A < t < -- < t*
2C c 8c < 4c*


0 < t < X/2C 0 < t* < 8*
8c*


ri'(x*,t*) n=
n=1,3,5,...


sin(nnx/2.) H
n (c* -1) [2t ) H2(t)


n'(x*,t*) has now been described for all t*. It can be seen that n*

can be described as a function of three dimensionless quantities

c*, X*, and t*.

Figure 2-17 shows profiles of n' on the shelf for various values

of t* and particular values of ,\ and c*.

Figure 2-18, perhaps of more interest, shows n'ax, the

maximum storm surge that occurs as a result of the storm passage as a

function of both c* and \*.

To determine the accuracy of both theory and calculations, an

independent approach was utilized to obtain equivalent dimensionless

surge heights. The one-dimensional tidal equation of motion and

equation of continuity in finite difference form, as developed in


sin(nrx/2g.)
n'(x*,t*) = n(c ) H2(t)+H2(tl)-2H2(t2)]
n= I,3,5,... *










0




0 0 z
0




O ()
co



rr

Q)
0
r- W
S\ / '

0












0-
LLL
,o,)

0


z D




00




*o Xo &


uIn
z


-O < L
d 0 .0
o 0 6

11 0




u0 0
0

soun AJDJliqj 't 'Iqb!lH binS i ~
























.5

4-

.3

2




o
.2 .4



FIGURE 2-18


DIMENSIONLESS
MAXIMA AS A
DiMENSIONLESS
VELOCITY


2 1.4 1.6 18 20



SURGE HEIGHT
FUNCTION OF
STORM WIDTH AND








Section 3-a, were used to obtain a numerical solution to this problem.

The equations were linearized (D h) to agree with the assumptions

in this section, and the bottom friction neglected. The "no flow"

boundary condition qx = 0 was prescribed at the coast with a boundary

condition r, = 0 assumed at the shelf break. Note that the qx = 0

condition is comparable to -= 0. Figure 2-19a depicts a shelf of

constant depth divided into 20 one-dimensional grid elements, as used

for the calculations. The "storm" or triangular shear stress dis-

tribution as illustrated in Figure 2-19b was moved over the shelf

waters at a constant speed, C. A slight approximation was used here.

Once the edge of a shear stress step encounters a new grid element,

the stress was assumed constant over that element for a time At = Ax/C,

or the time for the leading edge to traverse the grid element. The

numerical results were put into dimensionless form as in Equation (2-23)

and compared to the analytical results yielding Figure 2-20.




54



Notes:
(I) -r Lineorly Exlropoloted to Coast

(2) Boundary Conditions
(o) q =0 at Coast

(b) 77 0 in Grid (I)


- --- --- l.,!


h I,
I _


,Shear Stress Distribution


Moves at Constant
i/ \ Velocity C.
C z

-77 =. 0-


I i I I i I I
I I
-. I I I


20 19 18 17 16 15 14 13 12 II 10 9 8 7

20 Grid Elements
x


FIGURE 2-19a DEFINITION SKETCH

CALCULATIONS


To



'/2 To

'/6 r


I I
I I I I
I I I i
S I I

6 5 4 3 2 1=I


- [VIL--
Ax~


T- h-* CO


FOR NUMERICAL
FOR NUMERICAL


Idealized Distribution

-----.----------- I Discrete Approximation

/ I I I
"-------------

........ l
./ "I i I
/III' x


FIGURE 2-19b FINITE DIFFERENCE MODEL OF

TRIANGULAR SHEAR STRESS DISTRIBUTION


_ __


I



































.I f I I I I I I I I I
.1 .2 3 .4 5 .6 .7 8 9 1.0
r (Analytic)







FIGURE 2-20 COMPARISON OF NUMERICAL AND
ANALYTICAL SOLUTIONS












CHAPTER 3

THE NUMERICAL MODEL

a. One-dimensional Numerical Model


In most prototype situations it is unrealistic to make the simpli-

fying assumptions necessary for a successful analytical treatment of the

storm surge problem. To incorporate such features as an irregular coast-

line, irregular bathymetry, islands, and arbitrary wind stress patterns,

a numerical approach becomes necessary. The actual numerical approach or

algorithm is frequently referred to as a numerical model. Models may be

roughly categorized as one-, two-, and three-dimensional. One- and two-

dimensional models are presented in this chapter.

The one-dimensional model considers motion in only one direction,

e.g. along a perpendicular to the coastline. It is presented here for two

reasons: its simplicity allows physical insight into the problem and it

provides a check of the two-dimensional model, which can approach a

one-dimensional model as a limiting condition.

The development of the one-dimensional numerical model begins by

referring to Equations 2-6.


at 20(sing)qy -- gD an + ) (2-6a)
at P *x yx p ny b)
x. x


y + 2w(sing)qx = -D r gD n+ ( ) (2-6b)
t P 3y 6-y 0 n y b y







aqx q n
l+ + (2-6c)
ax ay at

Considering only motion in the x-direction leads to Equations (3-la) and

(3-1b)


S_ gD a -b ) (3-la)
at pX x p xI b

aqx +n
X +


xa +t- = 0 (3-lb)
ax at
pfqxlql
Using the quadratic friction law, Tb 8D Equations
x
(3-la) and (3-lb) can be written as


aq x D n n TnX fqx qx
- gD 3n + ~ (3-2a)
at p aX aX p 8D

aq
-- = 0 (3-2b)
Dx at

These are the equations to be used for the one-dimensional calculations.

The finite difference integration scheme used for the numerical

calculations is similar to that of Reid and Bodine [10] and Verma and

Dean [28]. Figure 3-1 provides a schematic diagram of a typical coastal

region showing the coordinate system and defining the variables in use.

The vertical coordinate, z, is oriented positively upward and referenced

to mean sea level, MISL; the horizontal coordinate x is oriented toward

shore with x = 0 occurring at the shelf break. In the computations n,

h, T, and p are considered uniform over the grid width LN,and centered

at the midpoint of the segment. The quantity, qx, denotes the flux in

volume per unit width through the plane representing the seaward bound-

ary of a grid element.







































C


I C


O



z
0



UI


O
z
-J







LL.
Y-

0
I -






0
z

-U
wo
0B


..........----- -- ---








For clarity the variables are placed into an indicial form where

qx(i+l,n) indicates the volume transport per unit width through the
plane section at the location x = i Lx and at time nt, i.e. flow into

the grid designated by index i+l. Similarly, rj(i,n+1/2) indicates the

value of n for the i- grid element. It is assumed to be centered at

x = (i+1/2)Ax and to be evaluated at time t = (n+1/2)At.

Following the scheme in Appendix A, the equations are placed in

finite difference form as Equations (3-3a) and (3-3b)

,n l) ,n) + t [D(i- n+ /2)+D(i,n+1/2)
q (i,n+1) [qx 1,n) + at { 2,x
x

{[n(i-l,n+1/2) r(i,n+l/2)] g + [p ((i-l,n+1/2) p (i,n+1/2)]}


+ at [T (i-l,n+1/2) + r (i,n+1/2)]] (3-3a)
2p r1
Sx x

where

Fx= 1 t f
F = 1+ 2 Iqx(i,n)I[D(i-l,n+1/2)+D(i,n+l/2)]-2


The continuity Equation (3-2b) provides, in finite difference form,

the equation for the free surface,


n(i,n+3/2) = n(i,n+1/2) + -' [qx(i,n+l) qx(i+l,n+l)] (3-3b)


To illustrate the calculation procedure assume that all values of

qx are available up to and including t = nAt. Beginning at the
seaward boundary, calculate qx's at time (n+l)at using Equations (3-3a).

Next, calculate n's at (n+3/2)at using Equation (3-3b), then advance the

pressure and shear stress to (n+3/2)At and calculate qx's at time








(n42)At, and so forth. In effect, averages over small segments are being

made in both space and time to approximate the true solution. This re-

sults in a type of leapfrog method, since qx's at time (n+l)At are used

to extend the solution from n(i,n+1/2) to n(i,n+3/2) and similarly

n and D at (n+3/2) are used to obtain qx(i,n+2). The process is con-

tinued until a desired time is reached, such as a surge maximum at the

coast.

The seaward boundary condition, the barometric tide, is applied

to grid (1,n), where the wind setup is zero and h(l,n) o. This con-

dition is written as n(l,n) = 1.15 [p. p (l,n)] where n is in feet and

p. is the undisturbed pressure in inches of mercury far from the storm.

The landward boundary condition of no flow requires qx(ic,n) = 0. The

initial conditions used will require all qx's = 0 and the grid to be in

barometric equilibrium. The only remaining data required for the calcu-

lations are the pressure and shear stress distributions as functions of

time. Figure 3-2 shows the movement of such a system at constant

velocity and with constant form. The required wind and pressure data are

extrapolated linearly from the storm model described in Section 2-c.

The model also allows for the "flooding" or "draining" of a grid

element. Consider the element (ic,n) in Figure 3-1, if n(ic-l,n)

<-h(ic) then the system is as shown, if n(ic-l,n) > h(ic) then the grid

element (ic) must be flooded or activated. In practice when activating

a grid element an instability can occur which results in switching be-

tween the activated and inactivated states. Also remember that at the

moment of grid activation D = 0 causing trouble with the calculations at

that point. To avoid these problems, the following steps are taken.





















































































































L." -.
-'Q~ ) (If

) 0. ) L
IJ)








Once activated the grid element is constrained to remain activated fo.r a

specified number of time steps. Five time steps has been found satis-

factory in the present model. In addition to the above, D is set equal

to .01 ft. in the just-flooded grid element prior to the calculations for

the next time step to prevent zero divide problems. This introduces a

small error into the system. Since flooding occurs infrequently with

respect to the number of time steps this error is considered acceptable.

To calculate a water velocity u from qx an average depth over ad-

jacent grids is taken and u(i,n) is approximated as

I
u(i,n) = 2qx(i,n) [D(i-l,n+1/2) + D(i,n+1/2)]- (3-4)

Since u(i,n) is not used in the calculations, the approximations in

Equation (3-4) are not considered significant. A more accurate value,

centered at time n, which must be calculated one time step after the

qx's is given by Equation (3-5).

u(i,n) = 4qx(i,n) [D(i-l,n-1/2) + D(i-l,n+1/2) + D(i,n-1/2)

1
+ D(i,n+1/2)] (3-5)

It is assumed that the hurricane or pressure and shear stress dis-

tributions approach normal to the coast. Numerical stability for the
nt -1/2
one-dimensional model requires that < (gD ) where D is the
Ax max max
maximum total depth expected anywhere in the system. Thus given a bottom

profile and, shear stress and pressure, as a function of time, the surge

can be calculated using this reasonably simple technique. This model is

particularly useful in quantifying effects of the numerical model that are

not inherent in the equations of motion but result as an effect of dis-




63


cretization. The model is also easily checked and thus provides a basis

for comparison with other models which are not easily checked.













b. Development in Two Dimensions


The development in two dimensions is analogous to that

dimension and therefore much of the detail will be omitted.

that the two-dimensional model considers motion in the x and

refer to Equation (2-6), and employing the quadratic bottom


_q_ an D P +qq
x 2mw(sin,)q + gD -n D ;P x f
-t y X P X p -CD7


aq dp ri fjqq
y + 2w(sinP)qx + gD r D '1 y
at x dy p y P 8D


dq aq
x + V. ; = 0
;X Dy at


in one

Recalling

y directions,

friction obtain



(3-6a)



(3-6b)




(3-6c)


These equations are now in forms which can be applied and where no

other approximations need be made. Expressing Equation (3-6a) in finite

difference form as described in Appendix B leads to the following for the

x-equation.

qx(i,j,n+l) = [qx(ijn) + 2w(sino)qy (i,j,n)&t


+ At { (i-I,j,n+1/2) + T (i,j,n+1/2)}
x x

+ {D(i-1,j,n+1/2) + D(i,j,n+1/2)) { [pn(i- ,j,n+1/2)

-P (i,j,n+1/2)] + g [n(i-1,j ,n+1/2) n(i,j,n+1/2)]}] (3-71)







where

F = 1 + { [qx(ij,n]2 + [qy(i,j,n)]2}1/2


[D(i-l,j,n+1/2) + D(i,j,n+1/2]-2 (3-7b)

Equations (3-7a), (3-7b), (3-8a), and (3-8b) are adapted using the

form given by Reid and Bodine [10], with changes in notation.

Following the same procedure for the y-equation yields


q(i,j,n+) [qy(i,j,n) 2w(sin) q- (i,j,n)At
y y


2pA { (i,j,n+1/2) + (i,j-1 ,n+1/2)}
y y

At 1
+ x {D(i,j,n+l/2) + D(i,j-l,n+l/2)} {- [pn(i,j-l,n+1/2)


pn(i,j,n+1/2)] + g[n(i,j-l,n+1/2) n(i,j,n+1/2)]}] (3-8a)

where
fAt 1f2
F = 1 + ~ {[qx(i,j,n)]2 + qy(i,j,n)]2} 1/2


[D(i,j,n+1/2) + D(i,j-l,n+1/2)]-2 (3-8b)

Here, Ax and Ay are the grid spacings as shown in Figure 3-3. Figure

3-3 also defines the variables in Equation (3-7a) and the sense of

operation of a grid element.

The values for the variables, n, p h, D, T are assumed to be

centered in each grid as shown in Figure 3-3. The values are, after

that, assumed uniform over the entire grid element for the calculation

at that time step.

Since qx(i ,j,n) and q (i,j,n) occur at different points it becomes




66










c cI
Lw






C.O
---- x ---------.i

or
3 c f c L


U5 )

---



0 I

c c2 LL
-


ow w
S oo
C "': I O

000

a EE
E5 .e 2
.- Q




0 "L

0 E0 O
W V3 C


L. 0
-, )<>- -J


. ....t Ll
a o o

Sr- a. .





"r')
S E Eu





-" ." 3"
.C Lo







necessary to obtain spatially averaged values q (i,j,n) and q (i,j,n)

for use in the Coriolis term and the friction terms Fx, F Thus,

q,,(i,j,n) and qx(i,ji,n) are considered to be at the same point and
J
qx(i,j,n) and q (ij,n) are considered to be at the same point. This
is illustrated more clearly in Figures 3-4a and 3-4b. q Ci,j,n)

is defined as


qy(ij,n) = [qy(i,j,n) + q (i-l,j,n) + q (i,j+l,n) + q (i-l,j+l,n)]

and qx is expressed by

-1
q (i,j,n) = [qx(i,j,n) + qx(i+l,j,n) + qx(i+l,j-l,n) + q (i,j-l,n)]


The equation of continuity, Equation (3-6c), becomes in finite

difference form,


n(i,j,n+3/2) = n(i,j,n+1/2) + A- [qx(i-l,j,n+l)-qx(i,j,n+l)]


+ [q (i,j-l,n+l) qy(i,j,n+l)] (3-9a)


The total depth D is defined by


D(i,j,c) = h(i,j) + n(ij,E) (3-9b)

where is an arbitrary time increment.

To carry out an actual computation, a grid network is set up

as illustrated in Figure 3-5.

The coordinate convention is the same as Figure 2-1 with x-

positive toward shore and y-positive to the left while facing the

shoreline. The boundary conditions are illustrated in Figure 3-5 and

are described as follows: 1) the fluid flux across any sea-land






















COMPONENTS OF y(i, j, n)


-t-
qx (+ j)

(i, j n)
qy(i, j, n)


qC(i, j, n)

( i-I, j)


qx(i+1, j-I)



q (i, j 1, n)
-----t----
qx(i, j-1, n)

(i-1, j-l)


FIGURE 3-4b


COMPONENTS OF (i, j, n)


FIGURE 3-4a




69



r7 Centered in Seaward and Lateral Boundary Rov's is
Equated to the Sum of Static Borometric Tide
and Astronomical Tide


I-
qs Normal to Land Boundaries
are Equated to Zero


FIGURE 3-5


SCHEMA OF GRID SYSTEM AND

BOUNDARY CONDITIONS EMPLOYED IN

TWO DIMENSIONAL CALCULATIONS


L.


S"0

o
-Ju
L5c
aI UJ







boundary is equal to zero; 2) the boundary condition at the seaward

boundary is given as n equal to the barometric tide plus the astro-

nomical tide at that point; 3) at the lateral boundaries qx = 0;

4) also at the lateral boundaries n is set equal to the barometric

tide plus the astronomical tide. Furthermore it will be assumed that

the lateral boundaries are far enough away from the region of interest

that the lateral boundary conditions have little effect. It should be

reemphasized that the lateral conditions are those conditions arrived

at empirically, which allow calculations without instabilities or other

problems with the model. Initially the system is considered quiescent

and in barometric equilibrium. This condition is of little effect

providing the storm begins sufficiently far offshore. In this model

the bathymetry and shoreline are arbitrary and the presence of islands

is allowed.

Currents in the coastal waters are of interest for several

reasons, among them: 1) to aid in making and evaluating future current

and surge calculations, 2) at the coastline hurricane-generated

currents are of interest as a cause of beach erosion, 3) hurricane-

generated currents may be instrumental in producing extreme forces on

ocean structures, and 4) currents are instrumental in dispersion and

pollution studies. Thus, it is of further interest to describe the

scheme used to calculate the water velocity at a specified grid

element. As in the one-dimensional case the total depth at time n-1/2

is used to calculate u(i,j,n) and v(i,j,n). In order to maintain a

consistent approach, an average in both the x and y directions is

taken. This leads to the following value of u(i,j,n) and for v(i,j,n)

as described in Figure 3-6.





















Note:
Shaded Areas
Weighted Twice



Figures Indicate Grid Elements
Involved in Calculating Weighted
Total Depth for Water Velocity
Calculations at Grid (i,j)


FIGURE 3-6


FILTERS USED TO CALCULATE WATER
VELOCITIES


( (i,) (j-1)



(i-l, j~ l) j-l 1- )


(i+l, j) (i+1, j-1)



(i,j) 4- (i, j-1)
q y(i, j,n)


(i- j) (i-I, j-I)




72


u(i,j,n) = 8qx(i,j,n) [2D(i-l,j,n)+2D(i,j,n)+D(i,j+l,n)+D(i-l,j-l,n)
+D(i ,j-1 ,n)+D(i-l ,j+l ,n)]-'

v(i,j,n) = 8q (i,j,n) [2D(i,j-1l,n)+2D(i,j,n)+D(i-1l,j,n)+D(i-1l,j-1l,n)
+D(i+1 ,j,n)+D(i+1 ,j-1 ,n)]-'

Adjacent to the land-sea boundary it is necessary to modify this

procedure because of the no flow boundary condition.














c. Effects of Discretization and Checks of the Numerical Model

In any model of this type it is desirable to obtain some

estimate of the effects of the grid size, time increments, and lateral

boundary conditions on the model output. A check of the programming

is also necessary. One way to accomplish this is with a series of two

or more programs which although of different concept or design, converge

to identical solutions for some limiting case. This does not neces-

sarily prove all are correct but provides confidence in the results.

A one-dimensional time dependent program should give results

similar to those of a two-dimensional time dependent program if the

radius scale of the hurricane is large enough with respect to some

length scale characteristic of the system. A comparison of this type

has been carried out for an extremely large storm, with the results

and storm characteristics as shown in Figure 3-7. One would, in

general, expect the one-dimensional program to over predict water

levels far from the coast since a cross-section of the storm is used

and this effect extends laterally to infinity, but in the two-dimen-

sional case the pressures and stresses die off with distance from the

center line. Near the coast the Coriolis effect and the two-dimensional

effects dominate. The one-dimensional plot shown in Figure 3-7a repre-

sents a cross section of maximum surge heights. The bottom was a

linear slope approximated by a finite difference form.




















10


Legend:
-- One Dmensional Model
Two Dimensional Model
No!es:
Surge Heights Plotted ore Maximum
Values Occuring During Passage of Storm


V=10 Knots f = 02
R=40 Nauticol Miles ,---
---
AP= 432 Inches of Mercury

Bottom Slope =.0005 .-

------------ -
-.. ... . .. I


100 90 80 70 60 50 40 30 20
Distance From Coast (Nautical Miles)


FIGURE 3-7a COMPARISON OF MAXIMUM SURGE HEIGHTS
CALCULATED BY THE ONE- AND TWO-

DIMENSIONAL MODELS


S80


S60
C()

* 40
5


)


/-




',


Legend:.
x-Component of Wind Velocity
-- Surface Pressure

Note:


30 >
U
(U

0

29 L





28
&


20-
Profile Token Through Point of Maximum Winds 0
R,V, and LP c"(e Given in Figure 3-70
0 L - I 1 '27 U)
-100 -80 -60 -40 -20 0 20 40 60 80 100
Distance, Along x-Coordinale, From Storm Center (Noutical Miles)
FIGURE 3-7b PROFILES OF STORM CHARACTERISTICS
USED IN COMPARISON


I


''


\
.\


,\
.\
.\


r
S\
'\
,\
'\
Ii

I0


I


120


100




















+120

+80



U')
" +40


rI

-, -40
-5
80


-120


Y- Component of
Wind Velocity


00 -80 -60 -40 -20/ 0 +20 +40 +60 +80 +100
Distance Along X-Coordinole
from Storm Center (Nautical
Miles)


Y-COMPONENT OF WIND VELOCITY


FIGURE 3-7c







As a further comparison, a time independent or steady state

model was developed as a further check of the computer modeling. It
aqx
is one dimensional and time independent in that is equal to zero.
L is
qx must then be equal to zero and the bottom friction term is

eliminated. The numerical technique consists of a Runge-Kutta

starter and Adam's Method predictor [29]. The stepwise bottom profile,

with slope .0005 as before, was included and flooding allowed, also

as before. The storm characteristics for the one-dimensional time

independent case are shown in Figure 3-7b, c, with the storm center

placed 10 nautical miles offshore. Note that although the storm is

stationary it has the characteristics of a storm traveling with

v = 10 knots. The numerical results and the steplike bottom profile

are illustrated by the dotted line in Figure 3-8.

In order to compare these results with those of the one-

dimensional time-dependent routine the time-dependent routine was

run with identical storm characteristics and then stopped when the

storm center was 10 nautical miles offshore. In this case the

bottom friction, f, was set equal to zero and as a result the solution

oscillated. The shaded boxes in Figure 3-8 represent the range of

oscillations. The stepwise approximation to the uniform slope,

.0005, was again used in this case.

How well does the 10 nautical mile step-like bottom profile

model a true linear slope? This question can be answered as follows.

The time independent model was run with all variables as before, with

the single exception of a uniformly sloping bottom. These results are

plotted as the solid line in Figure 3-8. In deep water the results are










t'-




L o 0
') 7
01)


U)ciJ
co 0 0
GJ ( -


-~0
t o o
Sc _c J

-2 U )0



c
I c I a
w1 1


SI I

\- \ O



UJ
O





I 0
C E W \ cr




w0 6
c : 0










o O
0 00










S i LUjo
I -.
--~

_I I Lii




c 0
L o I















(iOj) 3BJns uJ1S
0l~j 0EJn 0 0I







nearly identical. In the shallow region the differences become

noticeable. At the risk of over generalizing, this stepwise grid

network is permissible for values of D = h + n greater than 20 feet;

for values less than 20 feet the errors can become significant.

In the two-dimensional problem another "error-dimension"

exists, that of the effect of placement of lateral boundaries. This

is, in addition to, the effects of discretization. The lateral bound-

aries of the model are an artifice and if they are not a sufficient

distance from the zone of high wind stress then their proximity will

influence the calculated storm tides and velocities. In order to

evaluate this effect some basis for comparison must be available. A

two-dimensional analytical treatment of storm tides on continental

shelves is available from Dean and Pearce [21].

Consider the situation illustrated in Figure 3-9 in which a

"strip" wind field is uniform within a total width, w, and zero out-

side of this width. The Coriolis effect is omitted and the bottom

shear stress is expressed as proportional to the transport components

(i.e. linearized). Only shear stresses in the horizontal plane are

considered, i.e. adjacent water columns are considered to be laterally

uncoupled by stress terms; it may be anticipated that this will result

in a discontinuity in qx, the transport component in the x-direction.

The windstress in the region jyj < w/2 drives the water onshore and

results in a storm tide field which is a maximum at (x' E x/Z = 1.0,

y' = y/k = 0). The qx field in this region is a maximum at the

shelf break and decreases to zero at the coastline. It is clear, on an

intuitive basis, that the narrower the wind field relative to the

shelf width, the smaller the storm tide will be. For a very






















of Active
Stress


Elevation


FIGURE 3-9


IDEALIZED SHELF
WITH UNIFORM D


RESPONSE PROBLEM
EPTH


Shelf -
Break








broad wind field, the storm tide at the shelf should approach the value

corresponding to a one-dimensional condition. For wind fields which

are laterally limited, there will be a gradient in the storm tide

field such that the storm tide approaches zero as |y'| -+ The

surface gradient and wind stress fields will cause an onshore transport

in the wind stress region, an offshore transport outside of this

region and a y component of transport (" longshore transport") which

is directed outward from the region of onshore wind stress. Response

in terms of dimensionless n, qx and qy are presented in Figures

3-10a, b, c.

If the (artificial) lateral boundaries in the numerical model are

not a sufficient distance from the zone of high wind stress, the

proximity of the model boundaries will influence the calculated storm

tides and velocities. In order to evaluate this proximity effect,

numerical model computations were carried out for the shelf geometry

shown in Figure 3-11. The friction and total water depth terms in

the numerical model were linearized in order that the differences in

quantities calculated by the numerical model and the analytical

solution should be indicative only of artificial effects introduced

by the numerical model. The results of these calculations are de-

scribedin the following paragraphs.

The first set of calculations was carried out for a ratio of

shelf width to one-half the wind field breadth (c /(w/2)) of 1.83

and a grid size, Lx, which is variable. The lateral extent of the

shelf, incorporated into the model is denoted w* and the ratio of

wind field width to shelf extent, a*(i.e. a*=w/w*). The parameter







y
t

32










21


C.
o05



^ o\' ~~0~


IL-


0 0.2 0.4 0.6 0.8 1.0


FIGURE 3-10a


ISOLINES OF DIMENSIONLESS STORM
TIDE, ANALYTICAL SOLUTION


Noes
(I) See Figure 3-9 for
Definition Sketch

(2) Storm Tide Symmetric
About y'= 0

(3) Results Shown Apply
for a= -/(w/2) = 2.0


0
( L
tt,







Y
t

3- Notes:
(I) See Figure 3-9 for /
Definition Sketch
(2) q'x is Symmetric
About y'= 0 A
(3) Results Sho,vn Apply
For a = /(w/2) = 2.0
2C



I 7
/
Lq;









din d S l r e s s /
0 -
ILJ / / A ,


0 0.2 0.4 0.6 0.8 1.0


FIGURE 3-lOb ISOLINES OF DIMENSIONLESS TRANSPORT


IN x-DIRECTION, ANALYTICAL SOLUTION







I
t










2-



I
I
L5/


I


Notes:
(I) See Figure 3-9 for
Definition Sketch
(2) q' is Antisymmetricol
About y'= 0
(3) Results Shown Apply
for a = -e/(w/2) = 2.0




0
mn


0.2 0.4 0.6 0.8 1.0


FIGURE 3-10c ISOLINES


OF DIMlENSIONLESS


IN y-DIRECTION, ANALYTICAL


TRANSPORT
SOLUTION


&'J


<




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs