【二维泥沙数学模型及工程应用问题探讨】

【摘要】：从泥沙连续性方程和悬移质河床变形方程出发，分析恢复饱和系数的起因和悬移质含沙量垂线分布特点，给出了新的恢复饱和系数；对于动边界问题，通过阻力项的变换，给出稳定性好、精度高的迭代格式；对长江南京上元门河段河床冲淤变化进行了分析计算，结果表明本文所给泥沙数学模型和计算方法合理，能够反映河床冲淤变化规律。

【关键词】：恢复饱和系数；动边界；河床变形；泥沙数学模型

【作者简介】：赵明登(1963-)，男，武汉大学副教授。

　　1 引言
　　随着计算机存储容量和计算速度的飞速发展，河流泥沙数值模拟也得到长足的发展。一维泥沙数学模型虽然在理论上做了许多概化，但一维数学模型实践经验丰富，比较成熟，目前已能够解决许多长时期、长河段的河床冲淤问题。二维泥沙数学模型虽然在理论上较一维模型完善，又能够给出河床上任何部位的流速和冲淤变化情况，但由于实测资料和实践经验不足，目前还处于发展完善阶段。例如，二维糙率、挟沙力沿河宽分布问题、非均匀沙分组挟沙力及级配问题、恢复饱和系数和动边界问题等都有待进一步完善。本文从泥沙连续性方程和悬移质河床变形方程出发，分析恢复饱和系数的起因和悬移质含沙量沿垂线分布特点，给出了新的恢复饱和系数，并用长江南京上元门河段河床冲淤变化实测资料分析验证。

　　2 泥沙数学模型
　　2.1 水沙运动控制方程

　　2.2 悬移质河床变形方程
　　在河床表面上zb附近取一个微小棱柱体作为控制体研究，假设在床面上zb处悬移质时均含沙量为sb，在dt时间内，由于悬移质变化穿过床面zb的泥沙质量由两部分组成，一部分由于泥沙沉降引起，一部分由于泥沙紊动扩散引起，其大小为

4 模型应用
　　为了检验上述泥沙数学模型和恢复饱和系数合理性，对长江南京上元门河段河床冲淤变化进行了分析计算[5]。动量方程和泥沙连续性方程采用有限分析法离散，水流连续性方程和河床变形方程采用有限差分法离散[5]。
　　4.1 动边界处理方法
　　河道两岸的水边线位置是随着水位升降而变动的，这给数值计算带来了很大困难。河道非恒定流计算网格一般按高水位划分，当水位降低时，边滩露出水面，边滩上的网格处在水边线以外，从物理本质来说，此时边滩上水深、流速应为0，但在一般的迭代计算过程中，如果水深为0或很小时，常会出现计算机数值溢出或迭代发散。目前国内外许多学者常假设一个5～10cm的虚拟水深作近似处理，但由于虚拟水深较大，会引起水量、沙量和边滩冲淤量的计算误差。本文通过阻力项的变换，给出收敛性好、精度高的迭代格式，虚拟水深可以小到10-20cm，能够反映水沙运动的物理本质，新的迭代格式如下

　　4.2 非均匀沙处理方法
　　对于非均匀沙一般的处理方法是：将非均匀沙按粒径大小划分成若干个粒径组，将每个粒径组都看作是均匀沙，分别满足泥沙连续性方程和悬移质河床变形方程。对于非均匀沙模拟计算的关键问题是分组挟沙力计算问题。对于均匀沙来说，在泥沙粒径一定时，由于河床和悬移质中该粒径的含量都是百分之百，因此挟沙力只与水流强度有关。而对于非均匀沙来说，挟沙力不仅与泥沙粒径和水流强度有关，而且与河床中该粒径泥沙的含量多少有关，因为河床中含量越多，床面附近泥沙的紊动传递和悬富换机率越大。以某一粒径组泥沙在床沙中含量所占的百分比Pbn作为泥沙悬富换的机率，乘以该粒径组泥沙含量为百分之百时的可能挟沙力S′n*，即为该粒径组泥沙的实际挟沙力Sn*为

　　4.3 模型应用
　　模型验证资料为长江南京河段1997年10月至1998年10月实测水沙过程和1997年10月和1998年10月实测地形[5]。模型率定参数选取为：水流挟沙力公式中系数[6]Kn=0.03～0.07，m=1.0～2.0，恢复饱和系数中β=0.4～0.6。
　　图1分别为上元门附近河段1998年10月计算地形和实测地形比较图，图2为上元门断面计算地形和实测地形比较图。计算初始地形为1997年10月实测地形，1997年10月至1998年7月整个断面均为冲刷，1998年7月至1998年10月深槽有点回淤，深槽与边滩之间河床仍在冲刷。比较结果表明，经过一个大水中沙年的冲淤后，计算地形和实测地形符合良好，说明本文所给泥沙数学模型合理，能够反映河床冲淤变化规律。
　　利用本数学模型预测了京沪高速铁路南京长江沉管隧道工程对防洪航运及河床演变的影响[5]。计算结果表明，沉管工程阻水所产生的水位最大壅高为1cm，垂线平均流速最大变化为5cm/s，沉管附近河段左边河床受沉管控制不能冲刷，右边河床最大冲深0.5m，冲刷面积50m2左右。数学模型计算结果与河工模型试验结果相互印证，可以满足工程设计需要。

图1(a)上元门附近河段1998年10月计算地形图 　
　　　　　　　　　　　　　　　　　　　图1(b) 上元门附近河段1998年10月实测地形图

　　5 结论
　　从泥沙连续性方程和悬移质河床变形方程的推导过程可以看出，由于采用垂线平均含沙量代替河床底部含沙量而引入的系数α1、α2或α，从理论上讲，应该是大于1的，而实测资料率定的α小于1。产生这一矛盾的主要原因是，假定了非饱和悬移质垂线分布形式与饱和悬移质垂线分布形式相同，且α1=α2。本文通过分析系数的起因和悬移质含沙量垂线分布特点，给出了新的恢复饱和系数，可以避免上述矛盾且具有真正意义的调整恢复饱和程度。系数β取值可以由实测资料率定。

　　对于动边界问题，通过阻力项的变换，给出收敛性好、精度高的迭代格式，虚拟水深可以由通常的几厘米减小到10-20厘米；对长江南京上元门河段河床冲淤变化进行了分析计算，说明本文所给泥沙数学模型和计算方法合理，能够反映河床冲淤变化规律。

　　　　　　　　　　　　　图2 上元门断面地形比较图

参考文献

[1] 谢鉴衡。河流模拟，水利水电出版社，1988.6.
[2] Zhang Ruijin, and Xie Jianheng, Sedimentation Research in China, China Water & Power Press,1993.pp217-218.
[3] 杨国录。河流数学模型，海洋出版社，1993.10.
[4] 韦直林.付小平。论泥沙连续性方程的建立及相关的几个问题.武汉水利电力大学学报，1997.10，pp14-17.
[5] 赵明登。京沪高速铁路南京长江沉管隧道数值模拟研究报告，武汉水利电力大学，1999.6.
[6] 长江水利委员会长江下游水文水资源勘测局。南京河段泥沙基本特性分析。1999.1.
赵明登 李义天
（武汉大学河流工程系）

请教一下，迭代式（１４）中阻力项是否乘有有限分析系数cp

使用很方便哦,推荐一下,希望更多更好的全面资料上传分享,谢谢!楼主的精神值得我们大家学习。让我们共同提高水平

有没有软件啊 高手 东西很好 问题是 这些公式难道要手算出来吗？