Susanne Krüger
Laboratory of Hydraulics, Hydrology and Glaciology (VAW)
Swiss Federal Institute of Technology Zurich (ETHZ)
8092 Zurich
Switzerland
Tel: +41-1-632 4154, Fax: +41-1-632 1192,
E-mail: krueger@vaw.baug.ethz.ch
Nils Reidar B. Olsen
Department of Hydraulic and Environmental Engineering
The Norwegian University of Science and Technology
7491 Trondheim
Norway
Abstract:
Commonly, depth-averaged Navier-Stokes equations are used for shock-wave
computations. The classical approach is limited to regions where the hydrostatic
pressure assumption is valid. An extended approach taking into account
non-hydrostatic pressure distribution as well as momentum in the third direction
improves the simulation of the complex wave pattern. A comparison between
measurements, the classical and extended shallow water equations as well as with
the fully three-dimensional Navier-Stokes equations is undertaken for
computations of supercritical flow in a channel contraction considering Froude
numbers 4 and 6. The results show that the limitations of the classical approach
increase with the Froude number. The extended shallow water approach and the
fully three-dimensional Navier-Stokes equations yield both a good agreement with
the measurements.
Keywords: CFD, highly supercritical flow computations, channel-contraction, shock-waves
Supercritical flow
in channel contraction forms moving or standing shock waves for any change in
geometry, slope or friction. Commonly, shock waves are computed using the
depth-averaged Navier-Stokes equations, which yield the shallow water equations.
The classical approach assumes a uniform longitudinal velocity distribution,
vertical velocities are neglected such that a hydrostatic pressure distribution
results. An extension of the classical shallow water approach taking into
account non-hydrostatic pressure and momentum in the vertical direction [Krüger
et. al., 1998] shows significant improvements in the numerical simulation of
position and height of shock waves. This paper presents a comparison between
measurements, classical and extended computations as well as fully
three-dimensional simulations.
The classical shallow water equations [Vreugdenhil, 1994] are derived by depth- averaging the Navier-Stokes equations. Therefore, a uniform longitudinal velocity profile is assumed over depth. The vertical velocity components are neglected such that the pressure distribution is hydrostatic.
In contrast, the extended shallow water equations assume additional variables for the reproduction of the horizontal and vertical velocity profiles as well as for the pressure distribution [Steffler and Jin, 1993] (Figure 1). Using these distribution functions and depth-averaging the Navier-Stokes equations yields the following set of continuity and momentum equations:
(1)
(2)
(3)
(4)
with density
, depth-averaged turbulent normal stresses
, depth-averaged turbulent shear stresses
, bed shear
, bed elevation
and water-surface elevation
.The vertical bed and surface velocities are determined by the kinematic
boundary conditions
(5)
(6)
The distribution profiles of the extended shallow water equations allow the derivation of the first momentum such that the moment of the continuity equations and the three moment of momentum relations result
(7)
(8)
(9)
(10)
where
,
and
The depth-averaged stress terms
are neglected in this work. Setting
the additional variables equal to zero yields the classical shallow water
equations.
The finite element library FEMTOOL [Rutschmann, 1993] was used to solve the set of equations. Using a direct sparse matrix solver, both the classical and extended shallow water equations were solved simultaneously. The non-linearity was resolved by the Picard-iteration procedure. Numerical stabilization is obtained with a Galerkin/least-squares upwind scheme.
Solving the fully 3D Navier-Stokes equations, the k-e turbulence model was used [Rodi, 1980]. The pressure coupling was done by the SIMPLE method [Patankar, 1980]. A first-order upwind method was used for discretization of the convective terms. The water level was computed by using the water continuity deficit in the cells closest to the water surface. The SIMPLE method was not used in these cells, so water continuity was non-zero. The pressure in the top cells was instead computed by interpolation between the cells below and the zero level at the water surface.
The classical and extended shallow
water equations as well as the fully three-dimensional Navier-Stokes were
applied for shock-wave computations in a channel contraction. Two test cases
compare measurements and computations. Test case 1 considered a Froude number
, test case 2
. The measurements were conducted by Reinauer [1995]
in an experimental doctoral work. Figure 2 illustrates the contraction geometry
chosen for comparison, the solid area was used for test case 1, the contraction
was lengthened with the dashed area for test case 2.
Test Case 1: Froude number
The water depth at the upstream
boundary was constant
m. The upstream discharge
was set, corresponding to Froude
number
. For the shallow water computations, the extra terms were set to zero at the
inflow boundary corresponding to a hydrostatic pressure assumption. The
horizontal discharge variables were assumed to be tangential to the side walls.
Over the entire domain, the initial conditions correspond to the boundary
conditions on the inflow boundary. The problem was solved for unsteady flow
until the steady-state solution was obtained. Both the classical and the
extended shallow water equations were computed
using an unstructured grid containing 1240 nodes and 2261
elements. The fully 3D computations were
conducted with a structured mesh containing 100x21 nodes in the horizontal plane
and 11 nodes in the vertical direction.
Figures 3(a) and 3(b) compare the measurements with the computations of the classical and the extended shallow water equations as well as with the 3D equations along the side wall and along the axis. The classical shallow water computations show an underestimation of the wave heights and also large discrepancies in position. The measured water depth is only matched in regions, where the flow is undisturbed and a hydrostatic pressure is valid. Under-and overshooting occurs along the steep shock fronts, which shows that numerical diffusion is of minor relevance. The shock waves are computed in a plateau-like shape, whereas the experiments show a peak-like form.
Compared with the classical equations, the extended equations show a completely different shape of the shock waves. The extended equations present the peak-like form of the measured shock waves and match the wave peaks in position and height quite accurate.
The 3D computations present the same peak-like form of the shock waves as the extended shallow water computations and the measurements. Compared to the measurements, the computed wave along the side wall is too low and the computed wave front along the axis is less steep. Contrary to the extended shallow water computations, the fully 3D computations match the measured water depth downstream of the wave along the axis very accurate.
Both the extended solution and the 3D computations show differences with the measurements around the wave peak. This might be an effect of entrained air, which is always visible in hydraulic models and prototypes due to the interaction and reflection of the waves. In nature, the transition from water to air is visible as an air-water spray. Most numerical models are unable to reproduce this spray because of the mass conservation. Also the measurements are complicated in this regions causing inaccuracies.
Test Case 2: Froude number
The water depth at the upstream
boundary was set to
m as in Test case 1. The upstream
discharge
was defined so
that an upstream Froude number
resulted.
For the shallow water computations, the additional variables were set to zero at
the upstream boundary. Two sets of transformed boundary conditions were set for
the horizontal discharge variables to define the tangential flow along the side
walls. The initial conditions over the entire domain corresponded to the
upstream boundary conditions. The problem was solved with the transient
equations until the steady-state solution was obtained. For the classical and
extended shallow water equations, the resolution of the horizontal domain
consisted in an unstructured grid containing 1425 nodes and 2587
elements. The 3D computations were undertaken with a structured grid using
125x21 nodes in the horizontal plane and 11 nodes in the vertical direction.
Figures 4(a) and 4(b) compare the measurements with the computations of the classical and the extended shallow water equations as well as with the fully 3D computations. Along the side wall, the water depth was measured in several positions. Along the axis, the water depth was only fixed in one single point giving the wave maxima in height and position.
The classical equations again show steep shock fronts. The wave peak is extremely underestimated, the maximum along the axis is predicted with half of the measured wave. 3D effects increase with the Froude number and therefore cannot be reproduced with the classical shallow water approach.
The extended equations are able to reproduce these 3D effects. The measured wave is predicted in height and position accurate even if the wave maxima at the side wall is computed too high.
The 3D computations predict reliably the wave along the side wall but show a major underestimation of the wave along the axis. There seem to be a diffusion of the wave front. The finite volume approach computes the water-surface changes in the center of each cell. The grid line intersections are changed as an average of the water surface elevations in the four surrounding cells. The averaging procedure will cause the gradients to be smoothed. The finite element method computes the pressure in the grid cell corners, and the averaging process is not used. Using a finer grid resolution in the horizontal plane can reduce the diffusive effects of the averaging procedure. Preliminary tests with finer grids show instabilities and divergence.
Three numerical approaches are compared with measurements for two test cases with supercritical flows in channel contractions. The classical shallow water approach is based on a hydrostatic pressure assumption and negligible vertical velocities. Therefore, it is limited in regions where the water-surface slope is large or has high curvature. The 3D effects introduced with the shock-wave phenomena then cannot be reproduced.
The extended shallow water approach considers non-hydrostatic distribution and the velocities in all three directions. This approach significantly improves the simulation of shock waves. The waves are predicted in height and position quite accurate even for high Froude numbers. This model is able to consider these kind of 3D effects.
The
fully 3D computations enable the physically correct reproduction of the complex
flow pattern. However, they show too diffusive shock waves. The simulations
might be improved with a finer resolution of the horizontal plane.
Acknowledgement
The authors thank the Swiss National Foundation which founded the research project to numerically investigate shock waves in supercritical flows with the shallow water equations (SNF-Grant No. 20-45631.95).
References
Krüger, S., Bürgisser, M., Rutschmann, P. (1998). Advances in calculating supercritical flows in spillway contractions. Hydroinformatics’ 98, 163-170.
Patankar, S. V. (1980). Numerical heat transfer and fluid flow, McGraw-Hill Book Company, New York
Reinauer, R. (1995). Kanalkontraktionen bei schiessendem Abfluss und Stosswellenreduktion mit Diffraktoren. PhD thesis, ETH Zurich, Diss No. 11320.
Rodi, W. (1980). Turbulence modelling in hydraulic engineering, IAHR State-of-the-Art Monograph.
Rutschmann, P. (1993). FE solver with 4D finite elements in space and time. Proc. VII Int. Conf. On Finite Elements in Fluids, ed. K. Morgan, E. Onate, J. Periaux, J. Peraire, O.C. Zienkiewicz, Barcelona, 136-144.
Steffler, P.M., Jin, Y.-C. (1993). Depth averaged and moment equations for moderately free surface flow. Journal of Hydraulic Research, 31(1), 5-17.
Vreugdenhil, C.B. (1994).
Numerical methods for shallow-water flow. Water Science and Technology Library,
Vol. 13, Kluwer Academic Publ.: Dordrecht.
FIGURE
CAPTIONS

Fig. 1 Additional terms to extended shallow water equations

Fig. 2 Contraction geometry for Test case 1 (solid area) and Test case 2 (dashed area).

(a) along side wall (b) along side wall
Fig. 3 Water-Depth comparison for Test case 1

(a) along side wall
(b) along axis
Fig. 4 Water-Depth comparison for Test case 2