Applied Mathematics, 2010, 1, 104-117
doi:10.4236/am.2010.12014 Published Online July 2010 (http://www.SciRP.org/journal/am)
Copyright © 2010 SciRes. AM
Numerical Solution of the Rotating Shallow Water Flows
with Topography Using the Fractional Steps Method
Hossam S. Hassan1, Khaled T. Ramadan1, Sarwat N. Hanna2
1Department of Basic and Applied Science, Arab Academy for Science,
Technology and Maritime Transport, Alexandria, Egypt
2Department of Engineering Mathematics and Physics, Faculty of Engineering,
Alexandria University, Alexandria, Egypt
E-mail: hossams@aast.edu
Received April 16, 2010; revised June 2, 2010; accepted June 12, 2010
Abstract
The two-dimensional nonlinear shallow water equations in the presence of Coriolis force and bottom topog-
raphy are solved numerically using the fractional steps method. The fractional steps method consists of split-
ting the multi-dimensional matrix inversion problem into an equivalent one dimensional problem which is
successively integrated in every direction along the characteristics using the Riemann invariant associated
with the cubic spline interpolation. The height and the velocity field of the shallow water equations over ir-
regular bottom are discretized on a fixed Eulerian grid and time-stepped using the fractional steps method.
Effects of the Coriolis force and the bottom topography for particular initial flows on the velocity compo-
nents and the free surface elevation have been studied and the results are plotted.
Keywords: Shallow Water Equations, Fractional Steps Method, Riemann Invariants, Bottom Topography,
Cubic Spline Interpolation
1. Introduction
Shallow water equations form a set of hyperbolic partial
differential equations that describe the flow below the
pressure surface in the fluid, sometimes, but not neces-
sarily, a free surface. The equations are derived from
depth-integrating the Navier-Stokes equations, in the
case where the horizontal length scale is much greater
than the vertical length scale. They can be used to model
Rossby and Kelvin waves in the atmosphere, rivers,
lakes and oceans in a large domain as well as gravity
waves. The rotating shallow water equations including
topographic effects are a leading order model to study
coastal hydrodynamics on several scales including in-
termediate scale rotational waves and breaking waves on
beaches. Also, they are used with Coriolis forces in at-
mospheric and oceanic modeling, as a simplification of
the primitive equations of atmospheric flow.
Due to the nonlinearity of the model as well as the
complexity of the geometries encountered in real-life ap-
plications such as flow of pollutants, tsunamis, avalan-
ches, dam break, flooding, potential vorticity field…etc,
much effort has been made in recent years to develop
numerical methods to solve the equations approximately.
Bottom topography plays a major role in determining
the flow field in the oceans, rivers, shores, coastal sea
and so on. One of the most important applications of the
shallow water waves is the tsunami waves [1], usually
generated by underwater earthquakes which cause an
irregular topography of increasing or decreasing water
depth. In particular, the main problem in solving the
shallow water equations is the presence of the source
terms modeling the bottom topography and the Coriolis
forces included in the system.
The shallow water equations used in geophysical fluid
dynamics are based on the assumption H/L << 1, where
H and L are the characteristic values for the vertical and
horizontal length scales of motion respectively. These
equations are a two-dimensional hyperbolic system mod-
eling the depth and the depth-averaged horizontal veloci-
ties for an incompressible fluid.
Perhaps, rotation is the most important factor that dist-
inguishes geophysical f luid dyn amics f rom classical flu id
dynamics. If latitudinal varying Coriolis forces are in-
cluded in the shallow water equations, the resulting sys-
tem can support both Rossby and gravity waves. On the
other hand, by neglecting the Coriolis forces in the shal-
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
105
low water system, there is no Rossby-wave solution for
the system.
The Coriolis force is proportional in magnitude to the
flow speed and directed perpendicular to the direction of
the flow. It acts to the left of the flow in the southern
hemisphere and to the right in the northern hemisphere.
A somewhat inaccurate but helpful way to see why the
direction is different in th e two hemisph eres is related to
the principle of conservation of angular momentum. For
a given horizontal motion, the strongest horizontal de-
flection is at the poles and there is no horizontal deflec-
tion at the equator; for vertical motion the opposite is
true. The magnitude of the Coriolis force proportionally
depends upon the latitude and the wind speed. The direc-
tion of the Coriolis force always acts at right angles to
the direction of movement, which is to the right in the
Northern Hemisphere and to the left in the Southern
Hemisp h er e [ 2].
Many authors have used different numerical techniq-
ues to solve the shallow water equations such as finite
volume method, finite element method and fractional
steps method.
Lukácová-Medvid'ová et al. [3] pr esented a new well-
balanced finite volume method within the framework of
the finite volume evolution Galerkin (FVEG) schemes
for the shallow water equations with source terms mod-
eling the bottom top ography and the Coriolis forces.
Gallouët et al. [4] studied the computation of the shal-
low water equations with topography by finite volume
method, in a one-dimensional framework. In their paper,
they considered approximate Riemann solvers. Several
single step methods are derived from this formulation
and numerical results are compared with the fractional
steps method.
Dellar and Salmon [5] derived an extended set of sha-
llow water equations that describe a thin inviscid uid
layer above xed topography in a frame rotating about an
arbitrary axis.
Karelsky et al. [6] executed the generalization of clas-
sical shallow water theory to the case of flows over an
irregular bottom. They showed that the simple self-sim-
ilar solutions that are characteristic for the classical
problem exist only if the underlying surface has a uni-
form slope.
George [7] presented a class of augmented approxim-
ate Riemann solvers for the one-dimensional shallow
water equations in the presence of an irregular bottom,
neglecting the effect of Coriolis force. These methods
belong to the class of finite volume Godunov type
methods that use a set of propagation jump discontinui-
ties, or waves, to approximate the true Riemann solution .
Shoucri [8] applied the fractional steps technique for
the numerical solution of the shallow water equations
with flat bottom in the presence of the Coriolis force.
The method of fractional steps that he presented in his
paper has the great advantage of solving the shallow wa-
ter equations without the iterative steps involved in the
mul t i-dimensional interp olation , an d without the iteratio n
associated with the intermediate step of solving the
Helmholtz equation [9].
Abd-el-Malek and Helal [10] developed a mathemati-
cal simulation to determine the water velocity in the
Lake Mariut, taking into consideration its concentration
and the distribution of the temperature along it, by ap-
plying the fractional steps method for the numerical so-
lution of the shallow water equations.
Shoucri [11] applied the fractional steps technique for
the numerical solution of the shallow water equations to
study the evolution of the vorticity field. The method is
Eulerian [8], and the different variables are discretized
on a xe d grid.
Yohsuke et al. [12] presented two efficient explicit
schemes with no iterative process for the two-dimen-
sional shallow-water equations of a hydrostatic weather
forecast model. One is the directional-splitting frac-
tional-step method, which uses a treatment based on the
characteristics approach (Shoucri [8,11]).The other is the
interpolated differential operator (IDO) scheme (Aoki
[13]), which is one of the multimoment Eulerian schemes.
They compared the forecast geopotential heights ob-
tained from the fractional-step method and the IDO
scheme after 48 h for various resolutions with those of
the referenced scheme by Temperton and Staniforth [14].
Rotating shallow water equations including topogra-
phic terms are numerically dealt with the fractional steps
method. In most real applications there is variable bot-
tom topography that adds a source term to the shallow
water equations. There are several works, where both the
Coriolis forces as well as the bottom topography are
taken into account, see [3,15-17]. A standard and easy
way to deal with these source terms is to treat them in-
dependently by using the fractional steps method. It has
the great advantage of solving the equations without the
iterative steps involved in the multidimensional interpo-
lation problems.
In this work, we apply the fractional steps method to
solve the two-dimensional shallow water equations with
source terms (including the Coriolis force and bottom
topography) for different initial flows observed in the
real-life such as the tsunami propagation wave and the
dam break wave. The objective of the present work is to
simulate the influence of different profiles of the irregu-
lar bottom in case of neglecting and including the effect
of the Coriolis force on the velocity component in the
x-and y-directions, water depth and the free surface ele-
vation for different time. The results are illustrated
graphically for particular initial flows.
2. Mathematical Formulation of the Problem
The shallow water equations (Saint-Venant equations)
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
106
describe the free surface flow of incompressible water in
response to gravitational and rotational accelerations
(Coriolis accelerations), where the vertical depth of wa-
ter is much less than the horizontal wavelength of the
disturbance of the free surface [wave motion]. These
equations are often used as a mathematical model when
numerical methods for solving weather or climate pre-
diction problems are tested. Figure 1 illustrates the shal-
low water model where, “h” denotes the water height
above the non-flat bottom, “H” is the undisturbed free
surface level, ),( yx
denotes the bottom topog-
raphy, “h*” is the water height above the flat bottom and
Hh 
denotes the free surface elevation.
The Continuity equation and the momentum equations
for the two-dimensional shallow water equations system
taking into account the effects of topography and the
Earth’s rotation are formulated by Pedlosky [2] as
Continuity equation :
0
)()(
y
hv
x
hu
t
h (1)
Momentum equations:
vhf
x
hg
x
hg
y
vuh
x
uh
t
uh

)(
2
)()()(
2
2
, (2)
uhf
y
hg
y
hg
y
vh
x
vuh
t
vh

)(
2
)()()(
2
2
. (3)
The model is run over a rectangular domain centered
at the earth origin, such that the x-axis is taken eastward
and the y-axis is taken northward, with u and v the cor-
responding velocity components, respectively, g stands
for the gravitational constant,
sin2f is the Co-
riolis parameter, is the angular velocity of the earth
rotation,
is the geographical latitude of the earth ori-
gin coordinate and ),(ufvf represents the Coriolis
),( yx
h* h
H
Figure 1. Geometry of the shallow water model.
acceleration which is produced by the effect of rotation.
Let the geopotential height be
hg
(4)
Substitution (4) into (1)-(3), yields
0





y
v
y
v
x
u
x
u
t (5)
vf
x
g
xy
u
v
x
u
u
t
u


(6)
uf
y
g
yy
v
v
x
v
u
t
v


(7)
We concern ourselves to approximate the two dimen-
sional shallow water Equations (5)-(7) with source t erm s
modelling the bottom topography and the Coriolis
forces for different initial flows observed in the real-life
such as the tsunami propagation wave and the dam
break wave.
3. Solution of the Problem
Our solution is based on applying the fractional steps
method which was first proposed by Yanenko [18], to
the system of non-linear partial differential Equations
(5)-(7) by splitting the equations into two one-dimen-
sional problems that are solved alternately in x- and y-
directions [19].
3.1. Fractional Steps Method
Step 1. Solve for 2/t
in the x direction (without the
source terms):
Equations (5)-(7) will be
0



x
u
x
u
t, (8)
0

xx
u
u
t
u, (9)
0
x
v
u
t
v. (10)
Equations (8) and (9) can be written as
() 0
xx
RR
u
tx


 

(11)
where 2
x
Ru

are the Riemann invariants.
By applying the classical finite difference scheme to
(11), which is a first-order equation in time and space,
we get
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
107
(,,/2) (,,)
()/(2)
(,,)( , ,)0.
xx
xx
Rxyt tRxyt
utx
R xxytR xyt



 

 

 

(12)
The CFL stability condition is given by()u
1)2/(  xt . Let 2/)( tux , so, (12)
will be


,,/2
()/2,,
x
x
Rxyt t
Rxut yt

 (13)
The right-hand side of (13) is the value of the function
at time tnt  at the departure points of the charac-
teristics. The value of the function at 2/tt  at the
arrival grid points is obtained using cubic spline interpo-
lation from the values of the function at the grid points at
tnt  , (i.e. the value x
R after the time 2/t
from time t is the same value of x
R at that same time
but after the distance x is shifted by x).
Similarly, the solution of (10) for v at 2/tt
is
written as
,, /2/2,,,vxyttvx utyt  (14)
where 2/tux  , i.e. the solution of Equations
(8)-(10) can be written after time tn
)2/1( , n =
1,2,3…as the equality
),,())2/1(,,( tnyx
x
tnyx
x
v
R
v
R

where ()/2
x
xu t
  for the first row and
2/tuxx 
for the second row. Hence the CFL
stability condition will be satisfied automatically at any
time evolution tn  )2/1( .
So, The solutions of x
R, give the values of h and u
after time tn  )2/1( in the
x
direction. So, the
interpolated values in (13) and (14) are calculated using
the cubic spline interpolation, where no iteration is im-
plied in this calculation.
Step 2. Solve for 2/t in the ydirection (with-
out the source terms):
Equations (5)-(7) will be
0



y
v
y
v
t (15)
0
y
u
v
t
u (16)
0

yy
v
v
t
v (17)
Use the results obtained at the end of Step 1, to solve
(15)-(17) for 2/t
.
Equations (15) and (17) can be rewritten as
() 0
yy
RR
v
ty


 

, (18)
where 
2vR y are the Riemann invariants.
The solutions of (16) and (18) at 2/tt  are
,, /2,/2,,uxyttuxy vtt
 (19)
and

,,/2
,( )/2,,
y
y
Rxytt
Rxyv tt


(20)
i.e. after time tn
)2/1( we rewrite Equations (19)
and (20) as
),,())2/1(,,( tnyx
y
tnyx
y
u
R
u
R

where ()/2
y
yv t
  for the first row and
2/tvyy
for the second row.
Again, the solutions of y
R, give the values of h and
v after 2/t
in the y-direction.
Step 3. Solve for t
( for the source terms):
Equations (5)-(7) will be
vf
x
g
t
u

(21)
uf
y
g
t
v

(22)
Use the results obtained at the end of Step 2, to solve
(21) and (22) for t
.
Solutions of (21) and (22) are calculated at t
0
0
]),,(
),,([)/(
)1(,,
Vfntyxx
ntyxxtgU
tnyxu


(23)
0
0
]),,(
),,([)/(
)1(,,
Ufntyyx
ntyxytgV
tnyxv


(24)
where 0
U and 0
V
are the values of u and v at the be-
ginning of Step 3.
Step 4. Repeat Step 2:
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
108
The results obtained at the end of Step 3 are used to
solve for 2/t the equations in the y-direction as in
Step 2.
i.e. ))2/1(,,())1(,,( tnyx
y
tnyx
y
u
R
u
R


Step 5. Repeat Step 1:
The results obtained at the end of Step 4 are used to
solve for 2/t the equations in the x- direction as in
Step 1.
i.e. ))2/1(,,())1(,,( tnyx
x
tnyx
x
v
R
v
R


Shoucri [11] applied the fractional steps method for
the numerical solution of the shallow water equations
over flat bottom. He used the calculated variables (the
height and the velocity field) in studying the evolution of
the potential vorticity field.
Abd-el-Malek and Helal [10] developed a mathemati-
cal simulation to determine the water velocity in the La-
ke Mariut, taking into consideration its concentration and
the distribution of the temperature over it, by applying
the fractional steps method for the numerical solution of
the shallow water equations. The application they pre-
sented over flat bottom requires the variables at two
time-levels tand 2/t according to Strang method
[20] which is more accurate in time, since it has a sec-
ond-order accurate. They proved the convergence of the
fractional steps method and verified that the order of
convergen ce is of the first order.
In this problem, we use the same trend as done by
Abd-el-Malek and Helal [10] and Shoucri [11] in solving
numerically the two-dimensional nonlinear shallow wa-
ter equations in the presence of Coriolis force and bottom
topography using the fractional steps method. We study
the effects of the Coriolis force and the bottom topogra-
phy for particular initial flows on the velocity compo-
nents and the free surface elevation.
3.2. Cubic Spline Interpolation
To approximate the arbitrary functions on closed interv-
als by the aid of the polynomials, we used the most com-
mon piecewise polynomial approximation using cubic
polynomials between each successive pairs of nodes
which is called cubic spline interpolation. A general cu-
bic polynomial involves four constants, so, there is a
sufficient flexibility in the cubic spline procedure to en-
sure that the interpolant is not only continuously differ-
entiable on the interval, but also it has a continuous sec-
ond derivative on the interval, [21].
4. Results and Discussion
We apply the numerical scheme presented in Section 3 to
solve the problem at different values of time for different
initial flows describing the dam break wave and the tsu-
nami propagation wave over two main profiles of the
non-flat bottom in cases of neglecting and including the
Coriolis force. We represent the bottom with two differ-
ent shapes which are a hump and rising hill topography
respectively, as follows.
Hump topography: (25)
Rising up hill topography: (26)
The two forms of )(x
are defined for
100,0
y.
We present the bottom topographies considered by a hump
with maximum height 250 m and rising hill topogra-
phy up to 1850 m, Figur es 2 and 3. The problem is solved
for 2
/806.9 smg
and maximum depth 2
H km.
The scale is based with length 500x km, and width
100
y km. The horizontal grid cell length is
x
1
y km resulting in 100500 grid cells at time
step st1
. We assume the initial velocities (,ux
0)0,,()0,
yxvy .
The initial flows describing the tsunami propagation
wave and the dam break wave are illustrated as follows:

otherwise
xx
x
0
400350,]))375(
25
(cos1[125.0
)(
(25)

otherwise
x
xx
x
0
400,850.1
400200,)))400(
200
(cos1(925.0
)(
(26)
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
109
Tsunami propagation wave:
Tsunamis are long waves generated by submarine
earthquakes. Once they reach the open ocean and travel
through deep water tsunamis which have extremely small
amplitudes but travel fast. Tsunami propagation velocity
can be estimated by using the wave speed equation
CgH. For 2000 m water depth, the speed will be
about 504 km/hour. The most difficult phase of the dy-
namics of tsunami waves deals with their breaking as
they approach the shore. This phase depends greatly on
the bottom bathymetry. As a model of this initial tsun ami
displacement, we consider the wave presented by George
[7] which is given by (27).
Dam break wave:
The dam-break problem is an environmental problem
involving unsteady flows in waterways. The study of
flooding after the dam break is very important because of
the risk to life and property in the potentially inundated
area below the dam. The initial dam break model as-
sumed is given by (28).
Figures 4 and 5 represent the tsunami propag ation an d
dam break waves, respectively, as initial flows. In case
of the dam break problem, the water is assumed to be at
rest on both sides of the dam initially. At 0
t the
dam is suddenly destroyed, causing a shock wave (bore)
travelling downstream with 10
on x > 250 and a
rarefaction wave (depression wave) traveling upstream
with 10
on x < 250. The water pushing down from
above acts somewhat like a piston being pushed down-
stream with acceleration.
ζ
Figure 2. Profile of hump topography.
ζ
Figure 3. Profile of rising up hill topography.
η
Figure 4. Tsunami propagation initial wave.
η
Figure 5. Dam break initial wave.
2
100
10
*0.00235(100), 50150
(,,0)
H, otherwise
x
xeH x
hxy




(27)


500250)250(sgn)01.0(
2500)250(sgn)01.0(
)0,,(*xfor xH
xforxH
yxh (28)
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
110
From the geophysical fluid dynamics point of view, as
a first step toward understanding the role of bottom to-
pography, consider the flow in a periodic zonal channel
with solid boundaries to the North and the South with
idealized topography. Here 01.0f is chosen such
that the domain is resides in the northern hemisphere of
the earth. Sea surface velocities and water height over
the bottom topography at various time and different ini-
tial conditions are illustrated in the cases of neglecting
and including the Coriolis force
Equations (1)-(3) contain the most fundamental bal-
ances of shallow water flows, see [22,23]. The convec-
tive part on the left-hand side is a hyperbolic system of
conservation laws and the source term on the right-hand
side is due to gravitational acceleration and rotation of
the earth (Coriolis force). The steady state results from
a balance between the advection and decay processes
[24].
This suggests that we may have difficulties with a
fractional-step method in order to balance between the
advection terms and the source terms, where we first
solve the advection equation ignoring the reactions and
then solve the reaction equation ignoring the advection
[24]. Even if we start with the exact steady-state solution,
each of these steps can be expected to make a change in
the solution. In principle the two effects should exactly
cancel out, but numerically they typically will not, since
different numerical techniques are used in each step [24].
4.1. The Velocity Component in the x-Direction
Figures 6 and 7 illustrate the behavior of the velocity u
over the topo graphy )(x
represented by a hump given
by (25) for an initial flow represented by the dam break
wave given by (28), in case of neglecting and including
the Coriolis force at 400t, 1200 and 1800 sec res-
pectively. By neglecting the effect of earth rotation (i.e.
no Coriolis force), the coupling between Equations (2)
and (3) due to Coriolis force does no longer exists. Con-
sequently, it is expected that the velocity vector is in the
direction of wave propagation and that minor oscillation s
will appear in the results when the wave propagates
u
× 10
-4
t = 400 sec
u
× 10
-4
t = 1200 sec
u
× 10
-4
t = 1800 sec
Figure 6. Behavior of u with hump topography and dam
break initial flow without the effect of the Coriolis force at t
= 400, 1200 and 1800 sec.
over topography. This is in good agreement with the re-
sult shown in Figure 6. In fact, since the water at
x
250 km seems likely to acquire instantaneously, a ve-
locity different from zero, Figure 6.
u
× 10
-4
t = 400 sec
u
× 10
-4
t = 1200 sec
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
111
u
× 10
-4
t = 1800 sec
Figure 7. Behavior of u with hump topography and dam
break initial flow with the effect the Coriolis force at t = 400,
1200 and 1800 sec.
Figures 8 and 9 illustrate the behavior of the velocity
u over the topography )(x
representing a rising up
function given by (26) fo r an initial flow representing the
tsunami propagation wave given by (27), in case of ne-
glecting and includ ing the Coriolis force at t = 400, 1200
and 3000 sec respectively. As mentioned before, in case
of neglecting the Coriolis force, the velocity vector is in
the direction of the wave propagation and minor oscilla-
tions will appear in the results when the wave propagates
over the topography, Figure 8.
u
× 10
-4
t = 400 sec
u
× 10-4
t = 1200 sec
u
× 10
-3
t = 3000 sec
Figure 8. Behavior of u with rising up topography and tsu-
nami propagation initial flow without the effect of the Cori-
olis force at t = 400, 1200 and 3000 sec.
u
× 10
-4
t = 400 sec
u
× 10
-4
t = 1200 sec
u
× 10
-3
t = 3000 sec
Figure 9. Behavior of u with rising up topography and tsu-
nami propagation initial flow with the effect the Coriolis
force at t = 400, 1200 and 3000 sec.
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
112
4.2. The Velocity Component in the y-Direction
Figures 10 and 11 illustrate the behavior of the velocity
v over the topography )(x
representing a hump given
by (25) for an initial flow representing the dam break
wave give n by (28) and a rising up function given by
(26) for an initial flow representing the tsunami propa-
gation wave given by (27) with the effect of Coriolis
force.
4.3. The Water Height h*
Figures 12 and 13 illustrate the behavior of the water
height h* over the topography )(x
representing by a
hump given by (25) for an initial flow representing by
v
× 10-4
t = 400 sec
v
× 10
-3
t = 1200 sec
v
× 10
-3
t = 1800 sec
Figure 10. Behavior of v with hump topography and dam
break initial flow with the effect of the Coriolis force at t =
400, 1200 and 1800 sec.
the dam break wave given by (28) in case of neglecting
and including the Coriolis force at t = 400, 1200 and
1800 sec respectively at y = 0.
v
× 10
-4
t = 400 sec
v
× 10
-4
t = 1200 sec
v
× 10
-4
t = 3000 sec
Figure 11. Behavior of v with rising up topography and
tsunami propagation initial flow with the effect the Coriolis
force at t = 400, 1200 and 3000 sec.
h
*
t = 400 sec
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
113
h
*
t = 1200 sec
h
*
t = 1800 sec
ζ
Figure 12. Behavior of h* with hump topography and dam
break initial flow without the effect of the Coriolis force at t
= 400, 1200 and 1800 sec.
Figures 14 and 15 illustrate the behavior of the water
height h* over the topography )(x
representing a
rising up function giv en by (26) for an initial flow repre-
senting by the tsunami propagation wave given by (27)
in case of neglecting and including the Coriolis force at t
= 400, 1200 and 3000 sec respectively at y = 0.
h
*
t = 400 sec
h
*
t = 1200 sec
h
*
t = 1800 sec
ζ
Figure 13. Behavior of h* with hump topography and dam
break initial flow with the effect the Coriolis force at t = 400,
1200 and 1800 sec.
h
*
t = 400 sec
h
*
t = 1200 sec
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
114
h
*
t = 3000 sec
ζ
Figure 14. Behavior of h* with rising up topography and
tsunami propagation initial flow without the effect of the
Coriolis force at t = 400, 1200 and 3000 sec.
h
*
t = 400 sec
h
*
t = 1200 sec
h
*
t = 3000 sec
ζ
Figure 15. Behavior of h* with rising up topography and
tsunami propagation initial flow with the effect the Coriolis
force at t = 400, 1200 and 3000 sec.
Since the water at 250
x km seems likely to ac-
quire instantaneously, it is plausible that a shock would
be created instantly on the upstream s ide and a relatively
small propagated wave in the negative X-direction due to
bottom topography as shown in the circle in Figure 12.
By taking the Coriolis force ),( u
f
v
f
into con-
sideration, it becomes responsible for the oscillatory mo-
tion according to the solutions of Equations (23) and (24)
obtained in Step 3 in the fractional steps method. The
coupling between the velocity components u and v
causes the deflection of fluid parcels which are oscillat-
ing back and forth in the direction of wave motion and
causes gravity waves to disperse as shown in Figures 13
and 15. The magnitude of this deflection proved to be
independent of the absolute depth and depends only on
the slope of the bottom, as demonstrated in the momen-
tum Equations (2) and (3).
In case of the tsunami wave, the initial waves split into
two similar waves, one propagates in the positive x-di-
rection and the other in the negative x-direction. The hei-
ght above mean sea lev el of the two oppositely traveling
tsunamis is approximately half that of the original tsu-
nami, as shown in Figure 14. This ha ppened becau se the
potential energy that results from pushing water above
mean sea level is transferred to horizontal propagation of
the tsunami wave (kinetic energy) (i.e. the tsunami con-
verts potential energy into kinetic energy). It is well
known when the local tsunami travels over the continen -
tal topography, th at the wav e a mplitude in creases and the
wavelength and velocity decreases, which results in
steepening of the leading wave, see Figure 14.
In case of th e dam bre ak proble m, the ef fe ct of the Co-
riolis force on the water height at a sequence of times
after the breakage of a dam causes a wave travelling do-
wnstream and a wave travelling upstream as shown
in Figure 13. The rarefaction wave that developed in
Figure 13 is overtaken by dispersive wave which should
form shocks on both sid es of propagation.
In case of the tsunami waves, they propagate in coh-
erent wave packets, with little loss of amplitude over
very long distances as shown in Figure 15. As the water
depth decreases, the wave amplitude increases and the
wavelength and velocity decreases, resulting in steepen-
ing of the dispersed wave.
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
115
4.4. The Free Surface Elevation (η)
Figures 16 and 17 illustrate the effect of the Coriolis
force on the free surface elevation
when travelling
over the topography )(x
representing a hump given
by (25) and a rising up function given by (26) for an ini-
tial flow representing the dam break wave given by (28)
at 1200tsec.
The 500 kilometers long dam wall, which runs parallel
to the y-axis is 10 meters w ide and is centered at x = 250
kilometers. Along the boundaries, at x = 0 kilometers and
x = 500 kilometers and
is fixed at the upstream and
downstream water depth respectively. All other bounda-
ries are considered as reflective boundary conditions. In
case of neglecting the Coriolis force, a shock front al-
ways exists as shown in Figure 16, while the free surface
elevation
under the effect of the Coriolis force has no
shock travelling downstream and hence travelled in the
direction of propagation over the topography )(x
as
seen in Figure 17.
As seen in Figure 16, the jump at the gen erated shock
first decreases and then it increases as the shock ap-
proaches the bottom topography. It is observed that the
increasing in the shock wave when travelling over the
rising up hill topography is greater than when travelling
over hump topography.
It can be observed that the free surface elevation
in Figure 17 clearly differs from that in Figure 16 due to
the effect of the Coriolis force which is responsible for
the oscillatory motion in the direction of wave motion
which causes gravity waves to disperse.
Figures 18 and 19 illustrate the effect of the Coriolis
force on the free surface elevation
when travelling
over the topography )(x
representing a hump given
by (25) and a rising up function given by (26) for an ini-
tial flow representing the tsunami propagation wave
given by (27), at t = 2400 sec.
The effect of the Coriolis force on the tsunami propa-
gation wave is responsible for some part of the energy,
which is transmitted to the ocean with the seismic bottom
motions, to accumulate in the region of the disturbance.
This leads to a reduction of the barotropic wave energy
and tsunami amplitude. The direction of the tsunami ra-
diation varies and the energy flow transferred by the
waves is redistributed.
The effect of Coriolis force on transoceanic tsunami
with and without Coriolis terms shows differences in
wave height but not much difference in arrival time as
observed in Figure 18 and Figure 19.
η
(a)
η
(b)
Figure 16. Behavior of η for the dam break initial flow without the effect of the Coriolis force at t = 1200 sec over: (a) Hump
topography; (b) Rising up hill topography.
η
(a)
η
(b)
F igure 17. Behavior of η for the dam break initial flow with the effect of the Coriolis force at t = 1200 sec over: (a) Hump
topography; (b) Rising up hill topography.
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
116
η
(a)
η
(b)
Figure 18. Behavior of η for the tsunami propagation initial flow without the effect of the Coriolis force at t = 2400 sec over:
(a) Hump topography; (b) Rising up hill topography.
η
(a)
η
(b)
Figure 19. Behavior of η for the tsunami propagation initial flow with the effect of the Coriolis force at t = 2400 sec over:
(a) Hump topography; (b) Rising up hill topography.
Therefore the dispersion effects become more signifi-
cant as the wave energy is more spatially spread out and
scattered. The variability of the bathymetry is also quite
variable being an important parameter of wave disper-
sion as it controlled in the speed and the height of the
tsunami wave when there was a sudden change in the
depth of water as seen in Figure 18 and Figure 19.
As seen from Figures 16-19, the wave propagated
faster when travelling over the hump topography than
travelling over the increasing hill topography. This was
expected due to the wave speed equation HgC
which is proportional to the water depth and inversely
proportional to the bottom bathymetry.
5. Conclusions
The fractional steps method for the numerical solution of
the shallow water equations is applied to study the evo-
lution of the height and the velocity field of the flow un-
der the effect of Coriolis force and bottom topography.
The method consists of splitting the equations and suc-
cessively integrating in the x-and y-directions along the
characteristics using the Riemann invariants, associated
with the cubic spline interpolation. In this work, we ap-
plied the fractional steps method to solve the two-dim-
ensional shallow water equations with source terms (in-
cluding the Coriolis force and bottom topography) for
two different initial flows namely dam break wave and
tsunami wave. The presence of the Coriolis force in the
shallow water equations causes the deflection of fluid
parcels in the direction of wave motion and causes grav-
ity waves to disperse. As water depth decreases due to
bottom topography, the wave amplitude increases, the
wavelength and wave speed decreases resulting in steep-
ening of the wave. Th e effect of the Coriolis force is res-
ponsible for the oscillatory motion in the direction of
wave motion which causes gravity waves to disperse.
The overall performance of the fractional steps method
in solving the shallow water equations with source terms
is particularly attractive, simple, efficient and highly ac-
curate as our results verified the reality about the nature
of the dam break problem and the tsunami propagation
wave. In future, we shall apply the fractional steps
method to the shallow water equations with the source
term including Coriolis force and a movable topography
as it appears in oceanographic modeling as well as in the
river flow modeling is underway and comp are the results
analytically using Lie-group method, [25-29].
6. Acknowledgements
The authors would like to express their sincere thank the
H. S. HASSAN ET AL.
Copyright © 2010 SciRes. AM
117
reviewers for suggesting certain changes in the original
manuscript, for their valuab le co mments wh ich improv ed
the paper and for their great interest in that work.
7. References
[1] D. G. Dritschel, L. M. Polvani and A. R. Mohebalhojeh,
“The Contour-Advective Semi-Lagrangian Algorithm for
the Shallow Water Equations,” Monthly Weat her Review,
Vol. 127, No. 7, 1999, pp. 1551-1564.
[2] J. Pedlosky, “Geophysical Fluid Dynamics,” Springer,
New York, 1987.
[3] M. Lukácová-Medvid'ová, S. Noelle and M. Kraft, “Well-
Balanced Finite Volume Evolution Galerkin Methods for
the Shallow Water Equations,” Journal of Computational
Physics, Vol. 221, No. 1, 2007, pp. 122-147.
[4] T. Gallouët, J. M. Hérard and N. Seguin, “Some Ap-
proximate Godunov Schemes to Compute Shallow-Water
Equations with Topography,” Computers & Fluids, Vol.
32, No. 4, 2003, pp. 479-513.
[5] P. J. Dellar and R. Salmon, “Shallow Water Equations
with a Complete Coriolis Force and Topography,” Phys-
ics of Fluids, Vol. 17, No. 10, 2005, pp. 106601-106619.
[6] K. V. Karelsky, V. V.Papkov, A. S. Petrosyan and D. V.
Tsygankov, “Particular Solution of the Shallow-Water
Equations over a Non-Flat Surface,” Physics Letters A,
Vol. 271, No. 5-6, 2000, pp. 341-348.
[7] D. L. George, “Augmented Riemann Solvers for the Shal-
low Water Equations over Variable Topography with Steady
States and Inundation,” Journal of Computational Phys-
ics, Vol. 227, No. 6, 2008, pp. 3089-3113.
[8] M. Shoucri, “The Application of a Fractional Steps Method
for the Numerical Solution of the Shallow Water Equa-
tions,” Computer Physics Communications, Vol. 164, No.
1-3, 2004, pp. 396-401.
[9] A. Stainiforth and C. Temperton, “Semi-Implicit Semi-
Lagrangian Integration Scheme for a Baratropic Finite-
Element Regional Model,” Monthly Weather Review, Vol.
114, No. 11, 1986, pp. 2078-2090.
[10] M. B. Abd-el-Malek and M. H. Helal, “Application of the
Fractional Steps Method for the Numerical Solution of
the Two-Dimensional Modeling of the Lake Mariut,” Ap-
plied Mathematical Modeling, Vol. 33, No. 2, 2009, pp.
822-834.
[11] M. Shoucri, “Numerical Solution of the Shallow Water
Equations with a Fractional Step Method,” Computer Phys-
ics Communications, Vol. 176, No. 1, 2007, pp. 23-32.
[12] I. Yohsuke, A. Takayuki and M. Shoucri, “Comparison of
Efficient Expli cit Schemes for Shallow-Water Equations-
Characteristics-Based Fractional-Step Method and Mul-
timoment Eulerian Scheme,” Journal of Applied Meteor-
ology and Climatology, Vol. 46, No. 3, 2007, pp. 388-395.
[13] A. Takayuki, “Interpolated Differential Operator (IDO)
Scheme for Solving Partial Differential Equations,” Com-
puter Physics Communications, Vol. 102, No. 1-3, 1997,
pp. 132-146.
[14] C. Temperton and A. Stainiforth, “An Efficient Two-
Time-Level Semi-Lagrangian Semi-Implicit Integration
Scheme,” Quarterly Journal of the Royal Meteorological
Society, Vol. 113, No. 477, 1987, pp. 1025-1039.
[15] E. Audusse, F. Bouchut, M. Bristeau, R. Klein and B. Per-
thame, “A Fast and Stable Well-Balanced Scheme with
Hydrostatic Reconstruction for Shallow Water Flows,”
SIAM Journal on Scientific Computing, Vol. 25, No. 6,
2004, pp. 2050-2065.
[16] M. E. Talibi and M. H. Tber, “On a Problem of Shallow
Water Type,” Electronic Journal of Differential Equa-
tions, Vol. 11, 2004, pp. 109-116.
[17] V. R. Ambati and O. Bokhove, “Space-Time Discon-
tinuous Galerkin Discretization of Rotating Shallow Wa-
ter Equations,” Journal of Computational Physics, Vol.
225, No. 2, 2007, pp. 1233-1261.
[18] N. N. Yanenko, “The Method of Fractional Steps, the
Solution of Problems of Mathematical Physics in Several
Variables,” Springer-Verlag, Berlin, 1971.
[19] D. R. Durran, “Numerical Methods for Wave Equations
in Geophysics Fluid Dynamics,” Springer-Verlag, New
York, 1999.
[20] G. Strang, “On the Construction and Comparison of Fi-
nite Difference Schemes,” Society for Industrial and Ap-
plied Mathematics, Journal for Numerical Analysis, Vol.
5, No. 3, 1968, pp. 506-517.
[21] R. L. Burden and J. D. Faires, “Numerical Analysis,”
PWS Publishing Company, Boston, 1993.
[22] S. Noelle, N. Pankratz, G. Puppo and J. R. Natvig,
“Well-Balanced Finite Volume Schemes of Arbitrary Or-
der of Accuracy for Shallow Water Flows,” Journal of
Computational Physics, Vol. 213, No. 2, 2006, pp. 474-499.
[23] S. Noelle, Y. Xing and C. Shu, “High-Order Well-Bal-
anced Finite Volume WENO Schemes for Sha llow Water
Equation with Moving Water,” Journal of Computational
Physics, Vol. 226, No. 1, 2007, pp. 29-58.
[24] R. LeVeque, “Finite Volume Methods for Hyperbolic Pro-
blems,” Cambridge University Press, Cambridge, 2004.
[25] Y. Z. Boutros, M. B. Abd-el-Malek, N. A. Badran and H.
S. Hassan, “Lie-Group Method for Unsteady Flows in a
Semi-infinite Expanding or Contracting Pipe with Injec-
tion or Suction through a Porous Wall,” Journal of Com-
putational and Applied Mathematics, Vol. 197, No. 2, 2006,
pp. 465-494.
[26] Y. Z. Boutros, M. B. Abd-el-Malek, N. A. Badran and H.
S. Hassan, “Lie-Group Method Solution for Two-dimen-
sional Viscous Flow between Slowly Expanding or Con-
tracting Walls with Weak Permeability,” Applied Mathe-
matical Modelling, Vol. 31, No. 6, 2007, pp. 1092-1108.
[27] Y. Z. Boutros, M. B. Abd-el-Malek, N. A. Badran and H.
S. Hassan, “Lie-Group Method of Solution for Steady
Two-Dimensional Boundary-Layer Stagnation-Point Flow
Towards a Heated Stretching Sheet Placed in a Porous
Medium,” Meccanica, Vol. 41, No. 6, 2007, pp. 681-691.
[28] M. B. Abd-el-Malek, N. A. Badran and H. S. Hassan,
“Lie-Group Method for Predicting Water Content for
Immiscible Flow of Two Fluids in a Porous Medium,”
Applied Mathematical Sciences, Vol. 1, No. 24, 2007, pp.
1169-1180.
[29] M. B. Abd-el-Malek and H. S. Hassan, “Symmetry
Analysis for Solving Problem of Rivlin-Ericksen Fluid of
Second Grade Subject to Suction,” submitted for publica-
tion.