A COUPLED DISCRETE SPECTRAL WAVE HINDCAST MODEL

By

LI-HWA LIN

A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL

OF THE UNIVERSITY OF FLORIDA IN

PARTIAL FULFILLMENT OF THE REQUIREMENTS

FOR THE DEGREE OF DOCTOR OF PHILOSOPHY

UNIVERSITY OF FLORIDA

1988

ACKNOWLEDGEMENTS

The author would like to express his sincere gratitude to his adviser, Dr. Hsiang

Wang, Department Chairman of Coastal and Oceanographic Engineering, for his

ideas, advice and support throughout this study, and also to coadviser Dr. Bent A.

Christensen, Professor of Civil Engineering, for his guidance and encouragement.

The author is also grateful to Dr. Ashish J. Mehta, Dr. Daniel P. Spangler, and

Dr. Louis H. Motz for their participation as members of his graduate committee.

The technical assistance and programing suggestion provided by Sidney Schofield

and Subarna Malakar are also appreciated.

The author wishes to thank the U.S. Army Corps of Engineers, Coastal En-

gineering Research Center, for furnishing the WIS model and for supporting the

Coastal Data Network program. The Kings Bay data were collected under the

sponsorship of U.S. Navy, Ship Research and Development Center.

Finally, the author would like to thank his wife Li-chu for her support, encour-

agement and patience.

TABLE OF CONTENTS

ACKNOWLEDGEMENTS ..................... .......... ii

LIST OF FIGURES ................................... v

LIST OF TABLES ........................... .... viii

ABSTRACT ...... ............... ............... ix

CHAPTERS

1 INTRODUCTION ................. .............. 1

1.1 Background .... ............. ..... .......... 1

1.2 Scope of W ork .............................. 4

2 REVIEW OF EXISTING MODELS AND THEIR PERFORMANCE 8

2.1 Introduction .................. ............... 8

2.2 Spectral Definitions ... .............. ............. 9

2.3 Spectral Shape Function Models ........... ......... 11

2.4 Deepwater Parametric Model ......... ........ 16

2.5 Deepwater Discrete Model ..... ............... . 21

2.6 Shallow Water Model .. ............... .......... 26

2.7 Wave Data and Comparisons ...................... 34

3 A UNIFIED COASTAL WAVE PREDICTION MODEL ......... 44

3.1 Introduction .. .. ... ........ .......... 44

3.2 UCWP Deepwater Submodel ...................... 45

3.2.1 Generation terms ........ ......... .......... 45

3.2.2 Dissipation Terms ................. ........ 62

3.3 UCWP Shallow Water Submodel .............. ....... 64

3.3.1 Wave Energy Growth in Shallow Water ............ 64

3.3.2 Wave Energy Dissipation in Shallow Water ......... 66

3.4 Wind Model .................... ............ 66

4 MODEL PERFORMANCE ............................ 69

4.1 -Introduction ............. ..... ............. 69

4.2 Model Performance in Hypothetical Cases ............... 70

4.3 Model Performance in Real Cases .... .............. .. 76

4.4 Wind-Wave Relationship in a Fully Developed Sea ......... 96

4.5 Comparison of Various Wave Energy Flux Elements .......... 102

4.6 Influence of c! in Shallow Water Wave Transformation Process 105

5 CONCLUSIONS ... ... .... ... ... ............... 115

5.1 Summary .... ........... ................. 115

5.2 Conclusions ................................ 119

5.3 Future Studies .......... ............ ........ 122

APPENDICES

A SUMMARY OF PARAMETERIZED PROGNOSTIC EQUATIONS ... 123

A.1 Prognostic Equations ..... .... .. ..... ........... 123

A.2 Elements of Di and Gi .... ..... ................ 123

B PROCEDURE IN SOLVING PARAMETERIZED Gin,1(f) ........ 125

C SUMMARY OF COMPUTED AND MEASURED NON-DIRECTIONAL

FREQUENCY SPECTRA AT THE BUOY #41006 LOCATION ... 128

C.1 Introduction .................. ............. 128

C.2 Summary of Figures .... ........... ............ 128

D PLOTS OF WAVE ENERGY FLUX ELEMENTS ............. 150

D.1 Summary of Figures ................... ......... 150

BIBLIOGRAPHY ........ ......................... 157

BIOGRAPHICAL SKETCH ..................... ...... 161

LIST OF FIGURES

2.1 Comparison of nondimensionalized shapes of the P-M spectrum,

the Kitaigorodskii spectrum, the JONSWAP spectrum and the

modified JONSWAP spectrum..................... 13

2.2 Plots of Hs(f,0) versus 00 .................... 15

2.3 Nonlinear transfer of energy (Hasselmann et al., 1973). ...... 20

2.4 Schematic illustration of growth of a spectral component, E(f),

at (a) fetch-limited case and at (b) duration-limited case . 30

2.5 Wave gage stations ...................... .. 35

2.6 (a) Weather map showing northeast winds generated by a high

pressure system (1/12/84) . . . . .. 36

2.6 (b) Weather map showing northeast winds generated by a high

pressure system (1/21/84) . . . . .. 37

2.7 Deepwater grid coordinate system . . . ..... 38

2.8 Comparison of NOAA and WIS wave spectra for January 1984

at the Buoy #41006 location ..................... 40

2.9 Comparisons of NOAA and WIS modal periods and significant

wave heights at the Buoy #41006 location. . . . 41

2.10 Comparison of NOAA and WIS wave energy scale parameters

at the Buoy #41006 location ..................... 42

3.1 Combined sources of nonlinear term and generation term in WIS

m odel .................................. 47

3.2 Plots of E(f.) versus f. ............. ......... 50

3.3 Plots of ay versus fo/f,2t. ................... 52

3.4 Plot of a versus f,/fo ......................... .. 53

3.5 Plot of 7 computed from Eq.(3.8) versus y from the least square

method............... .... ............... .. 55

3.6 Plots of aeb versus fm/fo ......................... 56

3.7 Simulations of spectral growth at duration-limited case; U10 =

30 knots. ............................... 61

3.8 Comparison of spectral decays between measured and computed

values .................................. 63

4.1 Comparison of Deepwater hindcast models in fetch-limited cases. 71

4.2 Comparison of Deepwater hindcast models in duration-limited

cases. ................... .. ........... 74

4.3 Shallow water grid system. . . . ... .. 76

4.4 Comparison of hindcasted and measured winds at the Buoy

#41006location ................... ........ 78

4.5 Comparison of NOAA and UCWP wave spectra for January 1984

at the Buoy #41006 location .......... .. ...... .. 80

4.6 Comparisons of NOAA and UCWP modal periods and signifi-

cant wave heights at the Buoy #41006 location. . ... 81

4.7 Comparison of NOAA and UCWP energy scale parameters at

the Buoy #41006 location. ................... .. 82

4.8 Comparison of CDN and UCWP wave spectra for January 1984

at the Kings Bay gage location...................... 83

4.9 Comparison of CDN and UCWP marginal directional spectra at

the Kings Bay gage location. ..................... 85

4.10 Comparisons of CDN and UCWP average wave direction, modal

periods and significant wave heights at the Kings Bay station.. 86

4.11 Tracks of Hurricanes Diana, Isidore, and Josephine (1984). 87

4.12 (a) Weather map showing northeast winds generated by the Hur-

ricane Diana (9/8/84).......... ......... ..... 88

4.12 (b) Weather map showing northeast winds generated by the Hur-

ricane Isidore (9/27/84) ............ ........... 89

4.12 (c) Weather map showing northeast winds generated by the Hur-

ricane Josephine (10/11/84) .................... 90

4.13 Comparisons of NOAA and UCWP modal periods and signifi-

cant wave heights at the Buoy #41006 location. . ... 92

4.14 (a) Comparisons of CDN and UCWP average wave direction,

modal periods and significant wave heights at the Kings Bay

gage location. ...... ....................... 93

4.14 (b) Comparisons of CDN and UCWP average wave direction,

modal periods and significant wave heights at the Marineland

gage location ................... ......... 94

4.14 (c) Comparisons of CDN and UCWP average wave direction,

modal periods and significant wave heights at the Cape Canaveral

Station. ................................. 95

4.15 Plots of (a) H. versus U10, and (b) Tm versus Ul0 at various fetch

and duration............................. 98

4.16 Plots of 27rUo/gT, versus 10'H,/gT, in (a) fetch-limited case

and in (b) duration-limited case. . . . . 101

4.17 Plots of 2rUio/gT, versus 10SH./gT, based on the NOAA buoy

#41006 1984-1985 data. ....................... 103

4.18 Comparisons of wave heights and modal periods computed in

September 1984 for cy=0.0, 0.005, 0.01, 0.02, and 0.04 .. 107

4.19 (a) Plot of individual shallow water effect sources for cy=0.0 109

4.19 (b) Plot of individual shallow water effect sources for c =0.005 110

4.19 (c) Plot of individual shallow water effect sources for cj=0.01 111

4.19 (d) Plot of individual shallow water effect sources for c1=0.02 112

4.19 (e) Plot of individual shallow water effect sources for cf=0.04 113

C.1 Summary of computed and measured frequency spectra. 129

D.1 Plot of deepwater energy fluxes in growing period; Case 1. 151

D.2 Plot of deepwater energy fluxes in decaying period; Case 2. ... 152

D.3 Plot of deepwater energy fluxes during small winds; Case 3. .153

D.4 Plot of individual shallow water effect sources during growing

period; Case 4. ............................ 154

D.5 Plot of individual shallow water effect sources during decaying

period; Case 5 ............................ 155

D.6 Plot of individual shallow water effect sources during decaying

period; Case 6. ............................ 156

LIST OF TABLES

4.1 Summary of Correlations between computed and measured val-

ues of H, and Tm at the buoy #41006 location. . .... 92

4.2 Summary of Correlations between computed and measured val-

ues of H, and T. at three shallow water stations. . ... 96

4.3 Summary of energy flux budget at the buoy #41006 location. 104

4.4 Summary of energy flux budget at the Marineland station. .105

4.5 Comparison of 6,,,.,,, 6,,., ,u)u E..,.,, and T,.,T,,m for

cf=0.0, 0.005, 0.01, 0.02, and 0.04. . . . ... 108

4.6 Summary of energy flux budget for cases with c1=0.0, 0.005 0.01,

0.02, and 0.04, respectively. ..................... 114

Abstract of Dissertation Presented to the Graduate School

of the University of Florida in Partial Fulfillment of the

Requirements for the Degree of Doctor of Philosophy

A COUPLED DISCRETE SPECTRAL WAVE HINDCAST MODEL

By

LI-HWA LIN

April 1988

Chairman: Dr. Bent A. Christensen

Cochairman: Dr. Hsiang Wang

Major Department: Civil Engineering

In late 1976 a numerical wave hindcast model was developed by the Wave In-

formation Study (WIS) at the U.S. Army Engineer Waterways Experiment Station

(WES) and it was employed to generate 20 years of wave conditions in 10-meter

water depth along the eastern seaboard of the United States. The model comprises

of deepwater and shallow water submodels. The significant nonlinear wave-wave

interaction process is included in the deepwater submodel based on the study of

parametric model. Recently, more detailed field wave information became available

in both offshore and nearshore regions, particularly shallow water directional wave

information reduced from data collected along the Florida Coast. Based upon those

data, the WIS model is reexamined.

A new wave model, named here as the Unified Coastal Wave Prediction (UCWP)

model, is developed. This model is of a discrete type and it links a wind field

model to a combined deepwater and shallow water wave prediction model. The

deepwater wave submodel follows the WIS model but with revised source terms

and dissipation terms. The shallow water wave submodel utilities the refraction

law and the equation of conservation of energy to solve for the shallow water waves.

The UCWP model represents a significant improvement over existing models,

at least for applications to the southeastern region of the United States. The im-

provement is particularly evident in wave period and wave direction. In addition to

the development of the new numerical model, a number of fundamental aspects of

air-sea interactions are also examined in light of the field data. New semi-empirical

formulas are proposed for a number of important parameters, including the JON-

SWAP spectral parameters, required to define wind wave spectra.

CHAPTER 1

INTRODUCTION

1.1 Background

Wave information is required for a variety of purposes, including design of

coastal and offshore structures, planning of marine operations, ship design, coastal

process studies, and coastal erosion prevention, etc. Because of a general lack of

measured wave data, particularly long-term data, wave hindcast, which is a tech-

nique in calculating the events of wind-generated water surface waves in the past, is

often the best means to estimate representative wave statistics and is almost always

the only means to estimate characteristics of the infrequently occurring extreme

waves.

In early days, methods utilized for hindcasting wind waves involved empirically

based nomographys and tables to evaluate wave height and wave period for promi-

nent waves. A well-known example is the Sverdrup-Munk-Bretschneider (SMB)

model which was established by matching some predescribed math functions to

measured wave heights and wave periods. The model was established in a period

during and after World War I and was used exclusively for either hindcast or fore-

cast of ocean waves before more modern wave models were developed later in the

1970s. With more physical processes governing the nature of waves being known

and concepts of wave energy spectra becoming widely accepted, spectral models

were then developed for wave hindcast, or forecast.

The spectral models are commonly categorized into two major branches: the

discrete model and the parametric model. With wind information as input, both

the discrete model and the parametric model generate wave energy as a function

2

of frequency and direction, propagate energy across oceanic region, and consider

energy spreading and dissipation. The difference between the two models is that

the discrete model evaluates wave spectra by summing all individual spectral com-

ponents corresponding to different physical processes, while the parametric model

constructs wave spectra by solving numerically spectral parameters from the prog-

nostic equations which represent the projections of wind-wave relationships in the

spectral parameter space.

The advantage of using parametric model to hindcast the wind wave spectra is

that its numerical solution procedure is simple. Computations in solving spectral

parameters in the parametric model are solely performed in the frequency domain.

The directional spreading of wave energy is considered, in the model by multiplying

the resulting frequency spectrum with an angular distribution function. The as-

sumed angular distribution function is a normalized bell-shaped, symmetrical func-

tion with respect to local wind direction. The 'disadvantage of using the parametric

model is that it does not include swell and, hence, cannot rigorously consider all the

physical processes involved in wave transformation. Furthermore, the parametric

model can not properly generate wave spectra in shallow water where bottom effects

become important. The discrete model in more flexible and can include swells just

as easily as that of the regular wind waves. It can also handle easily the shallow

water effects, such as refraction, shoaling, bottom friction, etc. However, in the

discrete model computations involved are usually tedious and time-consuming since

they have to be carried out in the combined frequency and direction domain.

Currently, there are a number of spectral models, both discrete and parametric

types, available for wave hindcasting. Most of them are for deepwater wave applica-

tion. These wave models are based upon similar theoretical considerations but their

numerical methods and computational requirements are vastly different. They all

contain empirical functions or coefficients derived from actual wave measurements

3

and are calibrated on limited data base. Due to the scarcity of the measured wave

data, none of the present models were tested against a wide range of wind and wave

conditions. Comparisons for extreme waves during severe storms and hurricanes

were few and comparisons with measured wave energy directional spectra were ex-

tremely rare. Accordingly, there is considerable uncertainty about the accuracy of

the existing wave models.

A few comparisons of different wave models were made with each other in the

past. Among them, the work by Resio and Vincent (1979) offered an interesting

comparison among the SMB model, Hasselmann's (1976) parametric model, a ver-

sion of the Pierson et al.(1966) model, Barnett's (1968) model, and a discrete model

by Resio and Vincent (1977, 1978). Comparisons of wave height versus fetch and

duration for selected constant wind speeds indicate significant differences among

these models and, thus, raise concern for their accuracy.

There are several problems which need to be addressed if wave models are to

improve their accuracy and to be applicable for a wide range of conditions in both

deep and shallow waters. First, source and sink terms should be continuously re-

fined with availability of more and better wave data. Second, wave models should

be calibrated with directional wave information and with a wide range of wind con-

ditions, i.e., large systems, small systems, extreme events, etc. Finally, there is a

need to acquire the accurate wind information since wave model relies on wind field

as input. The wind information can be either measured from the field or generated

by a wind model, based on the sea surface pressure distribution. Accurate specifi-

cation of wind fields through a numerical model is sometimes difficult. Therefore, it

is necessary to improve the wind model or to properly incorporate measured wind

data as input to the wave model.

4

1.2 Scope of Work

Recently, more detailed field wave information became available in both deepwa-

ter and shallow water regions. For instance, long term wave data have been collected

in deep water by several NOAA (National Oceanic and Atmospheric Adminstration)

maintained buoys in the Northwest Atlantic Ocean and in shallow water by twelve

CDN (Coastal Data Network) underwater stations along the Florida coast. Based

upon these data, a new coupled wave hindcast model is developed. This model,

referred to here as the Unified Coastal Wave Prediction (UCWP) model, is of a dis-

crete type. It links a wind field model to a combined deepwater and shallow water

wave prediction model. The deepwater submodel follows the same technique used

in WIS with revised source and dissipation terms. The shallow water submodel

utilizes the algorithm developed by Chen and Wang (1983) with modified source

term proposed by Vincent (1984). The model is unified in the sense that the deep-

water model and the shallow water model are based on the same set of fundamental

equations; both employ discrete finite different scheme with compatible source and

dissipation terms; the models are on the same level in detail and are expected to

yield same degree of accuracy.

A reexamination of two previously developed deepwater wave models, the WIS

model and the parametric model by Hasselmann et al.(1973), has revealed some

differences about the source functions utilized by the two models. First, the non-

conservative term of wave energy dissipation is completely neglected from the para-

metric model while it is included in the discrete model to account for energy loss

of low frequency components. Second, the conservative wave energy transfer, as

a result of nonlinear wave-wave interactions, is included in the parametric model

based on JONSWAP spectrum while it is included in the discrete model based on

Kitaigorodskii spectrum. Third, the wave growth due to wind input is considered in-

dependent to the nonlinear wave-wave interactions in the parametric model, while

5

it is assumed to depend on the nonlinear wave-wave interactions in the discrete

model. A well-known problem about the parametric model is its inability of treat-

ing swell components. The parametric is preferred by some due to its simplicity in

computations, while the discrete model is preferred by others because of its straight

forward physical interpretation.

Inspection of field wave data collected by NOAA buoys reveals that nonconser-

vative turbulent dissipation of wave energy induced by wave motion should not be

neglected in a wave model. This dissipation is incorporated in the UCWP model

to consume wave energy corresponding to an exponential decay rate. The non-

linear wave-wave interactions are included in the UCWP model based on a mod-

ified JONSWAP spectrum proposed by Donelan et al.(1985). The migration of

peak energy frequency caused by this wave-wave interaction mechanism is deter-

mined in the model based on the idealized prognostic equations derived by Hassel-

mann et al.(1976). The threshold spectrum which indicates a minimal condition for

the occurrence of nonlinear wave-wave interactions is empirically determined based

upon the field data.

The UCWP shallow water submodel is developed with the inclusion of most

of the known significant shallow water effects, such as refraction, shoaling, perco-

lation, bottom friction, and wave breaking due to depth limitation, etc. All these

shallow water effects have been documented in the past (Chen and Wang, 1983;

Shemdin et al., 1978). The turbulent dissipation of wave energy is also included in

the model. The nonlinear wave-wave interaction is, however, not included in the

model. The growth of waves due to wind input is evaluated in the model based on

a TMA spectrum, as suggested by Bouws et al.(1983). The five parameter TMA

spectrum has been correlated to wind and wave conditions by Vincent (1984).

During the autumn and winter seasons extratropical storms and hurricanes fre-

quently occur in the Northwest Atlantic Ocean. Waves generated by the storms and

6

hurricanes in general are the major cause of severe coastal erosion along the south-

east coast of the United States. The model's performance is, thus, tested against

such severe weather conditions. The results obtained from the UCWP deepwater

submodel are compared with measured wave data from NOAA buoys deployed in

the Northwest Atlantic Ocean. The results obtained from the UCWP shallow water

submodel are compared with measured wave data from CDN wave gages.

In summary, the detailed research included in this study can be briefly described

as follows:

(1) A number of standard spectral functions and angular spreading distribution

functions are reviewed. These are important in defining an appropriate general

directional wave spectrum.

(2) Two existing deepwater wave energy spectral models, the discrete model

by the Wave Study Information (WIS) group of Waterways Experiment Station

(Resio and Tracy, 1983) and the parametric model by Hasselmann et al.(1973), are

rederived to highlight the basic assumptions and fundamental principles.

(3) Deepwater wave generation mechanisms are reevaluated based upon the

currently acquired wave data.

(4) Equations that are either theoretically derived or empirically proposed as

shallow water source functions with respect to wave generation, propagation, and

dissipation, are examined and summarized. These are later utilized to build a

rigorous shallow water model.

(5) A new discrete model, named here as the Unified Coastal Wave Prediction

(UCWP) model, is developed for generating wave spectra in the deep and shallow

water. The model links a wind field model to a combined deepwater and shallow

water wave prediction model.

(6) Computer algorithm for the UCWP model is developed to facilitate numer-

ical computations.

7

(7) Results from the deepwater UCWP model are compared with the SMB

model, the WIS model, and Hasselmann's (1976) parametric model for fetch-limited

and duration-limited cases.

(8) Model outputs are comparedwith those measured by deepwater buoys and

shallow water wave gages.

(9) Relative importance of individual elements affecting the deep and shallow

water spectral transformations is analyzed.

(10) Different values of bottom friction coefficient are tested in the UCWP

shallow water model to investigate the sensitivity of model responses.

CHAPTER 2

REVIEW OF EXISTING MODELS AND THEIR PERFORMANCE

2.1 Introduction

The approach taken in the present study for ocean wave hindcasting along the

Florida coast is to develop a numerical model based on state-of-the-art deepwater

and shallow water wave prediction techniques. Currently, deepwater wave models

are either of discrete type which calculate wave energy spectra based on discretized

wave energy equation or parametric type which calculate a set of parameters that

in turn specifies the energy spectra. The existing shallow water wave models are

invariably of discrete type since the mechanisms involved in the wave transport in

the shallow water, such as wave refraction, shoaling, percolation, wave energy dis-

sipation, and shallow water wave-wave interactions, are too complicated to permit

parameterization.

The material presented in this chapter includes the definitions of directional

spectra as well as the unidirectional frequency spectra which play essential roles

in wave hindcast models, and a presentation of the typical theories and methods

applied to the current deepwater and shallow water wave models. Two widely known

deepwater wave models, the parametric spectral model by Hasselmann et al.(1973)

and the discrete spectral model by the Wave Information Study (WIS) group of

the Waterways Experiment Station (WES), are reviewed as a demonstration of

spectral methods. It is shown that the deepwater parametric spectral model in

general is meant for wind waves only and cannot account for the swells while the

deepwater discrete spectral model by WIS has no difficulty in tracing both. It is

further shown that the WIS spectral model has the same ability as the parametric

9

model to include the nonlinear wave-wave interactions in computing wave energy

transformation in the deep water; this is accomplished using the formula which

estimates the nonlinear wave energy transformation rate in the frequency domain

derived from the prognostic equations of parametric method.

The shallow water wave model is reviewed following the works done by Chen

and Wang (1983), Collins (1972), and Shemdin et al.(1978). Among the major

contributors of shallow water effects, wave refraction and shoaling are included in the

advection terms, and others are acting as a part of the forcing functions governing

the generation, dissipation, and nonlinear transfer of wave energy in shallow water

region.

The model performance is only shown here for the deepwater model of WIS.

The shallow water model which requires the output of deepwater model to specify

the boundary condition is not tested until the latter is modified. Its performance

will be discussed in Chapter 4.

2.2 Spectral Definitions

A directional wave spectrum is essentially a representation of the wave num-

ber spectrum of water surface waves that follow the linear wave theory dispersion

relationship between wave number and wave frequency. By linear superposition of

wave components, an irregular water surface can be approximated by taking the

stochastic integral over the wave number plane, a, which is taken here as identical

to the still water surface plane (Thomson, 1972).

(1, t)= Re[ da(k)Wei(=-1t)] = fj [da(k)ei('t) + da*()e-'i J] (2.1)

where rl is water surface elevation, t is time, s- = (zX,yi)-is the horizontal location

of interest in rectangular coordinates, a is radian frequency, k = (k.cos 0, k.sin 0) is

a vector wave number with magnitude k =1 kI and direction 0, da(k) indicates an

infinitesimal wave amplitude (complex number) which is the contribution of waves

in a wave number interval di, an asterisk appearing in superscript indicates the

10

complex conjugate, the operator Re[ ] indicates the real part of a complex number.

The wave amplitude component da(f) can be transformed into da(a, 0) since k can

be expressed by a direction 0 and a frequency a. To a first approximation, the wave

frequency satisfies the well-known linear wave theory dispersion equation

a' = gk-tanh(kh), k =Ik (2.2)

where g is the gravitational constant and h is the water depth. A small contribution

of the water surface waves representing a specific vector wave number is

d(, t) = (da()e4-*I1) + da*()e-'(I-)] (2.3)

Thus, the mean square value of dr can be evaluated, which has a form

[dn (,t)2 = 1 lim 1 / {Ida(Ck) 2e2( -Iat)+ 2da(c)da*(c)

4 T-o6 Tf-T1P

+[da*()]2e-i2(rl-t)}1d = da(da()da )

= | da(f) 1= da(o, 0) 1= E(O, 0)dad0 = dE(a, 0) (2.4)

2 2

where T is the time duration of interest, an overhead bar indicates the operation of

arithmetic mean; E(u, 0) is the directional frequency spectrum which is defined as

E(o, 0) = [dr(s, t)]2/dad0 = 1 Ida(o, 0) 12/dadO (2.5)

According to linear wave theory, the wave energy per unit water surface area per

unit water weight is equal to half the squared value of wave amplitude (Ippen, 1966).

Therefore, Eq.(2.4) specifies a quantity proportional to wave energy associated with

a particular frequency range, da, and directional interval, dO. However, Eq.(2.4)

may be regarded as a contribution from an infinitesimal amplitude wave with ampli-

tude da(a, 0) at the mid da and dO. Thus, E(a, 0) is interpreted as the distribution

of wave energy per unit water surface area per unit water weight per unit frequency

per unit direction. The corresponding one-dimensional frequency spectrum E(a) is

defined by integrating E(a, 0) over all 0.

E(a) = E(o, 0)dO (2.6)

The marginal directional spectrum is defined as

E(0) = E(ao,)do (2.7)

The total wave energy per unit water surface area per unit water weight is defined

by integrating E(a) over all a, or integrating E(0) over all 0, or integrating E(a, 0)

over all a and 0. That is,

Etota = E(a)da.= E(0)d = f E(, O)dad0 (2.8)

Since the spectral density defined in Eq.(2.5) or Eq.(2.6) has dimensions of rate

of change of wave energy density, the name 'power spectral density' is often used.

However, an alternate name 'energy spectral density' is more commonly used in

coastal engineering.

2.3 Spectral Shape Function Models

In the study of wave hindcast through spectral models, the standard spectral

shape functions are served as fundamental structure unit in building the wave en-

ergy spectral states in the sea. Three frequency spectral shape functions are well-

accepted in the past: (i) Pierson-Moskowitz (P-M) spectrum, (ii) Kitaigorodskii

spectrum, and (iii) JONSWAP spectrum. The Pierson-Moskowitz spectrum (1964),

EpM(f), is given by

eXp[.74(fo/f)4] apg(f

EpM(f) = (2 exp[--074(fo)] = ( exp[- (f/f)4] (2.9)

where f = a/21r is the wave frequency with units 1/time, a, is known as Phillips'

constant, fo = g/27rU1s is the frequency of the deepwater waves whose phase

speed is identical to the wind speed U19g. at a height 19.5 m above sea surface,

and f,, = 0.88fo is the wave frequency of the spectral peak. More often, the name

12

'spectral peak frequency' is referred to the frequency of spectral peak, f,, and the

wave period corresponding to this frequency is known as 'modal period'.

The Kitaigorodskii spectrum(1962), EK(f), is given by

f ag2fs(22r)-' exp[l (f//)f)], if f < f,

EKW = (2.10)

I ag2f-(2r)4, if f / f

which contains two free parameters, the energy scale factor a and the spectral

peak frequency f,.. The JONSWAP spectrum is derived from the Joint North Sea

Wave Project (Hasselmann et al. 1973) and constitutes a modification to the P-M

spectrum to provide for a much more sharply peaked spectrum. The JONSWAP

spectrum, Ej(f), is given by

1 f- f,,m

g2 5 exp[-( fm)2]

EJag( = exp[-(f f)'] 2 aL fm (2.11)

Ef)---(2r)4fs 4

with

= Ja, if f -f,,

i a, if f > fm

The function contains five free parameters, two scale parameters a and f,, and

three shape parameters -, oa and ab. The parameter 7 is the peak-enhancement

factor which is the ratio of the peak value of the spectrum to the peak value of the

corresponding P-M spectrum with the same values of a and f,,. The parameters

oa and ao yield the scales of left and right peak widths, respectively.

Figure 2.1 shows a comparison of spectral shapes, normalized by the spec-

tral peak energy of the P-M spectrum, of the P-M spectrum, the Kitaigorodskii

spectrum, and the mean JONSWAP spectrum, which has its spectral parameters

7, oa, and ab equal to 3.3, 0.07, and 0.09, respectively.. It is seen that the mean

JONSWAP spectrum and the Kitaigorodskii spectrum have similar shape and they

are much more sharply peaked than the P-M spectrum.

All three spectra shown above indicate an 'equilibrium range' on the rear face

of the spectral shape with the associated spectral energy densities being inversely

E(f)

EpM(fm)

1 2

f/fm

Figure 2.1: Comparison of nondimensionalized shapes of the P-M spectrum, the

Kitaigorodskii spectrum, the JONSWAP spectrum and the modified JONSWAP

spectrum.

14

proportional to frequency to the fifth power. However, more recent studies made by

Donelan et al.(1985) and Battjes et al.(1987) using both field and laboratory data

showed that an f-4 frequency dependence at the tail part of the wave spectrum

may be more adequate for the wind-generated waves than the f-5 spectral shape.

Based upon the assumption of the f-4 dependence on the rear face of frequency

spectrum, Donelan et al.(1985) suggested a modified JONSWAP spectrum

2 exp[-- f)]

EJM(f) = ag exp[-(fm/f)4 2 a fm (2.12)

(27r)4f41f

with all its free parameters a, f,, 7, and oa defined the same as in JONSWAP

spectrum. This modified JONSWAP spectrum is immediately applicable to any

discrete spectral method provided fm is known or defined in advance. This is

because all the empirical relationships between each of the free parameters a, -y,

and oa were established, based on field data, by Donelan et al.(1985) as functions of

fo/fm, A modified JONSWAP spectrum with its parameters and oab being equal

to the mean JONSWAP ones is also shown in Fig. 2.1. This spectrum has wider

width and higher peak than the mean JONSWAP spectrum, and it contains more

energy in the high frequency components than the mean JONSWAP spectrum, the

P-M spectrum and the Kitaigorodskii spectrum.

The directional dependence of the wave spectrum is usually modelled by a sym-

metrical directional distribution function which is also known as the 'angular distri-

bution function' or 'spreading function'. It is often convenient to define a normalized

spreading function as the ratio of (Cartwright, 1963)

E(f, 0)

H(f,) = E(f) (2.13)

where E(f,O) = Hs(f,O)E(f) indicates the directional spectrum, which has its

wave energy densities distributed in the (f, 0) domains according to Hs(f, 0).

A number of functional forms of Hs(f, 8) have been suggested but all with

similar shape, i.e., they are bell-shaped and symmetric about a central direction. It

1.5

S

s

20

15

1.0 -

10

Hs(f, ) 7

0.5 3

2

1

0.0

-90a -450 00 450 9dr

60-0

Figure 2.2: Plots of Hs(f, ) versus 0 0o.

is usually difficult to choose one in preference to another in practical applications.

A commonly used function, which is proportional to a power cosine function, is

given here (Sand, 1979).

S 22 2(S+l)os 25(O ), if 0 e001 /2

H (f,0)= r(2S+1) (2.14)

0, elsewhere

where S = S(f) is called the spreading parameter, which is always positive and a

function of frequency, 0o = 9o(f) is the central angle of the symmetric spreading

directions. Figure 2.2 shows several distributions of Hs(f, ) with the spreading

parameter ranging from 1 to 20. It is shown that Hs(f,0) associated with higher

spreading parameter can have very sharply peaked, narrow-band directional distri-

bution. When this Hs(f,0) is applied in a wave prediction model, the spreading

parameter S is commonly defined equal to 1 and 2, respectively, for all frequency

components of locally generated and nonlinearly transferred energies. It can be

16

seen in Fig. 2.2 that both spreading functions, Hs(f,0), corresponding to S = 1

and S = 2, respectively, have relatively wide distributions over a 1800 directional

domain.

2.4 Deepwater Parametric Model

In the absence of current in the ocean, the general equation governing the trans-

port of wave energy is given by (Hasselmann et al., 1973)

dE(f,0,S,t) aE(f,0, 8,t)+ (f, E(f, t) f, ,) (2.15)

dt a t + C,(f,a) = G(f,0,-, t) (2.15)

where E(f, O,, t) is the energy density of the spectral element with reference to

still water surface plane at horizontal location s and time t, with frequency f and

direction of propagation 0, i is the path of this spectral element in which di =

jdtje*', Cg(f, s) = aa/ak r = a(27rf)/aklI is the group velocity associated with this

spectral element, k is the magnitude of the associated wave number, and G(f, 0, ', t)

is the net source function representing the rate of energy transfer into or out of this

element. This equation is commonly referred to as 'wave energy balance equation',

or 'wave energy transport equation', or simply 'transport equation'. In the text to

follow, the notations of E(f, 0) and G(f, 0) are often used for simplicity.

The magnitude of wave group velocity Cg(f, ') can be derived from the disper-

sion relation given in Eq.(2.2) and expressed as

1 kh

C, (f,) = n C = ( + )C (2.16)

2 sinh(2kh)

where C = a/k = 27rf/k is phase speed of waves with frequency equal to f, h is

water depth at the location -, and n is the ratio of wave group velocity to wave

phase velocity. It is noted that both phase velocity C and wave number k are

implicit functions of f. The presence of the source function G(f, 0, 5, t) is thought

to include all the possible physical terms which are responsible for the local growth

and decay of water surface waves.

17

The deepwater wave is usually defined in coastal engineering by

tanh(kh) > 0.9962, or (kh) > %r (2.17)

Under the above condition, the wave properties are hardly influenced by the bottom.

Therefore, the shallow water effects, such as shoaling, refraction, wave dissipation

due to bottom friction, etc., are no longer important in the deepwater wave hindcast.

It is noted that in deep water, the dispersion relation is reduced to

o2 = (27rf ) = gk

and the deepwater wave group velocity, Co(f), can be expressed as

1 g g

Cgo(f) = 1 C= =

2 2a 47rf

which is not dependent on water depth and, hence, is not dependent on location.

The parametric spectral model was originally developed by Hasselmann et al.

at 1973. The model was later applied by many investigators without further mod-

ifications (Weare and Worthington, 1978; Gfinther et al., 1979). The presentation

of the parametric model in this section is mainly based on the work of Hassel-

mann et al.(1973).

To derive the equations solving the spectral parameters, the directional wave

spectrum used in the transport equation, Eq.(2.14), is assumed to have the following

expression:

Er(f )2 cos2(O 00), if 10 Ool_ _r/2

Ej(f,0) = E,(f)H,(f,e) = E, r (2.18)

0, elsewhere

where Hi(f,B), the spreading function, is a special case of Eq.(2.14) with S equal

to 1, and 0o coincides with the downwind direction. By first substituting Eq.(2.18)

into the transport equation, Eq.(2.15), and then intergrating the transport equation

over the direction domain, 0 < 0 < 27r, the direction-independent transport equation

is obtained as

dEj(f) E(f) 9E (f)

dt- + c( (f) =' G(f) (2.19)

with the directional averaged group velocity Cgo(f) equal to

',C o f) Ej,(f )dO 8 g

Co(f) -j ~ (2.20)

EJ (f) 37r 47rf

and the directional integrated source function G(f) equal to

G(f) = G(f, O)dO (2.21)

Note that symbols g and t, which stand for location and time, respectively, are

not shown in Eqs.(2.19), (2.20), and (2.21), and in the following derivations for

simplicity.

The source function G(f) is usually represented by three quantities, the energy

input from the atmosphere, Gi,(f), the dissipation of energy, Gdi,(f), and the

nonlinear energy transfer, G,,(f). That is,

G(f) = G,.(f) + Gdi,(f) + Gn(f) (2.22)

Since the spectral shape is assumed to be consistent with the JONSWAP spectrum

in the parametric model and is presumed in a state of equilibrium range, the dis-

sipation Gdi.(f) is in general dropped from Eq.(2.22). The rest two source terms

Gi,(f) and G,4(f) are thought to be equally important for the growth and transfer

of wave energy in the ocean. It was demonstrated by Resio and Vincent (1979) in

numerical model studies that the duration required to achieve fetch-limited condi-

tion in the case dominated by atmospheric input is of the same order as the case

dominated by the wave-wave interactions.

The growth rate Gn(f) is commonly approximated in terms of a linear func-

tion of one-dimensional frequency spectrum (Snyder and Cox, 1966). However, for

convenience the quantity Gi.(f) is often assumed to be proportional to E(f) in the

19

parametric model. For instance, the energy input source term utilized by Hassel-

mann et al.(1973) in their wave model is linearly proportional to one-dimensional

frequency spectrum and has the following expression:

G,.(f) = 5s f UEj(f) (2.23)

where s -0.0012 is the ratio of air density to sea water density, and Uo0 is the wind

speed at 10 m elevation above sea surface.

The general form of the nonlinear transfer source Ga(f) for a single peaked

spectrum is known to have a positive lobe appearing on the frontface of the spec-

trum, a negative second lobe appearing on the rearface of the spectrum closer to the

spectral peak, and a positive third lobe appearing on the rearface of the spectrum

away from the spectral peak (see Fig. 2.3). In wave hindcast, only the first two

lobes are important since the energy transfer into the third lobe is at such frequen-

cies that viscous and turbulent dissipations are assumed to keep the local energy

components small.

The nonlinear transfer function deduced by Hasselmann et al.(1973) has the

expression

g2as

Ga(f) = D1 ~t/(f/f,) (2.24)

where D1 is a proportional constant whose magnitude is on the order of 10-',

and k(f/f,,) is a non-dimensional function depending on the spectral shape. By

neglecting the energy transfer in the third lobe, where waves of high frequencies are

more likely governed by the equilibrium f-' law, Eq.(2.24) can be approximated by

the expression

( 2.22fmEjE(f), if < f,,

Gra(f) = -22a'fnE.(f), if f f(2.25)

With the source function given, the direction-independent transport equation,

Eq.(2.19), can be then projected onto the JONSWAP parameter space to yield

20

MEAN JONSHAP SPECTRUM

0.2

f= 0.3 (H) 3

a = 0.0095

S= 3.3

a, = 0.07

Ob = 0.09 2

E(f) 0.1 i 05. G,

m2

)Hz (mZ)

-1

I I I I I I I I

0.0 0.2 0.4 0.6 0.8 1.0

f (Hz)

Figure 2.3: Nonlinear transfer of energy (Hasselmann et al., 1973).

(Hasselmann et al., 1973)

abi iOb

At + Diji j = Gi, i,j = 1,2,3,4,5 (2.26)

where b1 = f,,, b = a, b3 = 7, b4 = ~,, and be = ob. The equations given in

Eq.(2.26) are known as the 'prognostic equations' in the parametric model. The

Dij matrix is a generalized propagation velocity, and Gi are the generalized source

functions. A list of the Dij matrix and G, functions so obtained, with Gdi,(f) equal

to zero, is given in Appendix A.

Thus, by solving Eqs.(2.26) simultaneously with given wind information, the five

free JONSWAP parameters can be obtained and the resulting JONSWAP spectrum

can be then expanded according to the H (f, 0) spreading function to generate a

directional distribution of wave energy densities.

The nature of the parametric approach is that only single-peaked spectra are

allowed in the hindcast, and as a result temporal and spatial variations of swells can

not be included in the model. To include the swell portion of the wave energy, it is

necessary to incorporate a discrete spectral method to accommodate the transport

of swells. Unlike the parametric method, which resorts to solving the free spectral

parameters of a predesignated spectral function, the discrete spectral method tracks

the spectral densities of both the preexisting and currently generated waves in the

wave field of interest. Accordingly, the spectra so obtained could have arbitrary

shapes including multiple peaks. The details of the discrete spectral method are

given in the next section.

2.5 Deepwater Discrete Model

The development of deepwater discrete spectral models is again based on the

wave energy transport equation given in Eq.(2.15). The method utilized in discrete

spectral models generates directly the directional distributions of wave energy ac-

cording to the local atmosphere energy input and then transfers them advectively in

the horizontal domain according to the governing equation. Since the fundamental

22

structure of a discrete spectral model does not change from one model to another,

the solution procedure used in a discrete spectral model will be described in this

section through the presentation of the discrete model developed by Wave Informa-

tion Study (WIS) group at the Waterways Experiment Station (WES), U.S. Army

Corps of Engineers.

The WIS model was originally created to hindcast the deepwater waves in the

Atlantic Ocean between 1956 and 1975. The documentation of the model can be

found in the series of WIS reports of Vol. 8 (Jensen, 1983a) and 12 (Resio and

Tracy, 1983). The version of the WIS model presented here is the one modified by

the same WIS group in 1984 by studying the deepwater wave hindcast in the Pacific

Ocean.

The WIS model solves Eq.(2.15) numerically based on the following discretized

form:

Ek+1(f ,0) = Ek(f,9) C(f) AEk') At + G(f ,0) At (2.27)

where the superscript indicates the time level. The second term on the RHS is

energy advection and the third term is the net energy source, and these two terms

are uncoupled.

The advective transfer of wave energy is computed in the model using the

method developed by Hasselmann et al.(1968). This method divides the entire

direction domain, 0 < 0 < 27r, uniformly into sixteen directional bands, each hav-

ing an angular span of 27r/16 radian, and wave energy contained in each band is

advected independently without coupling with others. To express the method in a

numerical scheme, the advective energy E(f, 0),d may be written explicitly as

Ek+1(f, ) = Ek(f, ) Co(f) AE(f, ) At

The source term is further divided into three components:

G(f, ) At = [G,, (f, ) + G ,,(f,0) + Gk(f, )] At (2.28)

23

where Gi,, Gdi,, and G,. represent atmospheric energy input, wave energy dissi-

pation and nonlinear wave-wave interaction, respectively. The Gi. term, which is

computed at the current time level, is further split into two parts: (i) Gi,,1(f,0),

the rate of growth of the prevailing waves and (ii) Gn,s(f, 0), the rate of generation

of new wind waves in the downwind direction.

In WIS model, the two atmospheric input terms are formulated as follows:

g4.4 Atfa cos4(8 ), if f : fm

G,,At = (27r)4f5 31r

0, elsewhere

(2.29)

with the mean propagation direction, 7, defined as

-tan{ ff sin OE(f, 0) dfdO

8 = tan { (2.30)

ff cosOE(f,0) dfde

where the scale parameters a and fm are to be determined later.

ag 2 fm4 2

S exp ( ) cos2(O 8,), if f, 0.25 (Hz), f < f,,

f,m 0.2 (Hz)

Gilzat = ( ag2 2

G2,At =2 2 cos2(O 0,), if f, > 0.25 (Hz), f > fn,

(2,)4f5 r. -

f, > 0.2 (Hz)

0, elsewhere

(2.31)

where 0, denotes downwind direction and f, denotes a peak energy saturation

frequency which is not necessarily equal to the peak spectral frequency and is defined

as

E(f.) = max {E(f) > E,(f)}

24

E,(f) is the one-dimensional saturation spectrum defined as

cag2

S(2.)4f5 if f < fo = g/2rUlo

E.(f) =

(2) 5 (1 (folfm)4] + EpM(f)(fo/fm),, elsewhere

(27r)4fr

(2.32)

where fo is the wind-wave resonant frequency to which the corresponding wave

phase speed is identical to the wind speed at 10 m elevation, U1o, and EpM denotes

the P-M spectrum.

There are two scale parameters a and f, to be determined first at time level

k+1. The scale parameter fm is obtained by solving the prognostic equations de-

rived in the (a, f,) plane as described by Hasselmann et al.(1976) in their para-

metric model. By allowing variations of only two free parameters a and f,, in the

JONSWAP parametric space, the a- and f,- component prognostic equations can

be derived as

aa aa a8fm fmUI0413 5 a3f

+ 0.47 C, (f) T. + 0.2 (f) = 0.005 ( ) f -5 a (2.33)

and

af, af foa

+ ,,(f) 0.07Co(f) a = -0.54 a2f (2.34)

If the changing wind speed can be expressed as a power law of time,

U10 = Uo( )

Uo

with Uo a reference wind speed, f. and a can be solved from Eqs.(2.33) and (2.34)

as

fmUIo = A(t)(3/7)(1) (2.35)

g Uo

and

a = B(fmo)2/3s (2.36)

g

25

with constants A and B given by (Hasselmann et al., 1976)

1 + 1.51 q

To obtain f, at time level k+1, the value at t in Eq.(2.35) is replaced by t + At,

which results in

U~ f A( + )(s/7)(q-1) (2.37)

g Uo Uo

Eliminating gt/Uo from Eqs.(2.35) and (2.37) yields

U /k+1 = ((fk)(7/3)(q-1) + gAt( -(/s)(-1) (s/) l) (2.38)

U7o Uo Ag"

For constant wind speed within the duration At, the above equation reduces to

fk+1 = {(fk)-7/3 + o)4/At}-3/7 (2.39)

where a,j is the nonlinear wave energy transport coefficient defined by Hassel-

mann et al.(1976) and suggested as equal to 0.00138. In WIS model, the fk+1

value is replaced by fo if the spectral peak frequency is less than fo. Also, instead

of using Eq.(2.36) to compute a, WIS model employs the following relationship

suggested by Resio and Vincent (1977):

a = 0.0489 ( ) -0.2 (2.40)

where U. is the surface frictional velocity obtained by assuming logarithmic over

water wind profile

U, 1 z

U=1 In z- (2.41)

U. 0.4 zo

where U, is wind speed at a height z above the sea surface, and zo is the equivalent

roughness at water surface given by

0.1525 U2

zo = + 0.0144-* 0.00371, (in cgs units) (2.42)

/ 9

26

Although Eq.(2.36) is not used in the WIS model, it can be linked to an

empirical relationship (g2Etoti/Uo)(f,,/Uo/g)4/a=constant suggested by Hassel-

mann et al.(1976) to yield a power law a ~ (g2Etoj/IU/o)-0'2 which is similar to

Eq.(2.40).

The transfer of wave energy due to a nonlinear wave-wave interaction source

G,~(f, 0) is formulated based upon Kitaigorodskii spectrum (1962) and a spreading

function with the central angle parallel to the mean wave direction 0 (Resio, 1981).

g2ct3 8

0.0023 At t(f/fm) cos'(0 0), if 10 0\ < r/2

GnAt = 3 (2.43)

0, elsewhere

where

f exp[1 (f,n/f)4], if f < f,

I-(M/ /f)5, if f > fm

The dissipation includes two parts: the energy loss of wind waves and the energy

loss of swell. The energy loss of wind waves is accounted for by restricting the growth

to saturated spectral condition with f-' law on the downslope side. The loss of wave

energy in swell is adopted from Resio (1981), in terms of one-dimensional frequency

spectrum

-4 At E( f E(f), if f,, < f < 0.7 f,

G At = ) (2.44)

0, elsewhere

where f,, indicates peak spectral frequency of swell.

2.6 Shallow Water Model

As deepwater waves propagate into area of shallow Water, they will sense the

bottom and wave energy will be redistributed and dissipated in the direction and

frequency domains due to refraction, shoaling, bottom friction, percolation, etc.

Among these shallow water effects, only refraction and shoaling can be handled

fairly readily in many situations, others require further theoretical and experimental

27

studies. Shallow water waves can also grow if winds are present and large enough

to add energy into waves as a result of air-sea interaction.

In shallow water the wave energy transport equation can be expressed as follows:

(Chen and Wang, 1983)

aE(f, ) aE(f, ) C,(f,8) cos 8 aE(f, 8) C, (f, ) sin

at + + = G(f ) (2.45)

at Qzx ay

where 0 is defined counterclockwise from the positive x axis. Here, unlike the

deepwater case, both C, and 0 are no longer independent at position (actually

water depth). In order to track energy transport, changes in wave direction due to

wave refraction must be considered simultaneously. This requires solving additional

equation governing the change of wave number vector. According to the refraction

law, the curl of wave number is identical to zero, or

(ak cos 0) (ak sin 0)

Vxk= ( a =0 (2.46)

Qy Qz

Therefore, if the source function G(f, 8) is known in advance, the directional

spectrum E(f, 8) in shallow water can be obtained by solving both Eqs.(2.45) and

(2.46). It should be noted that neither current nor tide is considered to influence

the wave spectral transformation governed by Eqs.(2.45) and (2.46). Although it is

possible to include the current and tide effects in solving the wave spectra in the

shallow water (Chen and Wang, 1983), the influences of currents and tides to the

wave spectral transformation in the shallow water will be exclusive in the present

study for the following two reasons. First, the numerical model which includes the

current and tide effects in estimating shallow water wave spectra is in general dealing

with a large matrix of equations such that execution of the model is not efficient

in the aspect of computational time. Second, the insufficient information about

interactions among waves, currents, and tides, according to the present knowledge,

may cause problem in defining the source terms associated with both current and

tide influences.

28

The source function G(f, 0) may still be represented by the summation of three

contributions, namely, generation, dissipation, and nonlinear wave-wave interaction,

although the formulas are expected to be different from the deepwater ones due to

the influence of the bottom.

G(f,0) = Gi,(f,0) + Gdi,(f,0) + G,(f,0)

Explicit formulas for individual source terms are well documented in the past (Chen

and Wang, 1983; Collins, 1972; Shemdin et al., 1978) and are briefly described in

the remaining of this section.

The generation term is usually defined by a combined Phillips' and Miles' mech-

anism (Barnett, 1968; Chen and Wang, 1983), i.e.,

G,,(f,o) = a(f,o) + PE(f,0) (2.47)

where a and / are linear and exponential growth coefficients, respectively. Contri-

butions of the first term (Phillips' mechanism) is small but essential in their model

to initiate wave growth. The coefficient / is proposed by Barnett (1968)

U10 U10

C/ = (2.48)

0, elsewhere

Under this growth, wave energy is increasing faster in the high frequency compo-

nents than the low frequency components. Therefore, saturation condition needs to

be invoked to exclude the excessive energy by means of wave breaking. In shallow

water region, the saturation condition has been suggested by Bouws et al.(1983)

to be the product of the deepwater saturated spectrum, such as the one shown in

Eq.(2.34), and an energy transfer function

(k,h) k3(f,oo) C,(f, o) 1 k2(f, oo) (2.49)

k3(f, h) C,(f, h) 2n k2(f, h)

The growth of a spectral component according to the combined Phillips' and

Miles' mechanism can be studied for the hypothesized fetch-limited and duration-

limited cases. For a steady state, one-dimensional, and frictionless situation, the

/

29

growth of a spectral frequency component in the downwind direction is governed

by the following ordinary differential equation

C,(f)dE(f) = a + 3E(f)

By specifying E(f) = 0 at z = 0 as the boundary condition, the solution for spectral

component E(f) is

E(f) = [exp( ) -1]

For small fetch, E(f) approximates (a/C,)z, which corresponds to the Phillips'

linear growth mechanism. For large fetch, E(f) approximates (a/f) exp(3z/C,),

which corresponds to the Miles' exponential growth mechanism. In the real situa-

tion the exponential growth is limited to a saturation condition and an overshoot

growth occurs between the exponential growth region and the saturation region. A

typical growth curve for one spectral component in the fetch-limited case appears

in Fig. 2.4(a).

For a duration-dependent, one-dimensional, and frictionless situation, the growth

of a spectral component can be determined from the following ordinary differential

equation

dE(f) = a + E(f)

By specifing E(f) = 0 at t = 0 as the initial condition, the general solution for E(f)

is

E(f)= ~[exp(Ot) 1]

For small duration, E(f) approximates at, which corresponds to the Phillips' mech-

anism. For large duration, E(f) approximates (a/f) exp(pt), which corresponds to

the Miles' mechanism. In the real situation the exponential growth is again limited

E(f)

EXPONENTIAL GROWTH

LINEAR

GROWTH

-E(f) = (a/C,)z

o fetch, z

(b)

OVERSHOOT SATURATION

PERIOD PERIOD

(I)

E(f)

EXPONENTIAL GROWTH

LINEAR

GROWTH

S- NE(f) = at

0

0 duration, t

Figure 2.4: Schematic illustration of growth of a spectral component, E(f), at (a)

fetch-limited case and at (b) duration-limited case.

31

to a saturation condition and an overshoot growth occurs between the exponen-

tial growth and the saturation periods. A typical growth curve for one spectral

component in the duration-limited case appears in Fig. 2.4(b).

The rate of dissipation in the shallow water, Gdi,(f,0), is in general thought

to come from five sources: bottom friction, percolation, whitecaps, breaking in

maintaining the equilibrium spectrum, and breaking due to the depth limitation.

The dissipation rate attributed to bottom friction can be determined by the

work done by the bottom stress, i.e.,

Gd,,b(f, 8) dfdO = drbdub(f, 8) /(pg)

where ub is horizontal velocity component at sea bed, r7 is the bottom shear stress,

and p, is the water density. By expressing the bottom shear stress based on the

empirical law

rT = -p, CUs(f,80) I|u(f,0)|

where c, is bottom friction coefficient, and expressing dub(f, 8) in terms of dr(f, 8)

according to linear wave theory,

dub(f, ) -- ckh dr(f,) ,

a cosh kh

the rate Gdi,,b(f, 0) may be approximated by the equation (Collins, 1972)

g k- cc (Us)

Gd(,) o cosh2 khE(f, 0) (2.50)

where (ub) denotes the root mean square value of horizontal wave orbital velocity

component at bed and can be computed from the expression

(ub) --= ub(t)]2 = j Idu (f,O)]2 = gk g(k2

c cosh2 kh [dr(f 8)]

Sc oh2 cosh kh

SE(f'dO df (2.51)

Joo -r o'2 cosh'2 kh

32

The rate of dissipation of wave energy in the shallow water due to effect of a

porous sandy bottom, Gd i,p(f, 0), can be estimated by the work done by the action

of wave-induced dynamic pressure at the surface of the bed to cause the fluid motion

in the porous sand layer.

Gda,,(f, ) dfdO = dp-dw- /(p,g) (2.52)

where pb and wb, respectively, stand for the dynamic pressure and vertical water

partial velocity at the surface of the bed. For porous sand layer of depth d this

rate of dissipation derived by Shemdin et al.(1978) is

--tanh(\/\i/\ kd)

Gdi,p(f, ) = -kX2 cs 'kh E(f, ) (2.53)

where A, and A2 are the horizontal and vertical percolation coefficients, respectively,

and they have the dimension of (length)/(time).

The dissipation of wave energy caused by damping due to kinematic viscosity,

ik, and turbulent eddy viscosity, 1t, according to a theoretical analysis made by

Lamb (1932), can be shown as

Gd,,, (f, 0) = -4 (v, + Lt)k2E(f, 0) (2.54)

The kinematic viscosity of sea water mik, which is approximately equal to 0.018

cm2/sec at 320F, is generally several orders of magnitude smaller than the turbulent

eddy viscosity vt. By neglecting the dissipation caused by the kinematic viscous

effect in Eq.(2.54), and assuming the magnitude of turbulent viscosity induced in

the water waves is linearly proportional to the product of wave phase speed and

wave length, the resulting equation may be expressed in the form

Gd,.,((f, 8) = -c, fE(f, 8) (2.55)

where c, is the constant of proportionality. This equation has been proposed by

Hasselmann (1974) to limit the growth of wave energy on the rearface of the spec-

I

33

According to the JONSWAP study, the dissipation rate defined in Eq.(2.55)

is generally found to be several orders of magnitude smaller than the generation

rate. However, actual dissipation due to whitecaps may not be small as indicated

by the JONSWAP study which utilized a hypothesized equilibrium spectral shape,

i.e., JONSWAP spectrum, to define the energy generation term that has included

the influence of dissipation due to whitecaps. Thus, without having a reasonable

estimate of dissipation due to whitecaps, the limit of wave growth in the shallow

water model is traditionally given by a saturation condition (Chen and Wang, 1983):

A(f) = B g2(2ir f)-'5 (k,h) (2.56)

where B is a coefficient to be of the order of 10-2, and Q(k, h) is defined in Eq.(2.49).

The breaking of shallow water waves due to depth limitation occurs as a conse-

quence of strong wave-bottom interaction in the nearshore region. It is this breaking

mechanism that produces the turbulent flow and enhances the significant movement

of bed materials responsible for either erosion or deposition of the bed. The estima-

tion of energy dissipation for breaking due to depth limitation is important since it

matters greatly the sediment transport in the coastal waters. However, an analyt-

ical estimation of the dissipation is not found and the determination of this value

must rely on empirical formula. An empirical criterion, which estimates the total

wave energy retained in waves after the waves break due to depth limitation, often

used in the shallow water model is (Chen and Wang, 1983)

(Etot,,)1/2 = 2 tanh(kbh) (2.57)

where the subscript b indicates the breaking condition. The value kb is usually

taken as the wave number corresponding to modal period although it is in general

not clear what value of kb should be used in the equation.

34

2.7 Wave Data and Comparisons

In 1977, the Department of Coastal and Oceanographic Engineering at Univer-

sity of Florida (UF) initiated a shallow water wave data collection program around

the State of Florida. By 1985, twelve wave stations have been established (See

Fig. 2.5). In the offshore region, there are two NOAA maintained wave buoys with

ID No. 41002 and 41006 as shown.

Of the twelve shallow water stations, data from Kings Bay, Marineland and

Cape Caneveral are eventually being used here for comparison purposes. At Kings

Bay, two PUV underwater packages as manufactured by Sea Data Inc. (Model

635-12) were deployed at 15 m and 18 m water depth. These packages are self-

contained with the ability to measure two axes water particle velocities U and V,

and the pressure information, P. Directional wave information is derived form

the P-U-V data based on the method proposed by Longuet-Higgins et al.(1963).

Detailed data handling and analysis procedures are given in the User's Guide by

UF (1984). The gages at Marineland and Cape Caneveral are shore-connected

underwater pressure transducers at 10 m and 8 m water depths respectively. Only

wave height information can be retrieved from them.

The deepwater buoy stations are located at 29.30N, 77.30W (#41006) and at

32.30N, 75.30W (#41002). The water depths of these two stations are 970 m and

3200 m, respectively. These buoys also collect 5 m-level winds every hour, based on

an 8.5 minute record.

One of the intentions of the field data measurement program is to provide wave

climate information along the south eastern seaboard of the United States. To

accomplish this goal, a numerical wave prediction model is required to fill in the

data gaps and to extend data extrapolation.

The WIS model was first tested against the measured deepwater wave data

at NOAA buoy #41006. A segment of the data covering the complete month of

Figure 2.5: Wave gage stations.

500N

350N

S\ 00N

O t-T AN IC OC A 250N

80OW 75'W 70W 65'W 60'W 550W 500W

Figure 2.6: (a) Weather map showing northeast winds generated by a high pressure

system (1/12/84).

January, 1984 was used for testing purposes for two reasons: (i) During this month,

directional wave data at Kings Bay was available for the entire month and (ii)

January is usually the month of high waves dominated by large scale high pressure

systems. Figure 2.6(a) and 2.6(b) show two examples of typical large scale high

pressure systems occurring in January. In these examples, the northeast winds

developed from the high pressure systems will generate waves that propagate to the

southeastern coast of the United States.

The grid system covers an orthogonal-spherical area with the east and the west

boundaries coinciding with 70.0W and 82.50W, respectively, and with the north

and south boundaries coinciding with 42.50N and 22.5N, respectively. The grids are

all parallel to the longitudes and latitudes such that each grid is 2.5 degrees apart

I -

5 OON

120Z 1100 //8N

4 50N

40 H 4 OON

36 35ON

20

28 4 02 30aN

1006

12

-20

80 W 75'W 70OW 65WO 600W 550W 500W

Figure 2.6: (b) Weather map showing northeast winds generated by a high pressure

system (1/21/84).

' 1

S 75' ~ 70W 65 NW 60W 55 5

12001 5/1/7

450ON

8 35N

20 30ON

.20 ATLAN IC CC 25ON

80OW 75'W 70'W 65*W 60'N 55OW 50W

Figure 2.7: Deepwater grid coordinate system.

from the next and a total of 6x9 grid points are found in the system (Fig. 2.7).

The wind information served as the input for the model was obtained by analyzing

surface pressures from daily weather maps. The theories of the geostrophic winds

and planet boundary layers were used in the analysis (see Section 3.4 wind model).

The 10 m-level winds were generated at the grid system for forty days starting on

December 21, 1983. The computed winds were then blended with the ones measured

by the NOAA buoys.

The radiative boundary condition is utilized in the wave model. However, two

modifications were made to the boundary condition. The first one was made to

the boundary where the grid system is connected with land. It simply regulates

I -

39

that no wave energy is transported to sea from land. The second one was made

to the boundary where the grid system is connected with the sea. It regulates

that the amount of wave energy advected from a grid point at the boundary to

the grid points inside the boundary is the same as the wave energy transported

from the area outside the boundary to the grid point at the boundary. Therefore,

a quasi-stationary wave energy transport process is assumed at the seaboard grid

boundaries.

The wave model was run at time increment equal to two hours for a total of

forty-one days from December 21, 1983 to January 31, 1984. The first ten days were

for initialization and the results were not used for comparison. Figure 2.8 shows the

comparison of the time series of the wave energy spectral evolution process of the

entire month. The correlations of the measured and computed modal period and

significant wave height are shown in Fig. 2.9. The significant wave height, which is

commonly denoted as H113 or H,, is defined as the average height of the one-third

highest waves. When the heights of individual waves in a deepwater wave record

are ranked from the highest to lowest, the frequency of occurrence of waves above

any given value is given to a close approximation by the cumulative form of the

Rayleigh distribution. Therefore, the significant wave height can be determined

approximately from the following equation:

H, = 4.0 (Ett)1/2

The correlation coefficient, 6, given in Fig. 2.9 is defined as follows:

N N N

6 = E[(Xc),(x~ ),-( .),(x) )]/{[l[(Xc),-(X)j' f[(x,),-(x.)m)2}1/ (2.58)

i=1 i=1 i=1 -

where X, is a computed X, Xm is a measured X, N is the sample size, and the

overhead bar indicates the operation of arithmetic mean.

The measured and hindcasted energy scale parameter, fH,./g, are plotted in

Fig. 2.10 for waves with H, >3 m, 3 m> H, > 2 m, and H, <2 m.

I -

m2

E(f) (),

(NOAA)

E(f) (2),

(WIS)

60

50

30 -

20 *

10

0,

0.2

0.1 0.0

f (Hz)

60

50

0

30 *

20

0 (,

0.0

f (Hz)

JAN. (1984)

1

10 JAN. (1984)

Figure 2.8: Comparison of NOAA and WIS wave spectra for January 1984 at the

Buoy #41006 location.

41

20

10 5---10----15----20----25--------

15 6 =-0.22 COMPUTED OATA(RIS)

TN

o *** *,-. -

(SEC) "

5

1 5 10 15 20 25 30

JANUARY (1984)

10

STATION: NORA BUOT *41006 HERSUREO DATA(NORA)

8 6 =0.49 COMPUTED DATA (MIS)

S6

(M) .4 *3

2I

1 S 10 15 20 25 30

JANUARY (1984)

Figure 2.9: Comparisons of NOAA and WIS modal periods and significant wave

heights at the Buoy #41006 location.

0.12

0.10 l

0.08 I-

10 *f ,

(WIS)

(WIS)

0.06 -

0.01

0.02 -

0.0

0.0

+ 2n>Hs 6 = -0.33

A 3m>Hs>2n 6 = -0.55

S Hsa3M 6 = 0.17

0

0

a^6 O

0.02

0.04

0.06

0.08

10. (NOAA)

g

Figure 2.10: Comparison of NOAA and WIS wave energy scale parameters at the

Buoy #41006 location.

I ~ ~. -

43

Clearly, the model performs better for hindcasting wave height ( with a cor-

relation of 0.5) than for wave frequency (with a negative correlation coefficient).

The model performance also improves progressively for higher wave environment

(Fig. 2.9) and this is particularly the case if the wind wave dominates. The model

fares the worst on hindcasting the swell component. As can be seen in Fig. 2.8,

both growth and decay of low frequency components hindcasted by the model are

clearly inadequate.

I r I

CHAPTER 3

A UNIFIED COASTAL WAVE PREDICTION MODEL

3.1 Introduction

It is shown in the last chapter that the parametric model and the discrete model

are two alternative methods to hindcast waves in deep water. Growth rates of wave

height predicted by the parametric model of Hasselmann et al.(1976) and the WIS

model were consistent in the case which waves are growing either along a fetch or

through time (Resio, 1981; Resio and Tracy, 1983). However, the comparison of the

hindcast based upon either model to the measured wave heights from buoy in the

open ocean often shows discernible deviations (Ginther et al., 1979; Corson and

Resio, 1981). The comparison of the modal periods computed from either model

to the buoy data also shows certain degree of discrepancy (Ginther et al., 1979).

Thus, it seems necessary to modify the existing models in order to improve the

accuracy of the hindcasted waves compared to the measured data.

A new unified wave model, named here as the Unified Coastal Wave Predic-

tion (UCWP) model, was developed in the present study. The model is of discrete

type and is comprised of deepwater and shallow water submodels. The deepwater

submodel is a modified WIS model. The modifications are mainly in the source

functions by incorporating some recent development in air-sea interaction. The

model is calibrated by utilizing the field wind and wavedata collected by NOAA

buoys and it is suitable for hindcasting ocean waves generated by either high pres-

sure or low pressure weather system. The shallow water submodel is based upon

the numerical algorithm developed by Chen and Wang (1983) and it is suitable for

hindcasting coastal water waves.

I-. -

45

In the following sections of this chapter, both the UCWP deepwater and shallow

water submodels will be introduced. The deepwater submodel is to be described in

Section 3.2 and the shallow water submodel is in Section 3.3. A wind model which

prepares the wind data for the generation of water surface waves in the UCWP

model is described in Section 3.4.

3.2 UCWP Deepwater Submodel

The deepwater submodel is a spectral model of the discrete type. This model

is more flexible in selecting source functions than a parametric model as it is not

always an easy task to introduce new source functions in the JONSWAP parameter

space.

The governing equations for the UCWP deepwater model are the same as the

WIS model described earlier, i.e., Eqs.(2.27) and (2.28). The advection term in

model is computed in the same manner as in the WIS model. The source terms

including both generation and dissipation in the UCWP deepwater model are mod-

ified from WIS model.

3.2.1 Generation terms

Recent studies made by Mitsuyasu and Honda (1982), and Mitsuyasu (1985)

indicate that the mechanism corresponding to wave growth due to wind input shall

be considered separately for the prevailing gravity waves and the locally generated

wind waves. Therefore, the utilization of Eq.(2.23) alone to estimate wave growth

does not seem to be a proper approach in the parametric model. On the other hand,

the WIS model which utilizes two separate mechanisms expressed by Eqs.(2.29) and

(2.31) to estimate, respectively, the growth of prevailing waves and new wind waves

seems to be more consistent with the nature of wave growth.

In the WIS model, the wave generation comprises of three terms growth of

prevailing waves, generation of new wind waves and nonlinear transfer of energy

among wave components. The growth terms are computed at time level k+1 and the

I -

46

transfer term is computed at time level k. The first and the last terms Eqs.(2.29)

and (2.43), respectively are of a similar nature, responsible for the evolution of

existing wave components. Their combined effect is shown in Fig. 3.1. Here, the

nonlinear energy transfer is seen to be the dominant mechanism and has an abrupt

cutoff at f = f,i. The growth of the down-slope face of the spectrum is shown to

be a minor process. At peak spectral frequency, f = f,, approximately 80 percent

of the energy growth is being transferred to lower frequencies. By comparison with

field data as discussed in the previous section, the predicted spectral shape is found

to be too sharp peaked due to the abrupt cutoff and the low frequency components

grow too fast because of the proposed nonlinear transfer process. A new growth

term is proposed.

Since we have shown that the wave-wave interaction source term behaves like

an apparent wind source term, it can be lumped into the wind source term for

the growth of the prevailing waves. Based upon the common practice that this

source term should be related to a stable spectral shape to maintain similarity,

we confine our search to established spectral functions by testing them against field

data collected during steady winds when the prevailing wave growth is the dominant

mechanism. After various trails, a new source function is proposed here:

2 8 r7r

c.,a fEM(f) cos'(0 ), if lO l< 3 and 10e 3 1

Gi,( f, 0) =37 2 2

0, elsewhere

(3.1)

where c, is a non-dimensional coefficient, regarded by some investigators to be a

function of density ratio of air to water, and EJM(f) is the one-dimensional modified

JONSWAP spectrum proposed recently by Donelan et al.(1985):

1 fI,)

ag2 exp[- ( )

EjM(f) = (2,r)4f4f exp[-(f,/f)4] 2 r f (3.2)

with all its free parameters a, f,n, 7, and aab defined the same as in the JONSWAP

I I

1.6

fm = 0.1(Hz)

1.4

a = 0.0075

-- Gi.,n + Gw, (WIS)

1.2 -- G,.,n, (UCWP)

106 G

1.0 I

a g/(27r)4f.'

0.8 -

I \

0.6 \

0.4f

I

I

0I

O.I \

j I \ \

/; \\

0.0 0.05 0.1 0.15 0.2 0.25 0.3

S(Hz)

Figure 3.1: Combined sources of nonlinear term and generation term in WIS model.

0.0 0.05 o.1 0.15 0.2 0.25 0.3

f (Hz)

Figure 3.1: Combined sources of nonlinear term and generation term in WIS model.

I_ _

48

spectrum. Substituting Eq.(3.2) into Eq.(3.1) and absorbing c, into a, we have

a g2 exp[--( f)2] 8

(2)fsf exp[-(f,,/f)4 2 abfm cos(O ),

Gin,,(f,e) = 1ifle-I and l,-1|<

2 -2

0, elsewhere

(3.3)

The peak frequency f, is determined from the following equation:

fM+ = {(f)-7 + a( hl cos(r- 0,)),/At}-' (3.4)

which is similar to Eq.(2.39) with the exceptions that the wind component parallel

to the mean wave propagation is used here instead of the absolute strength and

that the spectral peak frequency fk given in Eq.(2.39) is now replaced by a newly

defined peak energy transfer frequency, f,,. This peak energy transfer frequency is

defined as

E(f,m) = max {E(f) > Ec(f)}

where Ec(f) is a threshold spectrum below which no energy transfer will occur from

high frequency components to low frequency components. To use f", instead of fm

in Eq.(3.4) physically limits the energy transfer to those components exceeding the

threshold condition.

At present, there is no established guideline to determine the threshold spec-

trum, E,(f). However, based upon physical reasoning, we observe here that E'(f)

should, at least, satisfy the following criteria: (i) The threshold spectrum should

fall below the saturation equilibriumm) spectrum as energy transfer is expected to

take place before the spectral energy reaches saturation. (ii) The threshold spec-

trum should be a monotonically increasing function with decreasing frequency and

should approach equalibrium spectrum at low frequencies where energy transfers

become more and more difficult.

49

To propose a suitable Ec(f), two years measured data 1984, 1985 from

NOAA buoys #41002 and #41006 were examined. Since we are dealing with the

mechanism on the growth of prevailing waves, only those data corresponding to

increasing in both wave energy and modal period were used. Figure 3.2 shows the

plots of E(f,) versus f, from the data of buoys #41002 and #41006. If the mean

values of the data points are connected, they form a skewed bell-shaped curve with

a peak value, f,, occurring at frequency between 0.09 and 0.11 Hz. This indicates

that energy transfer experiences strong resistance when f,, becomes less than this

peak value.

Based upon the above observations, we propose here the following functional

form for Ee(f):

E,(f) = (2r)4f5 p- f)2] (3.5)

which is the P-M equilibrium spectrum modified by a multiplier. Here c,~ is an

empirical scale factor which is chosen as equal to 0.183 by matching Eq.(3.5) with the

data mean at designated f,=0.1 Hz. Both Eq.(3.5) and the original P-M equilibrium

spectrum are shown in Fig. 3.2, with the data. As can be seen, E,(f) has the basic

characteristics as discussed with the most active energy transfer process taking place

in the energy containing range.

It is noted here that plots of spectral peak energy versus spectral peak frequency

in the wave growth cases in Fig. 3.2 show that the data are scattered in the equi-

librium range indicated by the P-M spectrum. Although the data scattering could

be exaggerated by the presence of swells, which generally appear as long waves,

there is no indication that a simple curve, such as the P-M equilibrium spectrum,

can represent the relationship. The data scattering suggests active nonlinear wave-

wave interactions. That is, energy transfer from higher frequency components to

lower frequency ones is a nonlinear process and it is dependent upon multi-variables

including wind speed, wave phase speed, peak energy frequency, etc.

I II

E

E(fm) (-)

II.

U. .

0.00

m 2

E(fm,) (

(H)

0.05 0.10

0.5

fm (Hz)

0.20

0.25

0.30

FIELD DATA COLLECTED BT

SNORA BUOT *o1006 (1984-19851

-- PROPOSED THRESHOLD SPECTRUM

--- SPECTRAL EQUILIBRIUM RANGE

INDICATED BY P-H SPECTRUM

- AVERAGED SPECTRAL FORM

OBTRINEO FROM DATA SAMPLE

III l. soI..n

i

0.00

S

i

0

*I

0.05 0.10 0.15

fm (Hz)

0.20 0.25

0.30

Figure 3.2: Plots of E(fm) versus f,.

I

a FIELD DATA COLLECTED BT

O A BUOY *41002 11984-198!

S --- PROPOSED THRESHOLD SPECTRUM

----SPECTRAL EQUILIBRIUM RANGE

INDICATED TB P-H SPECTRUM

-- AVERAGED SPECTRAL FORM

OBTAINED FRH O DATA SAMPLE

*.

0

0

\I

*

a *

Jrto ~u^..

I

+

U m --

7

I -

5)

H

_- -- -

51

The nonlinear transfer coefficient, a,w, in Eq.(3.4) has been suggested by Has-

selmann (1976) to assume a constant value of 0.00138 under the condition of steady

wind within a duration At. However, the values of a,a evaluated from the present

NOAA buoy data appear to be varying instead of constant. Figure 3.3 shows aa

values from buoys #41002 and #41006 plotted against the non-dimensional param-

eter fo/(f,2At). It is clearly seen there exists a lower bound of the plotted data

which are approximated by the dashlines in the figure. Above the dashlines, most

data points are clustered in a narrow strip shown in the thickline box. Inside the

box, the separation of data points into two large clusters is caused by the shortage

in having the one-hour time interval, At, in computing aj in Fig. 3.3. Outside

the box, the data shows considerable scatter, which is due to the contamination of

swells moving into the region (data with dual spectral peaks). If a curve is fitted

to the upper bound of the clustered data and using the condition that the value of

a,, approaches the Hasselmann's constant of 0.00138 when fo/(f/At) -- 0 (steady

wind over long period), we arrive at

aa = 26[0.0195 + ( )3/5 ]5/2 (3.6)

This is drawn in Fig. 3.3 as a solidline curve which approaches asymptotically to the

dashline, the lower bound of the plotted data, as fo/(fAt) approaches infinity.

Equation (3.6) has been used in the new model as a unique parametric formula in

accounting the migration rate of spectral peak in the frequency domain.

The same data set excluding those exhibit the swell property is then employed

to determine the parametric formulas for the parameters a, , o, and ab. This is

accomplished by minimizing the square errors over the ensemble of measured data,

AE(f)/At, and the computed ones from the proposed Eq.(3.3). The values of a

are plotted against f,,/fo in Fig. 3.4. The data are seen more scattering when

f,r/fo is less than 0.5, which indicates that the generation of local wind waves

becomes important in this situation to effect the nonlinear transfer process. Since

100

10

at 0.1

0.01

0.001

0.0001

0.0001

100

10 -

0.01

0.001

0.0001

0.0001

0.001 0.01 0.1 1

fo

f2At

0.001 0.01 0.1 1

fo

f.At

Figure 3.3: Plots of an versus fo/f/At.

+ 2n>Hs

A 3n>Hsz2n

0 Hsk3m

+

-- PRESENT STUDY

--- HASSELHRNN ET AL. (1976)

........ DONELAN ET AL. (1985)

J. I

A

+ .+

*

++

A +

A + ..

0.005

0.002

0.001 1

0.001

1I I

0.01 0.1 1 1

Figure 3.4: Plot of a versus f,/fo.

0.05

0.02 -

0.01

54

the values of AE(f)/At from the data are contaminated by the advected wave

energy into the region, the best-fit curve of the lower bound values is believed to be

more representative of the growth rate. This curve can be approximated by

a = 0.0122 ( )/ (3.7)

or equivalently,

a = 0.0122 (U0 )7/6

Cm.

where C, = g/2rf,, is the wave celerity at f,,. The ratio fm/fo is named the

resonant parameter, which has been found important to estimate the nonlinear

transfers of energy in previously developed models,' for instance, it appeared in

Eqs.(2.36) and (2.37), respectively, in the parametric and WIS models. Also shown

in Fig. 3.4 are the empirical equations proposed for a by Hasselmann et al. (1976)

and Donelan et al.(1985). The equation proposed by Hasselmann et al. (1976) is

obtained by defining a constant wind speed during At in Eq.(2.36), which is derived

in the parametric model by allowing the variations of only a and f,.

a = 0.0091 (fM)2/

fo

The empirical equation for a proposed by Donelan et al.(1983) is

a = 0.006 (fM)o0.5

which is obtained based upon data collected in the Great Lakes. It should be noted

here that when the above two equations were compared with Eq.(3.7) in Fig. 3.4

they were multiplied by a factor (4.4)1/3 which is the constant absorbed in a in

Eq.(3.7). Both empirical equations proposed by Hasselmann et al. (1976) and

Donelan et al. (1985), as shown in Fig. 3.4, appeared to be less representative for

the lower bound values. The equation proposed by Hasselmann et al. (1976) is

more representative for the mean values. In the present study, Eq.(3.7) is adopted.

~

200 + 2n>Hs +

A 3n>Hs~2n ++4 *

100 0 Hs33n + +

+ ++

50 -A -4. +

+ + A +&++ +

++

20 + +

A + + +A

y, Eq.(3.8) + +4. +

10 + ++

S ^++ +

4 *!+ +

S A A + +

2 ,+ "++

1 I III I I I

1 2 5 10 20 50 100 200

-, (l.s.m.)

Figure 3.5: Plot of y computed from Eq.(3.8) versus y from the least square method.

The y values solved from the least square method can be approximated by the

following equation (see Fig. 3.5):

1= exp[2 i]/2.5( )3 (3.8)

Jo fo

The above equation gives a minimum '7 F 2.38 when f,,/fo =1.5. Although y is

not bounded for small f,,n/fo, the spectral peak value of Gin,, which is proportional

to a3y exp[2f,,/foJ(fm,/fo)1/2, is bounded for small f,,/fo.

Values a. and ab obtained from the field data are scattered (see Fig. 3.6). The

arithmetic mean values of oa and ab are 0.12 and 0.08, respectively, and they are

in the same order of magnitude as given by the JONSWAP spectrum. The values

56

0.6

,

A + 2w>Hs

0.5 A 3H>Hsa2M

HS3uN

0.q

0.4 -

+ +

S0.3 ++

4A .

S + + 0 +

+- +

0.2 + +

t +^. +" A -

+ AAH

AA. A:;,+-

# + *'t wAA+ a + -+ +

V4 &+

0.0 A I

0.0 0.5 1.0 1.5 2.0

fm

fo

0.6

+ 2Nw>Hs

0.5 A 3N >H;sa2m

+ 4 HSz3N

0.4i

Fo 0.3 +s

A +

0.2 A

+ +

0.1. + A + +

+ +

0 .0 +

0.0 0.5 1.0 1.5 2.0

f..

fo

Figure 3.6: Plots of ai versus fmfo.

I

57

of oa and as control the spectral width of the forward face and down-slope side of

the spectrum, respectively. Since at any instant, the energy contained in the wave

system, which is manifested by the magnitude and shape of the energy spectrum,

is the cumulative difference between the gain and the loss. The measured a, and ab

could also be viewed as parameters characterizing the net energy levels on each side

of the spectrum. Since on the down-slope side (high frequency components), the

loss of energy due to dissipation and transfer is expected to be very active during

the growth period, the measured ab values evidently are too small to represent the

gross gain. Based upon our estimate, a value of approximately four times as large

should be used, i.e., a, should be adjusted to 0.32. On the forward side, on the

other hand, the loss due to dissipation is minimal. Therefore, the arithmetic mean

value of measured oa should represent good estimate of energy gain due to energy

transfer and no adjustment is made.

The source function Gi,,1 defined in Eq.13.3) is also plotted in Fig. 3.1 for

comparison. The energy growth on the down-slope face of the spectrum is seen to

be much stronger than that in the WIS model, and the unnatural abrupt cutoff is

smoothened out.

In the calibrations of parameters a and 7, the past practice is to assume wind

direction has to be parallel to the mean wave direction. However, in the real situa-

tion, winds and waves may have different directions, and the wind-wave resonance

can be greatly reduced when the wind direction deviates significantly from the

wave propagation direction. Therefore, this difference between the wind and wave

directions needs to be properly incorporated into the resonance condition. This is

accomplished in the model by redefining fo = g/2trU1o cos( 09) in Eqs.(3.7) and

(3.8), which are used to estimate a and -y, respectively.

The rate of growth of newly generated waves in WIS model is treated in such a

way that the total energy growth during time step At is independent of the length

58

of At. This requires either the wind to be sufficiently strong or the duration to

be sufficiently long, so that in every time step the growth of local waves reaches

its upper limit. The WIS model was found to have a tendency of overpredicting

wave growth under mild wind condition (see Fig. 2.8) and this source term might

be partially responsible.

We observe here that the growth of local waves due to winds at the current time

step can take place in one of the following three forms:

a. If the shift in wind direction from the previous time step is sufficiently large,

say over a certain threshold value, e,, new fresh waves will be generated in the

direction of new wind.

b. If the shift in wind direction from the previous time step is small and the

wind direction is different from the prevailing wave direction, wind waves generated

from the previous time step will continue to grow in the direction of the wind.

c. If the shift in wind direction from the previous time step is small and the

wind direction is the same as the prevailing wave direction, this source term will

vanish and growth will be limited to prevailing waves only.

Accordingly, the source terms for the generation of new wind waves, Gin,2, are

proposed as follows:

Case a. AO, > 0, In this case, it is generally recognized that two wave generation

mechanisms are at work, the Phillips' (1957) resonance mechanism and the Miles'

(1957, 1959a,b, 1962) instability mechanism. It is also recognized that for practical

application, the resonance mechanism which is only important at the very initial

stage can be neglected. Therefore, it is common to assume that the instability

mechanism dominates, or Gin,s oc E(f,0). We further observe that the growth at

a location will eventually reach duration limited condition, the following functional

form should be satisfied:

Gin,2(f, 0) = b[E.,(f, ) E(f, 8)]

(3.9)

59

where b is a growth coefficient and Eo is the duration unlimited spectrum, and is,

in general, fetch dependent. For the present model of deepwater application b is

chosen as equal to 0.1sf(U1o/C)' with s the density ratio of air to water and C the

wave celerity at frequency f, and E,, is the P-M spectrum with the scale parameter,

op, equal to 0.015, i.e.,

f 0.1s/f( ( )[E(f) 2 cos'(B 0,) E(f, 0)], if 0 O,1< r/2

Gin,2(L,0) = C

0, elsewhere

(3.10)

The above equation is applied in the present model to initiate wave generation at

high frequency (fo > 0.2 Hz) and whenever the shift in wind direction exceeds 0c.

At present, there is no known criterion on selecting 0,. For convenience, we have

chosen

S= ( rad/hr) At (3.11)

Case b. AOn < OB In this case, the growth mechanism should be similar to that of

the prevailing wave components as discussed earlier with the exception that now the

growth should be centered around the wind direction, 0w, instead of the prevailing

wind direction, B. Thus, the source function is assumed to have similar form as

Eq.(3.3) with 0 being substituted by 0,, i.e.,

a exp 4 -8 cos4(O O,), if J1 0.,< /2

G,,(0) =(2)'/fmo 317

0, elsewhere

(3.12)

The spectral peak frequency f, and scale parameter a are determined by the same

equations used in the prevailing wave case, i.e., Eqs.(3.4). and (3.7), or

fk-1 = ((f)-7/3 + aU(U cos( 0B))'4/3At)-/7 (3.13)

g

and

a = 0.0122 (f )7/6 (3.14)

ao = 0.012--0

i

60

with fo = g/27rUlo cos(AOw) being defined in Eq.(3.14) to account for the influence

of A#, to the resonance condition.

Case c. 0, = In this case, the wind direction is the same as the prevailing wave

direction, we have Gin, (f, 0) = 0.

Particular interest should be noted to the source function given in Eq.(3.10). If

this is the only source term for spectral evolution at the duration-limited case, the

transport equation can be solved analytically to yield the spectrum

E(f, ) = E((f){1 exp[-0.1sf (.O)t] } cos2 (0- ) (3.15)

which vanishes at t = 0, and is equal to Eo(f) (2/7r) cos2(0 0,) at t = oo. A

simulation of the growth of wave spectrum based on Eq.(3.15) under 30-knot winds

is shown in Fig. 3.7. The result is compared to Barnett's (1968) model, also shown in

Fig. 3.7, and the comparison shows a good agreement between the two simulations

although the spectra predicted based on Eq.(3.15) have wider and flatter shape than

these predicted from Barnett's model.

If E.(f) is eliminated from Eq.(3.10) by using Eq.(3.15), the source function

Gin,2(f, 0) can be rewritten as

exp[-O.l s/f(Uo/C)^ t ] Ulo

Gin,,(f, 0) = 0.1sf E,,s(f,0) l(exp [_.ls/ o/C)4t]( (3.16)

1 expj-0.1sf(UJo/C)' t] C

This equation indicates that the growth rate is equal to the product of a spectrum

and a forcing function which is dependent on wind strength, wave frequency, and

time. The relation described in Eq.(3.16) is different from Barnett's formula de-

scribed in Eq.(2.48) which states that the growth rate is equal to the product of

a spectrum and a forcing function which is dependent only on wind strength and

wave frequency. Since growth rate proposed by Barnett (1968) can not predict a

steady state spectrum and it has an artificial cutoff of energy growth as the ratio

Ulo/C is less than 0.9, it is more appropriate to use Eq.(3.10) rather than Barnett's

formula in predicting the local generation of wind waves.

2

Ef) (j)

E(f) (M)

Hz

20 -

10 -

0

0.0

0 L.

0.0

0.1 f (Hz) 0.2

0.1 0.2

f (Hz)

Figure 3.7: Simulations of spectral growth at duration-limited case; Uo0 = 30 knots.

3.2.2 Dissipation Terms

The dissipation of wave energy caused by turbulence has been studied by a num-

ber of investigators. In deep water, the functional form suggested by Hasselmann

(1974) is represented by

Gd,.(f,) = -c, fE(f, ) (3.17)

where c, is an dissipation coefficient and its value is usually obtained by utilizing

field wave data obtained during successively decreasing winds with steady direction.

The value selected in the present model is s/6. Substituting this value into the above

equation and integrating at the duration-limited case, we obtain

E,(f, 0) = E,,(f, ) e- t/6 (3.18)

Figure 3.8 offers a comparison of wave decays between computed values and field

data. Owing to this calibration practice, the above dissipation function is found to

be more suitable for large weather systems such as related to high pressure cells, in

which the wind direction is relatively stable.

For low pressure systems, such as hurricane winds, the wind direction could be

completely opposite on the forward and the trailing phases of the eye while the

wind strengths are usually on the same order of magnitude. With such a weather

system the process of wave energy dissipation is much more vigorous as wind in the

trailing phase often blows against the prevailing waves and, thus, promotes wave

breaking. For low pressure systems, the dissipation function may be determined

based on that proposed by Lamb (1932):

GdA(f, 0) = -c, fE(f, 0) Ido E(f, ) (3.19)

pCIdT7|

where p is the normal pressure disturbance and r is the water surface fluctuation. By

assuming dp is of the same order of magnitude as the downwind flux of the horizontal

momentum of the wind, d(paU,2) (Phillips, 1967; Longuet-Higgins, Cartwright, and

30

20

m2

E(f) ( )

10 -

0

0.0

30

20

m2

Elf) ( )

0 -

0.0

0.1 0.2

f (Hz)

0.1 0.2

f (Hz)

Figure 3.8: Comparison of spectral decays between measured and computed values.

64

Smith, 1963), Eq.(3.19) can be approximated by the following expression:

sU72

Gdi,(f, 0) = -c, f E(f, ) 0.02 sU1E(f, 0) (3.20)

which is a modification of that in high pressure systems by adding an additional term

-0.02[sUjo/CH,]E(f, ). If this is the only source term in the transport equation,

the wave energy evolution at the duration-limited condition can be found as

E {,(f, ) = ELd, (f, ) e- /a-^t (3.21)

where A = 0.02sU2o/CH,. The characteristics of e-"t are such that it varies from 0

to 1 depending upon the maturity of the sea and that smaller waves dissipate more

energy against wind than that of larger waves.

3.3 UCWP Shallow Water Submodel

The shallow water submodel is based upon the numerical algorithm developed

by Chen and Wang (1983) which is a nonstationary spectral transformation model of

discrete type. The governing equations utilized in the model are given in Eqs.(2.45)

and (2.46). In the present application, a number of source terms are modified.

3.3.1 Wave Energy Growth in Shallow Water

In Chen and Wang's (1983) model, the generation term is again defined by

Eq.(2.47). In order to be consistent with the deepwater model, we propose to use

the same functional form here for local wave generations as for deepwater case:

S0.1sf( )^[E.(f cos( ) E(f, 0), if O 0,.< r/2

Gi, = C r (3.22)

0, elsewhere

where E,(f) is the duration unlimited shallow water spectrum. The spectral evo-

lution based upon Equation (3.22) at a duration-limited case can be expressed as

Ek+(f,9) = Ek,(f,O)expl-0.1sf ( )4At]

+ E,(f){1 exp[-0.1sf_ ( )1At ] }cos2(O 0) (3.23)

C 7T

65

This equation is actually used in the UCWP shallow water model, and it is used to

compute the wave components for Ulo/C > 0.9 such that it will not contaminate

swells.

Currently, our knowledge on E,(f) is limited. In the present model we adopted

the TMA spectral form proposed by Bouws et al.(1983):

1 f-fm

ag2 5 exp- -( )21

E,(f) = e [-- (f/f) 2 oa fm (k, h) (3.24)

(27r)4f5 4

where (k,h) is a transfer function defined earlier in Eq.(2.49):

(k, h) = kfoo) 0)

ks(f, h) C,(f, h)

The spectral parameter a is estimated in the model by the following formula:

a = 0.0156 (9 km)0.49 (3.25)

where km, is the corresponding wave number at peak spectral frequency, fi, which in

turn is related to the resonant frequency as f, = 0.88 fo. Equation (3.25) describes

a best-fit curve of the upper bound values of a computed, based on field data, by

Vincent (1984). The spectral parameters 7, and aab are also adopted from Vincent

(1984):

S= 2.47(U k)o.39 (3.26)

(0.07 if f < fm,

Oa = (3.27)

1 0.09 if f >fm.

However, since a significant amount of wave energy is continuously being dissipated

in the low frequency range under the influence of shallow water, the energy densities

in this range may not be sufficiently large to maintain the spectral peak near the

resonant frequency. Vincent (1984) suggested that the peak spectral frequency

should have a lower bound

S2 (E)/2 (3.28)

where A lies between 0.8 and 1.0.

66

3.3.2 Wave Energy Dissipation in Shallow Water

The dissipation of wave energy in the shallow water may be attributed to the

following effects: (i) flow turbulence, (ii) bottom friction, (iii) bottom percolation,

(iv) breaking due to depth limitation and (v) breaking due to energy saturation of

spectral components.

The dissipation mechanism due to flow turbulence is essentially the same as

in deep water, or Eqs.(3.17) and (3.18) still apply. The bottom friction effect is

estimated from Eqs.(2.50) and (2.51). The dissipation of wave energy due to perco-

lation of porous sand bottom is estimated from Eq.(2.53). However, this quantity

is found to be very small and can be neglected. The water depth limited condition

is given by Eq.(2.57).

The limiting of spectral growth due to energy saturation is based on Eq.(2.56),

or

A(f)= B g(2rf)-b5(k, h)

where B is chosen equal to 0.073 in this study.

3.4 Wind Model

Wind is an important factor for generating water surface waves. The wind

information required as input to a wave model is the time series of the overwater

wind vector fields. The wind information can be obtained either from the field

measurement or from a wind model. Since measurement of the winds over an entire

water surface region is almost impossible, numerical generation of the wind field is

generally the only means to prepare the wind information as input to a wave model.

The numerical model utilized to generate the wind information here follows the

method described by Clarke and Hess (1974). It generates wind field based upon

the input of sea surface pressure distribution. The procedures are summarized 'as

follows:

67

(i) The geostrophic gradient winds are produced directly from the pressure gradient

information by the following equation (Byers, 1974)

u2 1 ap

uf + = (3.29)

r p, an

where u=l|it is the magnitude of wind speed, r and n are, respectively, the radius of

curvature and the directional normal of the isobars, p. is air density, f= 20 sin 0 is

the local component of the planetary vorticity normal to the earth's surface, or the

Coriolis parameter with n the angular velocity of earth rotation and 4 the latitude,

and p is the atmospheric pressure. If the geostrophic gradient wind is to be solved

in the Cartesian coordinates, (z, y), 8p/an can be obtained from the expression

pp (Op p 211/

= (P)2 + (aP)2 1/2

Qn ax+ ay

and the radius of curvature, r, is obtained from

8ao apo a4o as a9oPy a9o asn

r = 1( ), with =+ = cos O, + sin

at a a t t ay at ax ay

where t is the distance along the path where air flows, and 0, is the downwind-

direction tangent to the isobar. It should be noted that in the Northern Hemisphere

a positive r indicates the cyclones and a negative r indicates the anticyclones. The

geostrophic gradient winds solved from Eq.(3.29) are then used to compute the

winds at a reference level above the sea surface.

(ii) The surface winds are then obtained by the corrections of frictional effects.

The wind speed, U,, and wind direction, 0,, are adjusted according to the following

equations:

U= In z (3.30)

U. 0.4 zo

and

0.4 u

w, O, = sin-'(Bo/ -- ) (3.31)

U.

68

where the frictional velocity, U,, and the frictional height, zo, are computed by

(Clarke and Hess, 1974)

(0 = (Bo)2 + [In( ) Ao]2 (3.32)

U.+ f

and Eq.2.42), or

0.1525 U2

Z + 0.0144- 0.00371, (in cgs units)

.* 9

where Ao and Bo are non-dimensional parameters, depending on horizontal temper-

ature gradient and vertical heat flux in the surface boundary layer. At adiabatic

condition Ao and Bo are contents and equal to 0.8 and 3.5, respectively. More in-

formation about the values Ao and Bo in the nonadiabatic condition can be found

in Resio et al.(1982), and Large and Pond (1982).

CHAPTER 4

MODEL PERFORMANCE

4.1 Introduction

The new UCWP model is to be tested here for hypothetical and real cases.

The hypothetical cases consider the wave generation at either fetch limitation or

duration limitation under stationary wind fields. Only the deepwater submodel

is tested for the hypothetical cases. The real cases consider the hindcastings of

deepwater waves in the Northwest Atlantic Ocean and shallow water waves along

the East Florida Coast. Both the deepwater and shallow water submodels are tested

in the real cases.

In Section 4.2 the UCWP deepwater submodel was tested against those of the

WIS model, the SMB model, and Hasselmann's (1976) parametric model. Such tests

were made for either fetch-limited or duration-limited condition at selected constant

wind speeds. In Section 4.3 the UCWP model was tested against the measured data

from NOAA buoy deployed in the Northwest Altantic Ocean and from CDN gages

installed along the East Florida Coast. In Section 4.4 the relationship between wind

and wave for a fully aroused sea is studied. In Section 4.5 the relative importance of

each factor affecting deep and shallow water spectral transformations is analyzed.

Finally, in Section 4.6 various values of bottom friction coefficient were tested in

the UCWP shallow water model to examine their effect on the wave transformation

process.

70

4.2 Model Performance in Hypothetical Cases

In this section, the UCWP model is compared with three standard wave hind-

casting models, the WIS model, the SMB model, and the Hasselmann's (1976)

parametric model, for fetch-limited and duration-limited deepwater wave genera-

tions. The deepwater WIS model and Hasselmann's parametric model have been

reviewed in Chapter 2. The SMB model is known as a semiempirical method for

deepwater wave prediction; it was first developed by Sverdrup and Munk at 1947

and then revised by Bretschneider at 1952 and 1958 with additional empirical data.

The wave properties used for comparisons are the significant wave height and modal

period.

The UCWP deepwater submodel was tested first in fetch-limited cases with

constant wind speeds equal to 20, 30, 40, and 50 knots, respectively. The results of

significant wave height H, and modal period T, are plotted against fetch at 20, 50,

100, 150, 200, 250, 300, 350, 400, 450, 500, and 550 kilometers in Fig. 4.1.

According to the SMB model the significant wave height, H,, and modal period,

T, along a fetch F are estimated from the following formulas, respectively:

H.= 0.283 tanh[0.0125(i-)3/ ] (4.1)

T 7.54 tanh 0.077(gF)1/4 (4.2)

10o U120

where Ua, the adjusted wind speed, is defined as a power function of U10 (Shore

Protection Manual, 1984):

U. = 0.71 U (in meter-second units) (4.3)

According to Hasselmann's parametric model the significant waves and modal

periods, along a fetch F, are estimated from the following equations:

=g H 0.0016 (g F)/ (4.4)

U'2 U.'2

8-

C -

UJ

g

DC

4J

t U

0

UJ

6-

g -

o L

U)

tr-

"

72

g 0.286 (F (4.5)

Uo0 U120

In the classical SMB and parametric models, the adjusted wind speed U4 is

equal to the wind speed U10. However, it is found later that a correction to Ua is

necessary since both the classical SMB model and the Hasselmann's model tend to

underestimate the wave heights. The utilization of U4 according to Eq.(4.3) in the

both models is suggested by the Corps of Engineers, U.S. Army, and U4 obtained

from the equation will be greater than U10 if U10 is greater than 4 meters per second.

Therefore, by using the adjusted wind speed U4 as defined in Eq.(4.3), the SMB

and parametric models will predict higher wave heights than the classical ones if

U0o is greater than 4 meters per second.

The significant wave heights and modal periods obtained from the WIS model,

the SMB model, and Hasselmann's parametric model are compared to the UCWP

model in Figure 4.1. At low and moderate wind speeds (20, 30 and 40 knots)

and large fetch (greater than 400 km), the significant wave heights obtained from

the UCWP model are quite consistent with these obtained from the WIS model.

At high wind speeds (40 and 50 knots) and small fetch (less than 400 km), the

significant wave heights obtained from the UCWP model are the highest among the

four models, while the significant heights obtained from the WIS model and the

parametric model agree with each others. As fetch increases, the significant wave

heights obtained from the UCWP model, WIS model, and SMB model tend to reach

to a height limit due to energy saturation. The parametric model, however, is not

bounded.

The modal periods obtained from the UCWP model are consistent with those

obtained from the SMB model at low and moderate wind speeds (20 and 30 knots).

At high wind speeds (40 and 50 knots), the modal periods obtained from the UCWP

model are much smaller than the other three models. For high wind speeds (40

and 50 knots) and small fetch (less than 400 km), the modal periods computed

73

from the WIS model and the parametric model agree with each other. At large

fetch, the modal periods obtained from the UCWP model, the WIS model, and the

SMB model again tend to reach an upper limit due to saturated conditions. The

parametric model does not have an upper bound on modal period.

The UCWP deepwater submodel was compared next in the duration-limited

cases with the WIS model, the SMB model, and Hasselmann's model at the constant

wind speed of 20, 30, 40, and 50 knots, respectively. The comparisons are shown in

Figure 4.2.

The significant wave height and modal period, respectively, for a duration t are

estimated in the SMB model from the following formulas:

g H.

U = 0.283 tanh[ 0.00087( )5/8] (4.6)

Tm= 7.54 tanh[0.016( -)3/81] (4.7)

Uo0 U10

where the adjusted wind speed, U., is defined in Eq.(4.3).

The significant waves and modal periods for a duration t, are estimated in the

Hasselmann's parametric model from the following equations:

g H = 0.000082 ()5/7 (4.8)

g (gt)s/7 (4.9)

"T = 0.06 (ZF) (4.9)

It is seen in Figure 4.2 that at short duration (less than 20 hrs) consistent sig-

nificant wave heights are found between the UCWP and SMB models, and between

the WIS and Hasselmann's models. At high wind speeds (40 and 50 knots) and

short duration (less than 20 hrs), the significant wave heights obtained from the

UCWP and SMB models are found to be much greater than the WIS and para-

metric models. As a matter of fact, the significant wave heights obtained from the

UCWP model are the highest among the four models at short duration. At long

duration (greater than 20 hrs), both the SMB and parametric models predict much

S4J

- n

a a a

=

a

I-

az

"

c

aC

I

ao

75

higher waves than the UCWP and WIS models. At long duration, wave heights

generated from the UCWP and SMB models have upper limits. The significant

heights generated from the UCWP model show 'overshoot' patterns in Fig. 4.2. In

other words, waves become over saturated before settle to a stable saturation con-

dition. The largest 'overshoot' 0.2 m associated at 30 knots winds; it is still very

small (only 3 %) compared with a 6.6 m height at saturation. The 'overshoot' of

wave height is overall small but it is essential to gain enough energy for establishing

the saturation condition in the UCWP model. The 'overshoot' does not occur in

the WIS model nor in the SMB and parametric models. None of the four models

shows overshoot behavior in the fetch-limited cases.

The significant wave heights predicted from the WIS model at 40 knots wind

show a peculiar periodic, irregular damping pattern when duration becomes longer

than 40 hrs. This damped wave height which is unlikely to occurr in the real

situation ranges from 10.6 m to 13.4 m and repeats its pattern every 12 hrs. It

appears that it is caused by an improper choice of saturation spectrum E,(f), or

Eq.(2.32), which function is to limit the growth of waves in the WIS model.

The modal periods obtained from the four models do not differ appreciably to

each others for small duration (less than 10 hrs). For moderate duration (between

10 hrs and 40 hrs), the modal periods obtained from the WIS model agree with the

SMB model at low wind speed (20 and 30 knots) but become closer to that predicted

by the parametric model at high wind speed (40 and 50 knots). For long duration

(greater than 40 hrs), consistent modal periods are found between the WIS and

SMB models. At long duration, the modal periods obtained from the parametric

model are the largest among the four models. At long duration and high winds (40

and 50 knots), the modal periods obtained from the UCWP model are found to be

much less than the other three models.

/ l 0N

SsDEEP WATER NAVE /

INFORMATION STRTION 1

^ /

Q iSNALLON WATER A VE

__ _1002 320N

KINGS B iT

-5 STbT ON\

STRTIIN l

1006 300N

~,, ATLANTIC

-- j OCERN

CAPE ANAVERRL /

STATI N 5

W 280N

FLORDR D 26N

t 260N

82 W 800W 780W 76OW 74 0

Figure 4.3: Shallow water grid system.

4.3 Model Performance in Real Cases

The UCWP model was first tested against the bench mark case of the period

December 21, 1983, to January, 1984 to offer comparisons with the WIS model

and field data. For the deepwater model, the test conditions including grid size,

time step, etc., are identical to that described in the WIS model testing as shown

in Section 2.7. For the shallow water model the grid system is given in Fig. 4.3.

The grid system covers an area 300 kilometers long from mid Georgia to the

northern edge of Bahama Islands, and about 100 kilometers wide from shoreline

to approximately 600 m contour line offshore.

The wind information was generated by the wind model with daily synoptic

77

surface pressure as input. Since the surface chart is given daily at 7:00 a.m., E.S.T.

only, linear extrapolations are used for winds at other times. The computed wind

directions and speeds for January 1984 at grid point 30.00N and 77.5"W are shown

in Fig. 4.4 together with the measured winds from NOAA buoy #41006, located

at 29.3N and 77.3W. The correlation, 6, given in Fig. 4.4 has been defined in

Eq.(2.58). The root mean square difference in wind direction, ORMSE, given in

Fig. 4.4 is defined as follows:

N

ORMss = { [(e C), (0,,m)i}) 2/N (4.10)

i=i

where 8,,e is the computed wind direction, 8,., is the measured wind direction, N

is the sample size, and the overhead bar indicates the operation of arithmetic mean.

The wind strength exhibits a periodical behavior with a periodicity of roughly

equal to 7 days. The wind direction associated with high winds are found to be

more stable and to have better agreement with the predicted values. The small

winds are less stable, thus, less predictable.

The measured wind speeds are seen in Fig. 4.4 to be slightly scattered around

the computed ones within 10 knots. This is because the measured winds are

reported in an every hour schedule while the computed winds are from the wind

model at an every 24 hrs interval. The fluctuations of either magnitude or direction

of hourly recorded winds, as seen in Fig. 4.4, may be caused by abrupt changes of

temperature between the air and sea water, and between the stratified air layers

above the sea surface. For simplicity the temperature effect is not include in the

wind model which assumes a neutral condition.

The wind information used as input to the wave model are the computed winds

adjusted by the measured values at the buoy. That is, the wind speed at any grid

point is the computed value multiplied by a factor which is the ratio of measured

and computed ones. However, this factor is limited in the range between 0.5 and

1.5 to limit the influence of abnormal data values.

1 5 10 15 20 25 30

JANUART (1984)

Uto

(KNOT)

Figure 4.4:

location.

1 5 10 15 20 25 30

JANUARY (1984)

Comparison of hindcasted and measured winds at the Buoy #41006

79

Figure 4.5 shows the comparison of the time series of the wave energy spectral

evolution process of the entire month. Figure 4.6 shows the comparisons between

the measured and computed modal periods and significant wave heights at the buoy

location. Figure 4.7 shows the comparison between the measured and computed

energy scale parameters at the buoy location. The model performs well for this

data set. The good agreement, however, should be viewed as partially artificial as

the model was calibrated on this data set. The model outputs in the ensuing cases

are data independent.

In the shallow water model, only one empirical value, the bottom friction coef-

ficient, needs to be defined as input. Shemdin et al.(1978) suggested values ranging

from 0.006 to 0.05 along the east Florida coast. A medium value 0.02 was chosen

for the model. Directional wave spectra from deepwater model are used as input

boundary conditions. The wind field is spatially linearly extrapolated from the ad-

justed output from the wind model as discussed earlier. All the input conditions are

upgraded at 2 hrs interval but the time step in the shallow water model is chosen

as 24 minutes to insure numerical stability. The dissipation of wave energy due to

percolation is not considered in the shallow water model since it was found to be a

minor effect to the wave transformation process in the water region near the East

Florida Coast (Shemdin et al., 1978).

For the case of January, 1984, the hindcast results from the shallow water model

are compared with the directional wave information collected at the offshore Kings

Bay Station (St.#5 in Fig. 2.4). Figure 4.8 shows the comparison of the time series

of the one-dimensional spectra. By comparing Fig. 4.8 with Fig. 4.5, the measured

deepwater and shallow water wave spectra, we noticed the significant reduction in

wave energy and the reduction is particularly severe in low frequency components.

This point will be further discussed later. The time series of the marginal directional

spectra, which are obtained by integrating the two-dimensional spectra over the

2

E(f) (f),

(NOAA)

m

E(f) (M),

(UCWP)

U.Z

0 JAN. 11984)

S(Hz0.0

f'(Hz)

Figure 4.5: Comparison of NOAA and UCWP wave spectra for January 1984 at

the Buoy #41006 location.

60

50

40 -

30 *

20 -

10 -

JAN. (1984)

I

f (Hz)

81

20

STATION NOAR BUOY q*1006 HERSURED ORTR(NORA)

15 6 =0.38 COMPUTED ORTR(UCHP)

T 10 **

(SEC)

.1 5 10 15 20 25 30

JRNUART (1981)

10

STATION: NORR BUOY *41006 ERSURED ORATA(NOAR)

8 6 = 0.85 COMPUTER ORTA(UCHP)

H3 6

(MI)

2

1 5 10 15 20 25 30

JANUARY (1984)

Figure 4.6: Comparisons of NOAA and UCWP modal periods and significant wave

heights at the Buoy #41006 location.

0.12

0.10 F

0.08 -

10* ,

gCW

(UCWP)

0.06

0.04 I

0.0

0.I

+ 2n>Hs 6 = 0.4l

A 3N>HsH2n 6 = 0.43

S Hs33n 6 = 0.67

A +

*A +++

--1-

+ -9. A +

+ +

4A A A. A

A

A sp O A+ 4

A AA 0A AA

I + + A

+ + +.. ++ A

0

0.02

0.04

0.06

0.08

10 (NOAA)

g

Figure 4.7: Comparison of NOAA and UCWP energy scale parameters at the Buoy

#41006 location.

0.02

2

E(f) (),

(CDN)

E(f) (M),

(UCWP)

Ak

JAN. (198U)

0.2.

0.2

f (Hz)

vt0 JAN. (1984)

7

f (Hz)

Figure 4.8: Comparison of CDN and UCWP wave spectra for January 1984 at the

Kings Bay gage location.

!0

r%70

84

frequency domain, are shown and compared in Fig. 4.9. Finally, the measured and

computed mean wave directions, modal periods, and significant wave heights are

summarized in Fig. 4.10.

In the case illustrated above, the winds were predominantly induced by high

pressure systems which are large scale and relatively stationary. In the following

case, wave hindcasting was performed in September and October, 1984 to assess

the model performance for low pressure system and other coastal locations.

In September and October, 1984, three hurricanes occurred in the Mid Atlantic

Ocean (see Fig. 4.11). The three hurricanes are Hurricane Diana, Isidore, and

Josephine. The cyclones found during these hurricanes are good examples of low

pressure system found in September and October; they are small scale and relatively

nonstationary. Synoptic charts showing low pressure systems of Hurricane Diana,

Isidore, and Josephine are given in Figs. 4.12(a), 4.12(b), and 4.12(c), respectively.

All three cylcones generate strong northeast winds and, hence, large waves heading

to the southeastern coast of the United States.

Hurricane Diana came near the Florida east coast in the morning of Septem-

ber 8, 1984. After looped off the northeast Florida coast for almost 36 hours, it

moved north toward the coast of North Carolina. The highest wave recorded at

the Marineland CDN station during the visit of Hurricane Diana reached 3 m in

terms of the significant wave height. Hurricane Isidore came ashore near Miami

Beach, Florida, during the morning of September 27, 1984. It then moved north-

ward across central Florida and returned to the Atlantic Ocean near the border

between Florida and Georgia. The highest wave caused by the Hurricane winds

was recorded at the CDN wave gauge station at Cape Canaveral and reached 3.2

m in terms of the significant waves. Hurricane Josephine was developed in the Mid

Atlantic Ocean on October 10, 1984. The path of the hurricane is nearly parallel

to the east coast of North America. The closest distance between the hurricane's

0.3r

E(O) ( ),

rad

(CDN)

m2

E(O) ( 2),

(UCWP)

(UCWP)

0.2h-

5

0.1

0. 0

E

15 JAN. 18

10 JAN. (198U)

0.3r

25

15

JAN. (198L)

SNE

Figure 4.9: Comparison of CDN and UCWP marginal directional spectra at the

Kings Bay gage location.

%0

S WN

1,

1 5

TN

(SEC)

15

JANUARY (1984)

1 5 10 15 20 25 30

JANUARY (1984)

1 5 10 15 20 25 30

JANUARY (1984)

Figure 4.10: Comparisons of CDN and UCWP average wave direction, modal peri-

ods and significant wave heights at the Kings Bay station.

KINGS BRY(FL.) STATION

.ers= 16.7*

* MEASURED ORTA(CON)

- COMPUTED DATA (UCHP)

111111111111t 11111 11II i4i II

500N

45"N

4 50N

'400N

350N

300N

25N

80 W 75W 70 W 65"W 60W 550W 500W

Figure 4.11: Tracks of Hurricanes Diana, Isidore, and Josephine (1984).

Figure 4.12: (a) Weather map showing northeast winds generated

Diana (9/8/84).

S OON

50'N

l5'N

350N

24

300N

A 25N

50 N

by the Hurricane

500N

'4 5N

350N

300N

S1SIDORE

SATLAN IC OC A 250N

12

800W 75'W 70'W 65'W 60'W 55W 50aW

Figure 4.12: (b) Weather map showing northeast winds generated by the Hurricane

Isidore (9/27/84).

40ON

350N

2PHIN 300 N

u100N

oATLAN OC A 25'N

80OW 75'W 70'W 650W 60'W 550W 50W

Figure 4.12: (c) Weather map showing northeast winds generated by the Hurricane

Josephine (10/11/84).