CMS-Flow Numerical Methods: Hydrodynamics

Momentum and Continuity

The governing equations are solved in a segregated manner in which each governing equation is linearized and solved separately in a sequential manner within an iteration loop in order to obtain a converged solution. Coupling between the velocity (momentum) and water level (continuity) is achieved with the SIMPLEC algorithm (van Doormal and Raithby 1984). The main difficulty in solving the momentum equations is that the water level is not known a priori and must be calculated as part of the solution. The solution algorithm procedure is described below. First, an initial pressure $p^{*}=\rho g{\overline {n}}^{*}$ is estimated based on the previous time step value. Then, the momentum equations are solved for the corresponding velocity

${\frac {\partial (hV_{i}^{*})}{\partial t}}+{\frac {\partial (hV_{j}^{*}V_{i}^{*})}{\partial x_{j}}}={\frac {\partial }{\partial x_{j}}}\left(v_{t}h{\frac {\partial V_{i}^{*}}{\partial x_{j}}}\right)-{\frac {h}{\rho }}{\frac {\partial p^{*}}{\partial x_{i}}}+S_{i}^{*}$

(1)

where $S_{i}^{*}$ includes all the remaining terms. The inter-cell velocities are calculated with a Rhie and Chow (1983) type interpolation method

$V_{f}^{*}=(H_{\perp }^{*})_{f}-\left({\frac {\Delta A_{P}}{a_{P}}}\right)_{f}\left({\frac {h}{\rho }}\bigtriangledown _{\perp }p^{*}\right)_{f}$

(2)

where $(H_{\perp }^{*})_{f}=V_{\perp }^{*}+{\frac {\Delta A_{p}}{a_{p}}}{\frac {h}{\rho }}\bigtriangledown _{\perp }p^{*}$ and $\left({\frac {h}{\rho }}\bigtriangledown _{\perp }p^{*}\right)_{f}$ denotes the linear interpolation operator. This method avoids the checkerboard oscillations associated with the collocated grid. It is noted that this approach is slightly different from Lai (2010) and others and was found to be significantly more stable.

Next, the velocity $V^{'}$ and pressure corrections $p^{'}$ are defined such that both the momentum and continuity equations are satisfied

$V_{i}^{n+1}=V_{i}^{*}+V_{i}^{'},p^{n+1}=p^{*}+p^{'}$

(3a,b)

Subtracting the initial velocity equation from the momentum equation leads to velocity correction equation

${\frac {\partial (hV_{i}^{'})}{\partial t}}+{\frac {\partial (hV_{i}^{'}V_{j}^{'})}{\partial x_{j}}}={\frac {\partial }{\partial x_{j}}}\left(v_{t}h{\frac {\partial V_{i}^{'}}{\partial x_{j}}}\right)-{\frac {h}{\rho }}{\frac {\partial p^{'}}{\partial x_{i}}}$

(4)

This equation can be written in the discretized form as

$a_{p}V_{i,P}^{'}=\sum a_{N}V_{i,N}^{'}-{\frac {h}{\rho }}\bigtriangledown _{i}p_{P}^{'}\Delta A_{P}$

(5)

In the SIMPLEC algorithm, the velocity correction is assumed to vary smoothly so that $\sum a_{N}V_{i,N}^{'}$ may be approximated as $V_{i,P}^{'}\sum a_{N}$ , which leads to the velocity correction equation

$V_{i}^{'}=-G\bigtriangledown _{i}p^{'}$(3-32)

where

$G=(h\Delta A_{p}/\rho )/\left(a_{p}-\sum a_{N}\right)$ . Using $\partial h/\partial t=\partial p/(\rho g\partial t)$

and substituting $V_{i}^{n+1}=V_{i}^{*}+V_{i}^{'}$ in the continuity equation yields the semi-discrete water level correction equation

${\frac {(1+0.5\theta )(p^{*}+p^{'})-(1+\theta )p^{n}+0.5\theta p^{n-1}}{\rho g\Delta t}}={\frac {\partial }{\partial x_{j}}}\left(hG{\frac {\partial p^{'}}{\partial x_{j}}}\right)-{\frac {\partial (hV_{j}^{*})}{\partial x_{j}}}$(3-33)

Note that at convergence, $p^{'}=0$ and the Equation 3-28 reduces to the continuity equation. Once the pressure correction equation is solved, the cell-centered water levels and current velocities are corrected. The cell face velocities are also corrected as $V_{f}^{n+1}=V_{f}^{*}+V_{f}^{'}$ , in which the velocity correction is given by $V_{f}^{'}=-G_{f}(\bigtriangledown _{\perp }p^{'})_{f}$ .

Summary of the SIMPLEC Algorithm:

1. Guess the water level and pressure field $p^{*}$

2. Solve the momentum equations (give Equation numbers here) to obtain $V_{i}^{*}$

3. Use the Rhie and Chow (1983) momentum to determine the velocities and fluxes at cell faces

4. Solve the pressure equation (give Equation number here) to obtain $p^{'}$

5. Use the correction equations to adjust the velocities and water levels

6. Treat the corrected water level and pressure field as a new guess, and repeat this procedure from Step 2 until convergence

Wetting and Drying

During numerical simulations of the surface water flows with sloped beaches, sand bars, and islands, the land-water interface changes with time. This means that it is possible for nodes at the land-water interface to be wet or dry throughout a given simulation. In CMS, a threshold water depth is used to judge drying and wetting. If the depth at the cell center is larger than the threshold value (i.e. a small value such as 0.02 m for field cases), then the node is considered to be wet. Similarly, if the depth at the cell center is smaller than the threshold value, then the node is considered to be dry. For the implicit solver, all of the wet and dry cells are included in the matrix solver. Dry cells are assigned with a zero velocity. Cell faces are classified as either open (if the two cells neighboring cells are wet) or closed (i.e. cell faces are not classified as wet or dry).