CMS-Flow:Hydro Eqs: Difference between revisions

From CIRPwiki
Jump to navigation Jump to search
Line 1: Line 1:
<big>
big
= Governing Equations =
= Governing Equations =
The depth-averaged 2-D continuity equation may be written as
The depth-averaged 2-D continuity equation may be written as
{{Equation|<math> \frac{\partial h}{\partial t}+\nabla \cdot (h\mathbf{U})= S </math>|2=1}}
{{Equation|math \frac{\partial h}{\partial t}+\frac{\partial (h{{U}_{j}})}{\partial {{x}_{j}}}=0 /math|2=1}}


where <math>h</math> is the total water depth <math>h=\zeta+\eta</math>, <math>\eta</math> is the water surface elevation, <math>\zeta</math> is the still water depth, <math> \mathbf{U}=\left( {U_1},{U_2},{0} \right) </math> is the depth-averaged current velocity, <math> S </math> is a source term due to precipitation and evaporation, and <math>\nabla</math> is the divergence operator.
where mathh/math is the total water depth mathh=\zeta+\eta/math, math\eta/math is the water surface elevation, math\zeta/math is the still water depth, mathU_i/math is the depth-averaged Langrangian current velocity defined as math U_i=U_i^E+U_i^S/math, where mathU_i^E/math is the phase- and depth-averaged current velocity (i.e. Eulerian velocity), and mathU_i^S/math is the depth-averaged Stokes velocity (Phillips 1977)
 
{{Equation|mathU_{i}^{S}=\frac{{{E}_{w}}}{\rho hc}{{w}_{i}} /math|2=2}}
   
   
The momentum equation can be written as
The momentum equation can be written as
{{Equation| <math> \frac{\partial (h\mathbf{U})}{\partial t}+\nabla \cdot (h\mathbf{U} \mathbf{U})+{{f}_{c}}\mathbf{k}\times h\mathbf{U}=-gh\nabla \eta -\frac{h}{{{\rho }_{0}}}\nabla {{p}_{a}}+\nabla \cdot \left( {{\nu }_{t}}h\nabla \mathbf{U} \right)+\frac{1}{\rho }\left( {{\tau }_{w}}+{{\tau }_{S}}-{{\tau }_{b}} \right) </math>|2=2}}   
{{Equation| math \frac{\partial (h{{U}_{i}})}{\partial t}+\frac{\partial (h{{U}_{i}}{{U}_{j}})}{\partial {{x}_{j}}}-{{\varepsilon }_{ij}}{{f}_{c}}h{{U}_{j}}=-gh\frac{\partial \eta }{\partial {{x}_{i}}}-\frac{h}{{{\rho }_{0}}}\frac{\partial {{p}_{a}}}{\partial {{x}_{i}}}+\frac{\partial }{\partial {{x}_{j}}}\left( {{\nu }_{t}}h\frac{\partial {{U}_{i}}}{\partial {{x}_{j}}} \right)+\frac{1}{\rho }\left( \tau _{i}^{w}+\tau _{i}^{s}-\tau _{i}^{b} \right) /math|2=2}}   
    
    
where <math>g</math> is the gravitational constant, <math>f_c</math> is the Coriolis parameter, <math>\mathbf{k} = (0,0,1) </math>, <math>p_a</math> is the atmospheric pressure, <math>\rho_0</math> is a reference water density, <math>\nu_t</math> is the turbulent eddy viscosity, <math> \tau_{w} </math> is the wind stress, <math> \tau_{S} </math> is the wave stress, and <math>\tau_{b}</math> is the combined wave-current mean bed shear stress.
where mathg/math is the gravitational constant, mathf_c/math is the Coriolis parameter, mathp_a/math is the atmospheric pressure, math\rho_0/math is a reference water density, math\nu_t/math is the turbulent eddy viscosity, math \tau_{w} /math is the wind stress, math \tau_{S} /math is the wave stress, and math\tau_{b}/math is the combined wave-current mean bed shear stress. math\varepsilon_{ij}/math is the permutation parameter equal to 1 for mathi,j/math = 1,2, -1 for mathi,j/math = 2,1; and 0 for mathi=j/math.


= Numerical Methods =
= Numerical Methods =
== General Transport Equation: Discretization ==
== General Transport Equation: Discretization ==
All of the governing equations may be written in general form
All of the governing equations may be written in general form
{{Equation| <math>\underbrace{\frac{\partial (h\phi )}{\partial t}}_{\text{Temporal Term}}+\underbrace{\nabla \cdot (h\mathbf{U}\phi )}_{\text{Advection Term}}=\underbrace{\nabla \cdot \left( \Gamma h\nabla \phi  \right)}_{\text{Diffusion Term}}+\underbrace{S}_{\text{Source Term}} </math>|2=3}}   
{{Equation| math\underbrace{\frac{\partial (h\phi )}{\partial t}}_{\text{Temporal Term}}+\underbrace{\nabla \cdot (h\mathbf{U}\phi )}_{\text{Advection Term}}=\underbrace{\nabla \cdot \left( \Gamma h\nabla \phi  \right)}_{\text{Diffusion Term}}+\underbrace{S}_{\text{Source Term}} /math|2=3}}   
where <math> \phi </math> is a general scalar, <math> t </math> is time, <math> h </math> is the total water depth, <math> \mathbf{U} </math> is the depth averaged current velocity, <math> \Gamma </math> is the diffusion coefficient for <math> \phi </math>, <math> \nabla =({{\nabla }_{1}},{{\nabla }_{2}}) </math> is the gradient operator, and <math> S </math> includes all other terms. Note that in the case of the continuity and momentum equations <math> \phi </math> is equal to 1 and <math> U_i </math> respectively.  
where math \phi /math is a general scalar, math t /math is time, math h /math is the total water depth, math \mathbf{U} /math is the depth averaged current velocity, math \Gamma /math is the diffusion coefficient for math \phi /math, math \nabla =({{\nabla }_{1}},{{\nabla }_{2}}) /math is the gradient operator, and math S /math includes all other terms. Note that in the case of the continuity and momentum equations math \phi /math is equal to 1 and math U_i /math respectively.  
=== Temporal Term ===
=== Temporal Term ===
The temporal term of the momentum equations is discretized using a first order implicit Euler scheme
The temporal term of the momentum equations is discretized using a first order implicit Euler scheme
{{Equation| <math> \int\limits_{A}{\frac{\partial (h\phi )}{\partial t}}\text{d}A=\frac{\partial }{\partial t}\int\limits_{A}{h\phi \text{d}A}=\frac{{{h}^{n+1}}\phi _{{}}^{n+1}-{{h}^{n}}\phi _{{}}^{n}}{\Delta t}\Delta A </math>|2=3}}   
{{Equation| math \int\limits_{A}{\frac{\partial (h\phi )}{\partial t}}\text{d}A=\frac{\partial }{\partial t}\int\limits_{A}{h\phi \text{d}A}=\frac{{{h}^{n+1}}\phi _{{}}^{n+1}-{{h}^{n}}\phi _{{}}^{n}}{\Delta t}\Delta A /math|2=3}}   
where <math> \Delta A </math> is the cell area, and <math> \Delta t </math> is the hydrodynamic time step.
where math \Delta A /math is the cell area, and math \Delta t /math is the hydrodynamic time step.


=== Advection Term ===
=== Advection Term ===
Line 25: Line 27:
where  is the outward unit normal on cell face f,  is the cell face length and  is the total water depth linearly interpolated to the cell face. Here the overbar indicates a cell face interpolation operator described in the following section. For Cartesian grids the cell face unit vector is always aligned with one of the Cartesian coordinates which simplifies the calculation. Defining the cell face normal velocity as  the above equation simplifies to   
where  is the outward unit normal on cell face f,  is the cell face length and  is the total water depth linearly interpolated to the cell face. Here the overbar indicates a cell face interpolation operator described in the following section. For Cartesian grids the cell face unit vector is always aligned with one of the Cartesian coordinates which simplifies the calculation. Defining the cell face normal velocity as  the above equation simplifies to   


{{Equation| <math> \int\limits_{A}{\nabla \cdot (h\mathbf{U}\phi )}\text{d}A=\oint\limits_{L}{h\phi \left( \mathbf{U}\cdot \mathbf{n} \right)}\text{d}L=\sum\limits_{f}^{{}}{{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{U}_{i}} \right)}_{f}}{{{\tilde{\phi }}}_{f}}} </math>|2=4}}   
{{Equation| math \int\limits_{A}{\nabla \cdot (h\mathbf{U}\phi )}\text{d}A=\oint\limits_{L}{h\phi \left( \mathbf{U}\cdot \mathbf{n} \right)}\text{d}L=\sum\limits_{f}^{{}}{{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{U}_{i}} \right)}_{f}}{{{\tilde{\phi }}}_{f}}} /math|2=4}}   
    
    
where <math> \mathbf{n}={{\hat{n}}_{i}}=({{\hat{n}}_{1}},{{\hat{n}}_{2}}) </math> is the outward unit normal on cell face f, <math> \Delta {{l}_{f}} </math> is the cell face length and <math> {{\bar{h}}_{f}} </math> is the total water depth linearly interpolated to the cell face. Here the overbar indicates a cell face interpolation operator described in the following section. For Cartesian grids the cell face unit vector is always aligned with one of the Cartesian coordinates which simplifies the calculation. Defining the cell face normal velocity as <math> {{U}_{f}}={{U}_{i}}\in f\bot i </math> the above equation simplifies to  
where math \mathbf{n}={{\hat{n}}_{i}}=({{\hat{n}}_{1}},{{\hat{n}}_{2}}) /math is the outward unit normal on cell face f, math \Delta {{l}_{f}} /math is the cell face length and math {{\bar{h}}_{f}} /math is the total water depth linearly interpolated to the cell face. Here the overbar indicates a cell face interpolation operator described in the following section. For Cartesian grids the cell face unit vector is always aligned with one of the Cartesian coordinates which simplifies the calculation. Defining the cell face normal velocity as math {{U}_{f}}={{U}_{i}}\in f\bot i /math the above equation simplifies to  
{{Equation| <math> \sum\limits_{f}^{{}}{{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{U}_{i}} \right)}_{f}}{{{\tilde{\phi }}}_{f}}}=\sum\limits_{f}^{{}}{{{n}_{f}}{{F}_{f}}{{{\tilde{\phi }}}_{f}}} </math>|2=5}}
{{Equation| math \sum\limits_{f}^{{}}{{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{U}_{i}} \right)}_{f}}{{{\tilde{\phi }}}_{f}}}=\sum\limits_{f}^{{}}{{{n}_{f}}{{F}_{f}}{{{\tilde{\phi }}}_{f}}} /math|2=5}}


where <math> {{F}_{f}}={{\bar{h}}_{f}}\Delta {{l}_{f}}{{U}_{f}} </math>, <math> {{n}_{f}}={{n}_{\bot }}={{\left( {{{\hat{e}}}_{i}}{{{\hat{n}}}_{i}} \right)}_{f}} </math>, with <math> {{\hat{e}}^{i}}=({{\hat{e}}_{1}},{{\hat{e}}_{2}}) </math> being the basis vector. <math> n_f </math> is equal to 1 for West and South faces and equal to -1 for North and East cell faces. Lastly, <math> \tilde{\phi }_{f}^{{}} </math> is the advective value of <math> \phi </math> on cell face f, and is calculated using either the Hybrid, Exponential, HLPA (Zhu 1991) schemes. The cell face velocities <math> U_f </math> are calculated using the momentum interpolation method of Rhie and Chow (1983) described in the subsequent section.  The advection value is calculated as <math> {{\tilde{\phi }}_{f}}=\tilde{\phi }_{f}^{L(\exp )}+  \tilde{\phi }_{f}^{H(\text{imp})}-\tilde{\phi }_{f}^{L(\text{imp})} </math>, where the superscripts <math>L</math> and <math>H</math> indicate low and high order approximations and the superscripts <math>(exp)</math> and <math>(imp)</math> indicate either explicit and implicit treatment. The explicit term is solved directly while the implicit term is implemented through a deferred correction in which the terms are approximated using the values from the previous iteration step.
where math {{F}_{f}}={{\bar{h}}_{f}}\Delta {{l}_{f}}{{U}_{f}} /math, math {{n}_{f}}={{n}_{\bot }}={{\left( {{{\hat{e}}}_{i}}{{{\hat{n}}}_{i}} \right)}_{f}} /math, with math {{\hat{e}}^{i}}=({{\hat{e}}_{1}},{{\hat{e}}_{2}}) /math being the basis vector. math n_f /math is equal to 1 for West and South faces and equal to -1 for North and East cell faces. Lastly, math \tilde{\phi }_{f}^{{}} /math is the advective value of math \phi /math on cell face f, and is calculated using either the Hybrid, Exponential, HLPA (Zhu 1991) schemes. The cell face velocities math U_f /math are calculated using the momentum interpolation method of Rhie and Chow (1983) described in the subsequent section.  The advection value is calculated as math {{\tilde{\phi }}_{f}}=\tilde{\phi }_{f}^{L(\exp )}+  \tilde{\phi }_{f}^{H(\text{imp})}-\tilde{\phi }_{f}^{L(\text{imp})} /math, where the superscripts mathL/math and mathH/math indicate low and high order approximations and the superscripts math(exp)/math and math(imp)/math indicate either explicit and implicit treatment. The explicit term is solved directly while the implicit term is implemented through a deferred correction in which the terms are approximated using the values from the previous iteration step.


=== Cell-face interpolation operator ===
=== Cell-face interpolation operator ===
The general formula for estimating the cell-face value of <math> \tilde{\phi }_{f}^{{}} </math> is given by  
The general formula for estimating the cell-face value of math \tilde{\phi }_{f}^{{}} /math  is given by  
{{Equation| <math> {{\bar{\phi }}_{f}}={{L}_{\bot }}{{\phi }_{N}}+(1-{{L}_{\bot }}){{\phi }_{P}}+{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{N}}{{L}_{\bot }}({{x}_{\parallel ,O}}-{{x}_{\parallel ,N}})+{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{P}}(1-{{L}_{\bot }})({{x}_{\parallel ,O}}-{{x}_{\parallel ,P}}) </math>|2=6}}
{{Equation| math {{\bar{\phi }}_{f}}={{L}_{\bot }}{{\phi }_{N}}+(1-{{L}_{\bot }}){{\phi }_{P}}+{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{N}}{{L}_{\bot }}({{x}_{\parallel ,O}}-{{x}_{\parallel ,N}})+{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{P}}(1-{{L}_{\bot }})({{x}_{\parallel ,O}}-{{x}_{\parallel ,P}}) /math|2=6}}
where <math> {{L}_{\bot }} </math> is a linear interpolation factor given by <math> {{L}_{\bot }}=\Delta {{x}_{\bot ,P}}/(\Delta {{x}_{\bot ,P}}+\Delta {{x}_{\bot ,N}}) </math> and <math> {{\nabla }_{\parallel }} </math> is the gradient operator in the direction parallel to face f. By definition <math> \parallel \,=2\left| {{{\hat{n}}}_{1}} \right|+1\left| {{{\hat{n}}}_{2}} \right| </math>. Note that for neighboring cells without any refinement <math> {{x}_{\parallel ,O}}-{{x}_{\parallel ,P}} </math> and <math>{{x}_{\parallel ,O}}-{{x}_{\parallel ,N}} </math> are zero and thus the above equation is consistent with non-refined cell faces.
where math {{L}_{\bot }} /math is a linear interpolation factor given by math {{L}_{\bot }}=\Delta {{x}_{\bot ,P}}/(\Delta {{x}_{\bot ,P}}+\Delta {{x}_{\bot ,N}}) /math and math {{\nabla }_{\parallel }} /math is the gradient operator in the direction parallel to face f. By definition math \parallel \,=2\left| {{{\hat{n}}}_{1}} \right|+1\left| {{{\hat{n}}}_{2}} \right| /math. Note that for neighboring cells without any refinement math {{x}_{\parallel ,O}}-{{x}_{\parallel ,P}} /math and math{{x}_{\parallel ,O}}-{{x}_{\parallel ,N}} /math are zero and thus the above equation is consistent with non-refined cell faces.


=== Diffusion term ===
=== Diffusion term ===
The diffusion term is discretized in general form using the divergence theorem  
The diffusion term is discretized in general form using the divergence theorem  
{{Equation| <math> \int\limits_{A}{\nabla \cdot \left( \Gamma h\nabla \phi  \right)}\text{d}A=\oint\limits_{S}{\Gamma h\left( \nabla \phi \cdot \mathbf{n} \right)}\text{d}S=\sum\limits_{f}^{{}}{\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{\nabla }_{i}}\phi  \right)}_{f}}} </math> |2=7}}
{{Equation| math \int\limits_{A}{\nabla \cdot \left( \Gamma h\nabla \phi  \right)}\text{d}A=\oint\limits_{S}{\Gamma h\left( \nabla \phi \cdot \mathbf{n} \right)}\text{d}S=\sum\limits_{f}^{{}}{\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{\nabla }_{i}}\phi  \right)}_{f}}} /math |2=7}}


The discritization of the cell-face gradient is described in the next section. On a Cartesian grid the above expression may be further simplified as
The discritization of the cell-face gradient is described in the next section. On a Cartesian grid the above expression may be further simplified as


{{Equation| <math> \sum\limits_{f}^{{}}{{{n}_{f}}\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{\nabla }_{\bot }}\phi  \right)}_{f}}}=\sum\limits_{f}^{{}}{{{D}_{f}}\left[ {{\phi }_{N}}-{{\phi }_{P}}+{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{N}}\left( {{x}_{\parallel ,O}}-{{x}_{\parallel ,N}} \right)-{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{P}}\left( {{x}_{\parallel ,O}}-{{x}_{\parallel ,P}} \right) \right]} </math> |2=8}}
{{Equation| math \sum\limits_{f}^{{}}{{{n}_{f}}\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}{{\left( {{\nabla }_{\bot }}\phi  \right)}_{f}}}=\sum\limits_{f}^{{}}{{{D}_{f}}\left[ {{\phi }_{N}}-{{\phi }_{P}}+{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{N}}\left( {{x}_{\parallel ,O}}-{{x}_{\parallel ,N}} \right)-{{\left( {{\nabla }_{\parallel }}\phi  \right)}_{P}}\left( {{x}_{\parallel ,O}}-{{x}_{\parallel ,P}} \right) \right]} /math |2=8}}
where <math> {{\nabla }_{\bot }}\phi </math> is gradient in the direction perpendicular to the cell face and  
where math {{\nabla }_{\bot }}\phi /math is gradient in the direction perpendicular to the cell face and  
<math> {{D}_{f}}=\frac{\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}}{\left| \delta {{x}_{\bot }} \right|} </math>.
math {{D}_{f}}=\frac{\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}}{\left| \delta {{x}_{\bot }} \right|} /math.
   
   
=== Cell-centered node-based gradient operator ===
=== Cell-centered node-based gradient operator ===
Line 78: Line 80:
----
----
= Variable Index =
= Variable Index =
{| border="1"
{| border=1
! Symbol !! Description !! Units
! Symbol !! Description !! Units
|-
|-
| <math> t </math> || Time || sec
| math t /math || Time || sec
|-
|-
| <math> h </math> ||  Total water depth <math> h = \zeta + \eta </math> || m
| math h /math ||  Total water depth math h = \zeta + \eta /math || m
|-
|-
| <math> \zeta </math> ||  Still water depth || m
| math \zeta /math ||  Still water depth || m
|-
|-
| <math> \eta </math> ||  Water surface elevation with respect to the still water elevation || m
| math \eta /math ||  Water surface elevation with respect to the still water elevation || m
|-
|-
| <math> U_j </math> || Current velocity in the jth direction || m/sec
| math U_j /math || Current velocity in the jth direction || m/sec
|-
|-
| <math> S </math> || Sum of Precipitation and evaporation per unit area || m/sec
| math S /math || Sum of Precipitation and evaporation per unit area || m/sec
|-
|-
| <math> g </math> || Gravitational constant || m/sec<sup>2</sup>
| math g /math || Gravitational constant || m/secsup2/sup
|-
|-
| <math> \rho </math> || Water density || kg/m<sup>3</sup>
| math \rho /math || Water density || kg/msup3/sup
|-
|-
| <math> p_a  </math> || Atmospheric pressure || Pa
| math p_a  /math || Atmospheric pressure || Pa
|-
|-
| <math> \nu_t  </math> || Turbulent eddy viscosity || m<sup>2</sup>/sec
| math \nu_t  /math || Turbulent eddy viscosity || msup2/sup/sec
|}
|}


----
----
</big>
/big
[[CMS#Documentation_Portal | Documentation Portal]]
[[CMS#Documentation_Portal | Documentation Portal]]

Revision as of 17:47, 25 June 2011

big

Governing Equations

The depth-averaged 2-D continuity equation may be written as

  {{{1}}} (1)

where mathh/math is the total water depth mathh=\zeta+\eta/math, math\eta/math is the water surface elevation, math\zeta/math is the still water depth, mathU_i/math is the depth-averaged Langrangian current velocity defined as math U_i=U_i^E+U_i^S/math, where mathU_i^E/math is the phase- and depth-averaged current velocity (i.e. Eulerian velocity), and mathU_i^S/math is the depth-averaged Stokes velocity (Phillips 1977)

  {{{1}}} (2)

The momentum equation can be written as

  {{{1}}} (2)

where mathg/math is the gravitational constant, mathf_c/math is the Coriolis parameter, mathp_a/math is the atmospheric pressure, math\rho_0/math is a reference water density, math\nu_t/math is the turbulent eddy viscosity, math \tau_{w} /math is the wind stress, math \tau_{S} /math is the wave stress, and math\tau_{b}/math is the combined wave-current mean bed shear stress. math\varepsilon_{ij}/math is the permutation parameter equal to 1 for mathi,j/math = 1,2, -1 for mathi,j/math = 2,1; and 0 for mathi=j/math.

Numerical Methods

General Transport Equation: Discretization

All of the governing equations may be written in general form

  math\underbrace{\frac{\partial (h\phi )}{\partial t ({{{2}}})

_{\text{Temporal Term}}+\underbrace{\nabla \cdot (h\mathbf{U}\phi )}_{\text{Advection Term}}=\underbrace{\nabla \cdot \left( \Gamma h\nabla \phi \right)}_{\text{Diffusion Term}}+\underbrace{S}_{\text{Source Term}} /math|2=3}}

where math \phi /math is a general scalar, math t /math is time, math h /math is the total water depth, math \mathbf{U} /math is the depth averaged current velocity, math \Gamma /math is the diffusion coefficient for math \phi /math, math \nabla =({{\nabla }_{1}},{{\nabla }_{2}}) /math is the gradient operator, and math S /math includes all other terms. Note that in the case of the continuity and momentum equations math \phi /math is equal to 1 and math U_i /math respectively.

Temporal Term

The temporal term of the momentum equations is discretized using a first order implicit Euler scheme

  math \int\limits_{A}{\frac{\partial (h\phi )}{\partial t ({{{2}}})

\text{d}A=\frac{\partial }{\partial t}\int\limits_{A}{h\phi \text{d}A}=\frac{{{h}^{n+1}}\phi _{{}}^{n+1}-{{h}^{n}}\phi _{{}}^{n}}{\Delta t}\Delta A /math|2=3}}

where math \Delta A /math is the cell area, and math \Delta t /math is the hydrodynamic time step.

Advection Term

The advection scheme obtained using the divergence theorem as where is the outward unit normal on cell face f, is the cell face length and is the total water depth linearly interpolated to the cell face. Here the overbar indicates a cell face interpolation operator described in the following section. For Cartesian grids the cell face unit vector is always aligned with one of the Cartesian coordinates which simplifies the calculation. Defining the cell face normal velocity as the above equation simplifies to

  {{{1}}} ({{{2}}})

\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{U}_{i}} \right)}_{f}}{{{\tilde{\phi }}}_{f}}} /math|2=4}}

where math \mathbf{n}={{\hat{n}}_{i}}=({{\hat{n}}_{1}},{{\hat{n}}_{2}}) /math is the outward unit normal on cell face f, math \Delta {{l}_{f}} /math is the cell face length and math {{\bar{h}}_{f}} /math is the total water depth linearly interpolated to the cell face. Here the overbar indicates a cell face interpolation operator described in the following section. For Cartesian grids the cell face unit vector is always aligned with one of the Cartesian coordinates which simplifies the calculation. Defining the cell face normal velocity as math {{U}_{f}}={{U}_{i}}\in f\bot i /math the above equation simplifies to

  math \sum\limits_{f}^{{}}{{{{\bar{h}}}_{f ({{{2}}})

\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{U}_{i}} \right)}_{f}}{{{\tilde{\phi }}}_{f}}}=\sum\limits_{f}^{{}}{{{n}_{f}}{{F}_{f}}{{{\tilde{\phi }}}_{f}}} /math|2=5}}

where math {{F}_{f}}={{\bar{h}}_{f}}\Delta {{l}_{f}}{{U}_{f}} /math, math {{n}_{f}}={{n}_{\bot }}={{\left( {{{\hat{e}}}_{i}}{{{\hat{n}}}_{i}} \right)}_{f}} /math, with math {{\hat{e}}^{i}}=({{\hat{e}}_{1}},{{\hat{e}}_{2}}) /math being the basis vector. math n_f /math is equal to 1 for West and South faces and equal to -1 for North and East cell faces. Lastly, math \tilde{\phi }_{f}^{{}} /math is the advective value of math \phi /math on cell face f, and is calculated using either the Hybrid, Exponential, HLPA (Zhu 1991) schemes. The cell face velocities math U_f /math are calculated using the momentum interpolation method of Rhie and Chow (1983) described in the subsequent section. The advection value is calculated as math {{\tilde{\phi }}_{f}}=\tilde{\phi }_{f}^{L(\exp )}+ \tilde{\phi }_{f}^{H(\text{imp})}-\tilde{\phi }_{f}^{L(\text{imp})} /math, where the superscripts mathL/math and mathH/math indicate low and high order approximations and the superscripts math(exp)/math and math(imp)/math indicate either explicit and implicit treatment. The explicit term is solved directly while the implicit term is implemented through a deferred correction in which the terms are approximated using the values from the previous iteration step.

Cell-face interpolation operator

The general formula for estimating the cell-face value of math \tilde{\phi }_{f}^{{}} /math is given by

  math {{\bar{\phi }}_{f ({{{2}}})

={{L}_{\bot }}{{\phi }_{N}}+(1-{{L}_{\bot }}){{\phi }_{P}}+{{\left( {{\nabla }_{\parallel }}\phi \right)}_{N}}{{L}_{\bot }}({{x}_{\parallel ,O}}-{{x}_{\parallel ,N}})+{{\left( {{\nabla }_{\parallel }}\phi \right)}_{P}}(1-{{L}_{\bot }})({{x}_{\parallel ,O}}-{{x}_{\parallel ,P}}) /math|2=6}}

where math {{L}_{\bot }} /math is a linear interpolation factor given by math {{L}_{\bot }}=\Delta {{x}_{\bot ,P}}/(\Delta {{x}_{\bot ,P}}+\Delta {{x}_{\bot ,N}}) /math and math {{\nabla }_{\parallel }} /math is the gradient operator in the direction parallel to face f. By definition math \parallel \,=2\left| {{{\hat{n}}}_{1}} \right|+1\left| {{{\hat{n}}}_{2}} \right| /math. Note that for neighboring cells without any refinement math {{x}_{\parallel ,O}}-{{x}_{\parallel ,P}} /math and math{{x}_{\parallel ,O}}-{{x}_{\parallel ,N}} /math are zero and thus the above equation is consistent with non-refined cell faces.

Diffusion term

The diffusion term is discretized in general form using the divergence theorem

  {{{1}}} ({{{2}}})

\Delta {{l}_{f}}{{\left( {{{\hat{n}}}_{i}}{{\nabla }_{i}}\phi \right)}_{f}}} /math |2=7}}

The discritization of the cell-face gradient is described in the next section. On a Cartesian grid the above expression may be further simplified as

  math \sum\limits_{f}^{{}}{{{n}_{f}}\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f ({{{2}}})

\Delta {{l}_{f}}{{\left( {{\nabla }_{\bot }}\phi \right)}_{f}}}=\sum\limits_{f}^{{}}{{{D}_{f}}\left[ {{\phi }_{N}}-{{\phi }_{P}}+{{\left( {{\nabla }_{\parallel }}\phi \right)}_{N}}\left( {{x}_{\parallel ,O}}-{{x}_{\parallel ,N}} \right)-{{\left( {{\nabla }_{\parallel }}\phi \right)}_{P}}\left( {{x}_{\parallel ,O}}-{{x}_{\parallel ,P}} \right) \right]} /math |2=8}}

where math {{\nabla }_{\bot }}\phi /math is gradient in the direction perpendicular to the cell face and math {{D}_{f}}=\frac{\bar{\Gamma }_{f}^{{}}{{{\bar{h}}}_{f}}\Delta {{l}_{f}}}{\left| \delta {{x}_{\bot }} \right|} /math.

Cell-centered node-based gradient operator

Cell-centered face-based gradient operator

Cell-face gradient operator

Source terms

Hydrodynamic Solver

Wetting and drying

In the numerical simulation of the surface water flows with sloped beaches, sand bars and islands, the water edges change with time, with part of the nodes being possibly wet or dry. In the present model, a threshold flow depth (a small value such as 0.02 m in field cases) is used to judge drying and wetting. If the flow depth on a node is larger than the threshold value, this node is considered to be wet, and if the flow depth is lower than the threshold value, this node is dry. Because a fully implicit solver is used in the present model, all the wet and dry nodes participate in the solution. Dry nodes are assigned a zero velocity. On the water edges between the dry and wet nodes, the wall-function approach is applied.

References

  • Buttolph, A. M., Reed, C. W., Kraus, N. C., Ono, N., Larson, M., Camenen, B., Hanson, H.,Wamsley, T., and Zundel, A. K. (2006). “Two-dimensional depth-averaged circulation model CMS-M2D: Version 3.0, Report 2: Sediment transport and morphology change,” Tech. Rep. ERDC/CHL TR-06-9, U.S. Army Engineer Research and Development Center, Coastal and Hydraulic Engineering, Vicksburg, MS.
  • Ferziger, J. H., and Peric, M. (1997). “Computational Methods for Fluid Dynamics”, Springer-Verlag, Berlin/New York, 226 p.
  • Huynh-Thanh, S., and Temperville, A. (1991). “A numerical model of the rough turbulent boundary layer in combined wave and current interaction,” in Sand Transport in Rivers, Estuaries and the Sea, eds. R.L. Soulsby and R. Bettess, pp.93-100. Balkema, Rotterdam.
  • Rhie, T.M. and Chow, A. (1983). “Numerical study of the turbulent flow past an isolated airfoil with trailing-edge separation”. AIAA J., 21, 1525–1532.
  • Saad, Y., (1993). “A flexible inner-outer preconditioned GMRES algorithm,” SIAM Journal Scientific Computing, 14, 461–469.
  • Saad, Y., (1994). “ILUT: a dual threshold incomplete ILU factorization,” Numerical Linear Algebra with Applications, 1, 387-402.
  • Saad, Y. and Schultz, M.H., (1986). “GMRES: A generalized minimal residual algorithm for solving nonsymmetric linear systems,” SIAM Journal of Scientific and Statistical, Computing, 7, 856-869.
  • Soulsby, R.L. (1995). “Bed shear-stresses due to combined waves and currents,” in Advanced in Coastal Morphodynamics, ed M.J.F Stive, H.J. de Vriend, J. Fredsoe, L. Hamm, R.L. Soulsby, C. Teisson, and J.C. Winterwerp, Delft Hydraulics, Netherlands. 4-20 to 4-23 pp.
  • Wu, W. (2004). “Depth-averaged 2-D numerical modeling of unsteady flow and nonuniform sediment transport in open channels,” Journal of Hydraulic Engineering, ASCE, 135(10) 1013-1024.
  • Wu, W., Sánchez, A., and Mingliang, Z. (2011). “An implicit 2-D shallow water flow model on an unstructured quadtree rectangular grid,” Journal of Coastal Research, [In Press]
  • Wu, W., Sánchez, A., and Mingliang, Z. (2010). “An implicit 2-D depth-averaged finite-volume model of flow and sediment transport in coastal waters,” Proceeding of the International Conference on Coastal Engineering, [In Press]
  • Van Doormal, J.P. and Raithby, G.D., (1984). Enhancements of the SIMPLE method for predicting incompressible fluid flows. Num. Heat Transfer, 7, 147–163.
  • Zhu, J. (1991). “A low-diffusive and oscillation-free convection scheme,”Communications in Applied Numerical Methods, 7, 225-232.
  • Zwart, P. J., Raithby, G. D., Raw, M. J. (1998). “An integrated space-time finite volume method for moving boundary problems”, Numerical Heat Transfer, B34, 257.

Variable Index

Symbol Description Units
math t /math Time sec
math h /math Total water depth math h = \zeta + \eta /math m
math \zeta /math Still water depth m
math \eta /math Water surface elevation with respect to the still water elevation m
math U_j /math Current velocity in the jth direction m/sec
math S /math Sum of Precipitation and evaporation per unit area m/sec
math g /math Gravitational constant m/secsup2/sup
math \rho /math Water density kg/msup3/sup
math p_a /math Atmospheric pressure Pa
math \nu_t /math Turbulent eddy viscosity msup2/sup/sec

/big Documentation Portal

Retrieved from "https://cirpwiki.info/index.php?title=CMS-Flow:Hydro_Eqs&oldid=7136"