珠三角网河低水位变化与控制方程

对于河道水流、泥沙的模拟，如何布置网格，使之贴合曲折边界，并反映随着水位变化的边界形状，同时克服计算域长宽比悬殊的困难，是一个关键问题。如果选择矩形网格，为了顾及边界形状及河道横断面上的地形，不得不采用众多尺度很小的单元，这就大大增加了对计算机内存的要求，计算工作量也显著增加。本文通过引入边界贴体坐标［21］，采用贴体正交曲线网格系统来克服边界复杂及计算域长宽尺度悬殊所引起的困难。

5.2.1.1　坐标转换关系的基本方程

平面直角坐标系下的任意形状的区域，通过边界贴体坐标，可转化为新坐标系下的规则区域。新旧坐标采用：ξ=ξ（x，y），η=η（x，y）函数关系联系起来。假定变换关系满足Poisson方程与Dirichllet边界条件。

Poisson方程：

Dirichllet边界条件：

新旧坐标的变换关系为

其中，P，Q是与ξ，η，x，y有关的某一函数，反映了（ξ，η）平面上等值线在（x，y）平面上的疏密程度，适当选择P，Q函数，可使坐标变换为正交变换。根据水流势函数与流函数的性质及水流等势线与等流线的正交性，可导出生成正交曲线网格的转换方程［22］：

这里，α=。

5.2.1.2　平面二维贴体坐标系下的水流运动方程

水流连续方程：

ξ方向动量方程：

η方向动量方程：

式中　ξ、η——正交曲线坐标系中二个正交曲线坐标；

　　　u、v——沿ξ、η方向的流速；

　　　h——水深；

　　　H——水位；

　　　Cξ、Cη——正交曲线坐标系中的拉梅系数；

　　　σξξ、σξη、σηξ、σηη——紊动应力；

　　　νt——紊动黏性系数，一般情况下，νt=αu*h（α=0.5～1.0，u*为摩阻流速）。

对于不规则岸边、整治建筑物、桥墩作用引起的回流，可采用的k—ε紊流模型为

式中　k——紊动动能；

　　　ε——紊动动能耗散率。

正交曲线坐标系下紊动动能输运方程为

紊动动能耗散率输运方程为

紊动动能产生项为

Pkv、Pεv与摩阻流速u*间的关系为(www.xing528.com)

式中　Pkv、Pεv——因底床切应力所引起的紊动效应；

　　　Cμ、σk、σε、C1ε、C2ε——经验常数，采用Rodi（1984）［23］建议的值：Cμ=0.09，σk=1.0，σε=1.3，C1ε=1.44、C2ε=1.92，σs=1.0。

5.2.1.3　贴体坐标系下的悬沙不平衡输移方程

非均匀悬移质按其粒径大小可分成n0组，用SL表示第L组粒径含沙量，PSL表示此粒径悬沙含沙量所占的比值，则

针对非均匀悬移质中第L组粒径的含沙量，二维悬移质不平衡输沙基本方程为

式中　——第L组泥沙的挟沙能力；

　　　——第L组泥沙的挟沙能力级配；

　　　S*（ω）——挟沙能力；

　　　ωL——第L组泥沙的沉速；

　　　K0——挟沙能力系数；

　　　αL——第L组泥沙的含沙量恢复饱和系数。

5.2.1.4　贴体坐标系下的推移质不平衡输移方程

非均匀推移质按其粒径大小可分成nb组，推移质不平衡输移基本方程为［24］

式中　——第L组推移质的挟沙能力；

　　　——单宽推移质输沙率；

　　　SbL——床面推移层的含沙浓度；

　　　——第L组推移质泥沙的恢复饱和系数；

　　　ωbL——推移质的沉速。

5.2.1.5　床沙级配方程

床沙级配方程为

此式是将CARICHAR混合层一维模型［25］扩展到二维模型的。式中，PmL、PmL0分别表示当时时刻床面混合层内以及混合层以下初始床沙级配；左端第四项的物理意义为混合层下界面在冲刷过程中将不断下切河床以求得下层床沙对混合层的补给，进而保证混合层内有足够的颗粒被冲刷而不至于亏损，当混合层在冲刷过程中波及到原始河床时取ε1=0，否则取1。

5.2.1.6　河床变形方程

河床变形方程为

河床总冲淤厚度：