Wave-averaged-energy transport, transformation to roler energy and generation of wave forces

. Wave-averaged-energy transport,

transformation to

rol~r

energy and

generation of wave forces

report of implementation and testing. March 1997

M3007.30 W1009

-Chapter 1. 1.1 1.2 1.3 1.4 Chapter 2. 2.1 2.2 2.3 2.4 Appendix

Transport of wave-averaged short-wave energy Introduction

Analysis Implementation Tests

Roller energy transport Introduction

Analysis Implementation Tests

-Chapter 1.

Transport of wave-averaged short-wave energy

1.1. Introduction

In this text we will discuss the actions that were taken to solve the 2D surf beat equations which prescribe the forcing of waves through dissipation of short waves. The equations were solved by using and adapting the already available subroutines in TRISULA which are used there to solve the transport of pollutants. The variations in short wave energy which result from the computations can be used to determine the radiation stress tensor. The divergence of this tensor is the wave force which can be used to drive the flow in the TRISULA model. The advantage of this approach is that determining the wave forces and the driving of the flow is detemlined in one model.

1.2. Analysis

The following equation, which describes the transport and dissipation of short-wave averaged energy, is to be solved in order to find the time dependent energy distribution in the domain.

(1)

here Cg,x and Cg,y are the cartesian components of the group velocity. The sink term in the right

hand side of the equation can be detemlined by: ( _( ,J8Eol (pg)

)n)

D=2af _p l - e yh E _{0 (2)}

Of the three empirical coefficients in this equation y and n are related to the probability of breaking, a is related to the rate of dissipation of the breaking waves.

The other two parameters p and g represent the mass density of the water and the gravitational acceleration respectively.

By using the HISW A model the components of the propagation vector C g for the short-wave averaged energy can be approximated. Note that, although the boundary conditions which prescribe the incoming energy can be non-stationary, the wave direction 8 is approximated as fixed in time as this is assumed in the HISW A model.

After each time step the cartesian components of the radiation stress tensor can be determined as

S (=5 )=(n -2+n cos2_8)E

11 xx c 2 c 0

-where: n = Cg

₌

~(1

_{+ kh}

c C 2

propagate.

1 - tanh

2

(kh) ) and

e

is the direction in which the waves tanh (kh)

From these the wave forces components are determined as minus the divergence of the stress tensor:

_ ) _ aSll aS₁₂ F (-F - - - _ - _ _ 1 x

ax

ay

_ ) _ aS

12

aS

22 F (-F - - _ _ - _ _ 2 Y

ax

ay

1.3. Implementation

In curvilinear coordinates the equations transform to:

aEo

+ 1

(~(Cl

fGE

)+~(C2

fGE ))=-D

at

JG

G a~1 gV~22 0 a~2 gV~ll 0

11 22

where

c~

and

C:

are the physical components of the group velocity. Since the grid is orthogonal the contra- and co-variant physical components are identical.

Since the TRISULA grid is orthogonal, the non-trivial components of the metric tensor are given by:

The sink term becomes:

h

(li8

Eo

7 (

p g)

)n

were Z=

yh

If we assume that the changes in the depth are relatively small, the expression can be linearized as:

-In the TRlSULA routine DIFU the sink and source terms just before lable 849 be determined as:

sink: - dD = -af (2 + (-2 + nz)

e-

Z )

dE_o p

where h, which is present in z, is determined as:h = hOld= SO (nm) + DPS (nm)

D dn E - f -zE

sour: - + dE O,old-(X pnze old

°

The terms volnew and volold were replaced by gsqs(*).

In curvilinear coordinates the (physical contravariant) components of the wave force can be written as:

Here Sail are the cartesian components of the stress tensor. Using the properties for 2 dimensional transformations:

the components can be written as:

d~2 __ 1 dX2

ax1 - fG a~1 '

de _

1 ax1

ax 2 -

IG

a~1 '

-va;; [

aSl1 ax 2 ax2 aS12 ax2 ax1 dS12 ax1 ax2

G -

a~1 a~2 a~1

+

a~1 a~2 a~1

+

d~1 a~2 a~1

-aS22 ax1 ax1 dSl1 ( dx2

)2

OS/2 dX? dX1 OS22 ( dx 1 )2]

a~

1

a~

2 dP +

O~

2

a~

1 ) - 2

O~

2

a~

1

ap

+

d~

2

a~

1

The discretization of (the components of) these expressions becomes: for equation (3) which is to be discretized in a u-velocity point:

- 5

-(3)

GXi _{. .}

- - ~X~ -X~

G~2 m,n m,n-l

GXi

1(

i i i i )

a~l ~ 4 xm+l,n+ xm+1,n-l - xm-1,n- xm-l,n-l

for equation (4) which is to be discretized in a v-velocity point:

GXi i i

--~x -x 1

a~l m,n m- ,n

Several subroutines were introduced:

QKWCG

QKWBN WAVEU DIFUWE RADSTR

detennines the physical components of both C and C_g on the computational

grid.

extrapolates the compouenis of C and C_g to the boundary.

detennines the source and sink tenns of the transport equation (1)

detennines the transport of the wave energy; this is a modified version of the subroutine DIFU which is used for the transport of pollutant.

detennines the wave force components as the divergence of the radiation stress tensor.

1.4.

Tests using analytic solutions

In order to check the implementation in the TRISULA code, several test computations were made. As no analytic solutions are known for the transport equation with the nonlinear sink tenn, in some

of the tests use is made of a simplified right hand side which is linear in Eo.

test 1

As a first test we check the propagation in the x-direction of the energy Eo modelled on a

-rectangular grid. The propagation velocity C_gx is taken constant.

aEo aEo

- - + C - - = 0

_at

gx

ax

At the up-wave boundary we send in the signal:

(5)

The analytic solution of equation (1) is:

results of test

1

The computation was performed on a grid with dx =dy =2 m. The length of the area was 200 m the with was 20 m. At the locations x =20 dx, 60 dx, 100 dx a time series of the energy is determined. The initial condition for the energy is 2 on the entire grid.

The value of the period in the HISW A computation was chosen 2.5616 s on a constant depth of

10 m. The resulting group velocity is then C_g = 2

m/ s.

The direction of the group velocity was in

the length of the grid. The time step used in the TRISULA computation was IJ.. t == 1 . 2 s resulting

C_gll t _

in a CFL number of I1x -1.2.

The results of the computation are shown in figure 1. Here we see that it takes the signal 20, 60, and 100 s to arrive at the three stations. The amplitude of the signal has not decreased significantly while travelling over the 1 wavelength distance to station 1. This also applies to the 3 wavelengths to station 2. In the third station at a distance of 5 wavelengths from the boundary where the signal enters we find a reduction of the wave height.

-test 2

Here we solve the equation:

aEo aEo

- - + C - - = a E

_at

g

ax

0

which has the general solution:

If the signal at the up-wave boundary is given by: g ( t)

the total solution becomes:

-C_gln4

the choice a;

=

100 Ax forces the amplitude of the solution to be reduced with a factor .25 over

the distance of 100 f1x.

E (100 Ax t) ==Jc ( t - 100 f1x)

o ' 4 g C g

At x=50L\xwefind: Eo (sOflx, t) ==Jcg(t- Salix)

2 C

g

results of test 2

The computation was made on the same grid as used for test 1. The values of the time step and of Cg we also identical. The incoming signal is given in equation (5). At x=O, 50 dx and at 100 dx time series were generated. As shown in the analysis the maximum values of the periodic signals at the three stations should be 3.0, 1.5 and 0.75 respectively. The results of the computation is shown in figure 2. We find for the averages of the maximum values of the signals:

2.98 at x=O 1.47 at x=50 dx 0.71 at x= 100 dx

The 1800

phase shift of the periodic part of the signal at 50 dx is caused by the fact that the station

is located at 2.5 wavelengths from the boundary.

-test 3

The same computation as in test 2 but now on a curvilinear grid

results of test 3

The curvilinear grid used in this computation is shown in figure 3a. The results of the computation are shown in figure 3b.

-x

o til If) 0 1./)

l{) "T f'l r'1 {\J

A <J

Curvi-linear grid

test 4

At the boundary we send in the signal as used in test 1.

Since the E distribution is only dependent of x and t the wave force has only an x component:

F (x, t) - - (2 -- Cq - - ) -1 aEo _ - - (2 -Cq - - ) 1 2 1t cos ( 2 1t (x-C t) )

x C 2 ax C 2 2

a

.1x 2

a

Ax g

results of test 4

In Figures 4a and 4b the results of computations on the equidistant linear grid and on the curvilinear grid are shown. Notice that the results in Figure 4b are given along the curved lines where T) is

constant.

-,

<D N D

0 0

X~ ...

<38----Wave forces as function of x

at a rect-angular grid De 1ft Hydraul ics '" 'T 0 0 0 , on "' r-0 0 0 0 CJ H 23-'12 o o

j'"

1-",0 0 Fig ~" x

---II II 1/ <l>"'<l> </ I

==

I

=

I ill to If) '" '" 0 0 a X~ <J

Wave forces as funct ion

at a curvi- linea r- gr id Del f t Hyclr-aul ic

s

N '" " i l l 0 a 0 0 0 0 0 0 0 of x H 2342 " 0 a F I g o a (\) tn (\) 0 a

j"'

crF' 0 0 4b X

r

test 5

This test corresponds with the test "The effect of wave induced forces in a closed container. Comparison with an analytical solution" from the Test document DELFT2D-MOR.

A description of this test is given in the appendix.

The wave forces were read from a file in that test; here they will result from the energy input at the boundary.

At the left boundary we send in :

Eo(x=O,t) =A+pt

After the wave has propagated into the domain the solution for the energy will be:

E (x t)

=A-~

(x-C t)

0 ' C g

g

For the resulting wave forces we find: ( ) _ Cg 1)

P

F x t - ( 2 - - -

-x ' C 2 C

g

By choosing ~ = pgCg 10-3 _{we will get the wave forces which were used in the closed}

2 Cg

_l.

C 2

container test.

In that test the wave forces were provided from file, here they are determined by the program.

Since the group velocity will change significantly if the depth varies from 1 m to

f2

m in the container, we will take the depth at the left boundary to be 4 m and the required group velocity equal to Cg= 1 m/s. This leads to a period ofT= 1.280975 sand Cg/C=1I2 and kb=9.8 indicating that the water is deep. The group velocity is not sensitive to small changes in the water level

aC

( ah

g = -5 *10-7 s -1). For these parameters we find for the water level at the right boundary (at 500 m) the value h (x_{R )}

=/17 =4

.123m. The initial water level should be equal to 4.06186 m but given the above arguments a difference in water levels of 12.3 em should be found for all initial water levels close to 4 m.

results of test 5

The computation was made on an equidistant grid with LlX= 10 m. A time step of Llt=3 s was used. The group velocity determined by the program showed to be 1.01 mls instead of the intended value of 1 m/s. This resulted in a flatter profile of the energy distribution in space. The wave forces

derived from this distribution were as a result of this somewhat too small. The reSUlting stationary position of the free surface shows a difference in water level of 12.15 em. Figure 5 shows both the analytic and the numerical solution as a function of x.

-0 If) 0

~ u,

0 0

0 0 0

Water level elevation due to

using a recti- I inear grid in

Oe I f t Hyclrau I I CS ,"j -u ro - v wave '" 0 0 force a closed c H a 0 2342 o CO If) o o o o 1£)0 ~ 0 0 Fig 5

test

6

In this test we generate a free standing wave in the closed container by sending in an oscillating energy signal at the left boundary of the container. At this boundary two waves, both travelling to the right, will emerge, a bound wave travelling with velocity Cg and a free wave with velocity

.J(gh). At the boundary the sum of these waves will satisfy the boundary condition u=O. The free

wave will reflect at boundaries the energy signal will not. This leads to the generation by the energy signal of a free wave at the right boundary (at x=L) as well since there too the condition u=O needs to be met.

The parameters were chosen such that, based on linear theory, a gradually growing standing wave should emerge.

On the horizontal bottom in the closed container the shallow water equations can for low waves be approximated by:

all

+ h

au

= 0

at

ax

au

all _

1

as=.

+ g

-at

ax

ph

ax

(6) (7)

where Sxx is the only relevant component of the radiation stress tensor, u is the depth averaged horizontal velocity, h is the (constant) depth, g the gravitational acceleration and p is the mass density of the water.

S =(2

C

g

_-.!)E(X-C

t)

xx C 2 g

where we made the choice:

By substitution of l1b (x - C_g t) and u_b (x-C_g t) for the bound wave and requiring for the time averages fJ b = 0 and U b = 0 we find:

1 S -S

fJ (x-C t)

= -_

xx xx

b g P gh-C~

With the choice we made for the energy density we find:

The free wave, generated at the left boundary, is found from equation (6) and the condition u=O at this boundary, yielding:

C

Yj (x=O) = - - g -'11 (x=O).

f 1 r;:;-r::. • I b

yg11

-From this we fmd the free wave to be:

c

c

11 (x-.;gEt)=--g-n (-g-(x-/gEt».

fl

&Ii

".b

&E

For the free wave generated at the right boundary (x=L) we find, again using (6) and the condition u=O:

C

n (x=L) = - g -11 (x=L)

"I f2 .;gJi b

This wave becomes:

The free waves reflect at the boundaries x=O and x=L.

For the length of the container we chose L = 500 m. The depth was chosen to be h = 10m. For the

gravitational acceleration g=lO mls2 . C_g=2.5 m1s , To=100 s , C_g/C=0.5025119 , p=1000

kg/m3 and A=20000 J/m2 .

results of test 6

The computations were made on a grid with dx= 10 m and a time step of dt=3 s.

In the figures 6a, 6b and 6c the results of the computation are shown as time series at x = 10 m, x=240 m and x=490 m respectively. The analytical solution which holds for waves with a small wave height is given here as well. In the first of the figures we see the bound wave and the generated free wave. After 50 s the free wave generated at the left boundary is back at the left boundary after being reflected at the right boundary. As the free waves keep reflecting at the left and the right boundary and the length of the container is half the wavelength of the free waves, we find that the signal keeps growing. At the right boundary we see that it takes 50 s for the free wave generated at the left boundary to arrive at this boundary. After 200 s the bound wave has arrived here as well. At x=240 m the analytical solution shows no growth, here we only find the bound wave. The numerical solution deviates from the analytical solution after some time since the analytical solution does not satisfy the nonlinear equations solved in TRISULA. The analytical solution can be used to check the numerical solution only for a short time, as long as the wave height remains small.

-- -- -- -- analytica! - - - - - numerical 1.0 a 5 ·0 5 - I 0 o 200. ~oo.

Free SUI-fa c e elevation

Comparison of analytical and n ume ric a I solution

Del f

t

HycJl"dUI I cs I \ \ I \ I \ I \ I \ I \ I \ I 600. { I \ \ I ~ 1\ 1\ I II I \ \ I \ \ I \ \ I \

\ I \

i

J I 800. test 6 H 23-42 x 10 m Fig. 6a

- - - - analytical - - - - - numerical 1.0 O.S -0 5 - 1.0 o. 200. 400

Free surface elevation

Comparison of analytical and n ume ric a I 50 I uti 0 n

Delft Hydraul cs

600. 800.

test 6 X 240 m

- - - - analytical - - - - - numer icai 1.0 o 5 0.0 I \ I - 0 \ 1 o. 200. 100.

FI-ee surface elevation

Comparison of analytical and numerical solution Delft Hyclraul cs I \ I \ I \ I \1 600. I I I \ I I \ I _{\ I} \ I _{\ J} \ I _{\ I} \ I _{\ I} I 800. test 6 x = 490 m H 23"12 Fig 6c

Chapter 2.

_{Roller energy transport}

2.1

Introduction

In shallow regions modelling the reduction and transport of wave energy using Equation (1) results in too much reduction of the energy. This is caused by the fact that wave energy is in shallow regions transformed into roller energy and the reduction of the total energy is postponed to even less deep areas. The transformation from short-wave-averaged energy into roller energy is modeled through the introduction of a second equation which prescribes the transport of the roller energy with a velocity of the phase.

2.2

Analysis

The equations which are used to model the energy transport are:

This is a copy of equation (1).

The right hand side of this equation is given by Equation (2) in the first chapter. The equation for the roller energy is given by:

BEr

a

a " _

p

- + - ( 2 C E)+-(2C _{a t ax} ~ )-D---E

x r By Y r C r

(8)

(9)

Here the Cx and Cy are the cartesian components of the wave celerity. The value C in the right

hand side is the modulus of the wave celerity.

P

is an empirical parameter.

The boundary of the model has to be chosen outside the shallow region so that the boundary condition for the roller energy can be fixed to zero there.

Note that the reduction of energy per unit of time and per unit of surface which is given in the right hand side of equation (8) is the source term of equation (9) so the energy lost through breaking is added to the roller energy. Travelling to the coast the sink term in equation (9) wiII rapidly grow with the reducing depth as the celerity will reduce to zero at zero depth.

In order to determine the radiation stress tensor both the contributions of the wave energy and the roller energy have to be included. The cartesian components now become:

-The expressions for the cartesian components of the wave force again become:

as

as

F (=F)=-~--E 1 x

ax

ay

_ ) _ aS

12

aS

22 F (F -2 Y

ax

ay

2.3

bnplementation

The route followed to solve equation (9) is similar to that of solving equation (8) which is described in section 1.3. The prescribed components for the transport velocity now become equal to the components of the celerity. The source term, which is treated explicitly, is determined from the value of E. The sink term which is treated implicitly now contains the locally determined value for

C.

Two extra subroutines were introduced: ROLU

DIFURO

detemlines the source and sink terms of the transport equation (9)

determines the transport of the roller energy; this is a modified version of the subroutine DIFU which is used for the transport of pollutant.

2.4

Tests

In the following tests the sink term in equation (8) is simplified to D= -

a

Eo and the source term in equation (9) to A Eo' This effectively means that not all energy lost in the breaking process contributes to roller energy provided that A

<

IX •

On a horizontal bottom both group velocity and phase velocity become constant. For this simple case it is straightforward to fmd analytical solutions of the differential equations which describe the propagation of the organized wave energy and the roller energy.

Test

7, stationary solution

In this test we solve the stationary problem and choose A = IX :

a

ox(CgEo) +aEo =0 (10) and

a

p

Er -(CE)+

_ax

--('1.E =0 r C 0 (11)

where Eo represents the organized wave energy and Er the roller energy.

-As boundary conditions we use at x=O:

Eo ( 0) = EO and (12)

(13)

In equation (10) it is assumed that the wave height is larger than y h so that

a

becomes a constant.

From the first equation we find:

Multiplied by a this becomes the forcing term in equation (11). The solution of equation (11) now becomes: ( u

1

~_ a etc x a E (0) C (0) x

£

C2(T) Cg(-c)

-J

~du

o g l e du e oC2(u) c(x) 0 Cg(U)

In the case of a constant depth both C_g and C are constant and the solution becomes: a

- - x

Eo (x) =Eo (0) e cg ,

If

1

JL -

~

I'"

0 this expression can be replaced by Er (x) -<- 0.:

P

Eo (0) x 2 e

.~

x.

C2 _C

g C3

For the test we choose the parameters as follows:

(14)

(15)

(16)

h=S.O (m) ,g=10.0 (m/s2 ) , T=3..rc (s) =1.256637061 (s) from which we find:

5

C=2 .0 (m/ s) and Cg=l. 0 (m/ s) .

a=O.008 (S-l), P=0.128 (m/s 2) and Eo(O) =6.0 (J/m2).

The expressions for equations (15) and (16) now simplifY to:

4x 4x 16x

Eo (x) = 6

e

500 and Er (x) = e -500 - e 500

-result of test

7

In the TRlSULA computation use was made of an equidistant grid with ~ x =: Ll y == 10 (m) . Using

a time step of D.

t

= 6 (s) the stationary solution was determined.

In figure 7 and 8 we show the comparison of the numerical and analytical results.

6.0, 4.0 2.0 analytlca ilumer ica/

.O~--L---L-

_ _

L-~L-~

_ _

~

_ _

~-=~~~==~

o. 50. 100. 150. 200. 250. 300. 350. 400. 450 500.

Figure 7, EO stationary solution

.50 . 10 X analytical numericai OO/L---~----~---L--__ ~ ____ L -_ _ ~ _ _ _ _ - L _ _ _ _ -L ____ ~ __ ~ o. 50. 100 150. 200. 250. 300. 350. 400 450. 500. >::

Figure 8, Er stationary solution

-test

8,

Instationary problem

The differential equation which prescribes the development of the organized wave energy is here simplified as: aEo aEo - - +C - - =-aE

at

g

ax

0 (17) at x=o (18) Eo (x, t:) =0 at

t=

a

(19)

The equation which describes the roller energy is given by:

aEr aEr _

p

~

- - + C - - - - - E + J\ E

at

ax

C r 0 (20)

at x=o (21)

Er(x,t:) =0 at

t=o

(22)

To solve equations (17), (18) and (19) for x > 0 and t> 0, we introduce the following transformation:

1:=t:-~ , ~=~

C_g C_g

Equation (17) now transforms to:

oEo

-~ _=-aEo

The tilde indicates that the function is considered with respect to the new coordinates. We have:

Eo (-c J

0

=

Eo ( t:, x) .

The general solution of the equation is:

Transforming back to the originai coordinates yieids: a

- - x

Eo ( t: I x) = f l ( t: - ~) e Cg

C_g

In order to satisfy equations (18) and (19) we must have:

1: == t-

~

,

~

==

~

Equation (20) now transforms to:

- 28

aE

r _

P _

~­

- - - E +II.E

a~ C r 0

The tilde refers to the latest defined coordinate system. For the general solution we find:

E, ("

<>

+.

If, (

..

(1-

~) ~)

e

(~-" ,~), d~

+

f,

h)

1

e~'

after transforming back and applying the initial and boundary conditions (21) and (22) we find:

x

-Lx C (J!.-u-.f.)f!

Er(t,x)=)..e c2 Jf1(t-~+(1-~)!l)e C cg _dll

o C Cg This solution can also be written as:

Er (t, x)

p cg _-aC

C ( t - x)

C cg C

In our test the function f1 is prescribed as:

{

o

; t<o

f (t) =

1 A(l-e-

qt ) ;

t~O From equation (23) we now find:

o

Eo ( t, x)

(

-q(t-:!»)

-~x A 1 -e g e. g t -3.. C

J

f1 (V) e f:-~ _{, C} g

From (24) we find the solution for the roller energy:

- 29

-c-r -g

o

Er (t, x) = result of test 8 -q(t-~) e c ~e C (a ~-f!.. p~q) C+qc C2 _g t:O,~ C

With an interval of 96 s and the choice of parameters: A=l. 0 (J/m2), C=2. 0 (m/ s) ,

g=10.0 (m/82) , T=1.2566371 (s) and h=5.0 (m) both the result of numerical computations

and the analytical solution is shown in figures 9 to 16 .

-1.0.---, . B Eo

'\

.. -<;.,)2,-1 t::95 S ____ onol~tLC{]l __ ... numari..col

'\

°OL.--~5~O-. ~,~OO-:-.--,c!,c-O.--~?O"'O-. ~2~"~.~J~OO-:-.---,Jc!,~O.--,~0~0-. ~45"'O-. ---c-'~OO.

Figure 9, Eo at t=96 s Figure 11, Eo at t= 192 s X ",-0.02,-1 t=-192 5 >'-{).02 .. - 1 P-iJ.Ot • ,._2 '1-0 •01 ~ I _ _ _ anol~t\.col .. '" __ ... n'J""""rLcoL t=-288 5 >.-Q.Ol s:' ,-<l.Ot .. s·2 q· ... {I.OI.-1 _ _ _ or;al9l~col ... __ ... nU"'Bri..col .0 OL.---C5~O-. --c,-OOJ"C".--,c!,c-O=""'O"'O-. -o'~50-:-.--C3-00-:-0.----'Jc!'-oeO-. -4~O~0-. ~"~O-. ---c-':sOO.

Figure 13, Eo at t=288 s , 0 , , - - - ,

s\

Eo \

,

\\

'r

! l=384 5 ),·0.02 ;.-' p-ll.Ot .. 0;-2 q-O.O] ~-, ____ onoL!jlLcol ... nu ... ert.coL . 0 ~:-.----;5~O-. --;'~OO;;-.--c!15-:;-O.:c. ==C, OC;;O=. ~?~50;-. ---',"'0""0 .--3'"'"'0-. -':, 0""-. --:C"~O:-. -c!:500.

X Figure 15, Eo at t=384 s _ _ oniJl~ti.col '0 i~ . ____ nur,:;ri.col x Figure 10, Er at t=96 s . 3 0 , - - - , 50. 100. 150, 200 250. Figure 12, Er at t= 192 s Figure 14, El at t=288 s Figure 16, Er at t=384 s ,.-a.ro.-) t= 192 s l.-o.ru.~'" 300. 350. x 1-(1.01 .... _1 -100. -150. ,-(J.OI ~ ~..} q-<l,OItl _ _ onoL;.,ti.coL ... nU"'eri.caL x SOD,

The computation for which the results are shown here were made on an equidistant grid with grid

SIzes /1x= f1y

=

10m. For the time step we chose /1 t

=

3 s so the CFL numbers for the transport

of the organized wave energy Eo is 0.3 and for the roller energy it is 0.6.

On a curvilinear orthogonal grid with a length of 500 m and a width of 100 m the same computation was made. The time step was taken the same as in the former computation. The grid is shown in Figure 16.

Figure 16, Curvilinear computational grid

The results of the computation are again shown as a function of x. It should be noted however that the results were taken from a line of constant 1} , namely the sixth grid line.

I.O,---~ .B 6' Figure 17, EO at t=96 s.

'l

'0

'f\

x l:: 192 S ... -<).02:,.--1 p-o.Ot " ,,-J: q-O.OI.-· ... ...:1.02,·1 ).,-{i.oz ... •• jI'-{\.OI .. $_2 '1-0.0),,-1 2 4t \ \ _ _ onoL~ti.coL

I

---

n""erccol .0

~"

I " ,

Q. 50. 100. lJO. 200. 250. 3(10. 350 ~oo. -150. :sao.

y Figure 19, Eo at t= 192 s. - 32 -'1"'11.&t,-1

",

"f

_I co-<J.02.-1 t"""96 s p-o.OI _.-"' Er

I

tOf(\

j ,

D. ::10.

:~

100. ;5(1. 200. '250 301}. ]50. 400 150. 500 _____ nU!leri.O:JL x Figure 18, Er at t=96 s. 0 0 , -.20 Er

(~'\

..

,/

.~

i/

~

"I

O. 51). 100 150

"~

200. 2:>0. 300 350. ~~.OI-,-1 _ _ anol~t..lcol __ • __ nv ... ;;rt.col ~o-o. 150. !, 00. Figure 20, Er at t= 192 s.

·2 100. 150. 200. 250. Figure 21, Eo at t=288 s. l.:::384 5 .. .-Q.02,·1 1--0,1)2.,-1 1'I-<I.Ot • ~_2 ,,0.01,,-1 _ _ onol~l~col _____ nvmerl.col 300. 350. iDa. +-ISO. X SOD. 1.0,.r---______________________________ -,

'\

Eo

\\

.~

~Q. 100 150. 200. 250. 2

J

• J. Figure 23, Eo at t=384 s. t::-3B4 s ",-{j.02 ~-I \-{).O2: .. - 1 6-I).Ot • ... -l q-O.{)J",-l _ _ onolyl~col _____ (\I,H"Br i..coL 300. 350. 400 450 . X .500. .}iJ r - - - -__________ _ sao. X Figure 22, Er at t=288 s. . JO , - - - -___________________________ _ .. o{).(Q.-1 x Figure 24, Er at t=384 s. 33

APPENDIX

The effect of wave induced forces in a closed container.

Comparison with an analytical solution.

Contents

1 Introduction

2 Claim to be tested

3 Set-up

3.1 Analytical solution of the wave force effect in a closed container

3.2 Grid for flow computation

3.3 Control model

4 Results

5 Conclusions

-1

Introduction

In order to get an impression of the effect of wave forces on the behaviour of the system, an analytical solution was determined for the simple case where other wave effects are discarded. Comparison of the analytical solution to the result of computations gives insight in the accuracy of the system and how this is effected by grid sizes.

2

Claim to be tested

CLAIM HHHHHHHHm .m.m.m

The effect of wave induced forces in a closed container.

The system is capable to accurately determine the stationary free surface which is caused by a constant wave force. All other wave generated effects are discarded in this experiment so that an analytical solution is available for comparison.

-3

Set-up

3.1

Analytical solution of the wave force effect

in

a closed container

In a stationary situation the forces generated by the pressure gradient and those generated by waves are in equilibrium. The resulting velocity will be zero provided all other wave effects can be discarded.

In Figure 1 the horizontal forces are shown for the case where the wave force F is constant, we find: FL'"X

-c---~-- -

---l-~~pg~~:~

....

tpgh~X+6X)

x

X +.6. X

Figure 1 wave force effect

Integrating from XL to X yields:

The initial water level is found from:

H= l XJ!h(X)dx= 1

pg[(h2(x)+2F(X~X))%~h3(X)]

X ~X X ~x 3F L pg R L L R L ~ R L If we choose m, m and - 36

-we find for the initial water level:

H==~(2~

-1) ==1.218951416

m

3

and for the final water level at the right boundary:

h (x_{R )} ==/2

=

1.414213562 m

-3.2

Grid for flow computation

The values for XL and x_R were chosen as 0 and 500 m respectively. The grid size in x direction was dx = 10 m in the first computation and dx = 20 m in the second. In y direction the dimensions of the computational domain was 0 m

<

y

<

60 m with grid size dy=lO m and dy=20 m in the first and the second computation respectively.

3.3

Control model

In Figure 2 the tree structure of the control model is shown. Here we see that only one controller is needed to decide whether the stop criterion has been satisfied.

Controller 1 is set to perform another TRISULA computation if the difference in the x components of the velocities of the last two computations was larger than 1 % .

4

Results

The results of the computations and the analytical solution, h - H

, for the cases where dx= 10 m and dx=20 m were used, are shown in Figure 3. We see that the error is extremely small for both cases.

5

Conclusions

For the simple 1 dimensional case the numerical treatment of the Figure 2 Control-modei tree-constant wave force works very well. This suggests that in the structure

interior of the computational domain the wave force is handled

correctly. With respect to the boundaries this test shows that at closed boundaries no problems are introduced.

-(I) u ~ >. (I) c <{ 0 0 0 0 0 0 0 In 0 U> 0 U> eJ 0 0 CJ CJ Co CJ 0 0 CJ

Figure 3 h(x)-H, wave force effect

E E 0 0 0 0 C\J X x V U m

'"

U U (lJ Q) E E ~ ~ Z z

I

I

I

I

J I 0 0 0 0 In CJ N CJ D 0 o o u. 0 U> 'T 0 0 T 0 U> (T) 0 0 (T) 0 U> [,J o o ('J o on o o o L{) 0 0 In ['.J 0