UNCORRECTED PROOF
5. MODEL VERIFICATION1
CNM 1061
M.-C. HUNG ET AL.
5. MODEL VERIFICATION 1
5.1. Wind-induced flow
Wind-induced flow in a closed basin in which wind shear acting on the water surface results in a 3
drift current in the direction it blows and a bottom return flow in the opposite direction is chosen to verify the model’s accuracy and applicability. The analytical solution for steady laminar flow 5
under no-slip condition, given by Heaps[16], for constant vertical eddy viscosity is u(z)=
s
4v
z(3z −2H)
H (43)
7
where u(z) is the velocity at any depth z; H the water depth; v the vertical eddy viscosity;sthe wind stress acting on the surface; and the water density.
9
This numerical experiment is carried out to simulate a steady wind-induced flow in a uniform rectangular channel of length 400 m, width 160 m, and depth 2 m. The water initially at rest was 11
subject to a uniform wind stress, increasing linearly from 0 N/m2at t=0s to 0.5N/m2at t=10s, and maintaining 0.5N/m2at t>10s. The following parameters were adopted: x =20m, y =8m, 13
=1000kg/m3,v=0.015m2/s.
A constant viscosity leads to a parabolic velocity profile in the case of laminar flow. Effect of 15
number of layers was also carried out with four cases of 5, 10, 20, and 30 layers. Comparison of model predictions at the middle length of the basin to the analytical solution is shown in Figure 4.
17
The major improvement of predictions of velocity profile takes place between 5 and 10 layers, and no more significant difference takes place between predictions over 20 layers. Comparison of 19
the u profile at middle length using 20 layers to the analytical solution is as expected. The upper layer currents are in the wind direction while the lower layer currents are in the opposite direction.
21
It can be seen that the prediction has a slightly larger velocity. The predicted depth with zero
Figure 4. Effect of number of layers on velocity profile prediction.
23
UNCORRECTED PROOF
CNM 1061
APPROACH FOR SHALLOW WATER FREE-SURFACE FLOW COMPUTATION 15
Figure 5. Effect of number of layers on water depth prediction.
Figure 6. Variation of velocity in the wind-induced flow at middle length(x =200m).
velocity is at 13 of water depth under the free surface. It is consistent with the analytical solution.
1
The wind setup of the water surface has the same tendency as shown in Figure 5. With increase in the surface wind shear or decrease in the eddy viscosity, the maximum Reynolds number tested 3
for limited stable condition can be raised up to around 300.
The evolution of the wind-induced circulation could be classified into two steps: (a) surface 5
flow and water setup induced by wind shear and (b) bottom reverse flow induced by pres-sure difference. The transient process was also carried out to examine the physical mecha-7
Copyrightq 2007 John Wiley & Sons, Ltd. Commun. Numer. Meth. Engng (2007)
UNCORRECTED PROOF
16
CNM 1061
M.-C. HUNG ET AL.
Figure 7. Variation of water depth in the wind-induced flow.
nism. Figures 6 and 7 show the simulated evolution of velocity profile at middle length and 1
water surface setup varying with time with 20 layers, respectively. From Figure 6 it is obvious that the shear-induced surface flow in the beginning and the pressure-induced reverse flow at 3
middle length were successfully predicted by the proposed model. These simulation results include surface setup and pressure-induced circulation, two significant components of the physical 5
mechanism.
5.2. Effect of viscosity distribution on wind-induced flow 7
Owing to the distribution of vertical eddy viscosities existing in natural flows not being uniform, three more types of viscosity distribution were considered to further investigate the capability 9
of the proposed quadratic shape function for solving free-surface flows. The parabolic law of eddy viscosity distribution is the most frequently assumed type in natural flows [15]. It has the 11
maximum viscosity at middle depth and the minimum near the surface and the bottom. Figure 8 shows the simulated velocity profile with parabolic viscosity distribution compared to the one with 13
uniform viscosity. The larger velocity gradient is found near both surface and bottom due to its smaller viscosity. Contrarily, the smaller velocity gradient is found at middle depth due to its larger 15
viscosity. The simulated profile behaves in a well-consistent tendency with its physical mechanism.
To further examine the flexibility of the auxiliary shape function, the trapezoidal distribution of 17
viscosity was also conducted. The comparison of the velocity profiles due to change of the viscosity distribution is shown in Figure 9. For the case of trapezoidal type I with increased viscosity from 19
top to bottom, the simulated velocity profile is well adapted as its gradient decreases from top to bottom by the quadratic shape function. By preserving the condition of zero net flow, the point 21
of zero is risen a little compared with the uniform law. Figure 9 also shows that the quadratic shape function can also adapt the velocity profile to the trapezoidal type II of viscosity distribution 23
as well.
UNCORRECTED PROOF
CNM 1061
APPROACH FOR SHALLOW WATER FREE-SURFACE FLOW COMPUTATION 17
Figure 8. Velocity profiles corresponding to their own viscosities (uniform vs parabolic).
Figure 9. Velocity profiles corresponding to their own viscosities (uniform vs trapezoidal).
Copyrightq 2007 John Wiley & Sons, Ltd. Commun. Numer. Meth. Engng (2007)
UNCORRECTED PROOF
18
CNM 1061
M.-C. HUNG ET AL.
Figure 10. Three types of viscosity variation distribution employed in the study.
5.3. Examination of proposed eddy viscosity relation 1
Because there is no net flow in the wind-induced flow in a closed basin, the non-uniformity of viscosity distribution does not affect the mean water depth but the velocity profile. As far as 1D 3
or plane 2D is concerned, the normal depth is well known as a function of flow discharge and roughness in natural rivers. However, it is not sufficient in the 3D free-surface flow computations 5
because the water depth and velocity profile are also affected by the eddy viscosity. The following case is designed hypothetically to investigate the effect of viscosity variation on channel flow.
7
A steady flow is carried in a uniform rectangular channel with length 5000 m, width 100 m, slope 0.0005, Manning’s roughness n of 0.035, and discharge of 3.987 cm per unit width. The 9
normal depth determined by Manning’s equation is 3.0 m. The slip condition was used at the side walls but no-slip condition at the channel bed. The water body was divided into four layers of 11
equal thickness. The following parameters were specified:x =500m, y =10m, =1000kg/m3. Upstream, the lateral and vertical velocities v and w were set as zero, and the longitudinal 13
velocity u was equal to the inflow discharge divided by the instantaneous flow area. Downstream, zero longitudinal velocities gradient*u/*x =0, *v/*x =0 were applied. Three types of viscosity 15
variation distribution employed as shown in Figure 10 and the coefficient in Elder’s relation can be expressed as follows:
17
Type I: = v
H u∗=0.0304
UNCORRECTED PROOF
CNM 1061
APPROACH FOR SHALLOW WATER FREE-SURFACE FLOW COMPUTATION 19
Figure 11. Velocity profiles corresponding to the three types of viscosity variation distribution.
Type II: =0.0304
0.5+4 z H
1− z
H
Type III: = z H
1− z
H
, =0.41
where type I is calculated by Equation (42); type III is a parabolic distribution function proposed 1
by Elder[14]; and type II is an arbitrary function.
Figure 11 shows the velocity profiles corresponding to the three types of viscosity variation 3
distributions, where the shape of velocity profile apparently varies with the type of viscosity variation distribution employed. The shapes of velocity profiles for type I and type III are parabolic 5
and logarithmic, respectively, while the shape for type II is in between.
The water depth in 1D and 2D shallow water flow computation is determined by the flow 7
discharge, energy slope, and channel roughness. From Figure 11 one can see that the water depth is affected not only by these three factors but also by the viscosity variation distribution in the 9
3D free-surface flow computations. However, the flow regime is indeed turbulent because of their high roughness and irregular bed forms and the eddy viscosity is a function of flow field, 11
not a fluid property. The difficulty of the 3D channel flow simulation is that the eddy viscosity changes anytime and anywhere so that manually calibrating eddy viscosity in the model is almost 13
impossible. Before the turbulent model is considered in the proposed model, an estimation of eddy viscosity, related to Manning’s roughness n, water depth, and flow discharge is needed. The case 15
Copyrightq 2007 John Wiley & Sons, Ltd. Commun. Numer. Meth. Engng (2007)
UNCORRECTED PROOF
20
CNM 1061
M.-C. HUNG ET AL.
Table I. Comparison of computed water depth at mid-length from simulation and Manning’s equation.
n=0.03 n=0.035 n=0.04 Manning’s equation 2.74 3.00 3.25
Simulated 2.79 3.00 3.28
Figure 12. Simulation of water surface profiles with different n values.
designed hypothetically is the same as above, used to examine the eddy viscosity relation derived 1
previously in this study. is estimated to be about 0.0304 by Equation (42), and the corresponding eddy viscosity is about 0.011m2/s.
3
This numerical experiment with same conditions as the previous case is carried out to verify the relation between flow depth and Manning’s roughness n. Three different values 0.03, 0.035, 5
and 0.04 of Manning’s roughness n were considered in the numerical experiments. Table I shows the comparison of computed water depths with those calculated by Manning’s equation. Both the 7
water depth and velocity variation obviously demonstrate the response of the change of Manning’s roughness n. The normal depths are 2.74, 3.0, and 3.25 m for n=0.03, 0.035, and 0.04, respectively.
9
Computed water depths are a slightly larger than those from Manning’s equation. Comparison of the water depths under different n values is shown in Figure 12. The velocity profiles of inner 11
cross section at middle length for various Manning’s n are also compared in Figure 13. The larger the value of roughness considered, the larger the depth and smaller the velocity obtained in the 13
simulation.
5.4. Channel flow with wind shear 15
The final case demonstrated in the following is to solve the wind shear effect on the flow patterns.
All conditions were the same as the previous experiment with n=0.035 except that the wind 17
shears acting on the surface are 0.0, 2.0, and 5.0N/m2directing to the upstream. Figure 14 shows
UNCORRECTED PROOF
CNM 1061
APPROACH FOR SHALLOW WATER FREE-SURFACE FLOW COMPUTATION 21
Figure 13. Comparison of velocity profiles with different n values at x=2500m.
Figure 14. Simulation of water surface profile under different wind shears.
the comparison of wind-induced setups under different shear stresses. The water surface profile 1
behaves like the M2 curve under the action of wind shears because the rating curve was used as downstream boundary conditions. Increase of the water depth drives the flow against the additional 3
resistance on the surface. The maximum velocity takes place at one-fourth of depth under free surface due to wind shear resistance acting on the water surface. The converged water depth is 5
Copyrightq 2007 John Wiley & Sons, Ltd. Commun. Numer. Meth. Engng (2007)
UNCORRECTED PROOF
22
CNM 1061
M.-C. HUNG ET AL.
Figure 15. Comparison of velocity profiles under different wind shears at x=2500m.
3.0 m downstream and increases gradually to 3.17 m upstream under wind shear of 2.0N/m2due 1
to wind-induced setup. The velocity profiles were also compared in Figure 15 under different values of wind shears.
3
6. CONCLUSION
There were a few 3D free-surface flow models developed in the past decades but least applied 5
to field problem due to its heavy computational time consumptions. A novel interface continuity-based approach proposed in this study makes the layer approach applicable in the shallow water 7
free-surface flow computation. Use of a quadratic shape function guarantees that there is no discon-tinuity situation taking place between two adjacent layers further. The modeling procedure adopted 9
provides a simple and robust method for simulating 3D free-surface shallow water flows. The verification performed in the case of wind-induced flow in a closed rectangular basin by comparing 11
with the analytical solution shows a reasonable agreement in the simulated steady-state profiles.
The evolution of water surface slope and the velocity profiles at the middle length of the basin 13
were examined to investigate whether the model responding to the wind shear is consistent with the physical nature. Alternate shapes of velocity profile formed by their correspondent distribution 15
of vertical viscosity also show the flexibility of the proposed approach either in wind-induced flows or in channel flows. The variations of velocity profile due to change of viscosity distribution 17
and near bed velocity are also investigated. Further, a simplified relation between discharge, water depth, and eddy viscosity is derived as the closure of the governing equations. The proposed model 19
UNCORRECTED PROOF
CNM 1061
APPROACH FOR SHALLOW WATER FREE-SURFACE FLOW COMPUTATION 23 applicability for the 3D shallow free-surface water flow simulations was examined through two 1
hypothetical cases and shows convincing results. Further development of considering the turbulent model is needed for real-world application.
3
NOTATIONS
h layer thickness
n Manning’s roughness
p hydrostatic pressure u x-component of velocity
u∗ shear velocity
v y-component of velocity zB elevation of bottom bed
A coefficients of the algebraic equation for water columns
C arbitrary constant
F layer across flux
H water depth
K number of layers
S slope
Su, Sv source terms
U depth-averaged velocity
arbitrary constant for Elder’s relation
specific weight of water
dynamic eddy viscosity
kinematic eddy viscosity
water density
shear stress
s wind shear
Subscripts
e,s,w,n denote faces of the considered cell h horizontal direction
k denote the layer interface quantity between kth and(k −1)th layers k+1/2 denote the kth layer-averaged quantity
m dummy index used for integration along water depth
v vertical direction
E, S, W, N denote centroid of neighbors of the considered cell P denote centroid of the considered cell
5
ACKNOWLEDGEMENTS
Partial financial support of this study from National Science Council of Taiwan, Republic of China, through Contract NSC-89-2211-E-009-091 is greatly appreciated. The computer resources used in this study were provided by National Center for High-Performance Computing in Taiwan, Republic of China.
Copyrightq 2007 John Wiley & Sons, Ltd. Commun. Numer. Meth. Engng (2007)
UNCORRECTED PROOF
24
CNM 1061
M.-C. HUNG ET AL.
REFERENCES
1
1. Demuren AO, Rodi W. Side discharges into open channels: mathematical model. Journal of Hydraulic Engineering 1983; 109(12):1707–1722.
3
2. Hervouet JM, Van Haren L. Recent advances in numerical methods for fluid flows. In Floodplain Processes, Chapter 6, Anderson, Walling, Bates (eds). Wiley: New York, 1996.
5
3. Meselhe EA, Sotiropoulos F. Three-dimensional numerical model for open channel with free surface variations.
Journal of Hydraulic Research 2000; 38(2):115–121.
7
4. Hirt CW, Nicholls BD. Volume of fluid (VOF) method for dynamics of free boundaries. Journal of Computational Physics 1981; 39:201–221.
9
5. CFX. CFX-4.2: Solver. AEA Technology: U.K., 1997.
6. Faure JB, Buil N, Gay B. 3-D modeling for unsteady free-surface flow in open channel. Journal of Hydraulic
11
Research 2004; 42(3):263–272.
7. Li YS, Zhan JM. An efficient three-dimensional semi-implicit finite element scheme for simulation of free surface
13
flows. International Journal for Numerical Methods in Fluids 1993; 16:187–198.
8. Kim CK, Lee JS. A three-dimensional PC-based hydrodynamic model using an ADI scheme. Coastal Engineering
15
1994; 23:271–287.
9. Shankar NJ, Chan ES, Zhang QY. Three-dimensional numerical simulation for an open channel flow with a
17
constriction. Journal of Hydraulic Research 2001; 39(2):187–201.
10. Li CW, Gu J. 3D layered-integrated modeling of mass exchange in semi-enclosed water bodies. Journal of
19
Hydraulic Research 2001; 39(4):403–412.
11. Lin B, Falconer RA. Three-dimensional layer-integrated modeling of estuarine flows with flooding and drying.
21
Estuarine, Coastal and Shelf Science 1997; 44:737–751.
12. Tsai TL, Yang JC, Huang LH. An accurate integral-based scheme for advection–diffusion equation.
23
Communications in Numerical Methods in Engineering 2001; 17:701–713.
13. Delft Hydraulics. User Manual of Delft3D-FLOW—Simulation of Multi-dimensional Hydrodynamic Flows and
25
Transport Phenomena, Including Sediments. Delft Hydraulics, 2003.
14. Elder JW. The dispersion of marked fluid in turbulent shear flow. Journal of Fluid Mechanics 1959; 5:544–560.
27
15. Nezu I, Rodi W. Open-channel flow measurements with a laser Doppler anemometer. Journal of Hydraulic Engineering 1986; 112:332–355.
29
16. Heaps NS. Vertical structure of current in homogeneous stratified waters. Hydrodynamics of Lakes, Hutter K (ed.).
1984; 154–202.
31
17. Jobson HE, Sayre WE. Vertical transfer in open channel flow. Journal of Hydraulic Division 1970; 96:703–724.