珠江三角洲网河低水位变化的二维水沙数学模型

河流泥沙一维数学模型，只能给出水力、泥沙因子沿河道方向的变化，不能得到水力、泥沙因子沿河宽方向的变化。即使采用流管法可得到准二维的计算结果，李义天和吴伟明［12］的研究表明，冲淤沿河宽的分布对于平顺冲积性河段计算结果基本上能反映实际，但对于河宽变化较大的河段计算结果可靠性较差。分析其原因，主要是河宽变化剧烈的河段洪枯水主流位置和流速分布差别较大，且在突嘴、人工建筑物或江心洲尾部常有各种各样的副流（如回流等），而这些在一维数学模型中不能体现出来。因此，需建立平面二维数学模型对河道水流、悬移质含沙浓度、推移质输沙率平面分布及河床的平面变形进行比较细致地模拟。

5.1.2.1　基本方程（控制方程）

控制水沙运动的基本方程在沿水深按静水压力分布等的假定条件下，沿水深积分三维Reynolds方程和泥沙对流扩散方程得到

水流连续方程：

水流运动方程：

悬沙对流扩散方程：

式中　σs——经验常数，采用Rodi建议的值。

河床变形方程：

水流挟沙能力公式：

推移质输沙率公式：

式中　h——水深；

　　　u、v——垂线平均流速在x、y方向的分量；

　　　H——水位；

　　　νt——紊动黏性系数；

　　　n——满宁糙率系数；

　　　g——重力加速度；

　　　f——柯氏力系数。

5.1.2.2　方程组的闭合

式（5.1）～式（5.7）构成了描述水沙运动平面二维模型的控制方程，在采用式（5.1）～式（5.7）时，常忽略了水流运动的紊动效应，νt=0。这种做法对岸线比较平顺的河段而言是允许的，但对岸线不规则，有可能产生回流的河段而言，就不允许了［7］。这是因为回流的存在主要是以铅直面摩阻力的存在为前提的，主流区和回流的泥沙交换是通过紊动扩散作用实现的。

实际工程计算中，为了便于计算，常将νt取为常数，或采用经验公式，这样作法比不考虑νt要好一些，但往往回流的模拟误差较大。因此，准确模拟工程附近（如丁坝、港池开挖、挖槽）的流场及泥沙运动所引起的河床变形应采用完整的紊流模型（如K—ε模型，νt=Cμk2/ε）模拟水沙运动。(www.xing528.com)

紊动动能K输运方程：

紊动动能耗散率ε输运方程：

Pkv、Pεv与摩阻流速u*间的关系为

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

　　　Cμ、σk、σε、C1ε、C2ε——经验常数，采用Rodi建议的值；

源项Ph——紊动应力与平均流速梯度相互作用所产生的紊动能量，表达式为

5.1.2.3　阻力

水流运动方程中的阻力项通常用曼宁糙率系数n来表示，许多研究者认为，最可靠的办法是根据河流实测的流量及水面线资料反求n。这样做虽能较好反映出河流的客观情况，但无法预报出在水沙条件改变后或在工程修建后河槽中的阻力规律。对于冲积河流来说，n与床面形态（沙波）有关，而床面形态是随水流及输沙强度而变化的。

目前关于冲积河流阻力计算公式都是针对一维水流的。在复式河槽情况下，常将滩地与主槽分开，分别计算其糙率。在二维水流计算中如何计算阻力的研究很少。通常假定沿河宽不变即等于断面平均糙率。李义天曾经提出过一个糙率沿河宽变化的计算方法［17］，但问题远未获得解决（陈国祥、郁伟族，1989）［7］。

在建立模型时如何选用阻力公式问题，遵循的原则应是（Brownline，1983）［18］：①实验室和野外资料符合很好；②包含使用范围及误差分析；③易于在计算机模拟中使用；④能提供大范围独立变量解。因此，应收集大量实验室及野外资料尤其是工程前后水位、流速、河床形态、床沙组成等资料，采用统计分析方法建立不同河型河段阻力沿河宽方向变化的公式。

5.1.2.4　水流挟沙能力

近岸处垂线平均流速相对河道中间要小，由于河岸的影响，水流的紊动强度较大，使水流挟沙能力相对河道中间要大。

对于弯道及分汊河道，水流速度沿河宽分布很不均匀，导致水流挟沙能力沿河宽分布不均匀，对于不规则的岸边及工程附近常存在回流，其水流紊动及挟沙机理与通常河槽是不同的。

5.1.2.5　非均匀推移质输移及床沙交换计算

大多数平原河流河床变形主要是由悬移质运动引起的，推移质运动引起的河床变形仅占总量的5%～10%。因此，已有的二维全沙数学模型（窦国仁等，1987；周建军，1988；李义天，1989；陆永军等，1993b）［16，17，19，20］实际上大多验证了悬移质所引起的冲淤变化，由于缺乏野外及室内推移质运动引起的河床变形资料，对推移质运动的二维模型研究得不够。山区河流中河床多为卵石覆盖，当上游有砾石或粗沙通过某一河段时，这种砾石或粗沙推移质在卵石床面上输移，用一般公式计算推移质输沙率，由于卵石起动流速较大，水流速度常小于卵石起动流速而得到输沙率为零的错误结果，这就需要考虑有效床沙组成问题；其次是推移质与床沙总是在发生交换，交换的结果床沙级配发生变化，河床冲刷时床沙粗化，河床淤积时床沙细化，床沙交换引起的床沙级配变化机理及床沙交换层厚度怎样确定？港口航道工程中的丁坝及挖槽，常常由于推移质运动引起河床变形，推移质运动模拟问题在二维数模中尚未解决，同样在河工模型试验中也没有妥善解决。

5.1.2.6　计算方法

目前求解二维泥沙数学模型的数值方法很多，主要有差分法和有限元法。差分法主要有形式上比较简单的直接差分法，该方法简易明了，易于掌握和理解。鉴于时间和空间差分形式不同，出现了多种形式和差分格式，在时间处理上主要分为显式和隐式两种，有推导严谨的特征线法，该方法对原始方程进行处理，转化成沿特征方向的常微分方程再求解；以及交替方向法（ADI法），这一方法采用双时间层的方法将一个时间步分成两个半步长，在上半时步，对x方向的运动方程和连续方程采用隐格式，而y方向函数采用显式求解，在后半步反之；破开算子法是通过引进一个或若干中间变量，将偏微分方程中的时间微商破开成两个或更多的部分，得到相应的便于求解的若干个偏微分方程，且对每个偏微分方程选择不同的格式求解，具有一定的灵活性；控制体积法的基本思想是将计算区域划分成一系列控制体积，将待解的微分方程对每一个控制体积积分，便得出一组以网格上因变量为未知量的离散方程，由于对任意控制体积都体现准确的积分守恒，控制体积法的时间步长较ADI法及破开算子法要长得多，且计算稳定性好。

有限元法采用局部近似的低阶多项式作为试函数构成关于因变量结点值的代数方程。其优点是网格划分灵活，便于处理不规则边界，但其数学推演和编制程序比较复杂。传统的有限元法需要整体存储，因而占用内存大，计算费用昂贵。