Efficiënt berekenen van vertikale strukturen

I

I

I

I

I

I

l!~'

TU

Delft

I

·

Delft University of Technology

Department of Civil Engineering

Hydraulic andGeotechnicalEngineering Division

I

I

I

I

Efficiënt berekenen van vertikale strukturen

M.D.J.P. Bijvelds en G.S. Stelling rapport no. 10-96

I

I

In opdracht van RijkswaterstaatjR.I.K.Z.

I

I

Technische Universiteit Delft, Faculteit der Civiele Techniek, Sectie Vloeistofmechanica Postbus 5048, 2600 GA Delft.

,

t

_P dichtheid

I

Pa dichtheid van de atmosfeer

Po referentiedichtheid

T, Tb bodemschuifspanning Tij schuifspanningstensor

T6 windschuifspanning aan het vrije wateroppervlak

( wateroppervlakteuitwijking t.o.v. referentieniveau

Lijst van tabellen

I

1 Overzicht methoden voor het verkorten van de rekentijd van 3D-modellen.

2 Overzicht toepassingen en rekentijden voor verschillende numerieke methoden .

.

,

I

I

t

I

3

I

I

1

Inleiding

I

Bij het oplossen van vragen over de morfodynamische processen die zich voor de kust afspelen, is inzicht in de dynamiek van vooroever en brandingsruggen onontbeerlijk. Een instrument dat hierbij aan bij kan dragen is het ARGUS videosysteem dat in het kader van het KUST*2000 onderzoekprogramma bij Bergen Egmond zal worden opgesteld. Een vergelijkbaar videosys-teem is reeds in gebruik in Noordwijk. Ondanks de vorderingen die met dit sysvideosys-teem worden behaald, kent het systeem zijn beperkingen. Integratie met een golfvoortplantingsmodellijkt een geschikte benadering voor kustlangse uniforme banksystemen. Echter de morfologie van deze banken is minder kustlangs dan wordt aangenomen. Het onderkennen van muigeulen zal een uitbreiding van de interpretatietechniek naar drie dimensies vragen, waarbij naast het golfvoortplantingsmodel ook een 3D-stromingsmodel zal moeten worden ingezet.

I

Door de snelle ontwikkeling van de computertechnologie is het sinds enige tijd mogelijk de be-wegingsvergelijkingen voor stroming met een vrij oppervlak (zoals bijvoorbeeld zeeën en estu-aria) volledig drie-dimensionaal op te lossen. Echter, voor praktische situaties kan een volledig drie-dimensionaal model inefficiënt en duur zijn. Dit geldt zeker wanneer men geïnteresseerd is in de morfologische ontwikkeling van kustprofielen waarvoor geldt dat de karakteristieke tijdschalen een orde van grootte van maanden of zelfs jaren hebben. Versnelling van de bere-keningen is in deze gevallen vanuit praktisch oogpunt gezien een noodzaak. Een voorwaarde hierbij is echter wel dat informatie over de vertikale strukturen van de stroming hierbij niet verloren gaat.

I

t

In dit rapport zullen de mogelijkheden voor het efficiënt berekenen van vertikale strukturen in 3D-modellen besproken worden. Daarbij moet onderscheid gemaakt worden tussen het beperken van de rekentijd per tijdstap en het terugbrengen van de totale rekentijd door gro-tere tijdstappen te nemen. Voorts is het mogelijk de vereenvoudigingen te rangschikken naar mathematisch fysisch (mf) of numerieke (num) aard. In onderstaande tabel worden de te behandelen onderwerpen schematisch weergegeven. In hoofdstuk 5 zal een overzicht van de te verwachten rekentijden worden gegeven, gevolgd door de conclusies en aanbevelingen in hoofdstuk 6.

I

t

aard methode hoofdstuk

mf hydrostatische druk 2.1 mf vereenvoudigd k-E-model 2.2 rekenwerk per tijdstap num afhandeling advectie ~.2 num linearisatie continuiteitsvgl. 3.3 mf quasi-3D- modellering 4

grootte tijdstap num ADI 3.1

num AOI 3.1

Tabel 1: Overzicht methoden voor het verkorten van de rekentijd van 3D-modellen.

4

I

I

2

Modelvergelijkingen

I

De basisvergelijkingen voor de water beweging worden in hun meest elementaire vorm weer-gegeven door de Navier-Stokesvergelijking, die het behoud van impuls beschrijft, en de mas-sabehoudsvergelijking, zie bv. [1]. Een vereenvoudiging van deze vergelijkingen is mogelijk in de te beschouwen gebieden.

I

2.1

De ondiepwatervergelijkingen

In estuaria en kustgebieden wordt de stroming gekenmerkt doordat de karakteristieke ho-rizontale schalen van de stroming vele malen groter zijn dan de vertikale schalen. Onder deze omstandigheden is het mogelijk om de bewegingsvergelijking in de vertikale richting (de richting waarin de zwaartekracht werkt) te reduceren tot de hydrostatische drukaanname,

op

oz

=-pg (1)

I

Dit leidt tot een aanzienlijke vereenvoudiging van het stelsel vergelijkingen omdat de druk expliciet te berekenen is uit de waterstand. De bewegingsvergelijkingen in het horizontale vlak zijn hiermee te schrijven als:

OUi O(UiUj)

o(

9

1(

op

d 1 Orij

f

0

-+

+g-+-

-

z---+

i=

ot

OXj OXi Po z OXi PoOXj i=1,2 j=1,2,3 (2)

t

dichtheidsvariaties alleen in rekening worden gebracht als ze voorkomen in combinatie met dewaarbij gebruik is gemaakt van de Boussinesq-benadering wat betekent dat de invloeden van zwaartekrachtsversnelling; voor de andere termen wordt gebruik gemaakt van een referentie-dichtheid Po. Voorts geldt U}=U, U2=V, U3=W enX}=X, X2=Y, X3=Z. De vertikale snelheden

w worden berekend m.b.v de continuiteitsvergelijking welke in het geval van incompressibele stroming geschreven kan worden als:

I

iiu,-' =0

OXi (3)

Samen met de transport vergelijking,

oe

+

O(eui) _

...!_

(D

oe) =

0

ot

OXi OXi OXi (4)

I

levert dit het stelsel vergelijkingen dat opgelost moet worden.

I

2.2

Een vereenvoudigd k

-

é-model

De schuifspanningstensor, Tij, in vgl.(2) bestaat uit Reynolds-spanningen die ontstaan na

Reynolds-middeling van de bewegingsvergelijkingen (hierbij is er vanuit gegaan dat de schuif-spanningen tengevolge van de moleculaire viscositeit verwaarloosbaar zijn ten opzichte van de Reynolds-schuifspanningen). Om een gesloten stelsel vergelijkingen te krijgen moeten deze

5

,

I

termen gemodelleerd worden. Meestal is deze modellering gebaseerd op de zogenaamde Bous-sinesq hypothese die ervan uit gaat dat turbulente transporten evenredig zijn met gradiënten van de betreffende gemiddelde grootheden in de hoofdstroming,

I

(

OUj

OUi)

Tij

= -Vt

OXi

+

OXj

Voor de bepaling van de turbulente viscositeit, VI! wordt in de praktijk vaak gebruik gemaakt

van het k-E-model. Dit model heeft als voordeel ten opzichte van het mengweglengtemodel dat transport van turbulente grootheden in rekening wordt gebracht. De turbulente viscositeit wordt op grond van dimensieanalyse bepaald volgens:

(5)

(6)

I

De turbulente kinetische energie, k, en de dissipatie van turbulente kinetische energie, E,

worden bepaald met behulp van transport vergelijkingen voor deze grootheden. De volledige (gemodelleerde) vergelijkingen (zie bv. [6]) kunnen voor estuaria verder vereenvoudigd worden. Typische afmetingen voor de horizontale schalen in estuaria en kustwateren zijn van de orde van kilometers terwijl de vertikale schalen van de orde van meters is. Dit betekent dat in het algemeen het numerieke grid ook een een grote variatie van celgroottes in horizontale en vertikale richting zal hebben. Hierdoor zullen vertikale snelheidsgradiënten veel groter zijn dan die in horizontale richting. Dit maakt het mogelijk de vergelijkingen voor de turbulente grootheden k enE te vereenvoudigen zonder veel verlies aan nauwkeurigheid

Dk

= ~

(!2.

Ok)

+

p

+

!2.!!.

op _

e

Dt

OZ

(Tk

oz

(Tp

P oz

(7a)

I

DE = ~ (Vt OE)

+

CIE~P _ EE

Dt

OZ

(TE

OZ

k

waarbij

DjDt

=

oui/ot

+

Ujoui/OUj

en

P

=

Vt

[(oujoz)2

+

(ovjoz)2].

(7b)

I

In deze (niet-conservatieve) vergelijkingen zijn de horizontale diffusie en produktietermen verwaarloosd. Door deze grote roosterafstanden in de horizontale richting is het aannemelijk dat turbulentie met een karakteristieke lengteschaal van de waterdiepte en kleiner, verdwenen zal zijn in het volgende roosterpunt. Aan de hand van de relatie van Kolmogorov, kunnen we een tijdschaal voor de turbulentie schatten

(8)

I

waarbij U de karakteristieke snelheidsschaal en L de karakteristieke lengteschaal van de ma-crostructuur van de turbulentie zijn. Relatie (8) kunnen we interpreteren als: een turbulente eddie met energie rv

U2

verliest in een tijdschaal

T

rv

LjU

zijn energie. Nemen we

U

= 0.2

(mjs) en L = 10 (m) als karakteristieke schalen dan levert dit een leeftijd van de orde 1 mi-nuut. Voor advectie binnen het numerieke model geldt dat de het water in een tijd

ó.xju

de

6

I

I

I

volgende cel bereikt. Nemen we een celgrootte Llx = 500 (m) en u

=

1 (mis) dan geeft dit een reistijd van de orde 10 minuten. Dit rechtvaardigt het verwaarlozen van de advectieve termen in de vgl.(7a) en (7b). Omdat in 3D-modellering het berekenen van de advectieve termen een groot deel van de totale rekentijd in beslag neemt zal dit een significante verkorting van de rekentijd opleveren zonder dat de kwaliteit van de berekende resultaten ernstig aangetast zullen worden.

I

,

Voor het oplossen van de transportvergelijkingen voor k en

e

moeten randvoorwaarden wor-den voorgeschreven. Vlak bij de bodem heerst er een lokaal evenwicht tussen produktie en dissipatie van turbulente kinetische energie, wat tot de volgende Dirichlet randvoorwaarde leidt voor k en €:

I

(9a)

(9b)

I

De schuifspanningssnelheid bij de bodem, U.b, wordt bepaald onder de aanname dat vlak bij de bodem, in een dunne laag met bij benadering constante schuifspanningen, een logaritmisch snelheidsprofiel kan worden aangenomen:

I

U.b = ( )

In 1

+

~ZI

2zQ

(10)

waar het subscript 1 duidt op de eerste rekencel boven de bodem. Dit snelheidsprofiel levert een "no-slip"voorwaarde (u = v = 0) op z = Zo. Aan het vrije wateroppervlak worden, bij gebrek aan kennis over de effecten van het wateroppervlak op turbulentie, dezelfde rand-voorwaarden als vgl.(9a) en (9b) gebruikt. De schuifspanningssnelheid kan worden berekend als

I

₍₁₁₎

I

I

7

I

I

I

3

3D-modellering

I

Voor het oplossen van de ondiepwatervergelijkingen, wordt in numerieke modellen, zowel in twee als in drie dimensies, meestal gebruik gemaakt van de ADI methode. Deze methode maakt optimaal gebruik van het feit dat de beweging in het horizontale vlak niet beïnvloed wordt door hydrodynamische drukken maar alleen bepaald wordt door gradiënten in de wa-terstand en dichtheden. Een beschrijving van deze methode, en een variant hierop met kortere rekentijden, wordt in de volgende paragraaf gegeven. De mogelijke vereenvoudigingen voor de integratie in de ruimte wordt in de opeenvolgende paragrafen besproken.

I

3.1

Tijdsintegratie

Het volledig impliciet oplossen van het stelsel ondiepwatervergelijkingen iedere tijdstap is een kostbare zaak. Daarom wordt vaak ADI (Alternating Direction Implicit) factorisatle van een tijdstap toegepast voor de integratie van de ondiepwatervergelijkingen, zie appendix A. Hier-bij wordt een splitsing in richtingen toegepast; in de eerste halve tijdstap worden de termen in de x-richting impliciet geïntegreerd en die in de y-richting expliciet. In de tweede halve tijdstap gebeurt dit vice versa, expliciet in de x-richting en impliciet in de y-richting.

I

I

Tengevolge van de factorisatle komt het veelvuldig voor dat de tijdstap kleiner gekozen moet worden dan voor een goede beschrijving van de fysische processen noodzakelijk is. Deze tijd-stapbeperking is noodzakelijk om stromingen langs grillige geometrieën goed te beschrijven. Bij de ADI methode, die onvoorwaardelijk stabiel is, wordt de tijdstap dan ook volledig be-paald door de geometrie van het gebied (het zogenaamde AD I-effect , zie [7]) en niet door de tijdschalen van de fysische processen.

t

Deze tijdstapbeperking doet zich bij de AOI (Alternating Operator Implicit, zie appendix A) methode niet voor. De AOI-methode, die ook onvoorwaardelijk stabiel is, is gebaseerd op een splitsing in operatoren. De vergelijkingen worden opgesplitst in bijvoorbeeld in advectie-en viscositeitstermadvectie-en enerzijds en de termen die de golfvoortplanting beschrijven anderzijds. Ondanks het feit dat per tijdstap meer bewerkingen nodig zijn dan voor de ADI factorisatie, zal de AOI-methode minder rekentijd vereisen omdat de tijdstap groter kan worden gekozen. In de praktijk blijkt echter wel dat er voor zeer grote tijdstappen convergentieproblemen op-treden voor de iteratieve solvers die in de AOI methode worden toegepast. Dit geldt zowel voor de iteratieve solver voor de continuiteitsvergelijking als ook voor de iteratieve solver voor de bewegingsvergelijking. Deze laatste wordt echter ook toegepast bij de ADI-methode,

hetgeen voor deze methode ook een tijdstapbeperking betekent. Dit betekent dat er voor de AOI-methode ook een beperking is voor de grootte van de tijdstap. Deze is echter veel minder stringent dan voor de ADI-methode. Bij de ADI-methode kan men voor een model met grillige geometrieën en diepe geulen het Courantgetal, wat een maat is voor de tijdstap,

in het algemeen niet veel groter kiezen dan 5 á 10, terwijl het maximale Courantgetal voor de AOI-methode ongeveer 50 is. Naast het voordeel om grotere tijdstappen toe te passen bij de AOI-methode heeft het als voordeel, omdat de keuze van de tijdstap veel minder afhankelijk is van de geometrie, dat er minder inzicht van de gebruiker wordt vereist.

I

8

I

I

I

Het enige verschil tussen de beide methodes is de manier waarop de tijdsintegratie wordt uit-gevoerd. De ruimtelijke discretisaties kunnen identiek zijn. Per toepassing dient op basis van reken kosten en nauwkeurigheid afgewogen te worden welke methode de voorkeur verdient.

3.2

Discretisatie advectieve termen

I

I

.

Bij simulaties met drie-dimensionale modellen blijkt de discretisatie van de advectieve termen vaak een belangrijk deel van de rekentijd te vergen. Voor verscheidene toepassingen (bijvoor-beeld voor getijberekeningen) zijn de advectieve termen van ondergeschikt belang. Echter, er zijn situaties waarin de advectietermen een belangrijke rol spelen in de impulsvergelijking,

zoals bij sterk variërende bodemgeometrie en bij stroming rond obstakels. Het is daarom wenselijk een keuze te hebben uit verschillende discretisaties, waarin zowel een nauwkeurige (dure) optie als ook een minder nauwkeurige (goedkope) optie aanwezig is. Als discretisa-tie voor de minder nauwkeurige afhandeling kan een expliciete afhandeling in de tijd en een eerste orde upwind in de ruimte worden gebruikt. Hiermee vereenvoudigt de ADI-methode omdat de eerste stap hiervan, het impliciet berekenen van de v- of u-snelheden, niet meer nodig is. Door de vereenvoudigde C.q. expliciete afhandeling van de advectieve termen zal de tijdstap echter wel beperkt worden door de CFL-conditie. In situaties waarin advectie een onbelangrijke rol speelt is het dan ook denkbaar om deze termen geheel te verwaarlozen.

I

3.3

Linearisatie van de continuiteitsvergelijking

De waterstanden worden in ondiepwatermodellen berekend door de bewegingsvergelijkingen (2) en de continuiteitsvergelijking (3) te integreren over de diepte en gebruik te maken van de kinematische randvoorwaarde aan het vrije wateroppervlak en de bodem:

(12a)

I

U.llz=_h =0 (12b)

I

Door de geïntegreerde bewegingsvergelijking in de continuiteitsvergelijking in te vullen, ont

-staat een niet-lineair stelsel vergelijkingen. De oplossing van dit stelsel moet iteratief bepaald worden zoals bij de ADI en AOI methoden, wat extra rekentijd met zich meebrengt. Dit kan voorkomen worden door de vergelijkingen te lineariseren zonder iteraties uit te voeren; de totale waterdiepte Hn+l wordt daarvoor vervangen door de diepte H", wat de totale water-diepte in de voorgaande tijdstap is. Zolang de veranderingen van de totale waterwater-diepte tussen twee tijdstappen niet te groot is, hetgeen in het algemeen het geval is omdat anders de

fysi-sche processen niet nauwkeurig worden weergegeven, zal de fout die hierdoor geïntroduceerd

wordt gering zijn.

I

9

I

I

4

Quasi-3D-modellering

I

Het quasi-30-model concept is gebaseerd op het ontkoppelen van de horizontale en vertikale waterbeweging. In plaats van het volledig drie-dimensionaal oplossen van de vergelijkingen wordt de vertikale struktuur afgeleid van een berekende dieptegemiddelde stroming. De re-kentijd kan hierdoor t.O.V. een volledige 30-modellering sterk worden gereduceerd.

4.1

Dieptegemiddelde stroming

I

Door integratie van de ondiepwatervergelijkingen over de diepte worden de vergelijkingen voor de dieptegemiddelde (20H) berekening verkregen (zie bv. [9]). Tengevolge van deze integratie ontstaan echter wel extra termen. Integratie van de advectieve termen geeft dispersietermen zoals

8

j(

- (u -

u)(v -

v)dz

8y -h

(13)

I

waarbij de overstreping een dieptegemiddelde aangeeft. Omdat deze termen horizontale gradiënten van relatief kleine waarden aangeven, is het waarschijnlijk dat deze verwaarloosd kunnen worden in ondiep water, hetgeen in de meeste modellen gebeurt. Het meenemen van deze dispersietermen garandeert niet dat het model nauwkeuriger is; een dispersiecoëfficiënt moet worden voorgeschreven om de invloeden van vertikale middeling te verdisconteren. Voor relatief complexe stromingen zal dat onnauwkeurigheden met zich meebrengen. Bij toepas

-sing van quasi-30-modellering kan de dispersieterm voor de 20H-module eenvoudig worden berekend aan de hand van het vertikale snelheidsprofiel.

I

Het energieverlies wat optreedt door turbulente stroming over een bodem wordt in 20H-modellen vaak gemodelleerd met een kwadratische wrijvingswet,

(14)

I

waarbij u, de stroomsnelheid bij de bodem is. In 20H-modellen wordt hiervoor de dieptege-middelde snelheid genomen, hetgeen alleen juist is in geval van vertikaal uniforme stroming. Echter, de bodemschuifspanningen ontstaan door lokale gradiënten in de stroomsnelheden bij de bodem en moeten daaraan gerelateerd zijn. In stromingen met een niet-uniforme snel-heidsverdeling in de vertikale richting, zoals bijvoorbeeld bij windgedreven stroming of door Coriolis effecten, kunnen berekeningen van de bodemschuifspanning op basis van de diepte-gemiddelde stroomsnelheid sterk afwijken van de werkelijke stroming. In een quasi-30-model wordt de stroomsnelheid bij de bodem bepaald aan de hand van het nader te bepalen vertikale

snelheidsprofiel.

I

4.2

Vertikale struktuur

I

manier is het numeriek oplossen van de vertikale snelheidsprofielen door gebruik te makenEr zijn verschillende manieren om de vertikale struktuur van de stroming te verkrijgen. Een van een grid discretisatie, zoals voor het oplossen van de vergelijkingen in horizontale richting

10

I

I

I

wordt gedaan, zie bv. [5]. Om splitsen van de vergelijkingen in vertikale en horizontale richting (VHS) gemakkelijk toe te kunnen passen, worden bij de bepaling van de vertikale snelheden de horizontale diffusietermen verwaarloosd. Voor ondiep water zullen deze termen gewoonlijk klein zijn t.o. v. vertikale diffusietermen en het is daarom waarschijnlijk dat ze weinig invloed hebben op de vertikale verdeling van de stroming. Deze verwaarlozing kan aanleiding ge-ven tot smalle afwijkingen tussen lokale vertikale snelheden in de 2DH module en de vertikale struktuur module. Dit is vooral het geval in gebieden waar de horizontale uitwisselingstermen belangrijk zijn zoals bijvoorbeeld in een horizontale menglaag. Bovendien kunnen tengevolge van verschillende numerieke methoden in de twee modules afwijkingen ontstaan. De diep-tegemiddelde stroming uit de 2DH module kan gebruikt worden voor de correctie van het vertikale snelheidsprofiel zodat de dieptegemiddelde stroming in beide gevallen gelijk is. De vergelijkingen voor de bepaling van de vertikale struktuur worden dan:

I

i

=

1,2 j = 1,2,3 (15)

I

De vertikale snelheden, w, worden bepaald m.b.v. de continuiteitsvergelijking vgl.(3). Door het linkerlid van vgl.(15) impliciet en het rechterlid expliciet te berekenen (om toepassing VHS mogelijk te maken), ontstaat een tridiagonaal stelsel wat snel en eenvoudig met het Thomas-algoritme op te lossen is. De expliciete benadering van de barotrope druk term kan een beperking van de tijdstap betekenen. Ondanks het feit dat de advectietermen een kleine rol spelen in de impulsvergelijking voor een laag Froudegetal (hetgeen in de praktijk meestal zo is), worden deze toch meegenomen in beide modules. Deze termen hebben namelijk een cumulatief effect op de stroming, zie [5]. Het verwaarlozen van deze termen levert in 2D-berekeningen slechts een kleine besparing van rekentijd op. Nadat de vertikale profielen berekend zijn, wordt de bodemschuifspanning berekend aan de hand van deze profielen. De randvoorwaarde voor de schuifspanning bij het vrije wateroppervlak (tengevolge van wind) en de bodem (tengevolge van bodemwrijving) geeft:

I

I

r,

aul

- =

Vt(()-

.

p

a

z

z=( (16a) Tb

aul

- =

V

t

(-h)-p

a

z

z=-h (16b)

I

De turbulente viscositeit volgt uit het k-ê-model zoals beschreven in paragraaf 2.2. De

diep-tegemiddelde module en de vertikaal-profielmodule zijn aan elkaar gekoppeld middels de vrije wateroppervlakteuitwijking en de bodemschuifspanning. Soms is het mogelijk om de bereke-ning verder te versnellen door voor enkele tijdstappen de dieptegemiddelde stroming uit te rekenen zonder het vertikaal profiel te bepalen. Echter, een grotere tijdstap nemen in de snel-heidsprofielmodule is alleen toegestaan voor langzaam variërende stroming. Fysisch gezien

stelt deze module de vertikale uitwisseling van impuls via diffusie voor. Dit proces is gewoon-lijk een stuk trager dan de voortplanting van zwaartekrachtsgolven in de 2DH-berekening. Voor snel variërende stromingscondities moet voor een stabiele oplossing dezelfde tijdstap

11

I

I

worden genomen in beide modules omdat de bodemschuifspanning in de 2DH-module

expli-ciet berekend wordt.

I

De expliciete afhandeling van de bodemwrijving is in gevallen waar bodemwrijving een do-minante rol speelt niet mogelijk omdat hierdoor instabiliteiten in de berekening worden geïntroduceerd. Dit is bijvoorbeeld het geval in ondiepe gebieden. In deze gebieden is het noodzakelijk de bodemwrijving impliciet te berekenen. De discrete bewegingsvergelijking per laag in de x-richting, voor een volledig 3D-model luidt:

I

U~+l _ u~ (n+l _ (n+l IH/2,k+l/2 l+l/2,kH/2

+

iH i

+

tlt 9 tlx

tin 4>n+l _ tin 4>n+l

i+l/2,k+l i+l/2,k+l i+l/2,k i+l/2,k _ RHSn

tlz - i+l/2,k+l/2

(17)

waarbij het rechterlid alle expliciete benaderingen van de bewegingsvergelijking bevat. De snelheidsgradiënt 4>is gedefinieerd als:

I

(18)

I

Substitutie van vgl.(17) in vgl.(18) levert een tridiagonaal stelsel in</Jn+l. De snelheidsgradiënt bij de bodem kan worden verkregen door een enkelvoudige sweep van de double sweep me-thode [8]. Door deze betrekking in te vullen in de bewegingsvergelijking voor de eerste laag boven de bodem en deze vervolgens van de dieptegemiddelde bewegingsvergelijking af te trek-ken, wordt een lineaire relatie tussen de bodemschuifspanning

Ti~ii2,O

en de dieptegemiddelde snelheid

ur:tl1/2

verkregen, zie appendix B. Terugsubstitutie van de uitdrukking van de bo-demschuifspanning in de dieptegemiddelde bewegingsvergelijking levert een uitdrukking voor de dieptegemiddelde snelheid. Door de snelheid op deze manier te berekenen, wordt de keuze-mogelijkheid gecreëerd op een flexibele manier van numerieke methode te veranderen, zonder ingrijpende aanpassingen in de code aan te brengen. Zo wordt op een efficiënte en robuuste manier het drie-dimensionale karakter van de waterbeweging berekend, waarbij indien nodig, eenvoudig op een volledige 3-dimensionale berekening kan worden overgegaan.

I

I

Een tweede groep zijn de zgn. spectraal modellen (zie bv.

[4],

[10]) waarbij de (niet-lineaire) advectietermen in vgl.(15) verwaarloosd worden of worden gelineariseerd om de spectrale methode toe te kunnen passen. Dit zorgt ervoor dat de stelsels vergelijkingen voor de ho-rizontale en vertikale richting inconsistent zijn. De wateroppervlakteuitwijking, (, is naast de horizontale snelheden u en v de grootheid die de stelsels voor de twee richtingen koppelt. Deze grootheid wordt gebruikt voor het consistent maken van de stelsels; omdat u(z) en v(z) gebruikt worden voor het bepalen van de bodemschuifspanning is het een logische keuze om voor deze grootheden consistentie te behouden.

I

In deze spectraalmodellen worden de stroomsnelheden uitgedrukt als een lineaire combinatie 12

I

I

van een aantal orthogonale functies, de vormfuncties:

w(z

,

t)

=

L

Fk(t)fk(z)

+

E(t)e(z)

k

(19)

I

waarin

!k(Z)

de vormfuncties zijn en

H(t)

de wegingsfactoren. De tweede term in de rechter-zijde van vgl.(19), bestaand uit een vormfunctie

e(z)

en wegingsfactor

E(t),

zijn gelijk aan nul als er geen schuifspanning tengevolge van wind aan het vrije wateroppervlak is. Deze functies worden meestal afgeleid van de eigenfuncties van de vertikale diffusieoperator. De oplossing met een eindig aantal vormfuncties wordt in de basisvergelijkingen gesubstitueerd en dan geïntegreerd over de vertikale richting. De resulterende functies worden dan geminimaliseerd m.b.v. een Galerkin techniek. Hiermee wordt een stelsel van onafhankelijke vergelijkingen verkregen, waaruit de wegingsfactoren bepaald kunnen worden. De eigenfuncties hoeven maar een maal voor iedere situatie bepaald te worden, hetzij analytisch,voor simpele vertikale ver-delingen van de turbulente viscositeit, of numeriek voor praktisch relevante gevallen. De bodemschuifspanning volgt uit de randvoorwaarde voor de schuifspanning bij het vrije wa-teroppervlak en de bodem, vgl.(16a) en (16b): Een nadeel van het "vormfunctiemodelïs dat verondersteld wordt dat turbulente viscositeit constant is in de tijd. Bovendien worden in deze modellen de niet-lineaire advectietermen verwaarloosd.

I

I

4.3

Gelaagde stroming

I

De quasi-3D-modellering is ook toepasbaar op situaties waarin dichtheidsvariaties in de verti-kale richting een rol spelen. Echter, hiervoor moet de vierde term van vgl.(2), die de invloed van deze variaties beschrijft, wel worden meegenomen in de berekening van de vertikale struk-tuur van de waterbeweging.

I

I

I

13

I

I

5

Rekentijden

I

In dit hoofdstuk zal een overzicht worden gegeven voor de rekentijden die te verwachten zijn bij toepassing van de verschillende vereenvoudigingen binnen een 3D-model, zoals in de voorgaande hoofdstukken uiteengezet. Binnen het kader van deze studie is het echter niet mogelijk een exact overzicht te geven van rekentijden voor verschillende situaties met de diverse modellen. Veel van de getallen in tabel 2 zijn dan ook gebaseerd op ervaringen die met verschillende versies van TRISULA en TRIWAQ, numerieke modellen voor waterbeweging en

stof transport van het Waterloopkundig Laboratorium en Rijkswaterstaat, zijn opgedaan. Bovendien worden alleen methoden vermeld die de rekentijd per tijdstap verkorten; ADI of AOI methoden zijn algemeen toepasbaar. Als referentie voor de rekentijden isin onderstaande

tabel uitgegaan van een 2DH-model.

I

I

variabelen type problemen rekentijd

2DH u, v, (, Tb, Ta - berekeningen waarbij vertikaal 1

profiel van ondergeschikt belang

IS

Q3D u, v, (, Tb, Ta, U,V, w, k, e - als 2DH 2.5

- stroming over bodem met flauwe hellingen

- gestratificeerde stroming

3D1ou U, V, w, (, k, € - als Q3D 4

3Dcom U, V, w, (, k,

e

- als 3D1ou 16

U, V, w, (, p, k, e - als 3D1ou

- stroming over steile bodems

3Dnhs - stroming nabij lozings- en 32

innamepunten

- korte golven

I

Tabel 2: Overzicht toepassingen en rekentijden voor verschillende numerieke methoden. 1

I

5.1 Vereenvoudiging modelvergelijkingen

I

o

ndiepwatervergelij kingen

Het gebruik van de ondiepwatervergelijkingen in plaats vande volledige N

avier-Stokesvergelij-king heeft als voordeel dat de druk zelf niet berekend hoeft te worden; deze volgt direct uit de waterstanden. Voor het oplossen van de hydrodynamische druk, die belangrijk is op plaatsen waar de stroomlijnen sterk afwijken van het horizontale vlak, is veel extra rekentijd vereist. Voor een 2DV-model, gebaseerd op de Euler vergelijkingen, neemt de rekentijd ongeveer met

een factor 2 toe [2]. Voorstroming van kustwateren en estuaria is het voldoende nauwkeurig

I

13D1ou

=

advectie eerste orde upwind expliciet

3Dcom

=

advectie hogere ordeimpliciet

3Dnh•

=

niet-hydrostatische berekening van de druk (volledige Navier-Stokesvergelijking)

14

I

I

I

om de waterbeweging met de ondiepwatervergelijkingen te beschrijven. Lokaal kunnen grote bodemgradiënten voorkomen (bijvoorbeeld op de overgang van geulen naar platen) zodat de hydrostatische drukaanname niet meer juist is. Echter, de hierdoor geïntroduceerde fout heeft slechts een lokaal karakter; globaal zal de kwaliteit van de berekening hierdoor niet aangetast worden.

I

I

Vereenvoudigd k-ê-model

Het vereenvoudigde turbulentiemodel, dat terug is gebracht tot een lDV turbulentiemodel, vergt een minimale rekeninspanning. Door het toepassen van een double sweep methode in de vertikale richting, wordt de rekentijd per tijdstap bijna verwaarloosbaar ten opzichte van de andere rekenroutines. Op deze manier is een turbulentiemodel verkregen dat qua rekentijd kan concurreren met een mengweglengtemodel maar dat in meer complexe situaties, zoals bijvoorbeeld stroming over een drempel, betere resultaten geeft.

5.2

Numerieke implementatie

I

ADI en AOI methode

I

Door het gebruik van de ADI-methode wordt een grote tijdwinst verkregen t.o.v. het si-multaan oplossen van het hele stelsel vergelijkingen. Door afwisselend termen impliciet en expliciet te behandelen in iedere halve tijdstap, ontstaan stelsels vergelijkingen die relatief eenvoudig op te lossen zijn. Echter, in gebieden met sterk variërende geometrie zal de tijd-stap voldoende klein gekozen moeten worden. Bij de AOI methode is dit niet het geval. Ondanks het feit dat deze methode per tijdstap duurder is, zal de AOI methode daarom toch voor een simulatie (met dezelfde nauwkeurigheid) goedkoper zijn dan ADI. De AOI methode beperkt dan ook niet de rekentijd per tijdstap maar maakt het mogelijk dat de tijdstap zelf groter kan worden gekozen zonder verlies aan nauwkeurigheid. Het blijkt dat de AOI methode tot een Courantgetal van ongeveer 50 voor barotrope termen convergeert. De karakteristieke relatieve rekentijdverkorting die wordt verkregen door toepassing van AOI in plaats van ADI is ongeveer een factor 2.

I

Advectietermen en linearisatie

I

Voor de vereenvoudigde implementatie van de adveetleve termen (expliciet in de tijd, upwind in de ruimte) in combinatie met het lineariseren van de continuiteitsvergelijking zonder ite-ratie, levert naar verwachting een reductie van de rekentijd met een factor 4 op. Hiermee is deze methode slechts 4 maal duurder dan een 2DH-berekening. Door de grote maaswijdtes waarmee in de praktijk gerekend wordt, zal de expliciete afhandeling van de advectie in het al-gemeen geen sterke tijdstaprestrictie opleveren. Indien dit toch het geval is kan dit in situaties waarin advectie onbelangrijk is voorkomen worden door advectie volledig te verwaarlozen.

I

15

I

I

I

Quasi-3D- modellering

Quasi-30-modellering is een zeer efficiënte manier om 30-stromingen op een redelijke nauw-keurige manier te voorspellen. Voor een berekening van een windgedreven stroming in het Ysselmeer blijkt de benodigde CP U-tijd voor de quasi-30-modellering met 5 rekenpunten in de vertikale richting ongeveer 2.5 maal zo groot te zijn dan voor een 2DH-model, terwijl de resultaten veel beter overeenkomen met metingen [5].

I

I

I

I

I

I

I

16

I

I

I

6

Conclusies en aanbevelingen

I

Voor het voorspellen van de morfologische ontwikkeling van kust profielen zal gebruik wor-den gemaakt van een drie-dimensionaal stromingsmodel. Echter, doordat de karakteristieke tijdschalen van deze ontwikkelingen maanden of zelfs jaren zijn, is het noodzakelijk de bere-keningen te versnellen. Gezien het karakter van de waterbeweging voor de kust ligt het voor de hand om de ondiepwatervergelijkingen als basis van het numerieke model te gebruiken. Ondanks het feit dat in kustgebieden lokaal het gebruik van deze vergelijkingen niet altijd gerechtvaardigd is, zoals bij de overgang van geul naar plaat, zal dit globaal gezien geen ver-slechtering van de resultaten opleveren. Lokaal kunnen echter wel afwijkingen ten opzichte van de werkelijkheid ontstaan. Naast deze vereenvoudiging van de modelvergelijkingen is het mogelijk het k-€-turbulentiemodel terug te brengen tot een lOV-turbulentie model. De geringe diepte in kustgebieden in combinatie met grote roosterafstanden zorgt ervoor dat de leeftijd van de turbulente eddies korter is dan de reistijd van deze turbulentie naar de volgende rekencel. Verwaarlozing van de advectietermen in de transportvergelijkingen voor de turbu-lente grootheden is in deze gevallen dan ook geoorloofd. Hiermee wordt de rekeninspanning voor het turbulentiemodel teruggebracht tot een fractie van de totale rekentijd.

I

I

I

De numerieke afhandeling van de advectieve termen heeft een grote invloed op de totale re-kentijd die per tijdstap nodig is om het stelsel vergelijkingen op te lossen.

In

combinatie met het lineariseren van de continuiteitsvergelijking reduceert bij een expliciete eerste orde up-wind afhandeling van de advectietermen de rekentijd met ongeveer een factor vier. Bovendien blijkt uit de ervaring die is opgedaan met deze vereenvoudigde aanpak, dat de resultaten in veel gevallen sterk overeenkomen met de complexe discretisatie. Nadeel is echter wel dat de expliciete benadering een beperking van de tijdstap kan betekenen.

I

Een verdere reductie van de rekentijd is mogelijk door toepassing van een quasi-30-model waarbij de horizontale en vertikale waterbeweging ontkoppeld worden. Voor de berekening van het vertikale snelheidsprofiel bestaat daarbij de keuze tussen eindige differentie methoden en spectraalmodellen. Deze laatste hebben echter als nadeel dat de vergelijkingen geline-ariseerd moeten worden om de spectrale methoden toe te passen, hetgeen bij de eindige differentie benadering niet het geval is. Bovendien zijn spectraal methoden ingewikkelder. Toepassing van de eindige differentie methode voor de bepaling van het vertikale snelheids-profiel levert samen met vereenvoudigde turbulentiemodel en de impliciete afhandeling van de bodemwrijving een zeer efficiënt en robuust model. Ondanks de eenvoud van het model, is het mogelijk een groot scala van stromingen hiermee te berekenen, zonder dat de resultaten sterk afwijken van een volledig 3D-model. Nadeel van het model is echter wel dat door de expliciete benadering van de barotrope drukterm bij de bepaling van het vertikale snelheidsprofiel, de tijdstap beperkt kan worden door stabiliteitscondities.

I

I

Doordat de quasi-30-formuleringen voor de bodemschuifspanning en dieptegemiddelde snel-heid, zoals beschreven in appendix B, worden afgeleid van de volledige drie-dimensionale bewegingsvergelijkingen, ontstaat er op een flexibele manier een keuzevrijheid tussen het gebruik van het volledige 30-model en het quasi-3D-model. Het toepassen van deze

imple-17

I

I

mentatie lijkt dan ook de aangewezen methode voor de modellering van de water beweging.

In

veel gevallen zal de eenvoudige en daarmee goedkope quasi-3D-formulering voldoen voor het beschrijven van de waterbeweging.

In

meer complexe situaties, zoals voor stroming over zeer onregelmatige bodems, is het gebruik van een volledig 30-model nodig hetgeen de flexi-bele implementatie eenvoudig mogelijk maakt. Voor de eenvoudige advectieformulering zal hierdoor de rekentijd slechts met ongeveer 50% toenemen.

I

I

I

I

I

I

I

18

I

I

I

Referenties

I

[1] G.K. Batchelor. An introduetion to fluid dynamics. Press Syndicate of the University of

Cam bridge, Cambridge, 1967.

[2] M.D.J .P. Bijvelds. Towards numerical modelling of cooling water discharges and related phenomena. Master's thesis, Technische Universiteit Delft/Waterloopkundig Laborato-rium, 1995.

I

[3] V. Casulli and R.T. Cheng. A semi-implicit finite difference model for three-dimensional tidal circulation. In M.L. Spaulding, K. Bedford, A. Blumberg,

R.

Cheng, and C.

Swan-son, editors, Estuarine and coastal modeling, proc. f!1d conf., pages 620-631, New Vork,

1992. Ameriean Society of Civil Engineers.

[4] A.M. Davies. On extracting current profiles from vertieally integrated numerical modeis.

Coastal Eng., 11:445-477, 1987.

[5] X.Y. Jin. Quasi-three-dimensional numerical modelling of flow and dispersion in shallow water. PhD thesis, Technische Universiteit Delft, 1993.

I

[6] W. Rodi. Turbulence models and their applications in hydraulics - A state of the art review. International Association of Hydraulie Research, Delft, Nederland, 1980.

[7] G.S. Stelling. On the construction of computational methods for shallow water flow problems. Number 35. Rijkswaterstaat communications, Den Haag, Rijkswaterstaat,

1984.

I

[8] G.S. Stelling. Compact differencing for stratified free surface flow. In Advances in hydro -science and -engineering, volume 2, pages 378-386, 1995.

[9] C.B. Vreugdenhil. Numerical methods for shallow water flow. Kluwer Academie Press,

Dordrecht, Nederland, 1994.

I

[10] T.J. Zitman. Quasi three-dimensional current modelling based on a modified version of davie's shapefunction approach. Conto Shel] Res., 12(1):143-158, 1992.

I

I

19

I

I

A

ADI en AOI methode

I

Bij de ADI factorisatie wordt de berekening van een tijdstap in twee richtingen gesplitst. In de eerste halve tijdstap (stage) worden de termen in de x-richting impliciet geïntegreerd en die in de y-richting expliciet. In de tweede halve tijdstap gebeurt dit vice versa. In de continui-teitsvergelijking wordt één fluxterm impliciet genomen. Voor de bijbehorende snelheid wordt in de impulsvergelijking de drukterm ook impliciet berekend. Dit resulteert in een gekoppeld stelsel vergelijkingen in één richting. De horizontale snelheden van de opeenvolgende lagen zijn aan elkaar gekoppeld door de vertikale advectietermen en de vertikale viscositeitsterm. Over een hele tijdstap genomen is de methode tweede orde nauwkeurig in de tijd

I

Stage 1:De semi-discrete vergelijkingen voor de eerste halve tijdstap zien er als volgt uit:

I

(20a)

I

u* - un 8un 8un êu: H*,q 8C·,q+l ~~t

+

u· 8x

+

v· 8y

+

wn ÖZ - Jv·

+

H*,q+19 8x

+

82un 82un 82un - lJh 8x2 - lJh 8y2 - lJv öz2 =0 (20b) (20c)

I

De eerste stage begint met het berekenen van de v-snelheden op

t

=

~~t,

aangeduid met het superscript

*.

De horizontale snelheden u en w kunnen worden opgelost m.b.v. een Red-Black Jacobi-iteratie in combinatie met een double sweep methode. Door de bewegings-vergelijking voor u, te integreren over de diepte en in te vullen in de dieptegeïntegreerde continuiteitsvergelijking (met inachtneming van de kinematische randvoorwaarden aan het vrije wateroppervlak en bodem) ontstaat een tridiagonaal stelsel in ( wat snel opgelost wordt met het Thomas-algoritme. De drukgradiëntterm in bovenstaande vergelijkingen zijn geline-ariseerd door vermenigvuldiging met H*,q/ H·,q+l. Op deze manier ontstaat een lineair stelsel dat massabehoudend is. De geïntroduceerde fout hierdoor wordt geminimaliseerd door het stelsel iteratief op te lossen. De barocline term is in bovenstaande vergelijkingen weggelaten. Omdat deze expliciet berekend wordt is deze eenvoudig aan het model toe te voegen.

I

I

Stage 2:

20

I

I

De semi-discrete vergelijkingen voor de tweede tijdstap luiden:

un+1 - u* 8un+1 8un+! 8un+1 8(*

--:----

+

u*--

+

v*--

+

w*-- - Jv*

+

g-

+

~~t 8x 8y 8z 8x

82un+1 82un+1 82un+1 - Vh 8x2 - Vh 8y2 - Vu 8z2

=

0 (21a)

I

I

(21b) (n+!,q+l _ (*

+

8(H*

J~h

u*dz)

+

8(Hn+l,q+l

J~h

vn+1dz)

~~t

8x 8y (21c)

I

De berekening van de waterstanden en snelheden gaat in de tweede halve tijdstap analoog aan de oplosmethode voor de eerste halve tijdstap. De vertikale snelheden volgen uit de di-vergentievrijheid van de stroming.

I

De Alternating Operator Implicit (AOI) methode is gebaseerd op splitsing in operatoren, in tegenstelling tot een splitsing in richtingen zoals bij de ADI methode. Bij deze methode wordt een splitsing toegepast in advectie- en viscositeitstermen enerzijds en in termen die de golfvoortplanting beschrijven anderzijds. In de eerste halve tijdstap wordt de continuiteits-vergelijking volledig expliciet geïntegreerd en in de tweede halve tijdstap volledig impliciet. De vertikale snelheden volgen wederom uit de divergentievrijheid van de stroming. Evenals de ADI methode is de AOI methode tweede orde nauwkeurig in de tijd.

Stage 1:

I

(22a)

I

(22b)

I

(22c) 21

I

I

I

Een benadering voor u", die nodig is voor het oplossen van de bewegingsvergelijking voor u· wordt verkregen door de volgende vergelijking op te lossen:

I

I I 2 I

v· - vn 8v· 8(n

8

v·

1

+w

n-8

+9-

8

-vV-82

+

Jun

=0

2~t z y z

De horizontale advectie- en viscositeitstermen zijn weggelaten om de reken kosten voor het oplossen van vgl.(23 te beperken; er behoeft op deze manier slechts een tridiagonaal stelsel opgelost te worden en iteratie is dus niet nodig. Door een benadering voor v·' op het tussen-niveau te bereken, wordt de tijdsintegratie van o.a. de Coriolisterm meer impliciet, hetgeen mogelijke oscillaties (kleine instabiliteiten) in geval van expliciete afhandeling onderdrukt. Bij gebruik van hogere orde upwind discretisaties in de ruimte voor de adveetleve termen, zal een niet symmetrisch stelsel ontstaan. Deze stelsels kunnen opgelost worden m.b.v. een iteratieve Red-Black Jacobi techniek in de horizontale richting gecombineerd met een directe tridiagonale solver voor de vertikale richting.

(23)

I

I

Stage 2:

un+1 - u· au· s«: aun+1 Hn+l,q 8(n+l,q+l

n+l • •

J •

+

~~t

+

u 8x

+

v 8y

+

W

---a;- -

v -H-n-+-l-,q-+-=-l

9

ay

+

82u· 82u· 82un+1

- Vh 8x2 - Vh 8y2 - Vv 8z2

=

0 (24a)

I

(24b)

I

(n+l,q+l _ (* 8(Hn+l,q+l

J~h

un+1dz) 8(Hn+l,q+l

J~h

vn+1dz)

--~---+

+---~~---~~t

8x 8y

In beide tijdstappen worden, vanwege stabiliteitsredenen, de vertikale advectie- en diffusieter-men volledig impliciet geïntegreerd

[3].

De twee bewegingsvergelijkingen vormen samen met de dieptegeïntegreerde continuiteitsvergelijking een gekoppeld, niet lineair, stelsel. Evenals bij de ADI methode worden de impliciete druktermen daarom in beide impulsvergelijkingen met H·,q / H·,q+l vermenigvuldigd om een lineair stelsel te verkrijgen. Invullen van de twee bewegingsvergelijkingen in de continuiteitsvergelijking levert een pentadiagonaal stelsel op. Dit stelsel kan geschreven worden als

(24c)

I

A(n+l,q =B Matrix A is hierbij

(25)

I

A

=

I

+

( ~t)2

""2

9

H (~xx

+

~yy)

(26)

22

I

I

I

waarbij ~xx

+

~yy de gediscretiseerde Lapla.ceoperator is. Het stelsel vergelijkingen kan

m.b. v. een iteratieve solver worden opgelost. Dit levert echter wel een restrictie in de tijdstap vanwege convergentieproblemen die optreden voor grote tijdstappen. Een mogelijke oplossing voor dit probleem is door gebruik te maken van de geconjugeerde gradiënten methode (CG).

Deze methode is echter wel relatief duur omdat voor iedere tijdstap een matrix moet worden opgesteld.

I

I

I

I

I

I

I

23

I

I

I

B

Uitdrukking

voor de bodemschuifspanning

De bewegingsvergelijking per laag, hier alleen gegeven in de x-richting, luidt:

I

n+l _ n (n+l (n+l

uiH/2.kH/2 ui+l/2.kH/2

+

i+l - i

+

~t

9 ~x

vn </>n+l _ vn _A.n+l

iH/2,k+l i+l/2,k+l i+l/2,k'f'i+l/2,k _ RHSn

~Z - i+l/2,kH/2

waarbij het rechterlid alle expliciete benaderingen van de bewegingsvergelijking bevat. De snelheidsgradiënt </>is gedefinieerd als:

(27)

I

</>n+l = uni+l/2,kH/2+1 _ uni+l/2,k-l/2+1

i+l/2,k - ~z

Een uitdrukking voor de snelheidsgradiënt boven de bodem, wordt verkregen door vgl.(27) en vgl.(28) te combineren en een enkele sweep van de double sweep methode uit te voeren,

(28)

I

Hierbij volgen de coëfficiënten</>1= -a</>o

+

b a en b uit de voorwaartse sweep. De subschrift 0 en 1 duiden(29) respectievelijk op de bodem zelf en de eerste laag boven de bodem. Bij de bodem wordt een logaritmisch snelheidsprofiel aangenomen waarvoor geldt

I

(30)

Zo is de bodemruwheidslengte. De snelheid in de eerste laag boven de bodem wordt daarmee

u.

(ZI)

UI = -

In

-

=c</>o

/'i, Zo (31)

I

van vgl.(29) en vgl.(31) in de bewegingsvergelijking van de bodemlaag levert (de eerste index,waarbij

ZI

de coördinaat van het u-punt in de eerste laag is enc =

Zo

ln(zd

zo).

Substitutie i, is weggelaten m.u.v. de barotrope term)

( ~t

_c

₊

_{_vn a}

₊

_{_vn}

~t)

_{_A.n+l=_vn}

~t

_b

₊

_un

₊

_{g_ ((n+l _}

~t

_(!l+l)

~z 1 ~Z 0 '+'0 ~Z 1 0 ~x t+l t (32)

I

De vergelijking voor de dieptegemiddelde waterbeweging luidt

(33)

I

Door vgl. (32) van vgl. (33) af te trekken ontstaat een lineaire relatie tussen </>0(en daarmee de bodemschuifspanning) en

u.

Terugsubstitutie van de uitdrukking voor </>0in de dieptege-middelde bewegingsvergelijking (33) geeft een uitdrukking voor de dieptegedieptege-middelde snelheid

u.

24

I

I

I

I

I

I

I

I

I

I