Unsteady flow

using
Method of
Characteristics
 Definition and History
 Briefly called MOC
 It is a graphical procedure for solving pde.
 Developed by Monge in 1789 and first used
by Massau in 1889 and Craya in 1946 for
analyzing surges in open channels.
 The concept of MOC is important in
understanding the boundary conditions,
flow control and wave propagations.
 Basic Revision before
starting MOC
Solve: What Do You Mean?!!

du
 5s
ds

• 1st order Ode.

• Single independent variable, s.
• Single dependent variable, u.
 Basic Revision before
starting MOC
du means Get the unknown “u” for different s
Solve:  5s
ds within the domain

 du   [Link] s2
u  5  const.
2
Basic Revision before
starting MOC
We can easily solve ODEs
like the one before.
But… What about PDEs?
 Basic Revision before
starting MOC
U . y y
 0
x t
Means:
U  (U ) y
For every x and t U  g.  g.( S o  S e )
t x x

for ve
get u and y

l
(y, u)

So
river domain

t
(x, t)

X
 Method of Characteristics
• The method of characteristics (MOC) is an old technique used to
solve first order pde.
• The main idea of MOC is to discover lines within the solution
domain (called characterstic lines or characteristics) along which
the PDE degenerates into an ordinary differential equation (ODE).
• In many applications, along the characteristics, the solution is
constants.
 MOC (in Brief)
MOC discovers paths within the domain where
The PDE reduces to ODE.
These paths are called “Characteristics”.

Characteristic lines
t

X
 An Example for (MOC)
y To apply MOC,
s
I have to change this pde  ode
through out certain paths within the
• 1st domain.
order pde
• Linear equation x
• Independent variables (x,y)
dy b( x, y )

• Dependent variables u
dx a ( x, y )

If and only if

u
 MOC for St. Venant Eqs.
y  (U ) ( y ) y  (U ) ( y)
y U 0 ( y U )0
t x x t x x
U
t
U
 (U )
x
g
( y )
x
 g (S o  S f ) +
 U  (U )   ( y) g ( y) 
 t  (U   y )     t  (U  )   g (S o  S f )
 x    x 

dU d ( y) dx dx g
  g (S o  S f )  (U  y )  (U  )
dt dt dt dt 

g
2  or
g y
1, 2 ~ g
c 1, 2   or
y
 MOC for St. Venant Eqs.
y  (U ) ( y)
y U 0
c  gy t x x

U  (U ) ( y)
U g  g (S o  S f )
c 2  gy t x x

2 dc dy
c 
g dt dt

dU g d ( y ) Over the path dx

~  g (So  S f )  (U  C )
dt
dt c dt

dU g d ( y ) Over the path dx

~  g (So  S f )  (U  C )
dt c dt dt
 MOC for St. Venant Eqs.
y  (U ) ( y)
y U 0
t x x

U  (U ) ( y)
U g  g (S o  S f )
t x x

d (U  2c~ ) Over the path

 g (S o  S f ) dx
 (U  C )
dt dt

d (U  2c~ ) Over the path

 g (S o  S f ) dx
 (U  C )
dt dt
 Applications
Assume: A horizontal channel with a frictionless bed.
d (U  2c~ ) Over the path dx
0  (U  C )
dt dt

d (U  2c~ ) Over the path dx

0 dt
 (U  C )
dt

Over the path

(U  2c~ )  const.
dx
 (U  C )
dt

Over the path dx

(U  2c~ )  Const.  (U  C )
dt
 Applications
Over the path
(U  2c~ )  const.
dx
 (U  C )
dt

Over the path dx

(U  2c~ )  Const.  (U  C )
dt

P
t
C+ C-

to
R L

x
Applications (Neglecting Friction)
P
t
C+ C-

to
L R

UL+2cL=Up+2cp
Two equations in two
unknowns (Up, cp)
UR-2cR=Up-2cp
 Applications (Considering Friction)
t 8 Assume we know the
solution at points 3 and 5
6
We need to get the
7 5 unknowns at point 6
4 3
2
to 1
x
For C+
(v6+2c6)-(v5+2c5)=t.g(So-Sf5) After addition get:
v6=0.5(v5+v3)+(c5-c3)+t.g(So-Sfav)
For C- After subtracting get:
(v6-2c6)-(v3-2c3)=t.g(So-Sf3) c6=0.25(v5-v3)+0.5(c5+c3)
 What dose the last
example mean?!!

Conditions at point P are determined solely by the conditions

bounded within the two characteristics C1 and C2.

The conditions at point Q determine the

conditions within the area bounded by
the characteristic lines
 How to get C+, C- paths
t
tP P (xP, tP)
C+ C-
(xL, tL) (xR, tR)
L xP R x
For C+ For C-
dx dx Get
 (U  c)  (U  c)
dt dt
Xp, tp
xp xp P
P

 dx   (U  c)dt  dx   (U  c)dt
xR R
xL L
P P

x p  xL   (U  c)dt x p  xR   (U  c) dt
L R

x p  xL  (t p  t L )[ (U L  cL )  (U p  c p )]
1 x p  xR  (t p  t R )[
1
(U R  cR )  (U p  c p )]
2 2
 Example #1
 Consider a long channel Assumptions:
with the following flow
measurements at two gauging stations (at t=0).
Calculate: the location, -time
So=Se=0.0
and flow properties
at the intersection of the- characteristics issuing
from the two gauge [Link] straight lines for C+,C-

Station 1 Station 2
X (km) 7.1 8.25
Water depth (m) 2.2 2.45
Flow velocity .5 .29
(m/s)
 Example #1 (Solution)
t P

C+ C-

L R
7100 8250 x
Vp+2cp=VL+2cL Vp=0.5(VL+VR)+(cL –cR)
Vp-2cp=VR-2cR cp=0.5(cL+cR)+0.25(VL –VR)

cL  9.81* 2.2  4.646m / s Vp= 0.5(.5+.29)+(4.646 –4.902)=0.13814 m/s

cR  9.81* 2.45  4.902m / s cp= 0.5(4.646+4.902)+0.25(0.5 –0.29)=4.826 m/s
 Example #1 (Solution)
t P Vp= 0.13814 m/s
cp= 4.826 m/s
C+ C-
1 t 1
L Vp+cp Vp-cp R
7100 8250 x
(Vp+cp)t (Vp-cp)t
(0.13814+4.826)t (0.13814-4.826)t
= 4.96414t =-4.6878 t

= 4.96414t + 4.6878 t
= 9.6519 t

t = (8250-7100) / 9.6519
t = 119.14 sec,
Xp= 7100 +4.96519x 119.14 = 7691.6 m
 Example#2
d (U  2c~ )
Assumptions:
 Consider a long channel with the following flow  g (S o  S f )
dt
measurements at two gauging stations
- So=0 2c~t=0).
d (U (at )
Calculate: the location, time and flow properties g (S o  S f )
dtissuing
at the intersection of the characteristics
- n=.03
from the two gauge stations.
- Assume straight lines for C+,C-

Station 1 Station 2
X (km) 7.1 8.25
Water depth (m) 2.2 2.45
Flow velocity .5 .29
(m/s)
 Example #2 (Solution)
t P
V = (1/n) y(2/3) Sf(0.5)
C+ C- Sf= (V.n)2/y(4/3)

L R
7100 8250 x
d (U  2c~ )
dt
 g (S o  S f ) Vp=0.5(VL+VR)+(cL –cR) + gt[So-0.5(SfL+SfR)]
d (U  2c~ )
dt
 g (S o  S f ) cp=0.5(cL+cR)+0.25(VL –VR) - gt(SfL-SfR)
Vp= 0.5(.5+.29)+(4.646 –4.902)+9.81t[10-5 -.5((0.5x0.03)2/2.2(4/3)+(0.29x0.03)2/2.45(4/3))]
Vp = 0.0898 m/s
cp= 0.5(4.646+4.902)+0.25(0.5 –0.29) -9.81t[((0.5x0.03)2/2.2(4/3)-(0.29x0.03)2/2.45(4/3))]
cp =4.76 m/s
t = 120.78 sec,
yp = 2.31 m
Xp= 7685.85 m
 Example#3
 Consider a long channel Assumptions:
with the following flow
measurements at two gauging stations (at t=0).
Calculate: the location, -time
So=0and flow properties
at the intersection of the- characteristics issuing
from the two gauge stations.n=.03
- Consider curvature for C+,C-

Station 1 Station 2
X (km) 7.1 8.25
Water depth (m) 2.2 2.45
Flow velocity .5 .29
(m/s)
 Effect of Fn on B.C.
Sub Critical Critical Super Critical

Fn  1 Fn  1 Fn  1

V V V
 1  V  gy  1  V  gy  1  V  gy
gy gy gy

V  gy  ve V  gy  0 V  gy   ve
V  c  ve V c  0 V  c  ve

t P P P
t t
R L R R
 Interior and Boundaries (sub-
One boundary condition
critical)
One boundary condition
On the left side (is required) On the right side (is required)

t3

t2
C- C+
t1

Initial values (known)

 Interior and Boundaries
(super-critical)
TWO boundary conditions Boundary conditions On the
right side (is NOT required)
On the left side (are required)

t3

t2

t1

Initial values (known)

We Can Easily Notice that:
Super Critical Flow is controlled by u/s sections.
Sub Critical Flow is controlled by d/s sections.
More Examples
 New Example
Open-channel flow  By Subhash C. Jain, Subhash
Solution
 Example #4
 A pump is installed to pump water from a
channel intake to a WTP. If it is assumed
initially that the water is stationary
throughout the intake channel with a water
depth of 2m. Estimate the maximum
discharge that can be pumped from the
channel at the early beginning of pumping.

Plan view
 Example#5
 Water flows in an irrigation canal at steady
state (V=.9m/s, d=1.64m). The flume is
assumed smooth and horizontal. The flow
is controlled by a d/s gate. At t=0, the gate
is very slowly raised and the water depth
u/s of the gate decreases at a rate of
5cm/min. until the water depth becomes
0.85m.
 Calculate: the discharge per unit width at
the gate at t=10min.
Numerical Methods for Solving
St. Venant Equations

Method of Characteristics Finite Difference Method

 Finite Difference Method
t

t+1

t
i-1 i i+1
x
y yit 1  yit y yit1  yit1 U U it1  U it1
  
t i t x i 2x x i 2 x
 Finite Difference (Cont.)
Continuity equation :
 By using FTCS scheme, the continuity
equation can be discretized using the
finite difference technique (assuming
rectangular cross sections) as follow :

yin 1  yin 
2x

t n n

yi (U i 1  U in1 )  U in ( yin1  yin1 )
 Finite Difference (Cont.)
 8.1.2 Momentum equation :
 The can also be discretized using the FTCS
scheme as follow :
1  n 1 n 1 2 
U in 1   i  i  4in 1  
2 

in 1 
 R
hi 
n 1 4 / 3

2
n g t

  U in 
t
2x

U in (U in1  U in1 )  g ( yin1  yin1 )  gtS o 
 Example #6
A rectangular prismatic channel, made of poor quality
concrete, has a width of B=5m and a normal depth of 1.2m.
During a flood event, the discharge in the channel increases
linearly to Q=50m3/s in a time period of 20 min. and
decreases back linearly to its initial value during a time
period of 60 min. The channel has So of 0.001 and
Manning’s coefficient n=0.02 (metric). The length of the
channel reach to be studied is 3Km.
Using the explicit finite difference method, determine :
the hydrograph at two stations (x=1.5, 3.0 Km) and compare
them with the initial hydrograph.
The times at which Q=Qmax. at x=1.5, 3.0 Km.
The rating curve at x=L/4, L/2, 3L/4, L, where L=3000m