A COUPLED DISCRETE SPECTRAL WAVE HINDCAST MODEL
By
LIHWA 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 Lichu 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 Work ................... ........... 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 . .
3.3.2 Wave Energy Dissipation in Shallow Water .
3.4 W ind M odel ................. .......
4 MODEL PERFORMANCE ..................
4.1 Introduction ........................
4.2 Model Performance in Hypothetical Cases . . . .
4.3 Model Performance in Real Cases . . . . . .
4.4 WindWave Relationship in a Fully Developed Sea .
4.5 Comparison of Various Wave Energy Flux Elements
4.6 Influence of cf in Shallow Water Wave Transformation
5 CONCLUSIONS ........................
5.1 Sum m ary .................... .....
5.2 Conclusions . . .. . .. . . . . . . .
5.3 Future Studies
. . . . 64
. . . . 66
. . . . 66
. . . . 69
. . . . 69
. . . . 70
. . . . 76
. . . . 96
. . . . 102
Process . 105
. . . . 115
. . . . 115
. . . . 119
. . . . . . . . . . . . . ... 122
APPENDICES
A SUMMARY OF PARAMETERIZED PROGNOSTIC EQUATIONS ... 123
A.1 Prognostic Equations ................... ....... 123
A.2 Elements of Dii and Gi ......................... 123
B PROCEDURE IN SOLVING PARAMETERIZED Gin,i(f) ........ .125
C SUMMARY OF COMPUTED AND MEASURED NONDIRECTIONAL
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 PM spectrum,
the Kitaigorodskii spectrum, the JONSWAP spectrum and the
modified JONSWAP spectrum ................... 13
2.2 Plots of Hs(f,0) versus 0 Oo. .................... 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) fetchlimited case and at (b) durationlimited 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
model... .............................. 47
3.2 Plots of E(f,) versus f,,........................ 50
3.3 Plots of al versus fo/fAt. ................... 52
3.4 Plot of a versus fm/fo. ..................... 53
3.5 Plot of y computed from Eq.(3.8) versus y from the least square
method. . . . . . . . . . . . . . . .. 55
3.6 Plots of ab versus fm/fo . . ................. 56
3.7 Simulations of spectral growth at durationlimited 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 fetchlimited cases. 71
4.2 Comparison of Deepwater hindcast models in durationlimited
cases. .................... ............ 74
4.3 Shallow water grid system. . . . . . . . ..... . 76
4.4 Comparison of hindcasted and measured winds at the Buoy
#41006 location ......................... 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 Uo1, and (b) Tm versus U10 at various fetch
and duration . . . . . . . . . . . . . . . 98
4.16 Plots of 27rU0o/gT, versus 10OH,/gT in (a) fetchlimited case
and in (b) durationlimited case. . . . . . . . . . 101
4.17 Plots of 27rUo/gTm, versus 10OH,/gT, based on the NOAA buoy
#41006 19841985 data ........................103
4.18 Comparisons of wave heights and modal periods computed in
September 1984 for cf=0.0, 0.005, 0.01, 0.02, and 0.04 .. . 107
4.19 (a) Plot of individual shallow water effect sources for cf=0.0 . 109
4.20 (b) Plot of individual shallow water effect sources for c =0.005 110
4.21 (c) Plot of individual shallow water effect sources for cf=0.01 111
4.22 (d) Plot of individual shallow water effect sources for cf=0.02 112
4.23 (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 T, at the buoy #41006 location. . . . ... 92
4.2 Summary of Correlations between computed and measured val
ues of H, and Tm 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 6bH,cH, 6Tm,c Tm,m 1EH*,cH,, and ET,,cT,,,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 cf=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
LIHWA 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 10meter
water depth along the eastern seaboard of the United States. The model comprises
of deepwater and shallow water submodels. The significant nonlinear wavewave
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
airsea interactions are also examined in light of the field data. New semiempirical
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 longterm data, wave hindcast, which is a tech
nique in calculating the events of windgenerated 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 wellknown example is the SverdrupMunkBretschneider (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 II 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 windwave 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 bellshaped, 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 timeconsuming 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 wavewave 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 wavewave interactions in the parametric model, while
5
it is assumed to depend on the nonlinear wavewave interactions in the discrete
model. A wellknown 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 wavewave 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 wavewave 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 wavewave 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 wavewave 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 fetchlimited
and durationlimited cases.
(8) Model outputs are compared with 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 stateoftheart 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 wavewave 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 wavewave 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, s, which is taken here as identical
to the still water surface plane (Thomson, 1972).
77(i, t) = Re[i da(k)ei(r1)]= =Jf !(da()ei(i) + da(k)ei(n )] (2.1)
where r7 is water surface elevation, t is time, 9i = (xi, yi) is the horizontal location
of interest in rectangular coordinates, a is radian frequency, k = (k.cos 0, ksin 0) is
a vector wave number with magnitude k =Ik I and direction 0, da(k) indicates an
infinitesimal wave amplitude (complex number) which is the contribution of waves
in a wave number interval dk, 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(k) 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 wellknown linear wave theory dispersion equation
a2 = gk.tanh(kh), k = kl (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
dr7? (, t) = [da(k)e'( at') + da'(k)eei(k a)] (2.3)
2
Thus, the mean square value of d7r can be evaluated, which has a form
1 1 T/2
d7(,t)]2 = lim {[da(k)] 2'I.t) + 2da(k)da'(k)
4 T.oo T _2
1
+[da'(k)]2ei2(k,1t)}dt = da(k)da'(k)
2
1 1
= Ida(k) 2= da(a, 8) 2= E(a, 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(a, 0) is the directional frequency spectrum which is defined as
1
E(a, 8) = [di7(s,t)]2/dadO = Ida(a, 0) 12/dad (2.5)
2
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 onedimensional frequency spectrum E(a) is
defined by integrating E(a, 0) over all 0.
E(a) = E(a, 0)dO (2.6)
The marginal directional spectrum is defined as
E(0) = E(o, 0)da (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,
Etotai = E(o)da = E(0)d0 = f j E(, ,0)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) PiersonMoskowitz (PM) spectrum, (ii) Kitaigorodskii
spectrum, and (iii) JONSWAP spectrum. The PiersonMoskowitz spectrum (1964),
EpM(f), is given by
EPM(f) = (2r)f5 exp[0.74(fo/f)4] = (2) exp[ (f./ f)4] (2.9)
where f = a/27r is the wave frequency with units 1/time, cp is known as Phillips'
constant, fo = g/2rU19.5 is the frequency of the deepwater waves whose phase
speed is identical to the wind speed U19.5 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
Sag2f,5(2r)4 exp[1 (fm/f)4], if f < fm
EK(f) = (2.10)
ag2f5(2r)4, if f fm
which contains two free parameters, the energy scale factor a and the spectral
peak frequency fm. The JONSWAP spectrum is derived from the Joint North Sea
Wave Project (Hasselmann et al. 1973) and constitutes a modification to the PM
spectrum to provide for a much more sharply peaked spectrum. The JONSWAP
spectrum, Ej(f), is given by
1 ffm
a 92 5 exp[ ( f 2
Ej(f) = ( ep[(fm/f)4 2 ob (2.11)
with
{C f fa, if f < fm
S b, if f > fm
The function contains five free parameters, two scale parameters a and fm, and
three shape parameters 1, a, and ab. The parameter 7 is the peakenhancement
factor which is the ratio of the peak value of the spectrum to the peak value of the
corresponding PM spectrum with the same values of a and fm. The parameters
ao and ab 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 PM spectrum, of the PM spectrum, the Kitaigorodskii
spectrum, and the mean JONSWAP spectrum, which has its spectral parameters
f, a,, 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 PM 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
_ = Ej(fm) 
EpM(f.)
 PM SPECTRUM
 JONSWAP SPECTRUM
 KITAIGOROOSKII SPECTRUM
 MODIFIED JONSWNP SPECTRUM
y= 3.3
a, = 0.07
ab = 0.09
f /fm
Figure 2.1: Comparison of nondimensionalized shapes of the PM spectrum, the
Kitaigorodskii spectrum, the JONSWAP spectrum and the modified JONSWAP
spectrum.
E(f)
EpM(fm)
2H
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' frequency dependence at the tail part of the wave spectrum
may be more adequate for the windgenerated waves than the fs spectral shape.
Based upon the assumption of the f4 dependence on the rear face of frequency
spectrum, Donelan et al.(1985) suggested a modified JONSWAP spectrum
1 I fm
g2_ _( exp[( f f ]
EJM(f) = 2exp(fm/f)4] 7 2 ab fm (2.12)
(2 7r) f4 fm
with all its free parameters a, f,, y, and ab 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, ,
and ab were established, based on field data, by Donelan et al.(1985) as functions of
fo/fm. A modified JONSWAP spectrum with its parameters and ab 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
PM 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, 8)
Hs(f,) = (2.13)
E(f)
where E(f,0) = Hs(f, )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,0) have been suggested but all with
similar shape, i.e., they are bellshaped and symmetric about a central direction. It
1.5
S
20
15
1.0
Hs(f, 0) 7
0.5 3
2
0.0
900 450 00 1150 90
0 80
Figure 2.2: Plots of Hs(f,0) versus 0 00.
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).
1) os (0 o0), if l0 Oo < 7r/2
Hs (f, )= 7 (2S + 1) (2.14)
0, elsewhere
where S = S(f) is called the spreading parameter, which is always positive and a
function of frequency, 80 = 0o(f) is the central angle of the symmetric spreading
directions. Figure 2.2 shows several distributions of Hs(f,0) 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, narrowband 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,9, t) E(f, 0, ,t) 9E(f,0,0 ,t)
at= + C,(f, s) = G(f, 8, s, t) (2.15)
dt at ae
where E(f, ,s,t) is the energy density of the spectral element with reference to
still water surface plane at horizontal location and time t, with frequency f and
direction of propagation 0, e is the path of this spectral element in which dt =
Idele'I, Cg(f, s) = a/8k = a(27rf)/akIj 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, s) can be derived from the disper
sion relation given in Eq.(2.2) and expressed as
1 kh
C,(f, ) = n C =( + k )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 s, 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, 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) > 7r (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
a2 = (27 f)2 = gk
and the deepwater wave group velocity, Cgo(f), can be expressed as
C,(f) = g g
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; Giinther 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:
Ej(f)2 cos2( o0), if 10 Oo\< 7r/2
Ej(f,O) = Ej(f)Hi(f,O) = (2.18)
0, elsewhere
where Hi(f,O), the spreading function, is a special case of Eq.(2.14) with S equal
to 1, and 00 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 directionindependent transport equation
is obtained as
dEj (f) Ej(f) 8Ej (f)
dEj) t) + Co(f) Y = G(f) (2.19)
dt at ae
with the directional averaged group velocity Cgo(f) equal to
.f) ,fCgo(f)Ej(f,O)d 8 g
(f) 37 47(220)
and the directional integrated source function G(f) equal to
G(f)= JG(f, )dO (2.21)
Note that symbols s 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, Gdis(f), and the
nonlinear energy transfer, Gia(f). That is,
G(f) = Gi,(f) + Gdis(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 Gf(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 fetchlimited condi
tion in the case dominated by atmospheric input is of the same order as the case
dominated by the wavewave interactions.
The growth rate Gi,(f) is commonly approximated in terms of a linear func
tion of onedimensional frequency spectrum (Snyder and Cox, 1966). However, for
convenience the quantity G;,(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 onedimensional
frequency spectrum and has the following expression:
G,,(f) = 5s f E(f) (2.23)
where s 0.0012 is the ratio of air density to sea water density, and U10 is the wind
speed at 10 m elevation above sea surface.
The general form of the nonlinear transfer source GI(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
2 3
Ga(f) = D,  (f//fm) (2.24)
where D1 is a proportional constant whose magnitude is on the order of 103,
and O(f/fm) is a nondimensional 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 fs law, Eq.(2.24) can be approximated by
the expression
S2.2a2f ,Ej(f), if f < fm
G.Y(f) = (2.25)
1 2.2a 2fE(f), if f > f,
With the source function given, the directionindependent transport equation,
Eq.(2.19), can be then projected onto the JONSWAP parameter space to yield
0.2
f. = 0.3 (Hz) 3
a = 0.0095
y = 3.3
a, = 0.07
Ub = 0.09 2
E(f) 0. 105 Gn
m
(MH (m2)
1
0.0 /G 0
1
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 9bi
t + Ditj = Gi, ij = 1,2,3,4,5 (2.26)
where b1 = f,, b2 = a, b3 = y, b4 = a,, and b5 = ab. 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 Dii matrix and Gi 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 HI (f, 0) spreading function to generate a
directional distribution of wave energy densities.
The nature of the parametric approach is that only singlepeaked 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+(f,0) = Ekf) Cgo(f) AEk At + G(f, ) At (2.27)
Ae
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,0) = Ek(f,0) Cgo(f) (f At
The source term is further divided into three components:
The source term is further divided into three components:
G(f,0) At = [Gk+'(f,0) + G ,(f,0) + GkC(f,0)] At
(2.28)
23
where Gi,, Gad,, and Ga represent atmospheric energy input, wave energy dissi
pation and nonlinear wavewave interaction, respectively. The Gi, term, which is
computed at the current time level, is further split into two parts: (i) Gi.n,(f, O),
the rate of growth of the prevailing waves and (ii) G,,2 (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:
23 8
4.4 Atf, cos(0 ), if f fm
Gi,nAt = (27r)4f 37r
O, elsewhere
(2.29)
with the mean propagation direction, 0, defined as
= tan j f sin OE(f, 0) d/fd
0 = tan1{ } (2.30)
ff cos0E(f, 0) dfd0
where the scale parameters a and fm are to be determined later.
ag2 fm) 2
(274 exp[ ( ) cos2(0 0B), if f, 0.25 (Hz), f < fm,
(2) fm f 7r
fm > 0.2 (Hz)
Gi. ,2/t = ag2 2
(Gi),2A < g 2 cos(0 0), if f, > 0.25 (Hz), f > fm,
(27r) 'f I 7
fm > 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 onedimensional saturation spectrum defined as
ag2
ag if f < fo = g/27rUjo
(27r)4f5
E,(f) =
ag2
(27r)5 [1 (fo/fm)4] + EPM(f)(fo/fm)4, elsewhere
(2.32)
where fo is the windwave resonant frequency to which the corresponding wave
phase speed is identical to the wind speed at 10 m elevation, U10, and EpM denotes
the PM spectrum.
There are two scale parameters a and fm 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 fm component prognostic equations can
be derived as
a e a fa Di fmUo_
S+ 0.47 C(f) + 0.2 Cgo(f) = 0.005 ( ) afm 5 a3 f (2.33)
at ae fm ae g
and
afm afm fm a2
+ C7o(f) 0.07Co(f) a = 0.54 a2fm (2.34)
at ae a e
If the changing wind speed can be expressed as a power law of time,
gt
U10 = Uo( )
Uo
with Uo a reference wind speed, fm and a can be solved from Eqs.(2.33) and (2.34)
as
fmUo = A( t)(3/7)(q) (2.35)
g Uo
and
a = B( mUO)2/3 (2.36)
25
with constants A and B given by (Hasselmann et al., 1976)
A = 16.8 (1 + 1.51q)3/7, B = 0.031 (1 + 1.33 q1/2
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
UlkO+1 kf1 = A( gt + )(s/7)(q1) (2.37)
U~flg =A(got + 0 (2.37)
g Uo Uo
Eliminating g t/Uo from Eqs.(2.35) and (2.37) yields
Sfk+l = {(fk)(7/3)(q1) + ( (7/3)(q1)3/7)q1)
U U o Ag (2.38)
For constant wind speed within the duration At, the above equation reduces to
fk+1 = {(f_)7/3 + a, (1)4/3At}3/7 (2.39)
where a, 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 (gE ) 0.2 (2.40)
where U, is the surface frictional velocity obtained by assuming logarithmic over
water wind profile
Uz 1 z
= In (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)
U. g
26
Although Eq.(2.36) is not used in the WIS model, it can be linked to an
empirical relationship (g2Etota /U4o)(fUlo/g)4 /a=constant suggested by Hassel
mann et al.(1976) to yield a power law a ~~ (g2Etota,/Uo)02 which is similar to
Eq.(2.40).
The transfer of wave energy due to a nonlinear wavewave interaction source
Gj (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).
g23 8
0.0023 At i(f/fm) cos4(0 ), if 10 \ x7r/2
GonAt = f m 73 (2.43)
0, elsewhere
where
{ exp[1 (f,/f)4], if f < fm
W//m) =
(fm/f)s, 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 f5 law on the downslope side. The loss of wave
energy in swell is adopted from Resio (1981), in terms of onedimensional frequency
spectrum
E2(f" )
4 At f(27r)g E(f), if f,, < f < 0.7 f,
G8t = (27r)8' (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 airsea interaction.
In shallow water the wave energy transport equation can be expressed as follows:
(Chen and Wang, 1983)
aE(f, ) aE(f, )C,(f,0) cos0 aE(f, )C,(f, ) sin
+ + = G(f, 9) (2.45)
at az ay
where 0 is defined counterclockwise from the positive x axis. Here, unlike the
deepwater case, both Cg 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
vx (ak cos 0) (ak sin) = 0 (2.46)
Vxay a0 (2.46)
ay ax
Therefore, if the source function G(f, 0) is known in advance, the directional
spectrum E(f, 0) 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 wavewave 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,0) = a(f, ) + 3E(f,0) (2.47)
where a and )3 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 3 is proposed by Barnett (1968)
S5sf( 0.9), if > 0.9
= 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
S= k3(f, oo) Cg(f, oo) 1 k (f, oo)
(k, h) (2.49)
k3(f, h)Cg(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 fetchlimited and duration
limited cases. For a steady state, onedimensional, and frictionless situation, the
29
growth of a spectral frequency component in the downwind direction is governed
by the following ordinary differential equation
dx
Cg(f) E(f) = a + 0E(f)
By specifying E(f) = 0 at x = 0 as the boundary condition, the solution for spectral
component E(f) is
E(f)= [exp() 1]
For small fetch, E(f) approximates (a/Cg)x, which corresponds to the Phillips'
linear growth mechanism. For large fetch, E(f) approximates (a/P) exp({z/Cg),
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 fetchlimited case appears
in Fig. 2.4(a).
For a durationdependent, onedimensional, and frictionless situation, the growth
of a spectral component can be determined from the following ordinary differential
equation
E(f) = a + 3E(f)
By specifing E(f) = 0 at t = 0 as the initial condition, the general solution for E(f)
is
E(f) = a[exp(3t) 1]
For small duration, E(f) approximates at, which corresponds to the Phillips' mech
anism. For large duration, E(f) approximates (a/3) exp(pt), which corresponds to
the Miles' mechanism. In the real situation the exponential growth is again limited
E(f)
EXPONENTIAL GROWTH
LINEAR
GROWTH
S E(f) = (a/C,)z
0
fetch, x
(b)
OVERSHOOT SATURATION
PERIOD PERIOD
E(f)
EXPONENTIAL GROWTH
LINEAR
GROWTH
E(f) = at
0
o duration, t
Figure 2.4: Schematic illustration of growth of a spectral component, E(f), at (a)
fetchlimited case and at (b) durationlimited 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 durationlimited case appears in Fig. 2.4(b).
The rate of dissipation in the shallow water, Gdi,(f, ), 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.,
Gdis,b(f, ) dfdO = drbdub(f, ) /(pg)
where ub is horizontal velocity component at sea bed, rb is the bottom shear stress,
and p, is the water density. By expressing the bottom shear stress based on the
empirical law
b = pCf b(f,0) jUb(f,0)I
where cf is bottom friction coefficient, and expressing dub(f, 0) in terms of dr7(f, 0)
according to linear wave theory,
g k
dub (f,0) = gk (f, 0)
a cosh kh
the rate Gdi,,b(f,O) may be approximated by the equation (Collins, 1972)
Gda,,b(f, )  gCf Ub E(f, 0) (2.50)
a2 cosh2 kh
where (ub) denotes the root mean square value of horizontal wave orbital velocity
component at bed and can be computed from the expression
(ub)2 = [Ub(S,t)]2 = [db(f, 0)l2 = g2k [d(f,0)
= osh2 dsh (2.51)
7r 9 E( 0do df (2.51)
Jo "r O cosh 2kh
32
The rate of dissipation of wave energy in the shallow water due to effect of a
porous sandy bottom, Gdi,p(f, 0), can be estimated by the work done by the action
of waveinduced dynamic pressure at the surface of the bed to cause the fluid motion
in the porous sand layer.
Gdi,p(f, 0) dfdO = dpbdwb /(pg) (2.52)
where pa 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( A1/A kd)
Gd,,,p(f, 0) = kX 2 osh2 kh E(f, 0) (2.53)
cosh' kh
where Ax 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,
Vk, and turbulent eddy viscosity, vt, according to a theoretical analysis made by
Lamb (1932), can be shown as
Gdi,,.(f ) = 4 (vk + Ut)k2E(f, 0) (2.54)
The kinematic viscosity of sea water vk, 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
Gdis,v(f ) = c, f Ef, 0) (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
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(2r f)s 1(k,h) (2.56)
where B is a coefficient to be of the order of 102, and ((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 wavebottom 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)
(Etotal,b)l/ = 0.2 tanh(kh) (2.57)
kb
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
63512) 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 PUV data based on the method proposed by LonguetHiggins 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 shoreconnected
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 mlevel 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
system (1/12/84).32
4H 24 450N
U ST S
350N
4 000N
R T RN IC OC R 250N
80OW 75OW 700W 65OW 600W 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 orthogonalspherical area with the east and the west
boundaries coinciding with 70.00W and 82.50W, respectively, and with the north
and south boundaries coinciding with 42.50N and 22.50N, respectively. The grids are
all parallel to the longitudes and latitudes such that each grid is 2.5 degrees apart
500N
350N
80N 75N 70N 65N 60N 550N 50W
Figure 2.6: (b) Weather map showing northeast winds generated by a high pressure
system (1/21/84).
20
4 002
1006
A T IC OC A 250N
20
800W 750W 700W 65aW 60aW 550W 500W
Figure 2.6: (b) Weather map showing northeast winds generated by a high pressure
system (1/21/84).
1200Z 5/11/71
2L 450N
UN ST ES
20 300N
4 006
20 j RTLRN IC OC R 250N
80OW 75aW 700W 650W 60OW 550W 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 mlevel 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
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 quasistationary 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
fortyone 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 H1/3 or H,, is defined as the average height of the onethird
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 (Eotal)1/2
The correlation coefficient, 6, given in Fig. 2.9 is defined as follows:
N N N
6 = Z[(X)i(Xm)i(X)([/{E( X(X.m)(Xm)i]2}12 (2.58)
i=1 i=l i=l
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, fHI,/g, are plotted in
Fig. 2.10 for waves with H, >3 m, 3 m> I, > 2 m, and H, <2 m.
E(f) (f ),
Hz
(NOAA)
m
E(f) ( ),
(WI)
(WIS)
f (Hz)
4I
30 
10 20. (1984)
10 JAN. (1984)
0.2
 u.u
f (Hz)
Figure 2.8: Comparison of NOAA and WIS wave spectra for January 1984 at the
Buoy #41006 location.
50 
30
20 L
1
10 
JAN. 1984)
1
10 I
STATION: NOAA BUOY z41006 MEASURE DATA(NORA)
6 =0.22 COMPUTER OATAP(WIS)
  * *..* "=
** **" .... "',".
S I I I I I I I I I I I I I I I I
1 5 10 15 20 25 30
JANUlRT (1984)
1 5 10 15 20 25 30
JANUARY (1984)
Figure 2.9: Comparisons of NOAA a
heights at the Buoy #41006 location.
nd WIS modal periods and significant wave
20
15
(SEC0
(SECJ
10
8
0.12
0. 10
0.08
10 fH,
g
(WIS)
0.06
0.04 
0.02 
0.0
0.0
0.02 0.04 0.06 0.08
f2H,
10 f (NOAA)
g
Figure 2.10: Comparison of NOAA and WIS wave energy scale parameters at the
Buoy #41006 location.
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.
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 (Gunther 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 airsea interaction. The
model is calibrated by utilizing the field wind and wave data 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.
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
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,. The growth of the downslope 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 wavewave 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:
8 a dr
G (f) c,a' fEjM(f)  os4( O), if O and IB0 <
G.n,l(/, ) = 37r 2 2
0, elsewhere
(3.1)
where c, is a nondimensional coefficient, regarded by some investigators to be a
function of density ratio of air to water, and EJM(f) is the onedimensional modified
JONSWAP spectrum proposed recently by Donelan et al.(1985):
a g2 exp[ ( fm
EM(f) = (2r) fexp[(fm/f)4 2 a fm (3.2)
with all its free parameters a, f,, 7, and ab defined the same as in the JONSWAP
1.6
fm = 0.1(Hz)
1.'4
a = 0.0075
G., + Ga, (WIS)
1.2  Gi,,,, (UCWP)
1.0 
106 G
0.8 \
\
0.6 
S 1
I \
I \
0.2 
0.0 / ...
0.0 0.05 0.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.
48
spectrum. Substituting Eq.(3.2) into Eq.(3.1) and absorbing c, into a, we have
1 f fm
a3g2 exp[ ( )f ] 8
Sexp[(f,/f)4] 7 2 abfm Cos4( ),
(2r)4ffm 37r
Gi,1(f, 0) = if 10 and 10  <5
2 2
0, elsewhere
(3.3)
The peak frequency fm is determined from the following equation:
fk+1 = {(/)7/3 + a( COS( ( )4/3At_3/ (3.4)
g
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, fma. This peak energy transfer frequency is
defined as
E(f,,) = max {E(f) > E,(f)}
where E,(f) is a threshold spectrum below which no energy transfer will occur from
high frequency components to low frequency components. To use fm, 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 Ec(f)
should, at least, satisfy the following criteria: (i) The threshold spectrum should
fall below the saturation (equalibrium) 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 E,(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(fm) versus fm from the data of buoys #41002 and #41006. If the mean
values of the data points are connected, they form a skewed bellshaped curve with
a peak value, fp, occurring at frequency between 0.09 and 0.11 Hz. This indicates
that energy transfer experiences strong resistance when fm becomes less than this
peak value.
Based upon the above observations, we propose here the following functional
form for E,(f):
E,(f) = c, (2)4e expH 2] (3.5)
which is the PM equilibrium spectrum modified by a multiplier. Here cm, is an
empirical scale factor which is chosen as equal to 0.183 by matching Eq.(3.5) with the
data mean at designated fp=0.1 Hz. Both Eq.(3.5) and the original PM 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 PM 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 PM 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 multivariables
including wind speed, wave phase speed, peak energy frequency, etc.
E(
E Hz
100
80
60
40
20
0.00
m2
Hz
FIELD OATA COLLECTED BT
SNORA BUOY .41002 (198419851
S  PROPOSED THRESHOLD SPECTRUM
 SPECTRAL EQUILIBRIUM RANGE
INOICRTEO BY PM SPECTRUM
 RVERAGEO SPECTRAL FORH
OBTAINED FROM OATA SAMPLE
*
*
e
\ *
I U
\* *
' i
U a *..
0.OS
0.05
0.00
0.10
0.10
0.15
f, (Hz)
0.15
fm (Hz)
0.20
0.20
0.25
0.30
FIELD DATA COLLECTED BY
NORA BUOTY 41006 (198419851
 PROPOSED THRESHOLD SPECTRUM
SPECTRAL EQUILIBRIUM RANGE
INOICATEO BY PM SPECTRUM
 VERAGEO SPECTRAL FORM
OBTAINED FROM DATA SAMPLE
I
*
. ; *
I I . .
*tt=4i
0.20 0.25 0.30
Figure 3.2: Plots of E(f,) versus f,.
51
The nonlinear transfer coefficient, a., 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 al evaluated from the present
NOAA buoy data appear to be varying instead of constant. Figure 3.3 shows an.
values from buoys #41002 and #41006 plotted against the nondimensional param
eter fo/(f,Aft). 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 onehour time interval, At, in computing an, 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'/,t)  0 (steady
wind over long period), we arrive at
an = 26[0.0195 + ( fo )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/(f/,2 t) 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, y, 7a, 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 fm,/fo in Fig. 3.4. The data are seen more scattering when
fm,/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
1 
ani 0.1
0.01
0.001
0.0001
0.0001
0.001 0.01 0.1 1
fo
f A(
0.01
0.001
0.0001
0.0001
0.001 0.01 0.1 1
fo
f At
Figure 3.3: Plots of anL versus fo/fo At.
+ 2M>Hs
A 3M>Hsa2n
0 Hsa3M
+
0.001 I
0.001
 PRESENT STUOT
 HASSELMRNN ET AL.(1976)
........ NELAN ET AL.(1985)
+ 4 +
0 A
0.01 0.1 1 1
Figure 3.4: Plot of a versus f,/fo.
0.05
0.02
0.01 
0.005 1
0.002 
54
the values of AE(f)/At from the data are contaminated by the advected wave
energy into the region, the bestfit 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 ( )7/6 (3.7)
or equivalently,
a = 0.0122 ( U1 )7/6
Cms
where C,, = g/27rfm is the wave celerity at fm,. The ratio fs/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 fi.
a = 0.0091 (/ )2/3
The empirical equation for a proposed by Donelan et al.(1983) is
a = 0.006 ( 55
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)'/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 + 2M>Hs +
A 3M>Hsk2M ++, +
100 0 Hs23M + +
+ ++
++
50 A +
+ + +A+ ++++ +
20 +
y., Eq.(3.8) + +
10 20 1+ ++
+ ++
..+ +
2 ~ 1/2++
1 2 5 10 20 50 100 200
Figure 3.5: Plot of y computed from Eq.(3.8) versus y from the least square method.
following equation (see Fig. 3.5):
S= exp[ 2 i/2.5 ) (3.8)
The above equation gives a minimum y ; 2.38 when fm,/fo =1.5. Although yf is
not bounded for small fm,/fo, the spectral peak value of G,,i1, which is proportional
to a37 exp[2f,./fo]((f../fo)1/2, is bounded for small fm,/fo.
Values a, and ab obtained from the field data are scattered (see Fig. 3.6). The
arithmetic mean values of a, 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
0.5
0.4
0.3
0.0
0.
0.4
0.3
0.0
o0.
0 0
a + 2n>Hs
+
S A a A 3m>Hs22m
a
0 Hs>3M
A
a
+
A
+ AA 0
+ A
+++
+a + +
0 .
+ 9a+ ++ +
Z ++ + ++
A 0 It
0 0.5 1.0 1.5 2
+ 2n>Hs
a 3m>Hs22m
+ 0 Hs23m
+
+
A +
++ ++ A
+ 4 A + ++
+ +.J_4 t ga ++
++
a + *
+
A 0.5 I 1 I
3 0.5 1.0 1.5 2
Figure 3.6: Plots of ab versus fm/fo.
57
of a, and ab control the spectral width of the forward face and downslope 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 downslope 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., ab 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 a, should represent good estimate of energy gain due to energy
transfer and no adjustment is made.
The source function Gi,, defined in Eq.(3.3) is also plotted in Fig. 3.1 for
comparison. The energy growth on the downslope 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 windwave 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/2irU1o cos(0 0,) in Eqs.(3.7) and
(3.8), which are used to estimate a and 7, 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, 0,, 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, Gi,,2, are
proposed as follows:
Case a. AO, > Oc 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 Gi,.2 oc E(f,O). 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,O) = b[E,(f,O) E(f,O)]
(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(Ulo/C)4 with s the density ratio of air to water and C the
wave celerity at frequency f, and E, is the PM spectrum with the scale parameter,
ep, equal to 0.015, i.e.,
0.1 sf( ) [Eo(f) 2 cos( 0,) E(f,)], if 0 0, 7r/2
Gin,2 (f, 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 0c. For convenience, we have
chosen
O= ( rad/hr) At (3.11)
Case b. A0, < 0c 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, 0,, instead of the prevailing
wind direction, 0. Thus, the source function is assumed to have similar form as
Eq.(3.3) with 0 being substituted by 0,, i.e.,
3 g f 8
S9 exp[( ")] cos4(O 0,), if 0 0,5 xr/2
Gn,2(f 0) = (27)4 f3 f,o 3r
0, elsewhere
(3.12)
The spectral peak frequency fmo 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
f {(o)7/3 aa cos(0 10,)4/3At})3/7 (3.13)
g
and
co = 0.0122 ( f)7/6 (3.14)
fo
60
with fo = g/27rUo cos(AO,) being defined in Eq.(3.14) to account for the influence
of AO, to the resonance condition.
Case c. 0, = 0 In this case, the wind direction is the same as the prevailing wave
direction, we have Gin,2(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 durationlimited case, the
transport equation can be solved analytically to yield the spectrum
E(f, 0) = E,(f){1 exp[0.1sf( )4 t] 2 cos (0 _,) (3.15)
C 7rT
which vanishes at t = 0, and is equal to EZ(f) (2/r) cos2(0 0,) at t = oo. A
simulation of the growth of wave spectrum based on Eq.(3.15) under 30knot 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
Gi,,2(f,0) can be rewritten as
exp[0.1 sf(Uo/C) ] Uxo
Gi,,(f, 0) = 0.lsf Ein,2(f, 0) exp[O.sf(U/C)4 t I )4 (3.16)
1 exp[0.1sf(UIo/C)4t] 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
Uio/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
E(f) ( )
m2
E(f) ( H)
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 durationlimited case; U1o = 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, 0) = 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 durationlimited case, we obtain
Ek+(f,0) = Ei,(f ,O) e sfat/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):
Idp
Gdi(f, 0) = c. f E(f, 0) IdI E(f,0) (3.19)
PwCldrl\
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'o) (Phillips, 1967; LonguetHiggins, Cartwright, and
30
20 
20
E(f) (m)
Hz
10 
0
0.0
30 
20 
20
E(f) ()
Hz
10 
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:
sUV
Gdis(f, 0) = c, f E(f, 0) 0.02 E(f, 0) (3.20)
C H,
which is a modification of that in high pressure systems by adding an additional term
0.02[sU20/CH,}E(f,0). If this is the only source term in the transport equation,
the wave energy evolution at the durationlimited condition can be found as
Ek (f, 0) = Ek,(f,0) esfat/6AAt (3.21)
where A = 0.02sU0o/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:
0.1 sf( U)4 [E(f) 2cos2(00,) E(f,)], if 100,5 7r/2
G,, = C 7r (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 durationlimited case can be expressed as
Ek+1(f,0) = Ei.,(f,O)exp(O.lsf( C) At]
C 7r
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 fm2
2 5 expl ( l)
E.(f) = g exp[4(f,/f)1 2 cab fm ](k,h) (3.24)
where $(k, h) is a transfer function defined earlier in Eq.(2.49):
S(k,h) = k3(f,) C(f, oo)
k3(f, h) C,(f, h)
The spectral parameter a is estimated in the model by the following formula:
a = 0.0156 ( k)O0.49 (3.25)
where km is the corresponding wave number at peak spectral frequency, f,, which in
turn is related to the resonant frequency as fm = 0.88 fo. Equation (3.25) describes
a bestfit curve of the upper bound values of a computed, based on field data, by
Vincent (1984). The spectral parameters 7, and ab are also adopted from Vincent
(1984):
=2.47(U k 39, (3.26)
g
0.07 if f < fm,
Cab = (3.27)
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
fmn A ) /2 (3.28)
27rwhere A lies between h
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 g2(2r f)5(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=[ifu 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= 22 sin 0 is
the local component of the planetary vorticity normal to the earth's surface, or the
Coriolis parameter with 2 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, (x, y), 9p/an can be obtained from the expression
ap_ 8 app 2 ( ),
+n az ay
and the radius of curvature, r, is obtained from
a(p~ ap (p aax ap ay a9p 9p
r = 1/( ), with =  + y os 0, + sin 0
ae z ae axa 1yae a ay
where e is the distance along the path where air flows, and 09 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, ,,, are adjusted according to the following
equations:
Uz 1 z
= In (3.30)
U. 0.4 zo
and
60.4 u
0, 0, = sin'(Bo/ 0 ) (3.31)
U,
68
where the frictional velocity, U., and the frictional height, z0, are computed by
(Clarke and Hess, 1974)
0.4u U
( ) = (Bo)2 + [In( ) Ao]2 (3.32)
U. f zo
and Eq.(2.42), or
0.1525 U2
o= + 0.0144 0.00371, (in cgs units)
U. g
where Ao and Bo are nondimensional 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 fetchlimited or durationlimited 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 fetchlimited and durationlimited 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 fetchlimited cases with
constant wind speeds equal to 20, 30, 40, and 50 knots, respectively. The results of
significant wave height H, and modal period Tm 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,
Tm along a fetch F are estimated from the following formulas, respectively:
g H, = 0.283 tanh(0.0125(F)3/7] (4.1)
U., U.
gT = 7.54 tanh0.077(F)'/4] (4.2)
U10 U120
where U,, the adjusted wind speed, is defined as a power function of U1o (Shore
Protection Manual, 1984):
U, = 0.71 Uo4 (in metersecond 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:
SH,= 0.0016 (L)1/2 (4.4)
 N
U,
an
I .
In)
N
71
0
O
, 0
, n
1 i
.4 u
0 2
\:i \ a j.
O l, 
so
\ I \
M 'u4a'4
2
j 4
D' \
\ E 3 )
LLL.
U)
\\ C,
\\n
> a \M
\11 0
\ \
ED \ V 0 e
 \ C ) \ 
'S
\ i
~o
.4.' ~
(0 IO
.; z
gT 0.286 (F)'/3 (4.5)
Ul0 U210
In the classical SMB and parametric models, the adjusted wind speed Ua 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 Ua according to Eq.(4.3) in the
both models is suggested by the Corps of Engineers, U.S. Army, and Ua obtained
from the equation will be greater than Uo1 if U10 is greater than 4 meters per second.
Therefore, by using the adjusted wind speed Ua as defined in Eq.(4.3), the SMB
and parametric models will predict higher wave heights than the classical ones if
U10 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 durationlimited
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 = 0.283 tanh( 0.00087( t)5/8] (4.6)
U2 U.
gT. = 7.54 tanh[ 0.016( t)3/8] (4.7)
U1o Uo0
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 (t 5/7 (4.8)
U2 Uo
g T 0.06( t )3/7 (4.9)
U0o U1o
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
=r 0
C N
cu
N M ~ 
I L
a 3 N C 3
cc
3
0
0:
aC
"o
Cz
0
0 ~z
N N
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 fetchlimited 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.
E tIOEEP WATER WAVE / /
INFORMATION STATION ''
0 iSMHALLO WATEB MAVE
REPORTING STATION
il 002
1 132 N
KINGS T _,
*S SIAT ON
30ON
HARINE LAN
LI2 10006W
A/TLANTIC
.. _OCERN
CAPE ;ANAVEAAL
STATII 280N
FLOR IR oo200M ) 260N
82W 80 W 780W 76 0 74'W
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.50W are shown
in Fig. 4.4 together with the measured winds from NOAA buoy #41006, located
at 29.30N and 77.30W. 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
ORMSE = {E[(,c)i (0w,m)i]}J//N (4.10)
i=l
where 80,, is the computed wind direction, 0w,m 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.
78
STATION: NORA BUOT *41006 MEASURED ATA (NORA)
8~me=59.70 COMPUTED OATA
N ** .
SI I I I
S i A" \ .* ** / t
1 5 10 15 20 25 30
JANUARY (1984)
50
STATION: NOAR BUOY 41006 MEASURED OATA(NOAA)
6 = 0.67 COMPUTED DATA
40 
10 30
(KNOT)i
20'* i .
10 I * .
1 5 10 15 20 25 30
JANUARY (1984)
Figure 4.4: Comparison of hindcasted and measured winds at the Buoy #41006
location.
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 onedimensional 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 twodimensional spectra over the
60
so0 
i
30
I
20 
10
I
E(f) ( ),
(NOAA)
m"
E(f) ( ),
(UCWP)
U. 2
Pxr20
0 JAN. 1 15
0 JAN. (1980)
f (Hz)
J0 N.198)
0 JAN. (1984)
U. I 0.0
f (Hz)
Figure 4.5: Comparison of NOAA and UCWP wave spectra for January 1984 at
the Buoy #41006 location.
81
20
STATION: NOAR 8UOT =41006 MEASURED DATA (NOA)
15 6=0.38 COMPUTED DATA (UCHP)
T 1I0 10
(SEC)
0 I II I I I I I I I I I tI I I I l i I II I I I
1 5 10 15 20 25 30
JANUARY (1984)
10
STATION: NORA BUOY .41006 MEASURED DATA (NOAR
8 6 =0.85 COMPUTED OATA (UCWP)
Hs 6
2
1 5 10 15 20 25 30
JANUART (1984)
Figure 4.6: Comparisons of NOAA and UCWP modal periods and significant wave
heights at the Buoy #41006 location.
0.12
0.10 I
0.08 F
10 f H
(CW
(UCWP)
0.06 H
0.04
0.02 
0.0
0.0
+ 2n>Hs 6 = 0.44
A 3M>Hs>2M 6 = 0.43
O Hs>3M 6 = 0.67
O A O A
A
A +
+ +
A + +
. A
A A
A A ^. 0 A A,
+ + ++ +
A A A
.4. ~ A "
0.02
0.04
0.06
0.08
f2 H
10 , (NOAA)
9
Figure 4.7: Comparison of NOAA and UCWP energy scale parameters at the Buoy
#41006 location.
2
E(f) ( ),
(CDN)
E(f) ( ),
HzP)
(UCWP)
0.2
30
3
hA
JAN. (1984)
n 1
S0.0
f (Hz)
FI,
0.2
0.2
f (Hz)
Figure 4.8: Comparison of CDN and UCWP wave spectra for January 1984 at the
Kings Bay gage location.
25
10 JAN. (1984)
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) ( m),
rad
(CDN)
E() ( ),
rad
(UCWP)
0.2
0.1 I
r I 1'0
10 (1. I98q)1
0.3r
0. 1 I
JRN. 11984)
Figure 4.9: Comparison of CDN and UCWP marginal directional spectra at the
Kings Bay gage location.
10.0
E W N E
730'
t 7
w E
KINGS BAT(FL.) STATION
E ns = 16. 7
* MEASURED DATA(CON)
 COMPUTED DATA(UCWP)
TN
(SEC)
1 5 10 15 20 25 30
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.
50N
45N
'450N
400N
350N
30N
250N
800w 750W 700W 650W 600N 550W 500w
Figure 4.11: Tracks of Hurricanes Diana, Isidore, and Josephine (1984).
500N
120DZ 9/(/84
45"N
UN TLNIC OCSTA 250N
220
80"0 750W 700W 65'W 600W 550W 500W
Figure 4.12: (a) Weather map showing northeast winds generated by the Hurricane
Diana (9/8/84).
I 12doz 9/27/4) 24
3 50N
UN STA
28
350N
16
1 002
300N
a 1S IDORE
0A TLAN IC OC R 250N
12
800W 750W 700W 650W 600W 550W 500W
Figure 4.12: (b) Weather map showing northeast winds generated by the Hurricane
Isidore (9/27/84).
