Modelling free surface incompressible flow around thin vertical wall incorporating hydrodynamic pressure on sigma-coordinate

H. Kim1 , D. H. Hwang2 , S. W. Baek3 , H. J. Yoo4 , C. Jang5 , J. Jin6

1, 2, 3Department of Civil and Environmental Engineering, Kookmin University, Seoul, Korea

4Geosystem Research Corporation, Busan, Korea

5Korean Intellectual Property Office, Daejeon, Korea

6Korea Institute of Ocean Science and Technology (KIOST, Ansan, Korea

1, 3Corresponding authors

Vibroengineering PROCEDIA, Vol. 4, 2014, p. 318-323.
Accepted 6 October 2014; published 3 November 2014

Copyright © 2014 JVE International Ltd. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Creative Commons License

Submerged thin walls are extreme case of submerged rectangular blocks, and could be used for many purposes at coastal zones or rivers, e.g. protecting tsunami. To understand flow characteristics including flow and pressure fields around a specific submerged thin wall a new algorithm in a numerical model, CST3D, is developed which includes computation of hydrodynamic pressure with solution procedure of a new divergence equation on moving sigma-coordinate. Sigma-coordinate is adequate for simulation of subcritical flow over mild-sloped beds. However sigma-coordinate is quite poor to simulate flows around sharp structures on the bed. Most sigma-coordinate approaches have adopted assumption of hydrostatic pressure. A few previous trials include hydrodynamic pressure by solving the Poisson equation on pressure. The other solution algorithm to treat hydrodynamic pressure was the divergence removing method, SOLA, developed by Los Alamos Lab., which is valid for fixed grid approaches only. The SOLA is modified here by including the flow unsteadiness within the divergence equation, so that accurate divergence is computed at every computational time step on sigma-coordinate. Computed flow field shows reasonable flow pattern including structure-scale vortices. Computed vorticity field reveals that the downstream vorticity is much stronger than the upstream one. The model was verified to laboratory experiments at a 2DV flume. Time-average flow vectors were measured by using one-dimensional electro-magnetic velocimeter. Computed horizontal velocity components of flow agree well with the measured ones within 6 cm/s error at 5 sections.

Keywords: submerged thin wall, CST3D, hydrodynamic pressure, sigma-coordinate, sola, vorticity.

1. Introduction

Submerged thin walls are designed to control river flood, or to protect tsunami at coasts. We need to know detailed flow pattern and force of fluid body acting on structure for proper design.

Fluid flow is unsteady, and three-dimensional from its nature. Coastal, estuary, or river flows usually have free surface, and are horizontal-velocity-component-dominant. Vertical distribution of water pressure in water column is close to hydrostatic pressure, which increases linearly from surface to bottom.

Most early numerical models for coastal or river flow were developed solving depth-average governing equations with assumption of hydrostatic pressure distribution [1-3]. The depth-average models were then expanded to three-dimensional form by adopting assumption of hydrodynamic pressure distribution to obtain three-dimensional information. Common vertical grids used for most three-dimensional flow models are z-grid, or sigma-grid, both of which have merits and demerits. Sigma-grid has a distinctive merit of simple free surface treatment for unsteady flows, and has been applied to variety of flows with mild-slope free-surface and bed level, e.g. tidal current flow, or density flow.

Methods to compute hydrodynamic pressure distribution in water column in models using z-grid have been developed and used widely [4-5], while methods to compute hydrodynamic pressure distribution in models using sigma-grid have not soundly been developed so far. Some sigma-grid models with hydrodynamic pressure distribution include Lin [6], Keilegavlen and Berntsen [7], and Zhu et al. [8].

Lin, Keilegavlen and Berntsen and Zhu et al. solved the highly implicit difference systems of equations from the Poisson equation regarding the hydrodynamic pressure driven from the continuity equation and two horizontal momentum equations. They used matrix solvers which require heavy computation. Furthermore Lin, Keilegavlen and Berntsen included high order terms generated from delicate non-orthogonality between horizontal and vertical grids, and included them in their numerical model, while Zhu et al. ignored the high order terms for simplicity. The non-orthogonality between horizontal coordinates and vertical sigma coordinate could be ignored, considering the initial adoption of sigma coordinate under the assumption of mild-slope of free surface and bottom morphology.

Hirt et al. [9] developed an explicit scheme called “Solution Algorithm” (SOLA) to solve systems of equations regarding velocity and hydrodynamic pressure fields for rectangular parallelepiped grid system on z-grid. SOLA starts from calculating fluid flow divergence field, and explicitly diffuses away unnecessary divergence by adjusting dynamic pressure and three velocity components surrounding each rectangular parallelepiped grid point in three-dimensional flows. SOLA’s procedure solving divergence corresponds to the solution of the Poisson equation regarding hydrodynamic pressure. Although the SOLA includes repetition to meet convergence criteria at each time step, it has a merit that the convergence level can be intentionally chosen with a maximum allowance of residual divergence. So far SOLA has not been applicable to sigma-grid, because the divergence was defined only on fixed z-grid, not on moving sigma-grid.

Here we modify Hirt et al.’s SOLA scheme by introducing divergence equation for incompressible fluid flow on moving sigma-grid.

2. Governing equations of new flow model

Governing equations to describe flows of incompressible fluid are the continuity equation, and three Reynolds’ momentum equations for turbulent flow of incompressible fluid. We introduce the Boussinesq approximation in buoyancy-driven flows.

New vertical variable, sigma (σ), replacing z, is defined as:

σ = z + d ζ ( t ) + d = z + d H ( t )   ,       0 σ 1 ,

where z is the vertical coordinate, ζ is the water level above a reference level, t is time, d is the water depth below the reference level, and H is the total water depth.

Continuity equation on sigma-coordinate is:

ζ t + H u x + H v y + w σ = 0 ,

where x, y are the two nearly horizontal coordinates linking constant sigma values; u, v are horizontal fluid velocity components in the x, y directions, respectively; w is the vertical fluid velocity component on the moving sigma grid.

Assume that non-orthogonality between two horizontal coordinates (x, y) and the vertical coordinate (sigma) is negligible. Then, momentum equations in the two horizontal directions become:

H u t + x H u u + y H v u + σ w u - f H v
            = - H x g ζ + p - d x - z H x p σ + σ H - 1 A v u σ + Q u   ,
H v t + x H u v + y H v v + σ w v + f H u
            = - H y g ζ + p - d y - z H y p σ + σ H - 1 A v v σ + Q v   ,

where f is the Coriolis coefficient, Qu, Qv are the momentum source or sink terms including horizontal dispersion, p is the pressure, and Av is the eddy viscosity in the vertical direction. The momentum equation in the vertical direction on sigma-coordinate reads:

w t + u w x + v w y + w H w σ = - 1 ρ H p σ - g ,
p = p d + p s ,
p s z = - ρ g ,
w = w * - w ^ ,

where w is the vertical fluid velocity component on sigma-grid, p is the pressure, pd is the hydrodynamic pressure, ps is the hydrostatic pressure, ρ is the fluid density, g is the gravitational acceleration, w* is the fluid vertical velocity on the fixed z coordinate, and w^ is the moving velocity of the sigma-grid.

The divergence-spreading procedure of the SOLA consists of two steps, i.e. computation of divergence step, and the adjustment step of pressure and fluid velocity components. First, the divergence on sigma-grid is newly defined as the following equation:

D = 1 H ζ t + H u x + H v y + w σ ,

which is based on the above continuity Eq. (2). This equation differs from the previous definition of the divergence on fixed grid in SOLA, that is:

D = u * x * + v * y * + w * z ,

where x*, y* are the two horizontal Cartesian coordinates, and u*, v* are the velocity components in x*, y* directions, respectively.

The above divergence, Eq. (6), is expressed in finite difference terms as follows:

D i , j , k = 1 H ζ i , j , k n + 1 - ζ i , j , k n t + ( u y z ) i + 1 2 , j , k - ( u y z ) i - 1 2 , j , k + ( v x z ) i , j + 1 2 , k
        - ( v x z ) i , j - 1 2 , k + ( w ^ x y ) i , j , k + 1 2 - ( w ^ x y ) i , j , k - 1 2 1 x y z ,

where subscripts i, j, k are defined at grid point centers, and ±1/2 describes grid border lines between two grid points. Theoretically divergence should be zero in the whole domain to satisfy the assumption of fluid incompressibility.

At the second step computed divergence field is reduced towards zero by adjusting the hydrodynamic pressure and the three velocity components with the following equations:

p i , j , k = - α 2 t D i , j , k 1 x 2 + 1 y 2 + 1 z 2   ,
u i + 1 2 , j , k = - t x p i + 1 , j , k - p i , j , k ,       u i - 1 2 , j , k = - t x p i , j , k - p i - 1 , j , k ,
v i , j + 1 2 , k = - t x p i , j + 1 , k - p i , j , k ,       v i , j - 1 2 , k = - t x p i , j , k - p i , j - 1 , k ,
w i , j , k + 1 2 = - t x p i , j , k + 1 - p i , j , k ,       w i , j , k - 1 2 = - t x p i , j , k - p i , j , k - 1 ,

where α is a parameter to control the adjusting speed of divergence (a number less than 0.14 for three-dimensional computation). The present modified SOLA scheme is adopted in a new flow and sediment transport model system, CST3D [10-11].

3. Verification of new scheme

Flows passing through a rapid bathymetric change could show the importance of the hydrodynamic pressure. As a validation step a non-uniform steady flow is chosen for comparison of the model results and laboratory experimental results. In this case the flow is considered two-dimensional vertical, ignoring micro-scale three-dimensionality of flow turbulence. A thin wall of height 30 cm blocks one-way water flow of depth about 50 cm in a channel of uniform width of 80 cm. A constant water flux was supplied from one end of the channel, and the flow arrived at a quasi-steady state after some time. Velocities at several heights along five sections in Fig. 1 were measured by using an anemometer, model LP1100 manufactured by Kenex.

Fig. 1. Measurement points of horizontal velocity components during experiment

 Measurement points of horizontal velocity components during experiment

Fig. 2. Measured and computed horizontal velocity components at two sections of experiment

 Measured and computed horizontal velocity components at two sections of experiment


 Measured and computed horizontal velocity components at two sections of experiment


Velocities were measured seven times at every point. Two large vortices were observed at both upstream and downstream, see Fig. 2.

Numerical computation was also carried out with several parameters chosen as: spatial increment in the x direction = 2.5 cm, spatial increment in the z direction = 2.5 cm, temporal increment = 0.001 s, total channel length = 4 m, and α= 0.14. Computed flow field is shown in Fig. 3. The vectors in the figure represent current direction and speed. The downstream vortex is outstanding, while the weak upstream vortex is seen on vorticity contours.

Fig. 3. Computed current field for experimental conditions: arrows are current vectors

 Computed current field for experimental conditions: arrows are current vectors

Computed and measured horizontal velocity components are compared in Fig. 2. The largest speed error of the computation is within 6 cm/s, which is satisfactory compared to the strongest current speed of 30 cm/s. The results demonstrates the overall behavior of the numerical model is reasonable.

4. Conclusions

Sigma-coordinate is useful for describing free-surface flows at coasts and rivers. Previous sigma-approaches have solved the Poisson equation on the hydrodynamic pressure, which involves heavy matrix-solving effort. Here computation of hydrodynamic pressure is included in sigma-coordinate approach with the divergence-dispersion method. Equation for divergence on sigma-coordinate is introduced by modifying the existing divergence-dispersion method, SOLA, which is valid for fixed z-grid only. The divergence is converted to a difference form, and adopted in a numerical model, CST3D.

The new scheme was validated against a laboratory experiment which is a non-uniform uni-directional flow in a channel. The model is believed to have reproduced experimental results well. The present modified SOLA could expand its usage for more general geometry or structure in the future.


This includes research results of “MIDAS” 2014, and “MEPS” 2013 supported by Ministry of Oceans and Fisheries (Korea).


  1. Tetra Tech Inc. The Environmental Fluid Dynamics Code: Theory and Computation Volume 1: Hydrodynamics and Mass Transport. US EPA Version 1.01, 2007. [Search CrossRef]
  2. Delft Deltares. Delft3d-Flow: Simulation of Multi-Dimensional Hydrodynamic Flows and Transport Phenomena Including Sediments: User Manual: Hydro-Morphodynamics: Version 3.15.33641, 2014. [Search CrossRef]
  3. DHI. Mike21 & Mike3 Flow Model FM: Hydrodynamic Module: Short Description. Horsholm, DHI, 2013. [Search CrossRef]
  4. Ho D., Boyes K., Donohoo S., Cooper B. Numerical flow analysis for spillways. Proceedings of the 43rd ANCOLD Conference, Hobart, Tasmania, 2003, p. 1-11. [Search CrossRef]
  5. ANSYS Inc 2012 ANSYS FLUENT. Theory Guide: Release 14.5. Canonsburg: ANSIS. [Search CrossRef]
  6. Lin P. A multiple-layer sigma-coordinate model for simulation of wave-structure interaction. Computer and Fluids, Vol. 35, 2006, p. 147-167. [Search CrossRef]
  7. Keilegavlen E., Berntsen J. Non-hydrostatic pressure in σ-coordinate ocean models. Ocean Modelling, Vol. 28, 2009, p. 240-249. [Search CrossRef]
  8. Zhu L., Chen Q., Wan X. Numerical modeling of nonlinear water waves with sigma coordinate and layer thickness optimization. Proc. XIX Int. Conference on Water Resources, University of Illinois at Urbana-Champaign, 2012. [Search CrossRef]
  9. Hirt C. W., Nichols B. D., Romero N. C. SOLA: A numerical solution algorithm for transient fluid flows. Los Alamos, Los Alamos Scientific Lab, 1975. [Search CrossRef]
  10. Jang C., Kim H., Lim N. A development of coupled wave-induced current modeling system and its application to the idealized shoreline with detached breakwater. Journal of Wetlands Research, Vol. 14, 2012, p. 439-455, (in Korean). [Search CrossRef]
  11. Kim H. Coastal sediment transport prediction model (CST3D): Version 2013-v1: User manual. Seoul, Saerom, 2013, (in Korean). [Search CrossRef]