Solution of the St Venant Equations / Shallow
Download
Report
Transcript Solution of the St Venant Equations / Shallow
Solution of the St Venant
Equations / Shallow-Water
equations of open channel flow
Dr Andrew Sleigh
School of Civil Engineering
University of Leeds, UK
www.efm.leeds.ac.uk/CIVE/UChile
St. Venant Equations
d
d u
g
u
g S o S f
x
x t
A
u
d
u
A
b
0
x
x
t
St Venant Assumptions of 1-D Flow
Flow is one-dimensional i.e. the velocity is uniform over the cross
section and the water level across the section is horizontal.
The streamline curvature is small and vertical accelerations are
negligible, hence pressure is hydrostatic.
The effects of boundary friction and turbulence can be accounted for
through simple resistance laws analogous to those for steady flow.
The average channel bed slope is small so that the cosine of the
angle it makes with the horizontal is approximately 1.
Cunge J A : Practical Aspect of Computational River Hydraulics
Control Volume
y
c
re
di
n
tio
A1
V1
f
of
w
lo
h1
free surface
h
2
A2
V2
z1
x1
s
z2
x
x2
bed
Continuity
CV in x,t plane
between cross sections x=x1 and x=x2
between times t=t1 and t=t2
conservation of mass
A A dx Q
x2
x1
t2
t2
t1
t1
x2
Qx1 dt 0
Momentum
Conservation of momentum
2
2
uA
uA
dx
u
A
u
Ax 2 dt
t2
t1
x1
x1
t1
x2
t2
g
t2
g
t2
g
t2
t1
t1
t1
I
1 x1
x2
x1
I1 x 2 dt
I 2 dx dt
AS
x2
x1
o
S f dx dt
Geometric Change Terms
I1
0
()
vertical change
in cross-section
h x
b
h
h-
d
hx x, d
change in width along the length of the
channel.
h x
I 2 hx
d
0
x hho
Integral / Differential Forms
Integral form do not require that any flow
variable is continuous
We will see later finite difference methods
based on this integral form.
Can derive differential form … but
Must assume variables are continuous and
diferentable
Replace variable with Taylor’s series
At 2
A
2 A t 2
At1
t 2
...
t
t 2
Differential form
-
A Q
0
t x
Q Q 2
gI1 gASo S f gI2
t x A
gI 1
h
gAx gI 2
x
x
Q
h
uQ gA So gASf 0
t x
x
In terms of Q(x,t) and h(x,t):
Where b = b(h), A=A(h)
Ah A h
h
b
t
h t
t
h 1 Q
0
t b x
Q Q 2
h
gA gAS f So 0
t x A
x
remember b = b(h), A=A(h)
Each of these forms are a set of non-linear differential equations which do not
have any analytical solution.
The only way to solve them is by numerical integration.
In term of u(x,t) and h(x,t):
using
A h A
Q
u
A
u
A u
A u
x
x
x
x
h
x
x hcons tant
h A u
h u A
u
0
t b x
x b x hcons tant
u
u
h
u
g
g S f S o 0
t
x
x
Characteristic Form
The St Venant equations may be written in a
quite different form know as the Charateristic
Form.
Writing the equations in this form enables
some properties and behaviour of the St
Venant equations become clearer.
It will also help identify some stability criteria
for numerical integration
will help with the definition of boundary
conditions.
Characteristic form
Consider a prismatic channel
h A u
h
u
0
t b x
x
u
u
h
u
g
g S f S o 0
t
x
x
Wave speed
Consider the speed, c, of a wave travelling in
the fluid.
c g
A
b
with respect to x and t gives:
c
h
2c g
x
x
2
c
c
u
2u c
0
t
x
x
2c
c
h
g
t
t
u
c
u
2c u
g S f S o 0
t
x
x
Combining
Adding equations (3) and (4) gives
v
v
c
c
v c 2 2v c g So S f
t
x
t
x
Subtracting equations (3) and (4) gives
v
v
c
c
v c 2 2v c g S o S f
t
x
t
x
Characteristics
Equations (5) and (6) can be rearranged to
give respectively
v c v 2c v 2c g So S f
x
t
v c v 2c v 2c g So S f
x
t
Total differential
For some function f of x and t
df sense.
f dx f
In a physical
dt
fx is
dt some
t property of the flow e.g.
If the variable
surface level,
if an observer is moving at velocity v,
the observer will see the surface level change
only with time relative to the observers' position.
The characteristic form of the St Venant
equations
If we take f = (v + 2c)
Total differential is
d v 2c v 2c dx v 2c
dt
x
dt
t
Compare with
v c v 2c v 2c g So S f
x
Clearly
and
dx
v c
dt
d v 2c
g S o S f
dt
t
Characteristic form of the St Venant
Equations
These pairs are known as the Characteristic
form of the St Venant Equations
d v 2c
g S o S f for
dt
dx
v c
dt
d v 2c
g S o S f
dt
dx
v c
dt
d v 2c
g S o S f
dt
for
for
dx
v c
dt
Meaning of the characteristics
The paths of these observers that we have talked about
can be represented by lines on this graph.
c2
dt
1
dx v c
t
dt
1
dx v c
8
c1
6
c4
7
c0
5
c5
4
3
2
Zone of quiet
1
x
Information paths
it takes time for information to travel
E.g. a flood wave at u/s end
The channel downstream will not receive the
flood for some time.
For how long?
The line C0 represents the velocity of flood
wave
Everything below C0 is zone of quite
Zones of Dependence and Influence
The idea that characteristics carry information
at a certain speed gives two important
concepts
c2
c1
t
c2
P
Domain of dependence of P
c1
Zone of influence of Q
Q
x
Stability
These zones imply a significance to
numerical methods and stability of any
solver
The numerical method must take only
information from within the domain of
dependence of P
this limits the size of time step
Calculation with characteristics
If we know the solution at points 3 and 5
Can determine the solution at point 6.
For C5-6 then
d v 2c
g i j
dt
1
dx v5 c5
v6 2c6 v5 2c5 t g i j
dt
For C3-6 then
d v 2c
g i j
dt
dt
1
dx v3 c3
v
6
2c6 v3 2c3 t g i j
Characteristics solution
Adding equations
v6
v5 v3
c5 c3 t g i j
2
Subtracting
c6
v5 v3 c5 c3
4
2
MOC on Rectangular Grid
t
P
t+t
t
t
W
L
O
R
E
x
x
P
t+t
t
t
W
L
O
x L
x
R
x R
x
E
Midstream discretisation
Away from boundaries
dx
v L c L
dt
x L
v O c O
t
P
t+t
dx
v P c P
dt
t
t
W
L
O
x L
x
x L
vO vW
v L vO
x
R
x R
E
x R
v O c O
t
x
x R
vO v E
v R vO
x
Solution
d v 2c
g So Sf
dt
vP 2cP vL 2cL t g So Sf L
d v 2c
g So Sf
dt
vP 2cP vR 2cR t g So Sf R
vP
vL vR
t g
So Sf L So Sf R
cL cR
2
2
cP
vL vR cL cR t g
So Sf L So Sf
4
2
4
R
Stability
Considering characteristics
C2
t+t
t
t
P
dt
1
dx = v+c
L
C1
dx
v c
dt
dt
1
dx = v-c
M
R
x
t
v c
x
x
t 0.9
v c max
Boundary conditions
Upstream boundary:
Downstream boundary:
forward characteristic
t+t
t
P
P
t+t
boundary
backward characteristic
boundary
t
t
t
O
R
E
W
L
O
x R
x
x L
x
Boundary Conditions
The second equation is a
boundary condition equation
Upstream depth boundary
c p ghp
hP = H(t)
vP 2cP vR 2cR t g So Sf R
vP vR 2cR cP t g So Sf R
Spuer-critical - mid
Right characteristic moves to left
P
t+t
t
t
W
L
R
O
E
x R
x L
x
x
solution method is exactly the same
Super-critical - upstream
No characteristics
P
boundary
L
R
t
O
E
Super-critical downstream
2 characteristics
No boundary condition
P
t+t
boundary
t
t
W
L
R
O
x R
x L
x
Finite Difference Schemes
Two basics classes
Implicit
Explicit
Commercial packages use implicit
Explict
for high accuracy (sometimes!)
Testing / understanding behaviour
Class examples!
Which Equations ?
Not always clear what equations a being used!
What are the shallow water equation?
We will look at schemes based on the integral
equations:
x1 At 2 At1 dx t1 Qx2 Qx1 dt 0
x2
t2
2
2
x1 uAt 2 uAt1 dx t1 u Ax1 u Ax 2 dt
x2
t2
g
t2
g
t2
g
t2
t1
t1
t1
I
1 x1
x2
x1
I1 x 2 dt
I 2 dx dt
AS
x2
x1
o
S f dx dt
Homogeneous Integral Equations
Without the gravity / frictions terms
A
x2
x1
t2
At1 dx Qx 2 Qx1 dt 0
t2
t1
2
2
Q
Q
dx
u
A
gI
u
A I1 x1 dt 0
1 x2
t1
x1 t 2
t1
x2
t2
Grid based
Consider the grid …
B
n+1
t
t
n
C
f
A
D
x
i-1
x
i
i+1
Integrate around the cell
Considering the cell ABCD,
Integral can be written in the general vector form :
fdx G f dt 0
ABCD
Q
f
A
Q
G f Q 2
gI1
A
Gridpoints / Variables
Variables are all known or will be calculated
at the grid points
xi represents the value of x at position i
tn represents the value of t at position n
Derivation approximates values:
Gxi , t Gin1 1 Gin
0 1
f x, tn fi n1 1 fi n
0 1
Substitute in approximations
And the equation become
f
xi 1
xi
n 1
i 1
t n1
tn
f
n 1
i 1
1 f i n 1 f i n1 1 f i n dx
G
n 1
i 1
1 f i
1 Gin1 Gin 1 1 Gin dt 0
n 1
f
n
i 1
1 f i
n
x
Gin11 1 Gin1 Gin1 1 Gin t 0
Difference equation
f
n 1
i 1
Divide through byΔx Δt
And can see that it is an approximation of
f G f
0
t
x
1 f i n1 f i n1 1 f i n Gin11 1 Gin1 Gin1 1 Gin
0
t
x
i.e. starting with integral form
discretisation is also valid for diferential form
Several schemes
This is a general discretisation scheme
Vary the parameters ψ and θ
Get a family of different finite difference schemes
Features are
They are implicit for values of ψ > 0. else explicit.
They link together only adjacent nodes.
Space interval can vary – no loss of accuracy.
They are first order, except for the special case of
ψ=θ=0.5 when they are second order.
Preissmann Scheme
ψ = 0.5 gives Preissmann 4-point scheme
Time derivative
f in11 f i n1 f in1 f i n
f
t
2t
Space derivative
G f Gin11 1 Gin1 Gin1 1 Gin
x
x
Equations become …
x n1
x n
Ai Ain11
Ai Ain1 t Qin11 Qin1 1 Qin1 Qin 0
2
2
x n1
x n
Qi1 Qin1
Qi1 Qin
2
2
n1
n1
n
n
2
2
Q2
Q2
Q
Q
t
gI1
gI1 1
gI1
gI1
A
A
A
A
i
1
i
i
1
i
I 2 AS f ASo n1 I 2 AS f ASo n1
j
j 1
gtx
n
n 1
2
1 I 2 AS f ASo j I 2 AS f ASo j 1
0
More common form
The terms I1 and I2 are often difficult
(expensive / time consuming)
(when originally attempted)
Usual form used in packages
h 1 Q
0
t b x
Q Q 2
h
gA gAS f So 0
t x A
x
Priessmann
Function represented by,
f
0.5 < θ ≤ 1
Continuity:
f
2
n 1
i 1
fi
n 1
1 n
f
n
f
i 1
i
2
h 1 Q
0
t b x
1 hin11 hin1 hin1 hin 2 Qin11 Qin1 1 Qin1 Qin
0
n1
n1
n
n
2
t
t x bi bi1 1 bi bi1
Priessmann
Momentum Equation
Q Q 2
h
gA gAS f So 0
t x A
x
1 Qin11 Qin1 Qin1 Qin
2
t
t
Q 2 n1 Q 2 n1 1 Q 2 n Q 2 n
x A i 1 A i
x A i 1 A i
1 n
n1
n 1
g Ai 1 Ai
Ai1 Ain E 0
2
2
Source term
So-Sf
hin11 hin1
hin1 hin
E
1
x
x
n
n
n n
Qin11 Qin11 Qin1 Qin1
Q
Q
Q
i 1 i 1
i Qi
1
2
2
n 2
n
K n1 2 K n1 2
K
K
i
i
i1
1 i1
2
2
2
Unknowns
There are the 4 unknowns
hnj1, hnj11, Qnj1, Qnj11,
Plus k, b, A which are functions of h and Q
n1
j
n1
j 1
n1
j
n1
j 1
h ,h ,Q ,Q ,
Need to linearise
“Linearise” equations
Ahnj1 Bhnj11 CQnj1 DQnj11 RHS
A, B, C, D and RHS are funtion of the
unknowns.
But not strongly
Boundary Conditions
In implicit scheme specify h or Q
Use characteristics to decide appropriately
Or relation between h and Q
2N unknowns (h, Q at each node)
2N – 2 equations from internal points
2 boundary equations
Junction
At a junction each chanel share the same
node
hjuntion 1 = hjuntion 2 = hjuntion 3 …..= hjuntion n = h
Continuity
Sum of inflow and outflow equal to zero
Iterative solution
Need to iterate updating coeficients
Cunge says …
iterates rapidly one or two iterations
Newton-Raphson methods used in packages
Stability
Formally unconditionally stable for all time
steps
0.5 < θ ≤ 1
Further away from Cr = 1 less accurate
Cr = 20 is common
Because of linearisation may fail for extreme
flows or those that are too far from original
assumptions
Explicit Schemes
Not used in simulations of real rivers
Time step limitations.
Important features of explicit schemes
they are simple to implement
allow experiment with weights, time-step and
space-step to understand behaviour of the
solution.
There are MANY schemes
Leap-Frog
Earliest scheme aplied to wave equations
Spatial derivative
f f
f
t
2t
n
i 1
gives
n
i 1
Temporal derivative
G f Gin1 Gin1
x
2x
n1
A
t
n1
n
i
f i n1 f i
Gin1 Gin1
2x
Qi
Stability
All explict scheme have time step limit
Courant confition (CFL)
dx
v c
dt
Cr < 1
x
t
v c
Lax Explict Scheme
Similar to leap from with weighting, α, in time
n
n
f
f
n 1
n
i 1
i 1
f i f i 1 1
2
f
t
t
For 0 ≤ α ≤ 1
Spatial derivative
G f Gin1 Gin1
x
2x
Lax Explict Scheme
Leads to function solution
n1
n
n
A
f
f
t
n1
n
n
n
i
i 1
i 1
f i n1 f i 1
Gi1 Gi1
2
2x
Qi
Boundary conditions must be applied using
method of characteristics
Some useful texts
Rather old- still basis of many commercial programs.
Cunge, J.A., Holley, F.M. and Verwey, A. (1980):
Practical Aspects of Computational River Hydraulics
Mahmood and Yevjevich (1975): Unsteady Flow in
Open Channels - Fort Collins, Colorado
Liggett & Cunge (1975)
Preissmann, A. (1960): Propogation des
Intumescenes dans les Canaue et Rivieres - 1st
Congress de l'Assoc. Francaise de Calcul, Grenoble