Front Cover
 Title Page
 Table of Contents
 List of Figures
 Theoretical review
 Governing equations
 Finite differencing
 Laboratory work
 Results and conclusions

Group Title: UFL/COEL (University of Florida. Coastal and Oceanographic Engineering Laboratory) ; 89/014
Title: Wave-current interaction over a submerged bar field
Full Citation
Permanent Link: http://ufdc.ufl.edu/UF00076140/00001
 Material Information
Title: Wave-current interaction over a submerged bar field
Series Title: UFLCOEL
Physical Description: ix, 102 leaves : ill. ; 28 cm.
Language: English
Creator: McSherry, Thomas Richard
University of Florida -- Coastal and Oceanographic Engineering Laboratory
Publication Date: 1989
Subject: Shore protection   ( lcsh )
Breakwaters   ( lcsh )
Ocean waves   ( lcsh )
Coastal and Oceanographic Engineering thesis M.S
Coastal and Oceanographic Engineering -- Dissertations, Academic -- UF
Genre: non-fiction   ( marcgt )
Thesis: Thesis (M.S.)--University of Florida, 1989.
Bibliography: Includes bibliographical references.
Statement of Responsibility: by Thomas Richard McSherry.
Funding: This publication is being made available as part of the report series written by the faculty, staff, and students of the Coastal and Oceanographic Program of the Department of Civil and Coastal Engineering.
 Record Information
Bibliographic ID: UF00076140
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved, Board of Trustees of the University of Florida
Resource Identifier: oclc - 20312025


This item has the following downloads:

UF00076140 ( PDF )

Table of Contents
    Front Cover
        Front Cover
    Title Page
        Title Page
    Table of Contents
        Table of Contents 1
        Table of Contents 2
    List of Figures
        List of Figures 1
        List of Figures 2
        List of Figures 3
        Abstract 1
        Abstract 2
        Page 1
        Page 2
        Page 3
        Page 4
    Theoretical review
        Page 5
        Page 6
        Page 7
        Page 8
        Page 9
        Page 10
        Page 11
    Governing equations
        Page 12
        Page 13
        Page 14
        Page 15
        Page 16
        Page 17
        Page 18
        Page 19
        Page 20
        Page 21
        Page 22
        Page 23
    Finite differencing
        Page 24
        Page 25
        Page 26
        Page 27
        Page 28
        Page 29
        Page 30
        Page 31
        Page 32
        Page 33
        Page 34
        Page 35
        Page 36
        Page 37
        Page 38
        Page 39
        Page 40
        Page 41
        Page 42
        Page 43
        Page 44
        Page 45
        Page 46
        Page 47
    Laboratory work
        Page 48
        Page 49
        Page 50
        Page 51
        Page 52
        Page 53
        Page 54
        Page 55
        Page 56
        Page 57
        Page 58
        Page 59
        Page 60
        Page 61
        Page 62
        Page 63
        Page 64
        Page 65
        Page 66
        Page 67
        Page 68
        Page 69
        Page 70
        Page 71
        Page 72
        Page 73
        Page 74
        Page 75
    Results and conclusions
        Page 76
        Page 77
        Page 78
        Page 79
        Page 80
        Page 81
        Page 82
        Page 83
        Page 84
        Page 85
        Page 86
        Page 87
        Page 88
        Page 89
        Page 90
        Page 91
        Page 92
        Page 93
        Page 94
        Page 95
        Page 96
        Page 97
        Page 98
        Page 99
        Page 100
        Page 101
Full Text




Thomas Richard McSherry










In the course of completing this thesis I have received abundant help in many forms

from a variety of people. The work required the extensive help of my Advisor, Dr. James

Kirby, and I am truly grateful for his assistance. I consider myself very fortunate to have

had the opportunity to study under Dr. Kirby. I also used the COEL extensively on two

separate occasions, and depended on the help of several people there. I thank Jim Joiner

and Sydney Schofield for their patient assistance, and the inputs from the entire group. A

great source of help at the COEL was in the way of mental release in the form of the daily

volleyball game, and although I jammed a few fingers in the less-than-benevolent net play,

I am grateful for the enthusiastic games.

I also thank Dr. Robert Dean for his help in every way, and I thank the secretaries who

have told me the logistical requirements for getting out from under the deadline hammer.

Also, thanks to the students that offered friendship and release, and provided help when I

really needed it. These people I wish the greatest of success.



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

ABSTRACT ..............


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

1.1 Problem Statement .......

1.2 Literature Review .......


2.1 Introduction ..........


3.1 Introduction ..........

3.2 Circulation Model ......

3.2.1 Governing Equations .

3.2.2 Radiation Stress Terms

3.3 Wave Model ..........

3.3.1 Governing Equation


4.1 Introduction ........ .

4.2 Circulation Model .......

4.2.1 Method of Solution .



Boundary Conditions .

Radiation Stress Subroutine

. . . . . . . .

4.3 Wave Model ...... ........ ........................ 35

4.3.1 Method of Solution .................... ........ ... .. 35

4.3.2 Finite Differencing .. ......... ................ .. 35

4.3.3 Boundary Conditions ......................... 37

4.3.4 Lateral Smoothing ..... ......... ........... .. 37

5 LABORATORY WORK ...... .... .... .............. 48

5.1 Introduction .. ... .... ...... ..... .. ... ... ... ...... 48

5.2 Equipment .... ...... ... .................. .. 48

5.2.1 Basin .............. ........... ........ 48

5.2.2 W avemaker .... ............ ............... 50

5.2.3 Beach and Cart System ........................ .. 50

5.2.4 Electronics .. . ........................ 52

5.2.5 Bar Field .................... .. .......... .. 54

5.3 Setup, Procedure and Data Analysis . . . . ... 55

5.3.1 Session One ....... ........... ............. .. 56

5.3.2 Session Two ... ... ... ... ... .. .. ... ... .... .. 62

6 RESULTS AND CONCLUSIONS .......................... 76

6.1 Introduction .... ........ ..... ....... ........... 76

6.2 Laboratory Work ................. .. ..... ........ .. 76

6.2.1 Reflection Results ..... ........ ............ .. 76

6.2.2 Bathymetry ...... ......... ..... ...... .... 86

6.3 Numerical Results .......... .. ................ .. 87

6.4 Conclusions and Status Report ........................ 92

BIBLIOGRAPHY ............... .................. 100

BIOGRAPHICAL SKETCH .. ............................ 102


2.1 Bar field in the presence of a depth-uniform current. From Kirby (1988). 6

2.2 Contour of Tr v. IFx v. w/Wro. ....................... 9

2.3 Surface projection of transmission Tr v. frequency and current. . 10

4.1 Flow chart for model MCSERRY. . . . ..... 25

4.2 Hyperbolic tangent startup function. . . . . ... 26

4.3 Grid for circulation model .......................... 28

4.4 Differencing of radiation stress gradients in x-sweep . . .... 30

4.5 Differencing of advective acceleration in x-sweep. . . ... 32

4.6 Wave calculated from parabolic equation on domain of width W. 39

4.7 Spectrum of A(y) without smoothing. Cutoff at nA = 9.051. . 41

4.8 Spectrum of y-direction driving force before smoothing. . ... 42

4.9 Spectrum of y-direction driving force after smoothing. . ... 43

4.10 IAI without smoothing. Bars are located on the bottom of this domain. 44

4.11 IAI after smoothing. Bars are located on the bottom of this domain... 45

4.12 y-direction forcing without smoothing. . . . ... 46

4.13 y-direction forcing after smoothing. . . . ........ 47

5.1 Plan view of large basin. ......................... .. 49

5.2 Plan view of cart track system. ........................... 51

5.3 Side view of profiling cart. ......................... 52

5.4 Sample profile. ....... ... ........ . . ... 53

5.5 Sample calibration curve for wave gage. .. . . . ... 54

5.6 Bar frame design. ............................... 55

5.7 Gage locations for two-gage reflection method. . . ... 57

5.8 Reconstructed bar profile for Session One. . . . ... 60

5.9 Plan view of barfield for Session One. . . . ..... 61

5.10 Gage locations for Session One. . . . ..... 62

5.11 Three gages inline with x-direction waves . . . ..... 63

5.12 Sample envelope, Session Two. Retrieved digitally from one of the three
on-offshore transects using the moving cart gage. . . .... 68

5.13 Sample envelope height versus X. . . . ..... 69

5.14 Sample Kr with offshore distance x from Session Two. . ... 69

5.15 Plot of Kr with transect location, Session Two. . . ... 70

5.16 Plan view of barfield, showing data collection locations. . ... 71

5.17 Sample time series of magnitude of current meter data, Session Two. 73

5.18 Profile for Session Two ............................ 74

5.19 Instrument setup for Session Two. . . . ..... 75

6.1 Plot of K, with wave period, Session One data and theoretical prediction
from 1-D model Kirby (1987) ......................... 78

6.2 Plot of Kr with period, data and theoretical for the case of no discernable
rip current. ....................... ............... 81

6.3 Plot of Kr with period for the case of a discernable shore rip current.
Data and theoretical .............................. 83

6.4 Two profiles, dotted line referring to profile used in the second and third
sets of testing, where a weak current formed. . . . ... 84

6.5 Plot of Kr with period for case of discernable shore rip current, moving
toward resonance ................... ........... 85

6.6 Bathymetry 1. Shoreward direction with increasing . .... 88

6.7 |A[ for initial amplitude of .03 meters and perfectly resonant waves. .. 89

6.8 IBI for initial amplitude of .03 meters and perfectly resonant waves,
without currents.......... ........... ......... 90

6.9 Steady state mean surface, ....................... 91

6.10 No-current wave model, incident amplitude of .02 meters. K, = 0.81. 93

6.11 No-current wave model, incident amplitude of .03 meters. Kr = 0.81. 94

6.12 No-current wave model, incident amplitude of .04 meters. Kr = 0.81. 95

6.13 No-current wave model, incident amplitude of .03 meters. K, = 0.81,
and bar half-length is 3 meters ....................... 96

6.14 No-current wave model, incident amplitude of .03 meters. Kr = 0.843,
and bar half-length is 2.25 meters. . . . ..... 97

6.15 No-current wave model, incident amplitude of .03 meters. Kr = 0.88,
and bar half-length is 1.75 meters. . . . ..... 98

6.16 No-current, tapered bars with .03 meter incident wave. Kr = .774... 99

Abstract of Thesis Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Master of Science




December 1989

Chairman: Dr. James T. Kirby
Major Department: Coastal and Oceanographic Engineering

People are moving to the coast in great numbers that seem to increase each year. Devel-

opment along the shoreline is therefore also on the uprise, and the public often assumes that

the shore environment is stable. Recognizing the appeal of the beach, and understanding

that this environment is anything but stable, the coastal engineer must attempt to curb the

natural destructive forces from the sea that would normally hinder man's use of the beach

and it's surroundings.

One present goal in this attempt is to achieve shore stabilization with a minimum

impact on the dynamic processes in the nearshore region. One candidate receiving much

recent investigation is the shore parallel low profile bar field. The bar field acts with the

incoming wave to actually reflect energy back out to sea. This action thus protects the

beach from some of the destructive wave energy that often erodes massive amounts of dune

material in heavy weather. If the concept could be perfected, the low profile bar field could

conceivably be a huge shield against any wave possible.

The way the bar field reflects energy is by having bar spacings that are one half the

incident wave wavelength. With this satisfied, the incident energy serves to form another

wave that propagates opposite to that of the incident wave. Knowing that the bar field

would be of finite longshore length, a low spot in energy would then exist shoreward of

the bar field. What would then be expected is to have a circulation begin by this energy

differential, much in the same manner as circulation into the lee of breakwater, or along

groins. It is also expected that the current would flow outward through the middle of the bar

field. Clearly, then, the current would alter the character of the incident wave, changing it's

wavelength. Then the resonant condition between the incident wave and the immoveable

bars would be upset. If the bar field is to considered for practical implementation, the

behaviour of the rip current driven by the energy differential will need to be investigated.

The most important aspect of the wave-driven current is how it responds to slight

shifts off resonant incident wave frequency. The main hypothesis of this paper is that

the bar-wave system seeks resonance, and therefore if the wave were to be initially out of

resonance a current would be required to tune the wave. By incrementally lowering the

initial wave frequency, and observing how the system responds and whether reflection is

enhanced or reduced, the hypothesis is proven or repudiated. Another factor in this scheme

is the available energy from the incident wave. If the waves are not big enough to create

sufficient energy differentials, then the current will never reach the strength required to

tune the waves that are out of resonance. The thesis thus looks at how current strength

responds to wave heights, and how reflection behaves during the process.

This paper studies results from laboratory data taken from a finite length bar field, and

also findings from a numerical model that incorporates coupled wave equations with the

momentum equations. The numerical scheme iterates between the two parts until a final

steady state is achieved. Mean values in surface elevation, lateral depth integrated currents

and wave amplitudes are solved for through a semi-implicit algorithm.


1.1 Problem Statement

The need for effective shore protection seems to be increasing with each year. The

worldwide dependence upon the coastline coupled with such chronic problems as rising sea

level has placed a responsibility upon coastal engineers to investigate truly effective methods

of shore stabilization. Effective means that it must reduce erosion, hold up under severe

conditions, be economically feasible, be environmentally sound, and not impede the regular

use of the coastline. Historically, structures meant to curb erosion have had some success,

but the side effects caused by the structures were oftentimes worse than the original problem.

Today, the main thought is to develop a device that works elegantly and efficiently. One

of the promising candidates is the low profile bar field that reflects incident wave energy

by use of a phenomenon now known in water wave theory as Bragg reflection, which was

named after the analogous process in crystallography.

By placing bars along the sea bottom parallel to the shore, and making their spacing

one half of the incident wavelength, a resonant condition occurs between these bars and

the surface waves. Specifically, a reflected wave is generated through the bottom boundary

condition. This reflected wave grows in amplitude with the number of bars placed until

theoretically complete reflection can occur if there were an infinite number of bars. The

findings to date have been very promising, with upwards of 80% reflection for waves travel-

ling over four sinusoidal bars. Previous laboratory data also confirm the expectations, and

these encouraging findings are fostering continued research into the low profile bar field.

These bars would necessarily be of finite longshore length, and if part of the incoming

energy were reflected, a low spot in energy would exist shoreward of the bar field, with

full amounts at the sides of the bars. This differential in the mean surface elevation along

the coast could possibly generate a large enough circulation to create a fairly strong rip

through the barfield and back offshore. What this means is that the incident wavefield

will be distorted over the barfield, and the resonant condition will be upset. An important

question that has been answered by Mei (1985), is how far off resonance the waves can be

before the bars quit reflecting energy. With the introduction of a current, however, the

question is whether reflection will be enhanced or reduced.

This thesis focuses on the wave-current interaction over the bar field, and what happens

when initial wave frequency is incrementally reduced. What the thesis tries to answer is

what happens to reflection given the current, and how the current changes as wave height

changes. Also, the thesis tries to prove that the wave-bar system seeks a resonant condition

by tuning off- resonant waves closer to resonance through the rip current. These questions

are answered by means of laboratory data and a numerical model. Some theoretical review

introduces the main content, and the problem is fully developed. Then the equations

necessary to the problem are presented, followed by a numerical solution to the equations.

After this, the laboratory work is covered, and the methods of data analysis are reviewed.

Finally, the results from the laboratory work and the numerical model are given, through

which the hypothesis will be proven or repudiated. Before all of this, however, the past

research that this thesis depends upon will be reviewed.

1.2 Literature Review

This paper is a direct extension of various avenues of previous research. The radiation

stress concept developed by Longuet-Higgins and Stewart (1964), is obviously important

in any study of waves and currents. These concepts were extended to terms involving

complex surface amplitude r' by Mei (1972). The expressions in Mei's book (1983) are used

here directly,and are the driving forces in the time averaged, depth integrated momentum

equations that govern the circulation and mean water surface.

The other terms appearing in the momentum equations have been included in the

model, and most of those formulations are based on previous work. The bottom friction

as developed by Longuet-Higgins (1970a) is represented in the model as a function of the

linear water wave maximum orbital velocity and the mean current U. In a companion

paper, Longuet-Higgins (1970b) looked at the transition of longshore current at the breaker

line, and developed an expression for the lateral mixing. This term also is included in the

present model.

The advective acceleration terms in the momentum equations are included in the model,

and are formulated as a four grid-point average of the current values, which is done by

Winer (1988). The momentum equations are solved by an Alternating Direction Implicit

(ADI) method which was introduced by Sheng and Butler (1982). The matrix set up for

this method can be seen also in Winer (1988). Although the matrix implicitly refers to

linear terms, the model introduces limitations in time step by including the nonlinear terms

explicitly. Iteration between wave model and circulation model was used by Noda et. al.

(1974), where the importance of wave current interaction in determining the final current

field was noted.

There has been much work on Bragg reflection in the last decade, as elegant and en-

vironmentally sound methods are sought for shore stabilization. A non-resonant theory

was developed by Davies and Heathershaw (1984), and partially successful comparisons to

laboratory data were made. The theory broke down at resonant frequencies, however, and

the reflection coefficients predicted were significantly larger than the data. Mei (1985) ex-

tended the theory to include resonant frequencies, and successfully predicted the reflections

seen in previous experimentation. Kirby (1986a) allowed for the presence of a mild slope in

the theory. The theory was further extended to include currents flowing over the bar field,

although on a flat bottom by Kirby (1988). Kirby (1986b) attempted to include currents

over a mild slope without the bar field. The work on flat topography noted that the pres-

ence of a current enhances reflection, and presented the background that this paper uses

as a starting point. The presence of a current allows incident waves of frequencies slightly

off resonance to tune with the bar field and actually achieve greater reflections than with

perfectly resonant waves. Transmission coefficients always drop in the presence of a current.

Thus with the joint action of a depressed energy zone driving a rip current that can alter the

incident wavelength, and transmission reduced against the current which further reduces

the low energy zone, the system may tend to want to go toward perfect tuning. The scheme

depends largely, however, on the incoming available energy to drive the current. Laboratory

data, will demonstrate cases where incident wave heights were too small to sufficiently create

a rip current strong enough to keep the reflective plateau.

Various methods are used to measure reflection from the laboratory data, and the

previous work that made these methods available will be cited in Chapter Five.


2.1 Introduction

To begin the task of answering some of the questions raised in the last chapter, specific

wave theory should be reviewed. It has been conjectured already that reflection will increase

when a current flows over the bar field against the waves.. The thesis will seek to prove this,

and the more important concept that the entire system will try to reach a resonant condition.

The mathematical basis for these hypotheses was developed by Kirby (1988).

Referring to the following figure 2.1 for the case of a current over a sinusoidal bar field

of length L, with bar wavenumber A = 27r/Lb, Kirby showed that if one considers incident

waves at or very near resonance, the ratio of transmission with a current to transmission

without a current is

SA(L)current cosh QL (2.1
T A(L)nocurrent coshQcL

The subscripts c refer to properties in the presence of a current. The bar field has length

L, where A(L) is the value of the transmitted wave amplitude at the back edge of the bars.

The lateral bar shape is described as Db. The other terms are defined as

Q = Q (2.2)

Qc = (2.3)

o gkDb (2.4)

Figure 2.1: Bar field in the presence of a depth-uniform current. From Kirby (1988).

( U 0)2 (2.5)
Q 2,,2 )

where C1 = Cgl + U and C2 = Cg2 U. Qc is defined later by (3.24). Note that when a

strong enough current opposes either wave, the wave can be stopped when Cg. + U = 0.

To simplify the expressions and to show what happens when the current strength changes,

consider a shallow water wave where the expressions reduce to the following;

ADb 2-l
Qc 8 -, F)-

where the Froude number in the x-direction is F, = U/vg-h. It appears that Qc has a

singularity as IFxI -- 1, which corresponds to the stopping Froude number. When IFx\ < 1,

from the above simplifications it is seen that Qc > Q, and therefore TA is less than unity.

What this means is that in the presence of a current, transmission for waves which are

perfectly resonant with the bars with an opposing current is less than the transmission

-~ x

---- U


for perfectly resonant waves without a current. The current is thus perceived to enhance

reflection, and it is interesting to note that this is independent of the direction of the current

along the x-axis.

Now for the shallow water limit consider resonant frequencies, where wr is in the pres-

ence of a current, and wro is without an imposed current. These are

Wr = )A(1-F (2.6)

1ro = IVghA (2.7)

which from the ratio wr/wro = 1 Fx describes a parabola. Refer to figure 5 in Kirby

(1988) for a picture of this ratio for various values of Ah, where Ah = 0 is also shown as the

shallow water limit. Now consider any frequency w shifted off the resonant frequency wr by

an amount 0, so that w = Wr + There is a cutoff condition for this shift, above which

for a given current, the resonance is lost between the wave and the bar field. This cutoff is

described by

S2 ; = ffg- C= h. (2.8)
cutoff = 64C2

The behaviour of the wave field depends'on the relation between Q and Qcutoff. Mei (1985)

and Kirby (1988) denoted four cases where the solutions to the coupled equations described

different wave fields, and these solutions for A(x) and B(x) can be put in terms of the

transmission coefficients, Tr = IA(Lb)/A(O)I. These expressions for the four cases are

case 1: 0 > Qcutoff

T, = 1 (2.9)
[(PC1)2cos2PLb + (S')2sin2PLb]/2

case 2: t = 2,cutoff

T,= (2.10)
[C2 + (fLb)211/2

case 3: Q < cutoff

T =QC (2.11)
[(QC1)2cosh2QLb + (QI)2sinh2QLb]1/2

case 4: Q = 0 (perfect tuning)

Tr = (2.12)

These relations may be evaluated close to shore using the approximations:

Ci = vgh(l + F)

C2 = Vh(-F)
gADb 1
8^ _( F.2

1- Fx
p2 = 2[ 02 2 f
r ~ \2gh(l F)

Q2 = -P2 when Q < ,cutoff

These expressions can be put into a program to plot the surface projection of transmission

Tr for any combination of w/wro and IFxI. Figure 2.2 shows the contour of the transmission

coefficient, and figure 2.3 is a surface projection plot of the same result.

Considering an incident frequency that is not far off wro, say 85% of it's value, then

by no-current resonant theory the transmission coefficient will be 0.8. The mean surface

elevation in the surfzone shoreward of the bar field would then be lower than surrounding

areas, and the differential could be enough to drive a current. By following this line of

constant wl/ro along the surface projection, clearly as the current increases the line heads


JR FROM 1.5B000E-1 TO 1.8101 CONTOUR INTERVAL OF 8.S1100E-11 PT13.31- 1.99934 LABELS SCALED BY 196B.
Figure 2.2: Contour of T, v. IFI v. /Wro.


Figure 2.3: Surface projection of transmission Tr v. frequency and current.

down the bank of the large trench toward the perfect resonance ratio Wr/wro, which again is

illustrated in figure 5 of Kirby (1988). As this minimum transmission coefficient is reached,

any more current would drive the transmission coefficient up the other side of the trench.

This will decrease the energy differential and reduce the current, allowing the system to fall

into the base of the trench again.

This assumes that there is enough available energy from the wave to drive the currents

necessary to reach the trench base, but clearly shows the process that the system seeks

an equilibrium value that is as close to the perfect resonance ratio as it can get. This

is the main point of this thesis, and Chapter Six will present results which are aimed at

investigating this hypothesis.


3.1 Introduction

This chapter addresses the governing equations that apply for the case of water waves

over mild topography, and the resulting mean currents and mean surface elevation. The

problem as stated must be simplified to be solved numerically, and therefore the equations

governing momentum balance in the domain between the mean values must be reduced to

two dimensions. Also, the wave equations that govern wave amplitude over the bathymetry

are presented as a parabolic equation. These implications inherently take away from the

full process, but without them the solution would be very difficult to achieve. First to be

presented in this chapter are the equations and terms that apply to the circulation. Then

the equations for the wave field will be presented.

3.2 Circulation Model

3.2.1 Governing Equations

The governing equations are presented in equations (3.1)-(3.2). They are the depth-

integrated and time averaged equations of motion, where time averaging is done over a wave

period. The derivation can be found in a number of references, most recently in Appendix

A of Winer (1988).

aU aU oU 0- 1 1
+ U- + V +g + b- s-
ft 9x 9y Dx pD pD

1 (dSId BS 1 drj
+ X +- - + 0 (3.1)
pD ozx dy p dy

OV OV dV Or 1 1
+ U + gV + + f y bi7s
9t 9x 9y 9y pD pD Y

+_1 + SYY) + o 0 (3.2)
pD \ Ox Oy px 0

The continuity equation completes the needed three equations for the three unknowns ?, U

and V.

+ (U D) + (V D) = (3.3)
Tt Ox By

These equations are simplified to the spatially two-dimensional case to render a realistically

solvable set of equations. The assumption does, however, take away from the complete

process which might include strong effects from vertically variable current fields. In these

equations, the terms are defined as

U = x component of mean current

V = y component of mean current

S= mean water surface elevation

p = vertically constant water density

ho = still water depth

D = total water depth = ho + 7

rl = lateral shear stress

Tbx = x component of bottom shear stress

Tby = y component of bottom shear stress

r,, = x component of surface shear stress

rsy = y component of surface shear stress


Sxx = x component of flux due to x-propagating wave

Sxy = y component of flux due to x-propagating wave

Syy = y component of flux due to y-propagating wave

Note that the equations will retain the nonlinear advective acceleration, gradients in the

radiation stresses, bottom friction, and lateral mixing. The numerical representation of

these terms are the same in this thesis as was presented in Winer (1988). Referring to that

dissertation, the bottom friction term is represented by equations (2.15) and (2.16), which

is a direct application of Longuet-Higgins (1970a). The lateral mixing term is in equations

(2.26) and (2.27), and also in Longuet-Higgins (1970b). The advective acceleration is rep-

resented in finite difference form in equations (4.24) and (4.25). The surface shear stress

is not included in the paper, but could be added easily enough. It was not included for

it remains out of the scope of the present problem dealing with energy differentials due to


3.2.2 Radiation Stress Terms

The original model by Winer (1988) did not address the reflected wave that is now

present. Since a reflected wave will be expected, the radiation stress terms must be altered

to include this wave to maintain a proper balance. Also, when considering waves on a

current, the incident and reflected wavenumbers will be different. If one looks at how the

waves interact with each other, there will be terms that will oscillate spatially at a spatial

frequency related to the difference of these wavenumbers. Therefore, deriving radiation

stresses would necessarily be involved. For demonstration of what types of terms are present

in a two-wave system at zero angle of propagation, consider the situation of zero current

and thus equal wavenumbers. To begin, the linear velocity potential for a two-wave system

is given by

(xytz,) ig coshk(ho + z) A(x,y) ei(kx-wt)
2 w cosh kho
ig coshk(h + z)B(, y) ei(-k-wt)+ +
2 w cosh kh,


and the instantaneous surface displacement 4' is given by

'(x,y,t) = 1 [A(x,y) eke + B(x,y)e-ikx] e-'t + *
,'(xY,t) = -


where denotes the complex conjugate. It has been assumed that the complex amplitudes

relate to a scaling parameter c defined as E = k lA, the steepness of the wave. The amplitudes

vary with x at O(e2), and in y at O(c). In keeping with the parabolic simplification to the

wave model, the wave angles are then assumed to be at or very near zero. The wave induced

velocities are then, by (3.6-8)

u'- 0 -

0 _
.-=^ =

w' =- z

g cosh k(h, + z) kA ei(kx-wt) i Axei(kx-wt)
2 w cosh kho I
+ kA*e-i(kx-wt) + i A e-i(kx-wt) kBei(-kx-wt)

iBxei(-k-t) kB*e-i(-kx-wt) + iBe-(-kx-w) (3.6)

g cosh k(ho + z) [i A, ei(kx-wt) + i A e-i(kx-wt)
- ~ e(2hkh ) i e(-x-)] (3.7)
- iBye i(-kx-wt) + i B* e-i(-kx-wt) (3.7)

i gk sinh k(ho + z) i(k-wt) +A* e-i(k-wt)
2 w cosh kho e-(-
B ei(-kx-wt) B* e-i(-kx-wt)

These wave induced velocities go in to form the radiation stress terms, along with the mean

surface elevation and wave surface displacement. The general expression for each of the

radiation stress terms are presented as a review from Mei (1983) as


s, = p i 2 dz + p g(- z) dz + (p /"-- d( dz -
J-ho J-ho -ho \ Jz X /

p /w'2 dz (h, + )2 + (77)2 (3.9)
J-ho 2 2

SXy = p u'v' dz (3.10)
3- ho

Syy = p v'2dz + pI g( z) dz + p d( dz -
.-ho -ho J-ho\ Jz O y
7 9 nP9
f w' dz P (ho + )2 + ) (3.11)
-ho 2 2

From these the derivation will now be completed for the two-wave system without currents.

First note that for any two complex number R and Q, that RQ* + R*Q = 2R[RQ*] and

RQ* R*Q = 2ij[RQ*]. Now consider each term in the radiation stress expressions


PfhoU'2 dz

p f 2dz = 2k 2|AI2 + 2k2BI2 1 4k2R [AB*e2ikx]

+ 21A 12 + 21Bx12 + 4 R[A B e2i"kx]

4 k[AA*] 4k [ABj e2ikr] 4 k [A B*e2ikx]

+ 4 k [B B*]} p (gz)2 dz (3.12)
2 w -h cosh kho

The integral in this expression has the solution,

S(cosh k(ho dz = 2 sinh 2kh, wCg
-ho cosh kho g Vk2 sinh 2kho gk


pf2ghog(i z)dz

p g(-- z) dz= -(h, + i7)2=p
-ho 2 2


P fh f'l7uW' d dz
Pjfho z ax --5--

p 7 o7- dC dz =
J-hoJz 1 OX

Pg { 2 IA2 + 2 IB|12 + 2 R[A Al] + 4 [A B e2ikx]

+ 2 Z[B BL] + 2 R[A BL e2ikx] 8k S[A Bd e2ik. ]

- 8 k A[Ax B* e2ikx] + 2 R[Axx B* e2ikx]

-82 [A ik sinh 2k(ho + ) d( dz
S8 -ho J sinh 2kho

The integral is evaluated as,

f f sinh 2k(ho + () 1
h ik-h + ( dC dz = [2kho coth 2kho 1]
-ho Jz sinh2kho 4k2



P f-h--2 dz


S 2dz = ( sinhk2kh ) { IA2 + IB12 + 2R[AB*e2ikx]} (3.17)

-P(ho + )2


Pg (/)2 = P.g { A2 + IB12 + 2,[AB*e2ikx]} (3.18)

Then writing the final expression for the radiation stress component Sxx,

4 (1 + sin h, k2A2 + k2 BI2 2k2 [AB*e2ik]

+ IA, 2 + iBx12 + 2R[AxxBe2ikx]
2 k S[A A*] 2 k .A[A B e2ikx] k Q[AB B* e2ikx] + 2 kQ [B B]}
+ -g(2khcoth2kh 1) { Ax12 + IB 12 + R[AA*x] + R[B B*]
+ [[A BZ e2ikx] 4 k S[A B e2ikx] 4 k b [Ax B* e2ikx]
+ R[AZ B* e2ik] 4 k2 R[A B* e2ikx] + 2 R[Ax B e2ik] }
sinh2kho) {IA2 + IB|2 + 2 R[AB*e2i]}

+ {A|12 + B|12 + 2R[AB*e2]ik} (3.19)

Now for the terms in the S,, expression, only the first applies.

P f hu' v' dz

and the expression for Sx follows as

sY = 4 1 sinh2k h) {[A A A] + [A, B; e2ik]

+ ~2[B Aye2ikx] + [Bx B] kQ [AA;]

k B[A Be2ikx] k [A B* e2ikx] + k3[B B*] (3.20)

For the Syy expression, which is similar to the Sxx term in the integral representations,
the final expression is

Syy = 1 + 2 k h ) {Ay2 + By2 + 2 [AB*e2ikx]
4k2 *sinh 2 k h,2
+ g(2khcoth2kh 1) {IAy\2 + IByI2 + 2l [AB*e2ik ] + R[A Ay]

+ ~[B B,] + R[AB*;e2ikx] + R[Ay B* e2ik]}
p sinh2k ho { A2 + IB12 + 2 R[AB*e2ik]}

+ {jAI2 + IB12 + 2 R[A B* e2ikx] (3.21)

The terms arising because of the two-wave system couple nicely when not considering
currents. After doing some preliminary testing of the model, it was found that the interac-
tion between the two waves in the radiation stress terms would give rise to small currents.
Since the only concern for the modelling effort was currents arising from differentials in
overall energy, these small disturbances did not matter. Thus it was decided to calculate
the wave field for a given current field and let the bar coupling extract the energy from

the incident wave, then set the reflected wave field to zero in the grid prior to calculating

the radiation stresses. This would not only remove the currents over the bar field due to

interplay of the partial wave, but also allow a simple form of the subroutine that included

the incident wavenumber. The equations are the same as (3.19)-(3.21), but contain only

incident A terms. The final expressions used in the model will be presented in Chapter


3.3 Wave Model

3.3.1 Governing Equation

The wave model emulates the wave equation presented in this section. It governs the

energy balance of a wave travelling over a mild bottom with fast undulations, and riding

on a current. The parent equations come from two sources, one being Kirby (1988) and

the other being a University of Florida Technical Report, Kirby (1986b). The first paper

developed coupled equations for waves on a current over a bar field on a flat bottom. The

other paper handled such a case, but on a mild slope. The forcing from Kirby (1988) is used

in the present thesis, and the main left hand side of Kirby (1986b) is also used to complete

the full equation. As a review, the coupled equations for the incident and reflected waves

from Kirby (1988) are the refraction approximation,

{TT, + (C., cosOl + Uo)Tx, + (Cg, sin0i + Vo)Ty,} = -c R (3.22)

2 {RT + (Cg cos02 + Uo)Rx, + (-C, sin2 + Vo)Ry} = fnT (3.23)

where the following terms are defined;

Qc = 'oc cos(01 + 2) -+ O1c,


o- = (3.25)
oc 4wcoshklhcoshk2h (3.25)

fl, = (4wa coshAh)-'I {AUo(12a1 lr2)Db

[gkik2 cos(81 + 2) + 02 ( 2 + (AUo)2)]Db} (3.26)
g g

a = AhF tanhAh
A = Barfield wavenumber =
Fx = Froude number =

Ii = kicosOi;

ai = intrinsic frequency = V/gki tanh kih

02 = 02
w = absolute frequency =
wave period
Db = Uniform bar amplitude

The amplitudes in these equations are surface displacement amplitudes divided by intrinsic

frequency, or T = A/la and R = B/a2. From the technical report, the parabolic coupled

equations are written as,

2ia(Cgi + Uo)T + 2oi(Cgi + Uo)(kl k1)T + 2iiTVTy

+ i {[l(Cg + Uo)]. + (aiV),} T + {(CCgl V2)(T),}

=coupling with R (3.27)

2iO2(C,2 U,)R, 2U2(Cg2 Uo)(k2 k2)R 2iC2VR,

+ i {[2(Cg2- Uo)] (02V)} R {(CCg2 V2)(R)y}

= coupling with T (3.28)

By dividing by 2i from (3.27-28), and multiplying w to the coupling from Kirby (1988), a

set of equations result for the case of a bar field on a mild slope with currents. These are

then solved using a Crank-Nicolson method for each wave until convergence is achieved. A

term w that handles the wave breaking energy dissipation is also included. This is added

to the incident wave only, since the reflected wave will not be expected to break anywhere.

This term is obtained from Dally et al. (1985). The concept is that if wave height exceeds

a criterion based on the local water depth, written as

2AI> (3.29)
then breaking commences, where K = 0.78. The term is then defined as

K [ 72h2
w = 2 1 A (3.30)

where 7 = 0.4 and K = 0.15 for the present paper. The full parabolic equations used in

the model are then

al(Cgi + U,)T. ioi(C,. + Uo)(ki k)T + olVT

+ {[al(C1 + Uo)]. + (aV),~ } T V)(T)y

+ -c7(Ci + Uo)T

= c R (3.31)

02(Cg2 Uo)R, + ia2 (Cg2 U,)(k2 k2)R 92VRY

+ [2(Cg2 U0)] (2V)y} R + f (CCg2 V2)(R)

= -w T (3.32)


This is recognized as the equations governing waves that are at perfect resonance with

the bars. If there were a slight shift Q, then the forcing would have a different form.


4.1 Introduction

The incorporation of the equations developed in Chapter Three involves finite differenc-

ing, decisions on initial and boundary conditions, and placement within a coherent algorithm

that performs efficiently. The latter requirement is best described by means of a flow chart.

So before moving to the details of the coding, an overall view of the model MCSHERRY

would be timely. Figure 4.1 provides the logic behind the main program.

After inputs of bathymetry and wave information, the model iterates between the wave

model WAVEMOD and circulation model CIRC. A hyperbolic tangent startup function is

multiplied to the initial offshore incident amplitude, which slowly increases it's value to the

final after a specified number of iterations. This is done so that the circulation model will

not be shocked by suddenly large forcing terms. The form of the startup function is

C(t) = 0.5[1 + tanh(t/cl c2)]

where the parameters cl and c2 are chosen to alter the offset and slope of the curve. For

the present case, cl = 3000 and c2 = 3 to render a curve as shown in figure 4.2.

The subroutine FLOOD allows the grid area to fluctuate depending in the mean surface

elevation. It is called at every iteration, and alters values at wetted or dried grids to conserve

mass flux in those grids. The subroutine CHECK simply checks for convergent values, and

sets a flag for the model to quit if the condition is satisfied at each grid point for U,V and

Figure 4.1: Flow chart for model MCSHERRY.






0 50 100 150 200
time (sec)

Figure 4.2: Hyperbolic tangent startup function.


4.2 Circulation Model

4.2.1 Method of Solution

The equations presented in equations (3.1)-(3.3) are solved using an implicit Alternating

Direction method that solves for the new U value and interim 7Y value in the x-sweep, and

the new 7 and new V in the y-sweep. This method applies to the matrix setup as suggested

by Sheng and Butler (1982). The detailed treatment of this method as used here can be

seen in Chapter 4 of Winer (1988).

The two unknowns in each sweep are applied at different places in the grid, where J

is at the grid centers and the mean velocities apply to the grid edges. Note that the first

offshore value is i7, at i = 1, and the first mean current U is at i = 2. The value of U and

V at i = 1 do not exist for the grid edge values, and are not used directly in the code. As

stated before, however, the wave model uses averages of the currents and for the i = 1 value

uses U, V at i = 2. The grid domain is shown in figure 4.3.

The x-sweep occurs in subroutine XCOEF, which sets up the coefficients, then passes

them to the double-sweep solver. At the lateral boundaries, the special subroutines XCOJ1

and XCOJN are called to handle the special coefficients when j = 1 and j = N. The

y-sweep is handled by subroutine YCOEF, except at the shoreline, where the variable grid

domain determined by FLOOD is deciphered by the subroutine YCVAR.

There is no time restriction for ADI method, but some of the terms are expressed

explicitly and added to the knowns of the implicit scheme, and restrictions are introduced.

The finite difference form of each of the terms in the momentum equations will now be


The radiation stress forcing is a function of the square of the amplitudes, thus being

nonlinear. Also, since the wave amplitudes are determined at grid centers in the wave

model, the radiation stresses also apply to grid centers.

The radiation stress gradients are put into the nth time level of the scheme, which

- offshore

1 4* 4


Figure 4.3: Grid for circulation model.


1* II


referring back to the last chapter for the x- and y-sweeps are

1 (OSxx Sy)
+ -Y
pD Ox Dy

pD Ox Dy

These terms are represented in finite difference form respectively by

1 SXXij SXXilj .SXYij+1 SXYij- + SXY + SXY j_
S'-"- "- .5 sXy,+ 2A -i-Ij-1

1 SYYij- SYYi, .SXY+j SXYi-_,j + SXYi+lj SXYl,j_
pD AY b 2AX

The symbols in the code for each variable are hopefully clear. The subroutine for the

radiation stresses is RADIAT, which is called just after the wave model and prior to the

subroutine CIRC.

Referring to figure 4.4, the terms are used to find the mean velocities, and are therefore

difference to straddle the known (i, j) velocity location. Note differencing is done over one

grid space, and is accurate to second order.

The bottom friction pFlUIU is nonlinear, but uses mean currents from unknown and

known time steps, and is thus represented implicitly. The terms are represented based on

the work of Longuet-Higgins (1970a) as

bx = pF luoblU + pF UIU

Tby = pF luoblV + pF IUIV

These are finite difference according to the sample x-component,

Figure 4.4: Differencing of radiation stress gradients in x-sweep.

= Unew 2fH +
x T sinhkh U2 + )

where the linear formulation for maximum orbital velocity has been used, T being the wave

period. Since the reflected wave is considered nonexistent when performing the radiation

stress subroutine, this also applies for this subroutine. Therefore, the linear formulation

of maximum orbital velocity for the one wave applies. The subroutine that calculates this

term is called SHEAR2, and is called just prior to CIRC.

The lateral mixing is expressed in the momentum equations as

p y y a U
P 9y Dy *Y 1y9

J-1 J J+1


I -- -- -- --

Il~l - -

D Orl a V\
p Ox Ox X x)

These are added to the knowns at the nth time level explicitly, and thus introduce a serious

timestep limitation. In the complete finite difference form they are expressed as

t AtEY U,j+1 2Ui, + Uij-1
Oy AY2

atr AtE +1j 2Uij + Ui-,j
Ox AX2

For the solution to remain in front of the fastest travelling error, the timestep must satisfy

< 1

These values are calculated in subroutine EXCOEF just prior to CIRC. Note that the

coefficient Ey is set to the highest value of Ex, which occurs at the breakerline.

The terms represented in the momentum equations by U and V2J in the x direction

are represented in finite difference form respectively as,

Uj Ui+l,j Ui-1,J
u'3 2AX


(Vi,j + Vi-,i + V-l J,+ + V1i-l,1+) U+l Ui-l,j
4 2AY

The same form applies to the y direction equations. A four point average is used for the

cross-derivative of the four surrounding velocities as shown in figure 4.5

Figure 4.5: Differencing of advective acceleration in x-sweep.

4.2.2 Boundary Conditions

The solution for the equations governing mean currents and mean surface elevation

require boundary conditions for each boundary in the grid. Winer (1988) used a fixed lid

condition at the offshore edge, with no flow into the shoreline, and free flow at the lateral

boundaries. Each condition had to be reviewed for the present case, and changes were made

where appropriate.

The original model constrained the mean elevation j to be zero at the offshore grid edge.

This was a fine assumption for the one wave system, but was found to cause instabilities

if a reflected wave occurred. The instabilities were due to a severe mean surface gradient

at the first two grid rows which gave rise to large velocities. For the standing wavefield, if

the grid edge happened to occur on the envelope antinode, the problem was accentuated.

Obviously a radiation condition was required to free that edge to allow disturbances within

J-1 j J+l


the domain to escape. The condition

+ V + 1 - = 0 (4.1)
ot Ox pD a x
was applied to the x-sweep of the model, where the interim 7 is solved for at the grid center

i = 1. After making the decision to not include the reflected wave, however, the original

fixed lid condition was again used and is what is specified now. If the model were to be

modified to include the reflected wave, an appropriate radiation condition would need to be

developed. Also, if a strong current were to be flowing across the offshore edge, an energy

set down condition might suffice.

The lateral conditions are reflective, that is the y-sweep sets V = 0 at j = 1 and j = N.

This is to remain consistent with the wave model.

There is no flow into or out of the boundary, so the x-sweep condition sets U = 0 here.

Also, there is no longshore current, so the y-sweep specifies that V = 0.

The model MCSHERRY begins by setting all unknowns to zero in subroutine INPUT.

Then, as the waveheights are ramped to their full value, the forcing slowly increases in the

circulation program. The model iterates until the difference between unknowns in successive

time steps is smaller that a prescribed tolerance.

4.2.3 Radiation Stress Subroutine

The finite differencing of the radiation stress components that were developed in Chap-

ter Three will now be presented. The calculation of these terms is performed in the sub-

routine RADIAT. This immediately follows the wave model and precedes the circulation

model. Since the wave model is only accurate to O(E2) in spatial derivatives, there are

terms based on derivatives in the radiation stresses that should be dropped as being too

small. Referring to the expressions for S.x, S.y and Syy from Chapter Three, the terms

that should be dropped are

Ox2 4

0 3
Ox Oy
The terms that remain are difference to remain at grid centers where the amplitudes are

calculated. As an example, the derivatives are expressed as

A Ai+l,j Ai-,l,j
A Aij+j Ai,j-1
Ay ,
A Aij+l 2Ai, + Ai,j-1
dyy -= y2 --

As was stated in Chapter Two, the reflected waves have been considered nonexistent, since

they do not physically affect the circulation landward of the bar field. To incorporate only

the differences in wave energy in the grid, and not wave-wave interaction effects, the final

expressions that exist in subroutine RADIAT therefore neglect terms consisting of reflected

wave B. The final forms used in the model are

+ pg [(A2)}] (4.2)(kh)'
SXYX 1+ (k2A 2) 2(AA)) (A A())
S4(k2)i sinh2(kho) i

Y pg ( 2(kh[ ) [( (A (l2)
4(l sinh2(kho)]
+ [(|A|2)i (4.2)

pg ( 2(kh))'
4(k2) sinh2(kho)

4 sinh2(kho)j


pg 2(ko)
4 sinh2(kh)
+ [(IA12i]

+ 8(k2) tanh 2(Ikh) A) [ 2+ R(A Ay,)] (4.4)

4.3 Wave Model

4.3.1 Method of Solution

The wave equations have been presented by (3.31) and (3.32). The finite differencing

follows a Crank-Nicolson method, where a tridiagonal matrix is formed for the three forward-

row amplitudes based on three known values at the present row. The solution requires only

an initial condition at the offshore edge, and lateral conditions to complete the matrix.

For the two waves, two loops make up the subroutine WAVEMOD. The first loop solves

for the incident wave A up to the shoreline. Then, the second loop starts at the shoreline

and solves for B out to the offshore edge. Due to the coupling through the bottom boundary
condition, the two sweeps must be iterated until a convergent solution is obtained.

4.3.2 Finite Differencing

Now the incident equation will be difference to serve as an example. Again, the

amplitudes T and R refer to amplitudes divided by intrinsic frequency. Using a Crank-

Nicolson method leads to the finite-difference approximation of the wave equation of

1 ([ai(Ci9 + IU)]+1 + ai(C91 + U)]} (Tj+1 J}
2AX j Tj
{[(Cgi + U)(ki k)]t+1 + [air(Cg, + U)(ki k)]}} {Tj+1 + Tj}

+ 1 [r(Cgi + U+1 r1(Ci + U)] T +1 + I'
+ A yI(' + 7(C [j+, ;_,
4AX '- I- J


+ 1Y {[V]i+l [rlV]-J + [rlV]}+l [ClV]_1I {Ti+1 + T}
16AY 31 )j ][+ I T;1]
2 {[(CC, V2)- vji (CCgl V2 +1 +1 [ +1
8Ay2 3 +1 j
-[(ccgl V2)+l + (CC,_ V2)+][+1 -_ T

+ [(c V2)j1 + (CCig V2)1[T1 T]

-[(CC, V2)i + (CC9 V2) -l][T T_-i]

+ 1 {[ia(Cg ) ] + U)wT] [al(C91 + U)wT]}}

fI{[wR]' + [QRi }
2 3

The reflected form is similar, with sign changes where appropriate. These are then put into

three coefficients for the i and i +1 rows, and put into a complex tridiagonal matrix double-
sweep solver. Since currents exists at grid edges, the values to place in the WAVEMOD

need to be averaged straddling the grid center. Thus Uj is really ((Uj+' + UJ), and likewise

for V. The grid domain is set so that the first calculated value of U is at i = 2, so then the

averaged value for velocities are set to U2 and Vj2 when marching from i = 1 to i = 2.

In WAVEMOD, a preamble defines the coded meanings of each variable. Then state-

ment functions are defined for the coefficients in the equations for the two waves. Before the

first loop, subroutine WAVNUM is called to determine the incident and reflected wavenum-

bers, average velocities, celebrities, group velocities, intrinsic frequencies, and forcing terms

Qoc and 1,c as defined by (3.25) and (3.26). After this is done, the iterations begin, with
a maximum of four iterations, but convergence is usually complete within 0.1% after three.

If the wave is perceived to meet the breaking criteria and start dissipating energy at some
grid point, then w is defined something other than zero, and the row is iterated again,

letting the energy out. The wave model uses depths of the uniform bottom, not including

the barfield. Therefore, shoaling over the bars is not included in the driving forces.

4.3.3 Boundary Conditions

The boundary conditions are that the lateral boundaries behave as walls. This requires

that no reflected wave disturbance reach the walls and bounce back into the domain, which

would be unrealistic.

4.3.4 Lateral Smoothing

When doing some preliminary testing on the wave model, the radiation stress gradients

that would be put into the circulation model would be plotted simply for observation.

It was found that the plots exhibited much irregularity in the driving forces, more than

was expected. In fact, the circulation model had stability problems, and the irregularities

were thought to be at fault. After some review, suggestion was made that the parabolic

approximation had altered the nature of wave properties in the transverse y-direction.

Physically, a wave can have a maximum wavelength for given frequency and depth. This

is uniform in all directions. But when assuming small incident angles, thus primarially x-

direction propogation, which is done to get from a hyperbolic equation to a feasibly solvable

parabolic equation, a modification is introduced. Consider a surface described as

7' = Aeikx (4.5)

where A is the complex wave amplitude. For a = IA|, the surface in the y- direction is

77 = aeiky (4.6)

If one considers this wave in the y-direction, the radiation stresses are proportional to

the wavenumbers, with the Sxy term being proportional to k1, and the term S,, being

proportional to k2. And then when considering the gradients of these terms as they appear

in the momentum equations, the term would be proportional to k3. Therefore, there

is a large sensitivity to the wavenumber in the transverse direction. By equating (4.5) and

(4.6), the amplitude A could be expressed as

A = aeie-ik = (y)e-ik (4.7)

Now consider a general parabolic approximation equation of form

2ikAx + Ayy = 0 (4.8)

and by making the substitution for A, the ODE is

Ay, + 2k2A = 0, (4.9)

which has the solution

A = aei\Vky (4.10)

Therefore, a limit has been placed on the maximum physically meaningful wavenumber in a

parabolic wave equation for the y-direction. Instead of the normal L, = 2, the new limit

is L, = '-- Any wavenumber in the transverse direction that is larger than this will have

severe effects on the driving forces in the transverse direction.

By showing this, consider a numerical grid with lateral width W, and some calculated

wave field signal along a transect in the y-direction, as shown in figure 4.6.

A(y) can be written as a Fourier series

A(y) = ~anein (4.11)
where A = J is the base spatial frequency. By taking the FFT of the signal, the spectrum

will show energy distribution over each component nA. For the parabolic model, however,

any energy that exists beyond frequency nA = VZ2k is not physically meaningful.

Now for the wavefields calculated by the numerical model WAVEMOD, a transect in y

was FFT'd, and was found to have a very small amount of energy beyond what was physical.

This is shown in figure 4.7, where the cutoff frequency component was nA = 9.051. Then a

y=O y



Figure 4.6: Wave calculated from parabolic equation on domain of width W.


transect of the driving forces ( + a ) were likewise decomposed spectrally, and a large

amount of energy existed beyond what was allowed. Shown in figure 4.8 for the y-direction,

this shows how a small amount of error in amplitudes can expressly alter what the circulation

model sees. The leakage of miniscule energy into high transverse wavenumbers in the

wavefield translated into major amounts of energy beyond the cutoff frequency component

for the driving forces. The process of the radiation stress gradients being proportional to

k3 shows the drawback when making a parabolic approximation. After taking the wave

field and doing the FFT row by row, and truncating the energy beyond nA = 9.051 and

inversing, a smoothed wave field was obtained. Figure 4.9 is the spectrum of the y-direction

driving after smoothing. By figures 4.10 and 4.11, the slight improvements are seen in the

wavefield, and in figures 4.12 through 4.13 the difference is remarkable. The variance in the

energy spectrum of the y-direction driving before the smoothing was calculated as 0.00804,

while after smoothing the total variance was 0.0030. In comparison to this 60% reduction

in variance, the wave field had an unsmoothed variance of 4.065 x 10-7 and a smoothed

variance of 3.80 x 10-7, about a 6% reduction.






5 10 15
frequency component

Figure 4.7: Spectrum of A(y) without smoothing. Cutoff at nA = 9.051.






20 40 60
frequency component

Figure 4.8: Spectrum of y-direction driving force before smoothing.





0 10 20 30
frequency component
Figure 4.9: Spectrum of y-direction driving force after smoothing.


Figure 4.10: IAI without smoothing. Bars are located on the bottom of this domain.

Figure 4.11: JA[ after smoothing. Bars are located on the bottom of this domain.

Figure 4.12: y-direction forcing without smoothing.

Figure 4.13: y-direction forcing after smoothing.


5.1 Introduction

The laboratory work was done to get a better idea of the process of Bragg reflection,

where effects not represented by the wave model were shown. The following chapter treats

the laboratory phase of this project specifically. Since the work was completed at two

separate times, the procedures and setup unique to each will be presented separately. The

results from both sessions, however, will be presented together, with comparisons to a

simple numerical model.

The work was done at the Coastal and Ocean Engineering Laboratory (COEL) of the

University of Florida. To begin the discussion, the basin and the other necessary equipment

will be described.

5.2 Equipment

5.2.1 Basin

The modelling basin is illustrated in figure 5.1. The internal walls are moveable, and

were positioned at right angles to the wave crests to minimize side wall reflection. The

basin can be filled to a depth of 60 centimeters, and allows for waves up to heights of

several centimeters. The main testing region lies inside the internal wave vanes, to minimize

intrusion of backscatter from the outer basin walls. A wavemaker sits at one end, with a

sloping beach at the other.

Figure 5.1: Plan view of large basin.

5.2.2 Wavemaker

The wavemaker has 80 paddles, each being 23.5 centimeters in width. Each paddle can

be adjusted for stroke and phase. The entire system has a frequency range of 0.5 Hz to 1.25

Hz. Since the beach was oriented at 10 degrees relative to the wavemaker, the paddles were

adjusted to generate an oblique wave at a 10 degree angle, which then arrives normally at

the beach. This is done by calculating the wavenumber for the wavefield over the desired

depth using a Newton-Raphson technique. Then the relation of paddle wavenumber kp to

the wavefield wavenumber vector kw is, by (5.1)

kp = k, sinO (5.1)

Then the paddle wavelength is Lp = -. For paddle width Wp, and paddle number np, the

phase of each paddle in degrees from 0 to 360 is, by (5.2)

p = L (5.2)

5.2.3 Beach and Cart System

The beach gently slopes to a maximum height of 60 centimeters above the floor. A

structure has been placed over the beach to aid in various types of testing. It allows a cart

to ride on tracks and be placed at any (x, y) coordinate in an 8 meter by 17 meter area on

the beach. For a plan view of this track system, see figure 5.2.

The track sits several centimeters above the greatest water elevation allowable in the

basin. The cart is positioned in the y-direction by the wire connected to the handcrank at

the side of the basin. The variable speed motor fitted onto the cross-member for the cart

controls the position of the cart in the x-direction.

The cart can carry wave gages, current meters, and whatever else is necessary in an

experiment. See figure 5.3 for a side view of the cart. Since profiling is important in the

conduct of moveable bed studies, the cart is fitted with a profiling system. The profiler arm

Figure 5.2: Plan view of cart track system.

Figure 5.3: Side view of profiling cart.

extends underneath the cart, dragging over the beach surface. A potentiometer produces

a signal which is converted to readings of arm angle, while another fitted to the trailing

wheel measures in voltage how many revolutions that wheel has made. With these two

instruments, the beach (x, z) position can be traced digitally. A sample profile is depicted

in figure 5.4.

5.2.4 Electronics

The computer for the basin is a PDP-11 microcomputer. It receives all basin data

which is amplified at a box next to the basin. A remote terminal at the basin runs all of the

programs. The main program is ATOD.FOR which takes multiple channels simultaneously

at any rate for any length of time. The PDP- 11 is located in an instrument room about

200 feet from the basin.

Wave gages were the main source of data. The capacitance gages were built at the


60 11 -TT I I I


N 20

0 j

0 200 400 600 800

x (cm)
Figure 5.4: Sample profile.

COEL. Figure 5.5 is a sample of the calibration points transposed on the calculated second
order curve. The response appears linear over the operational range. Calibration was done
daily, or if the temperature changed drastically. Along with the gages, a strip chart took
the signal from the cart gage in session one. It was used to trace wave envelopes from the
moving gage.

A Marsh-McBirney current meter took data for the second session. Since the calibration
was 10ft/sec = 1 Volt, and since the currents were small, an amplification of twenty
provided the necessary range to increase resolution. The data were quite erratic, and will
be addressed in a following section.
Lastly, a video camera was mounted above the beach. The camera was fitted with a re-
mote control to zoom, focus, point, and contrast. A taping device was included. Rhodamine
dye transport was recorded as a way to show current magnitude. The only drawback was
that the camera is black and white, and the dye does not show up as well as it could if the



-) 4



-2 I I I I I I I

1000 1500 2000 2500
Figure 5.5: Sample calibration curve for wave gage.

camera filmed in color.

5.2.5 Bar Field

The bars were made of wood and metal. The same bars served both sessions, and are
shown in figure 5.6 without the metal covering. They were sized according to the local
water depths at their position of placement. The planform dimensions were 3.048 meters
long and 0.318 meters wide, and each bar was made to be 40% of the local water depth.
Once the height h was determined by inspection, the bar radius in centimeters followed
equation (5.3).

(19.05)2 h2
(19.05)2 + h = R (5.3)


R 3.048

Figure 5.6: Bar frame design.

The full pattern was traced on paper, then transferred to the plywood frames and con-
structed. Since the main material was wood, the bars would float. They were anchored in
the first session using wires connected to the concrete floor. In the second session, the un-
dersides were filled with sand which was kept in place with a panel covering. Both methods
worked fine, but the first proved to be better in time. Another problem was wave orbital
velocities flowing under the bars, scouring large holes. Strips of sheet metal 4 inches wide
were extended off the bar edges as aprons sunk into the sand. This effectively stopped the

5.3 Setup, Procedure and Data Analysis

This next section discusses the actual testing methods. Since different methods were
used in each session, each will be discussed individually. The first session will be addressed,
with it's purposes and procedures, followed by the second session.


0 LL

h /\ \
, V \


5.3.1 Session One

The first session ran from September of 1987 to February of 1988. Since this was the

first work done in the lab on the currents attributed to Bragg reflection, the goals were

fairly loose. The purposes were specifically:

1. To measure Bragg resonant reflection.

2. To compare this reflection to theory.

3. To describe qualitatively the properties of the wave-induced current field.

4. To observe how reflection and current magnititude varied as wave period varied.

5. To observe sand movement during the experiment.

The first goal received most of the attention, with two methods being developed. The

first was to use two stationary gages, in line, and take data simultaneously. The data from

each then can be used to determine the reflection coefficient, and this method is from the

work of Goda, et. al. (1976). To review that method, for a wavefield which can be described

by a surface displacement

= a cos(kx wt + 61) + b cos(-kx wt + 62) (5.4)

the unknowns are a, b, 61, and 62. Four values for y from the data will therefore be needed.

Let two gages be located at x = 0 and x = 1, as in the figure 5.7.

If a data value is taken for 17 from each gage at times t = 0 and t = r, four values will

result. Each one can be expressed through (5.4) according to (5.5-8) as

71 = r(0,0) = a cos 61 + b cos 62 (5.5)

72 = (1, 0) = a cos(kl + 61) + b cos(-kl + 62) (5.6)

773 = 7(0, r) = a cos(61 wr) + b cos(62 wr) (5.7)

774 = 7(1, r) = a cos(61 + kl wr) + b cos(62 kl wr) (5.8)

Figure 5.7: Gage locations for two-gage reflection method.

Let the following relations hold;

c, = cos61 (5.9)

C2 = coS62 (5.10)

S1 = sin 6i (5.11)

s2 = sin b2 (5.12)

d = coskl (5.13)

e = sinkl (5.14)

f = cos wr (5.15)

g = sinwr (5.16)

Expanding out (5.5-8) and substituting (5.9)-(5.16) gives

?71 = a cl+ b 2 (5.17)

12 = (aci + b 2)d- (as b s2)e (5.18)

773 = (a c + bc2)f + (as + bs2)g (5.19)

7q4 = (aci -+bc2)f d + (bs2 asi)fe

+ (a s1 + b 82)g d + (a ci b c2)g e (5.20)

By substituting sequentially from (5.17) down to (5.20) and eliminating terms, the remaining

expressions are,

2acl = A = 171i+ 4 72 1fd (5.21)

2as = B = (5.22)
g e

2bc2=C = q1 4 713 712 + d (5.23)
2 b C2 = 7 ?i 1 --- ---- (5.23)

2b2 = D = 73 71f +2- d (5.24)
g e

From these expressions, the incident and reflected wave amplitudes are, by (5.25) and (5.26),

A2 + B2
a A2 B2 (5.25)

a C2 + 2
b C2 D2 (5.26)

The reflection coefficient is then Kr = b/a. The problem in this method is that if g or

e approach zero, the expressions for A, B, C and D blow up. Thus one must not choose

r = n( ) or 1 = n( ) where T = wave period, L = wavelength, n = integer.

To determine the accuracy of the method, a test was run where the reflection was from

the beach. To compare, a moving gage retrieved the partial standing wave envelope. It was

discovered that for various time increments r, the reflection coefficients were inconsistent.

This inconsistency was enough to scrap the two gage system, and use the envelopes as

Table 5.1: Bar heights and positions for Session One.
bar no. height (cm) x-position
1 7.62 420.
2 9.50 495.
3 11.20 570.
4 12.60 645.

the major source of reflection data. This required measuring by hand the minimum and

maximum envelope amplitudes, 7?min and r1max. Then the reflection coefficient is,

Kr = i]max r7min (5.27)
?7max + l7min

Envelope transects across the barfield were taken at four locations during the testing to get

a lateral average. The positions of these transects are presented later in this section. The

results are then presented and compared to the circulation model in the next session.

A camera was mounted above the testing site which filmed the movement of rhodamine

dye across the barfield at each wave period. The film is available for viewing through the


The tests were run at several wave frequencies around the designed resonant frequency.

Each run was made from a still water condition. With reflection and current observations

made at each run, a feel of the dependence was gained. The different runs are listed in the

next section for the periods and reflection coefficients.

Observations were made in regard to sand movement during the course of the testing.

These observations are listed as part of the results in Chapter Six.

The setup consists of where the bars were oriented. By using the equation for the bars

from the local depths, the bars were built according to Table 5.1.

Plots of the testing bathymetry were done as testing commensed, and are available as

of this writing, but the actual data that was written onto magnetic tape is not immediately




N o20

0 2 4 6 8
x (m)

Figure 5.8: Reconstructed bar profile for Session One.

available to the author. There is no other recourse than to reconstruct what the bathymetry
was at the time from the existing plots. This is the only handicap from losing the data,
since the envelopes are printed on strip chart paper. Figure 5.8 shows the lab bathymetry
with the bars placed. The depth was 55 centimeters, and the breakerline in the testing was
at about x = 3 meters. Figure 5.9 gives an overhead view of the barfield, with four lines in
the x-direction offshore which represent the transect lines where envelope data was taken.
The profile is available as a test bathymetry in the numerical model, although the bars need
to be sinusoidal.
Once the bars were placed, the wave gages were arranged. The barfield would create
a reflected wave that would propagate and diffract once beyond the bars. Two gages were
mounted on tripods in the center of the bar field longitudinally, and directly offshore of the
last bar. As a control, one gage was put off to the side, away from the influence of the
barfield. Another gage was on the moveable cart. Figure 5.10 shows a plan of the beach

1365 1300 1175 1110

16.0 16.0 ----16.0-

24.0 "0-4.0

32.0 0 __- 32.0 ---- 32.0

n f____

Figure 5.9: Plan view of barfield for Session One.


Figure 5.10: Gage locations for Session One.

region with the gage placements.

This concludes the discourse on session one. The data from the experiment will be

presented later, with the test run parameters, in Chapter Six. Now the second session will

be addressed with the procedures and setup.

5.3.2 Session Two

With the experience from the first session, the second experiment ran from September

of 1988 to December of the same year. The goals were more defined in this second session.

1. To repeat Bragg resonant reflection and a resultant current.

2. To extend the reflection to frequencies off the no-current resonant value.

3. To observe differences between moving away from resonance and shifting toward


Figure 5.11: Three gages inline with x-direction waves

4. To quantify current values.

In this session, two methods were again employed. The first was from the work of

Funke and Mansard (1980). Here, three in line wave gages took simultaneous data, and

the reflection coefficient came out of the spectra. The work covered waves travelling over a

flat bottom, with no current involved. Their results therefore needed to be revised for the

addition of current and variable topography. For the setup in figure 5.11, the algorithm was

modified in the following manner;

The wave profile for any probe p = 1,2,3 can be described as a Fourier series sum of

each component k by

N / 2rkt \
qp(t) = Ap,k sin + ap,k (5.28)


where Ap,k is the Fourier coefficient for frequency -, T is the time length of the data run,

ap,k is a random phase shift, and N is the desired number of components for the profile.
The Fourier coefficient is gotten from the transform of the time series at each probe p, rp(t),

which for any probe would be, by (5.29)

7 {lp(t)}= Bp,k = CI,kexp i 2(Xlk +iOk

+ CR,k exp {i 27(X1 2XR1 X1P) ik + O(k)}
I Lk
+ Yp,k exp{i(pp,k)} (5.29)

These equations are actually (1) and (10) in Funke and Mansard (1980). They point out

that the transforms represent functions of the complex amplitudes for the incident and

reflected waves. The three expressions for the FFT's of the three probes then lead to a

least squares determination of the incident and reflected amplitudes, Zi and ZR.

B1,k = ZI,k + ZR,k + ZN,1,k (5.30)

B2,k = Ks,I,2 ZI,k exp[qI,i,2] + Ks,R,2 ZR,k exp[--R,1,2] + ZN,2,k (5.31)

B3,k = Ks,I,3 ZI,k exp[~I,1,3] + Ks,R,3 ZR,k exp[-'Ra,,3] + ZN,3,k (5.32)

The terms in these three equations are defined as follows;

ZI,k = CI,kexp i27r +L iOk (5.33)
SLk .2X1+2XR1)
ZR,k = CR, exp {i 27(X1 Lk X +i(O, k)} (5.34)

ZN,p,k = Yp,kexp {iPp,k} (5.35)

The shoaling coefficients in the incident and reflected directions for each probe, and group

velocity Cg are,

Ks,I,2 = (5.36)

K,3 = (5.37)
V ^.2
sR,2 = (5.38)
V C9

Ks,,3 = Cg, (5.39)

The phase accumulations are obtained according to

PI,1,2 = j kidx (5.40)

,1,3 = kldx (5.41)

R,1,2 = k dx (5.42)

RR,1,3 = 3 kR dx (5.43)

These modifications reflect the different case of waves over a variable topography on a

current. Funke and Mansard (1980) addressed waves over a flat bottom where the incident

and reflected wavenumbers were equal, and wavelength Lk for each frequency was constant

in the domain. The reflective and incident wavenumbers k are found from the dispersion

relation that includes the current effect, using an iterative technique. Setting the ZN terms

to a small error parameter E, the least squares method is followed and the resulting equations

for unknowns ZI and ZR are

Zi {2 + 2K',I,2exp[2ii,1,2] + 2K',I,3exp[2i~P,i,3]}


+ ZR {2 -: /s,,,2AKs,R,2exp[i(II,1,2 R,1,2)] + 2Ks,I,3KS,,R,3exp[i(JI,1,3 T'R,1,3)]}

= 2B1 + 2Ks,I,2exp[i'@J,1,2]B2 + 2Ks,i,3exp[i',1 ,3]B45.45)

and then

ZI {2 + 2Ks,I,2Ks,R,2exp[i(PI,1,2 TR,1,2)] + 2K,I,,3Ks,,R,3exp[i('I,1,3 @R,1,3)]}

+ ZR 2 + 2I,R,2exp[-2iQR,1,2] + 2K ,R,3exp[-2iR,1,3]}

= 2B1 + 2Ks,R,2exp[-i'R,1,2]B2 + 2K',I,3 exp[-iQR,1,3]B3 (5.46)

These can be solved for the two complex amplitudes, Zi and ZR, and then Kr = ~IZR

The modifications were done numerically according to the following forms. For the

dispersion relation with opposing current,

w = -Uo k + V/gk tanh kh (5.47)

a Newton-Raphson iterative technique was used to solve for k, given depth, frequency, and

current Uo. It is assumed here that the incident wave angle 0 = 0, so that k cos 0 = k.

The simple integration scheme of the phase in space is performed by,

Wp = (kp + k~p-1)A -'-1 + Xp-1 (5.48)

where again p = the gage number, either 1, 2, or 3. Ax p-.p-1 is the interval size between

gages p and p 1. The intervals were on the order of 15 centimeters, so that the bottom

could be considered a plane sloping beach in these intervals, allowing the simple integration.

The group velocity used for the shoaling coefficient K, is, for each gage p,

Cgp P (1 + sinh2(h)p (5.49)

This method was tested for accuracy for reflection from the beach, being compared to

envelope data. It was found that the reflections differed between the two methods by as


much as 60%, and the development of the equations thus leaves something in error. Rather

than an in-depth search for an accurate form, it was again decided to use the envelope

method, as in Session One. Again, partial standing wave envelopes were retrieved from

a moving gage across the barfield, but were recorded digitally, and analyzed as per the

methods outlined below. The envelopes were taken at seven transects over the barfield. It

was found that the three center transects represented the highest consistent reflection for

a test. These locations are shown in the plan view in figure 5.16. After determining the

reflection from each transect, it was plotted with location in figure 5.15. As is apparent,

sections 850, 900, and 950 give a consistent representation of reflection. The method of

analyzation was to take each envelope, as shown in figure 5.12, and condense it into envelope

height versus x-position. A sample is shown in figure 5.13. The reflection coefficient was

calculated by taking a local maximum and immediately following minimum, using equation

(5.27). Thus a plot of Kr with x was obtained, as seen in figure 5.14. The maximum

reflection was taken from each envelope, and the three were averaged. This was done for

each wave period. The results are presented in the next chapter and show a repeatability of

the process, which is a good sign. Further similarity with theoretical output gives confidence

in the values themselves. From figure 5.14, there is a variability in measured reflection along

the bars, and the ending minimum should be disregarded since the cart may not have traced

the full wave envelope before stopping. The variability within the bar region, however, could

be attributed to focusing effects over the bar field.

This data was taken using program ATOD.FOR for a sampling frequency of 26.0485

IIz, and a run of 2048 points. The gage was the cart mounted gage, moved slowly through

the water over the barfield, insuring that the wave excursions stayed in the middle portion

of the gage wire.

The second priority was to get quantified results for current induced by Bragg reflection.

At each period, a line of current data was taken across the front edge of the barfield. The

locations were at the same locations used for the envelope measurements, except they were








400 600


x (cm)

Figure 5.12: Sample envelope, Session Two. Retrieved digitally from one of the three
on-offshore transects using the moving cart gage.















Figure 5.14:






I I1

x (cm)
Sample Kr with offshore distance x from Session Two.

400 600
x (cm)
gure 5.13: Sample envelope height versus X.




. I

I ...

' I I I


0.5 -





700 800 900

1000 1100

section number
Figure 5.15: Plot of K, with transect location, Session Two.

16.0 6 16.0 16.0

2 4.8 0_ ---10

Figure 5.16: Plan view of barfield, showing data collection locations.

stationary just inshore of the first bar. Figure 5.16 is a plan view of the session two barfield

with the locations used for current meter measurements and envelope data. A side view is

shown in figure 5.18, where the water depth was 50 centimeters, and the breakerline was at

x = 320 centimeters. The current meter was placed at 750, 800, 850, 900, 950, 1000, and


Data was taken at 26.0485 Hz for 512 points. The data proved very erratic from the two-

channel Marsh-McBirney current meter. The calibration for the instrument was 10 ft/sec =


1 Volt, and since the expected currents would be on the order of a few centimeters per

second, an amplification was done. An amplifier was built in house to magnify the signal

by a factor of twenty. The raw results from such a time series is given in figure 5.17, which

represents the velocity magnitude. Because of the nature of the signal, any trustworthiness

of discharge calculations would be in doubt. The irregular signal could be due to interference

from the metallic bars that were about 10 centimeters away from the probe. Still water

readings taken in a container away from the bars showed a fairly stable signal. Other reasons

could be excessive low frequency water displacements, but nothing like that was observed.

Coupled with a serious time limitation as of this writing, the choice was made to disregard

laboratory quantification of the rip current. A major goal of Session Two therefore is yet

to be realized, but the reflection results weakly imply the presence of a current, and the

dye experiments from Session One that are on film definitely show the rip current. Indeed,

observations made at the time were that the currents were very mild, and did not get past

the first bar.

Unlike session one, session two was a continuous test between wave period changes. The

intention was to extend the reflective peak as long as possible, thus keeping the induced cir-

culation from previous periods going through the enhanced reflection. Once it was observed

that circulation had dropped considerably, the runs began at rest from a high nonresonant

period, and moved at slow increments toward resonance. Thus a test was made to see if

reflection commenced upon return at the same place where it slacked when moving away

from resonance. As said, the current magnitudes were so small that it was difficult to see

a definite drop in reflection or current strength. Nevertheless, the results are presented in

Chapter Six and still offer some useful lessons for future study.

The bars were the same from session one, but some of the bars were buried slightly to

accommodate the new conditions. Table 5.2 lists the heights and positions of the bars in this







5 10 15
time (sec)

Figure 5.17: Sample time series of magnitude of current meter data, Session Two.

For an overhead view of the barfield, refer back to figure 5.16. A side view of the
barfield is shown in figure 5.18.

The instrument setup consisted of everything being put on the cart. The three gages
were positioned inline on the cart, along with the current meter. The overhead of the

Table 5.2: Session Two bar heights and positions
Sbar no. height (cm) x-position (cm)
1 6.25 436
2 7.00 500
3 11.00 565
4 12.50 631




N 20

0 200 400 600 800
x (cm)
Figure 5.18: Profile for Session Two.

instrument array is shown in figure 5.19. With the points made earlier, the three-gage
system was dropped in favor of the one moving gage to retrieve envelopes across the barfield
in the x-direction.
This concludes the description of Session Two. The next chapter describes the results
obtained from the laboratory work.

Figure 5.19: Instrument setup for Session Two.


6.1 Introduction

The results from the work done in the past year and a half are presented in this chapter.

The results are from laboratory data, where the methods that were used to obtain this data

are outlined in Chapter Five. The other source of results is the numerical modeling effort.

This consists of the wave- current model MCSIIERRY that as of this writing handles a

limited set of conditions. The preliminary findings could be construed as inconclusive, but

definitely suggest the need for more attention. The model should be running at full speed

in the next year (1989).

6.2 Laboratory Work

The data collected consisted of reflection data from wave gages, current meter data,

bathymetric data from the profiler, and some data from the camera. The results will be

presented in tabular form and pictorally where possible.

6.2.1 Reflection Results

The Session One reflection data is presented in table 6.1. Each run represents a Static

test, i.e., startup from a still basin. The averaged data over the middle transects is also

plotted in figure 6.1, transposed upon the theoretical prediction for a one dimensional

barfield without currents. This one-dimensional model is from Kirby (1987) that excludes

currents. This has been placed over the data to simply give an idea of the effect the

currents have on the data. The incident wave had a resonant period of 0.98 seconds, and

was 7 centimeters in height. The waves arrived normally onto the beach, with a breakerline

at x = 3 meters, referring to the profile in figure 5.8.



Table 6.1: K, and comments from Session One.
period (sec) 1110 1175 1300 1365 comments on rip
0.980 0.357 0.667 0.529 0.362 weak rip current
0.910 0.347 0.189 0.235 0.091 no current
0.963 0.172 0.579 0.368 0.238 no current
1.050 0.368 0.550 0.368 0.538 rip current
1.110 0.294 0.692 0.458 0.514 significant current.
1.180 0.333 0.529 0.300 0.176 no rip current
1.110 0.333 0.474 0.647 0.428 again strong rip repeated.
1.085 0.290 0.412 0.579 0.375 slightly weaker than T=1.ll.
1.020 0.474 0.438 0.400 0.405 weak rip, dye is stagnant.
1.000 0.333 0.474 0.538 0.333 weak current
0.942 0.273 0.444 0.350 0.170 no current

The predicted reflection from the one-dimensional model underestimates the reflection

for the given bathymetry, for the envelope data describes 50% reflection. The data were not

carried to sufficiently large periods so that the reflection clearly died away. From the films

of the current, the magnitude increases in strength until a period of 1.11 seconds where it is

quite strong, after which it dies off at 1.18 seconds. This would tend to suggest that while

the waves are moving increasingly off no-current resonance, the rip current increases to

keep the waves at resonance. This supports the hypothesis that the system seeks a resonant


From Session Two, the data were more abundant, but perhaps not as telling. Much

of the setup was the same as in Session One, except that the depths over the bar field

were less. This is due possibly to the state of the beach when testing began, in that the

slope was less than in Session One. To offset this the bars were buried into the sand about

10% of their heights, which was typically about a centimeter. Also, due to the testing that

had been going on before this work, some of the sand had been removed from the main

testing region, so water depths could not be as high. In Session One, the testing depth was

55 centimeters, and for Session Two the depth was 50 centimeters. To further offset this,




0.5 1.0 1.5 2.0
period (sec)

Figure 6.1: Plot of IK with wave period, Session One data and theoretical prediction from
1-D model Kirby (1987).


wave heights were decreased from 7 centimeters offshore to 5.5 centimeters offshore. This

seems to be the main value that has the greatest effect on the results, for a rip current never

developed as strongly as the one appearing in Session One. The reflections were high, but

it is conjectured that there was not enough available energy to drive the circulation.

The reflection data are presented in three tables, where the first set refers to reflection

with absolutely no discernible current throughout the course of the runs. The next set was

made after digging sand out from the region in between the bars and the land. This allowed

a weak current to form in some of the periods. The third set refers to runs that started at

low frequencies and slowly approached the no-current resonant case. In all the testing in

Session Two, the runs were continuous, without stopping the wavemaker.

The averaged data from the first set with no current are tabulated in table 6.2, and

illustrated in figure 6.2, which has the theoretical result transposed. Again, the theoretical

values are from the no-current, one- dimensional model of Kirby (1987). Three or less runs

were done at each period to get an idea of the repeatability of the experiment results.

Like the first type, the second type started at resonance and moved slowly outward. A

weak rip appeared after removing some sand from the foreshore. The tabulated results are

given in table 6.3, with the plot in figure 6.3. The two profiles from set one and two are

shown in figure 6.4.

From these plots, the data follows the theoretical fairly well, but the data drops away

at later periods. The sand did not move onto the first bar at all, whereas in Session One,

the first bar had to be repeatably uncovered.

The third type started off resonance, then moved toward resonance continuously until

reflection commensed and a weak current was excited. Table 6.4 gives the results, also

plotted in figure 6.5.

The difference between this test and the other is again not clearly evident, and may be

due to the lack of dominance of the current.

The one-dimensional algorithm from Kirby (1987) is over a variable mild topography.

Table 6.2: Kr. Reflections for case of no discernable rip current, Session Two.
Period (sec)__ Kr_
1.014 0.419 0.515 0.423
1.020 0.381 0.418
1.030 0.433 0.503
1.050 0.507 0.471
1.065 0.511 0.508
1.080 0.530 0.505
1.100 0.494 0.526 0.525
1.125 0.516 0.481
1.150 0.470 0.454 0.473
1.175 0.518 0.493
1.200 0.425 0.386
1.225 0.323 0.341
1.300 0.385 0.310
1.350 0.374 0.414
1.400 0.396 0.372
1.450 0.386 0.367
1.500 0.352 0.339

0.6 I I


1.0 1.5 2.0
0 0

period (sec)

Figure 6.2: Plot of K,- with period, data and theoretical for the case of no discernable rip
1.0 1.5 2.0

Table 6.3: Kr measurements from case of a discernable shore rip current. Session Two.
period (sec) i K,
1.014 0.425 0.464
1.030 0.391 0.442
1.050 0.521 0.536 0.539
1.080 0.466 0.510 0.483
1.100 0.483 0.503 0.514
1.125 0.482 0.521 0.482
1.135 0.447 0.486
1.150 0.532 0.520
1.160 0.461 0.490
1.200 0.402
1.250 0.349
1.300 0.349
1.350 0.231 0.289
1.400 0.215
1.450 0.257
1.500 0.259
1.550 0.323
1.600 0.310
1.650 0.272_


o o_

0.4 0



Figure 6.3: Plot of K, with period for the

1.5 2.0
period (sec)

case of a discernable shore rip current. Data and

1 I I 1

1 I I I I I



x (m)

Figure 6.4: Two profiles, dotted line referring to profile used in the second and third sets
of testing, where a weak current formed.

Table 6.4: Kr measurements for case of discernable shore rip current, moving toward reso-
period (sec) IKr
1.650 0.216
1.600 0.330
1.550 0.323
1.500 0.316
1.450 0.316
1.400 0.356
1.350 0.249
1.300 0.302
1.200 0.436 0.460
1.175 0.464 0.459
1.160 0.518 0.512
1.150 0.529 0.520
1.135 0.487 0.556







I1 I 1 1



0.6 I



1.0 1.5 2.0
period (sec)

Figure 6.5: Plot of K, with period for case of discernable shore rip current, moving toward


It again excludes current effects. The plots do not show the kind of reflection extension

that could occur because the currents were not as strong in the second session as the first.

Thus, only so much tuning was possible. Some notes about the comparisons of the first and

second sessions are worth noting to explain some problems with the execution.

1. The waves were slightly higher in the first session, thus involving more current

driving energy.

2. The anchoring system kept the bars very still in the first session, while there was

slight movement in the second.

3. Not enough time was allowed between period shifts in the second test to allow a

proper current maintainance.

The first reason seems the main cause of the decreased current generation. As is stated in

the theory, if a wave field does not have enough available energy, the rip current will not

be strong enough to achieve tuning of the nonresonant incident waves. Thus the efficiency

of the bar field is comparable to the case with no current, and affects only closely resonant

waves. For practical considerations, the most damaging waves will probably carry enough

energy to drive a sufficient current, and the bar field will be effective when most needed.

6.2.2 Bathymetry

This last section deals with the movement of sand during the experiments. There are

profile data on backup tape from Session One that is presently not available to the author.

Observations conclude that the current was sufficient to move sand over the barfield. In

fact, it was necessary to uncover the first bar periodically, which had been buried and

which lessened the effectiveness of the barfield. In the second Session, the current was too

weak to transport sand, and the profile shown in figure 5.18 was the same throughout the



6.3 Numerical Results

After much tinkering, some results are available from the numerical model. As of the

writing, however, the model requires another two months of work. Specifically, lateral

mixing was left out because the time step limitations would make the model very slow, too

slow for the given time allotment. Therefore there are circulation variations in the domain

that would probably not exist if the cells were allowed to influence on another. However, it

is thought that the overall pattern is similar for both cases.

The results presented now are from the no-current wave model coupled with the cir-

culation model. The process to achieve results was to run the wave model once and get

the final radiation stresses. These terms were then multiplied by the startup function and

inserted into the circulation model. The tests were run on a common bathymetry which is

shown in figure 6.6.

Without the wave-current equations to govern the wave field, the thesis cannot answer

the question about whether the system tends toward resonance through numerical results.

However, the no-current model can answer questions about how the magnitude of the current

is affected by wave height. Therefore, three wave inputs were chosen for the bathymetry,

where the period was the no- current resonant period. As a point of interest, the wave field

and mean surface elevation that are calculated by the wave model for an input of a .03

meter wave are shown in figures 6.7-6.9.

Then the following plots show how the current pattern varies as the amplitude is in-

creased. Velocity times total depth is plotted in the figures, which describes the discharge

across a vertical line. It is seen as expected that the rip current discharge increases as

wave height increases. The circular cell over the barfield ends is thought to be attributed

to a slow leaking of energy from the incident wave over the bar field. The action stated

in Chapters Three and Four of removing the reflected wave probably account for this lost

energy. The reflection coefficients remain the same for each wave height, being calculated

simply by the maximum reflected wave amplitude in the grid divided by the initial incident



-o. --


0 .3 i i 1 1 1 1 1 1 1 1 1 1 1 1
0 2 4 6 8 10
x (m)

Figure 6.6: Bathymetry 1. Shoreward direction with increasing x.


Figure 6.7: jA[ for initial amplitude of .03 meters and perfectly resonant waves.



Figure 6.8: IBI for initial amplitude of .03 meters and perfectly resonant waves, without

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

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