黄河河口平原多闸坝河道水流数学模型

陈学群，李福林，张瑞青，刘 健

(1.山东省水利科学研究院，山东 济南 250013； 2.中国市政工程华北设计研究总院，天津 300074)

黄河下游入海口黄泛平原区属黄河冲积平原区，隶属山东省东营市。由于自然条件的限制，河口平原区地势平坦，河流坡度较小，大多数为人工河流，河道水流受强人工调控作用明显，大量的坝、闸、泵等水工构筑物的存在破坏了河流的水沙条件和河床形态的相对平衡，使河道水流条件复杂多变，物理模型和实验手段对于揭示平原地区多闸坝河道水流变化、水沙运动及水质运移的规律问题似乎显得有些困难，而数值模型则在这方面表现出了较强的优势[1]。

平原河网不同于单一河流的特点在于多闸坝综合调控及河网错综复杂性，由此带来模型的数值离散和求解上的困难是多年来人们研究河网问题的一大难点。MIKE11河网水动力模拟软件在河口、河流、河网的水量模拟及闸坝运行调度处理方面具有较强优势，笔者利用该软件对黄河下游河口平原地区河网进行水流模拟计算，进一步揭示多闸坝调控对于河流水动力的影响，为深入研究及下一步的河道综合治理提供技术依据。

图1 广利河河网及闸坝分布示意图

1 研究区概况

研究对象为黄河河口平原地区广利河及其支流溢洪河、东营河、胜利干渠、老广蒲沟、五六干合排等7条河流，总长度186km。广利河为贯穿东营市中心城区的主要排水河道，是连接黄河与渤海的唯一河流，在供水、排水、防洪和改善生态环境等方面发挥着重要作用。从起端王营闸到入海前的防潮堤，高程从3.90m降为-1.57m，比降比较平缓，为0.016%。广利河属潮汐往复河流，水文条件复杂，受到海水上溯的影响比较严重，此外还受沿岸支流及农田灌溉引退水和市区段沿岸居民、企事业单位排放污水的影响。沿岸共有4条主要支流相接，主要支流河口均有水闸控制，其中主要闸门7座，橡胶坝1座。河流边界为：北起广利河、溢洪河共同的源头王营闸，南到广利河入海前的广利港，西到六干渠与黄河的节点，东到广利港，具体如图1所示。

2 数学模型

一维河网水动力模型是用数学方程模拟自然界明渠非恒定流在河道中的流动规律。任何复杂河网的水力数值计算问题，都可以归结为对描述单一河道的一维明渠非恒定流的Saint-Venani方程组的求解问题：

式中：x为距离；t为时间；A为过水断面面积；Q为流量；h为水位；q为旁侧入流单宽流量；n为河床糙率系数；R为水力半径；g为重力加速度。

Saint-Venani方程组在数学上属于一阶拟线性双曲型偏微分方程组，其解析求解是非常困难的，目前只能用数值离散的方法求其近似解。MIKE11利用有限差分法Abbott六点隐式格式进行计算。利用Abbott六点隐式格式离散上述控制方程组[2]，该离散格式在每一个网格节点水位和流量并不同时计算，而是按照顺序交替对水位和流量进行计算。连续性方程及动量方程经一系列变形整理后可简记为

3 模型概化

3.1 河网概化

河网概化要合理，既不能过于复杂，又需要最大限度地反映河流的实际情况。概化后的河网能基本反映河网的水力特性是河网概化的原则，也就是说概化后的河网在输水能力上必须与实际河网基本一致。在研究范围内，河道大多为人工河道，数量较少，主干河道明显，所以河道的概化比较简单。根据广利河实际情况，广利河本身不构成独立的水系，一方面它通过引黄干渠、闸、坝，同周围黄河及其他水系发生联系，受到黄河上游来水及灌溉退水的影响，另一方面广利河河水直接入海，容易受到渤海湾潮汐的影响，属黄泛平原感潮河网。因此，为反映引黄水对广利河水系的影响，上游边界概化为流量边界，其中广利河干流上游流量变化如图2所示。

图2 广利河干流上游流量变化

受渤海湾潮汐对广利河水系的影响，下游边界广利河口设置为潮汐水位条件，其年时间尺度变化趋势受上游来水影响较大，日时间尺度变化受到潮汐作用影响较大，下游资料采用交通运输部天津水运工程科学研究所在广利河口-1m水深处设置的潮位观测站资料[3]，在监测广利河河口潮位的同时进行了潮位、流速及流向等的监测，结果见图3、图4。

图3 下游广利河口日尺度(2009年6月)水位、流速及流向变化

图4 下游广利河口年尺度(2009—2010年)水位变化

3.2 断面设置

在计算断面的选取过程中需要考虑河势的变化和河道的走向等问题，当河段比较曲折时，断面选取应该在拐点处；对于坝、闸等水工构筑物，应该把其作为内部边界来处理。

天然河道断面不规则，MIKE11中的断面编辑器既可以处理不规则断面，也可以处理规则断面，该系统以矩形断面、梯形断面为主。断面编辑器将断面资料分成两种：①断面原始数据，可以把断面实测数据简化为折点高程、河面宽度、过水断面、河道边坡之后输入断面编辑器。②计算数据，根据断面原始数据计算出水位、断面面积、调蓄宽度、额外调蓄能力、水力半径及糙率系数。河道典型断面见图5。

图5 河道典型断面

3.3 闸坝水工构筑物概化

广利河水系属人工调控的多闸坝平原河网地区，几乎全面受到闸坝的控制。每个闸门都有其运行规则，考虑灌溉、防洪、景观及调水的需要，闸门的运行规则较复杂，对河网水体的流动起到决定性的作用。将闸门的运行规则输入到模型的河网编辑器中，模型就可以比较真实地反映闸门的运行对河网水动力模拟结果的影响。由于闸门调度资料缺乏，而且在日常运行中工作人员操作的随意性较大，所以只选择部分有调度资料的闸门进行模拟。其中明港闸位于溢洪河设计桩号46km+000m处，按100年一遇防洪，设计防洪流量300m3/s，防洪水位2.555m；5年一遇排涝，设计排涝流量128m3/s，排涝水位1.255m；50年一遇挡潮，设计挡潮水位3.630m。控制水位小于1.5m时，闭闸，上游各闸开启；控制水位大于2.2m时开闸放水，上游各闸关闭。秦家闸位于广利河上游，为广利河引黄河水控制闸门，开启后，可引黄河水进入广利河。其主要闸门情况及运行规则如表1所示。

表1 闸门情况及运行规则

4 参数选择

4.1 时间步长的选择

对于显示格式的计算，时间步长需要服从柯朗准则[4]。而MIKE11采用隐式格式，从理论上对时间步长没有限制。在实际模拟计算过程中，如果时间步长过大，模型运行计算后的结果容易过于坦化而失真；如果在模型中取的值太小，模拟结果就会因有些非线性的小扰动最后导致计算的失稳[5]。

在本次模拟中，出于对MIKE11软件模型计算稳定性的考虑及计算过程中对运行时间考虑，最后确定MIKE11水动力模型的时间步长取5～15min。

4.2 河道糙率

河道糙率与河道形态、河床粗糙情况、河道弯曲程度、植被成长情况、水位高低、河槽的冲击以及河道上人工构筑物等因素有关。一般情况下，对糙率的确定采用实测水位资料进行推算的方法，而对那些没有实测资料的河道采取类比相似河道糙率或用经验公式法来确定。对河网糙率参数进行率定时一般采取手工调试的方法，在河网规模不大的情况下可行，但如果河网规模比较大，这种方法的工作量将是巨大的，而且手工调试带有很大的随意性。为了克服上述缺点，韩龙喜等[6]从糙率的物理意义出发，根据水力特性和河网规模，将河网分成若干等级，同一级河网给与相同糙率，并采用最优化中的复合形法求解糙率。本次研究过程中，不同的河段按照河道的形状、水位等实际状况来确定糙率，在率定过程中再根据模型做出相应的调整。

5 模拟结果

模型利用2009—2010年期间对广利河干支流进行的调查及实测的资料，结合东营市水利局提供的河流水位、流量等资料进行验证。边界水位测点为王营闸、溢洪河源头、胜利干渠、五六干合排起点、老广蒲沟源头、供水公园、广利港。考虑了皇殿闸、秦家闸、五六干合排闸、明海闸、明港闸的调度运行过程。经过调算，广利河、溢洪河糙率在0.08～0.36之间，东营河、五六干渠等支流的河道糙率为0.10。模拟结果如图6所示。

由模型的率定结果可知，广利河干流的明海闸断面、广利港断面等测站的计算结果与实测值吻合较好，并且宏观上呈现出汛期水位较高、其他月份水位相对较低的情景，说明下游水位变化趋势受上游来水流量的变化趋势影响较大，水流状况受人工调控明显。同时东营河及溢洪河各个断面的水位计算结果与实测值吻合稍差，究其原因，河网闸坝较多，在实际调度过程中受不确定性因素的影响比较大，同时河流入海口受到潮汐作用的影响较大，但大部分断面的模拟误差仍能满足工程精度要求。

图6 广利河干支流典型监测断面模拟值与实测值对比

6 结 论

在充分考虑了闸坝的调度运行原则及潮汐作用影响的基础上，采用MIKE11模型对黄河河口平原河网概化后计算水位与实测水位进行了数值模拟，经过模型率定，广利河、溢洪河糙率在0.08～0.36之间，东营河、五六干渠等支流的河道糙率为0.10；同时，计算表明广利河河网下游水位变化趋势受上游来水流量的变化趋势影响较大，水流状况受人工调控明显。模拟结果较为合理，可以为下一步进行水量水质耦合模拟及河网库群与闸坝调度方案研究提供较为准确的水动力条件。

[1] 张瑞青.基于生态修复目标的黄泛平原河道水流水质数值模拟研究[D].济南：山东建筑大学，2011.

[2] Danish Hydraulic Institute(DHI).MIKEII：a modelling system for river sand channels user-guide manual[R].Copenhagen:DHI，2008.

[3] 交通运输部天津水运工程科学研究所.东营市滨海生态城防潮堤工程数学模型试验研究报告[R].天津：交通运输部天津水运工程科学研究所，2011.

[4] 周雪漪.计算水力学[M].北京：清华大学出版社，1994.

[5] 汪德.计算水力学理论与应用[M].南京：河海大学出版社，1989.

[6] 韩龙喜，金忠青.平原河网水量计算及参数率定：改进的湄公河模型[J].水动力学研究与进展，1991(6)：39-45.