专利名称：用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法
技术领域：
本发明涉及防洪减灾科学技术领域，特别是根据水量平衡方程式与动力方程式的 槽蓄增量关系求解圣维南方程组模拟洪水演进的方法。
背景技术：
我国是一个洪涝灾害频繁发生的国家，全国1/2的人口，1/3的耕地和主要大城市 都处于江河洪水位以下，受洪水威胁的地区工农业产值占全国2/3。为了达到减轻洪水灾害 以及利用水资源为国民经济与社会发展服务的目的，建立模拟洪水演进的数学模型，在河 道、渠道、分蓄洪区和水库合理确定防洪工程规模、洪水预报、科学调度洪水等方面发挥着 重要作用。描述河渠非恒定流运动规律的一阶拟线性双曲型偏微分方程，在1871年首先由 法国学者圣维南建立。在一般情况下，无法求出其普遍的解析解。目前，只能使用近似的计 算方法求解。近似的计算方法大致可分为两类1、水力学方法直接差分法就是其中的一种方法。直接差分法是用偏差商代替偏微商，将基本微 分方程离散化为差分方程，求在自变量域X t平面差分网格上各节点近似数值解的方法。 基本微分方程中的偏微商可用不同形式的偏差商即差分格式代替，从而得到不同的差分方 程。但从目前的技术手段而言，研究成果一般用于规划设计阶段，在流域水文预报模型中尚 未采用。2、水文学方法水文学方法如著名的“马斯京根法”，是将基本微分方程组分别简化为水量平衡方 程式和近似为蓄泄方程式，用差分形式合解水量平衡方程式与蓄泄方程式。缺点是①概化 相对较粗，难于详细分析河道内复杂水流情势以及主要特征；②难于提供河段内水位的时 空分布和动态模拟洪水在河段内的演进过程。
发明内容
本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法，是将圣维南 基本微分方程组在At时段简化为水量平衡方程式和动力方程式后，根据水量平衡方程式 和动力方程式可以分别计算河段槽蓄增量一是根据水量平衡方程式可计算槽蓄增量，二 是根据时段初、末的两个动力方程式也可计算槽蓄增量。由于两种方法计算的河段槽蓄增 量应相等，因此，这一发明就为水力学方法求解圣维南方程组模拟洪水演进提供了一种新 的、可靠的方法。本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法，于2009年11 月下旬已在湖北省水利水电勘测设计院进行的杜家台分蓄洪区可行性研究报告修编工作 中得到了初步应用，通过模拟洪水演进计算，结果发现没有必要兴建原报告的支洪道项目 (含出江口控制闸)，若实施，可节省投资2亿余元。
本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法，属于水力学方 法范畴。本方法既与水文学方法和现行的水力学方法有相似的地方，也有不同的地方。与水文学方法相似的地方是都涉及到了河段槽蓄内容；不同的是水文学方法如 “马斯京根法”是将动力方程式近似变为蓄泄方程式后与水量平衡方程式联解，而本方法是 直接将水量平衡方程式与动力方程式联解。与现行的水力学方法相比，相同之处是都涉及到了需要河道地形资料和依据水文 资料率定河段糙率；不同的是，本方法将基本微分方程在At时段简化为水量平衡方程式 和时段初、末的两个动力方程式后联解。具体方法如下现将圣维南基本微分方程组在△ t时段简化为水量平衡方程式(1)、动力方程式 (2)和动力方程式(3)
g±1 + g±2A^ = gT1 + gT2A + (ΡΓ 2 - fF 0-------------(1)
2 2 v 7
ΔΖ, = ^rL---------------------------⑵
ΔΖ2 = L---------------------------(3)
K22以上各式中=L为河段度长(m) ； Δ t为时段(s) ；Z为水位(m) ；Q为流量(m3/s) ；q 为区间流量(m3/s) ；W为河段槽蓄水量(m3) ； ΔΖ为河段上、下游计算断面的水位差(m) ； Q 为河段平均流量，AW= (W2-W1)；疋为河段平均流量模数；各量脚注1、2分别为时段初、 末的时刻。为了说明联解式(1)、(2)、(3)的方法，根据上述方程组绘制成河段槽蓄增量 △ W(水量平衡)和AW(动力)的图形。如图2所示，根据水量平衡方程式(1)绘制的图，槽蓄增量AW(水量平衡)是在 时段At内河段总入流量与总出流量的差值，即由入流过程Q±1 Q±2与出流过程Qti Qt2所包围的水体。如图3所示是根据动力方程式(2)、(3)绘制的，槽蓄增量AW(动力)是在时刻、 的水面线Z±1 Zti和时刻t2的水面线Z±2 Zt2所包围的水体，即由上、下游计算断面 分别在时刻^和t2的水位转换成过水断面积的差值平均后乘以河段长度得出的体积。当已知条件为Q±1、Z±1、QT1、ZT1、Qi2以及Qt2(或ZT2)时，可依据槽蓄增量 AW的关系求解Zi2和Qt2(或ΖΤ2)根据方程式(1)，如假定Qt2(或ΖΤ2)值，就可以计算出相应的河段槽蓄增量 Δ W (水量平衡)。根据方程式(2)、(3)，继续利用上述假定的Qt2 (或ΖΤ2)值，就可以计算出河段槽
5蓄增量AW(动力)。通过计算机程序按上述关系调节Qt2 (或ZT2)值，使两种方式计算的河段槽蓄增 量AW基本相等，此时的Ζ±2和Qt2 (或ΖΤ2)即为方程组的解。当计算断面个数为η时，则需要增加η-2个断面流量计算关系式求解方程组。如以上游起始断面的编号为i = 1，依次往下游的断面编号则为i = 2、i = 3-i =n_2、i = n_l、i = η。当断面编号1 < i < η时，根据水量平衡方程式(1)，如断面的流量和河段槽蓄增 量用二维数组表示，则断面i-ι与断面i之间、断面i与断面i+ι之间的水量平衡关系式分 别为Q(i-l，l)+Q(i-l，2)-Q(i，l)_Q(i，2) = 2 Δ W(i_l，i)/ Δ t---------------(4)Q(i，l)+Q(i，2)-Q(i+l，l)_Q(i+l，2) = 2 Δ W(i，i+1)/ Δ t---------------(5)联解方程式(4)、(5)Q(i，2) = (Aff(i, i+l)-Aff(i-l, i)) / Δ t+(Q (i_l，2)+Q (i+1，2))/2+(Q (i_l， l)+Q(i+l，l)-2Q(i，l))/2------------------------------(6)也可近似写成式(7)等形式Q(i，2) = (Aff(i, i+l)-Aff(i-l, i))/2/Δ t+(Q (i-1, 2)+Q (i+1, 2))/2------(7)式中二维数组流量Q括号内的第一项为断面编号，第二项为时段At的初、末代 号(初为1、末为2) ；二维数组槽蓄增量AW括号内的第一项为河段上游断面编号，第二项 为河段下游断面编号。因此，按水量平衡方程式计算的全河段累计槽蓄增量Σ AW(水量平衡)为Σ Δ W (水量平衡)=Δ t (Q (1,1) +Q (1, 2) -Q (η, 1) -Q (η, 2)) /2--------------(8)按动力方程式计算的全河段累计槽蓄增量Σ Δ W(动力)为
Σ 力)= ^W{i-\,i)-------------------------------(9)
图1为本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法假定糙 率系数为0. 031和0. 035时的1983年洪水杜家台闸计算与实测水位过程线对比图2为本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法根据水量平 衡方程式绘制的图3为本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法根据动力方 程式(2)、(3)绘制的图。
具体实施例方式以下对本发明的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法的具体 实施方式结合图1作进一步详细的说明，但本发明不局限于以下实施例。实施例一根据图1所述的用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法得出的1983年洪水杜家台闸计算与实测水位过程线对比图，分析杜家台闸至周帮河段1983年洪 水时的糙率系数杜家台闸至周帮河段(即洪道)长约21km，起止地点分别为杜家台闸和周帮，洪道 是引导汉江洪水进入分蓄洪区的通道。采用1998年实测洪道横断面28个，测量断面编号为1 28，测量精度为IV等水 准，断面高程为黄海基面。杜家台闸下和周帮的观测水位为冻结吴淞基面，换算关系为冻结吴淞高程=黄 海高程+1. 85m。洪水资料时间为1983年10月的8日7时至15日9时，水位观测时段为每 小时1次，共计171组数据。杜家台闸下断面的流量过程根据观测的闸上、下游水位按流量 公式计算而得。分析计算杜家台分蓄洪区洪道糙率系数的技术路线是首先应假定河段的糙率系数，如在表一中假定的糙率系数为0. 035 ；按照河段累计槽蓄增量Σ AW(水量平衡)和Σ AW(动力)应基本相等的原则 假定洪道下游末端周帮站的流量，如在表一中1983年10月13日11时假定断面编号为 28(周帮)的流量为2238. 14(m3/S)时，迭代计算累计河段槽蓄增量Σ AW(水量平衡)/ Δ t禾Π Σ Δ W (动力)/ Δ t均为-351· 39 (m3/s)，河段27至28的槽蓄增量Δ W (动力)/ Δ t 为-0. 49 (m3/s)。断面编号小于28 (周帮)大于1(杜家台闸下)的计算断面的流量按照式(7)迭代 计算，如在表一中1983年10月13日11时按式(7)计算断面编号27的流量为2238. 32 (m3/ s)；断面编号小于28的计算断面的水位按照式(3)迭代计算完成，在表一中1983年 10月13日11时按式(3)计算断面编号27的水位为28. 12m ；完成上述计算后，若在假定河段糙率系数的情况下，计算的杜家台间下水位过程 线与实测的水位过程线对比误差较小，则可认为该糙率系数即为所求；否则，应重新假定河 段糙率系数。本实施例假定河段糙率系数分别为0. 031和0. 035，从图1可以看到，当杜家台 闸下水位较高，糙率系数为0. 031时计算的杜家台闸下水位过程线与实测的水位过程线 对比误差较小，糙率系数为0. 035时的对比误差较大；当杜家台闸下水位较低，糙率系数 为0. 031时计算的杜家台闸下水位过程线与实测的水位过程线对比误差较大，糙率系数为 0. 035时的对比误差较小。洪道糙率系数随杜家台闸下水位涨高而减少的特点符合一般的 规律，可以采用。由于洪道的计算断面有28个，不便全部录入，因此，在表一中只列入了 4个计算断 面(即编号为1、26、27、28)，表一中反映的是假定糙率系数为0. 035情况下的部分计算结 果，节录计算的时间段为1983年10月13日的10时至21时，假定糙率系数为0. 031情况 下的计算结果也不便录入表一中。表一假定η = 0. 035时的洪水演进计算结果(节录)
权利要求
用槽蓄增量关系求解圣维方程组模拟洪水演进的方法，其特征是将圣维南基本微分方程组在Δt时段简化为水量平衡方程式和动力方程式后，根据水量平衡方程式(1)可计算河段槽蓄增量，根据时段初、末的动力方程式(2)、动力方程式(3)也可计算河段槽蓄增量。方程组如下 <mrow><mi>Δ</mi><msub> <mi>Z</mi> <mn>1</mn></msub><mo>=</mo><mfrac> <msup><msub> <mover><mi>Q</mi><mo>‾</mo> </mover> <mn>1</mn></msub><mn>2</mn> </msup> <msup><msub> <mover><mi>K</mi><mo>‾</mo> </mover> <mn>1</mn></msub><mn>2</mn> </msup></mfrac><mi>L</mi><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo></mrow> </mrow> <mrow><mi>Δ</mi><msub> <mi>Z</mi> <mn>2</mn></msub><mo>=</mo><mfrac> <msup><msub> <mover><mi>Q</mi><mo>‾</mo> </mover> <mn>2</mn></msub><mn>2</mn> </msup> <msup><msub> <mover><mi>K</mi><mo>‾</mo> </mover> <mn>2</mn></msub><mn>2</mn> </msup></mfrac><mi>L</mi><mo>-</mo><mo>-</mo><mo>-</mo><mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo></mrow> </mrow>以上各方程式中L为河段度长(m)；Δt为时段(s)；Z为水位(m)；Q为流量(m3/s)；W为河段槽蓄水量(m3)，ΔW＝(W2 W1)；ΔZ为河段上、下游计算断面的水位差(m)；为河段平均流量；为河段平均流量模数；各量脚注1、2分别为时段初、末的时刻。当已知条件为Q上1、Z上1、Q下1、Z下1、Q上2以及Q下2(或Z下2)时，可依据槽蓄增量ΔW的关系求解Z上2和Q下2(或Z下2)根据方程式(1)，如假定Q下2(或Z下2)值，就可以计算出相应的河段槽蓄增量ΔW(水量平衡)。根据方程式(2)、(3)，继续利用上述假定的Q下2(或Z下2)值，就可以计算出河段槽蓄增量ΔW(动力)。通过计算机程序按上述关系调节Q下2(或Z下2)值，使两种方式计算的河段槽蓄增量ΔW基本相等，此时的Z上2和Q下2(或Z下2)即为方程组的解。当计算断面个数为n时，则需要增加n 2个断面流量计算关系式求解方程组。如以上游起始断面的编号为i＝1，依次往下游的断面编号则为i＝2、i＝3…i＝n 2、i＝n 1、i＝n。当断面编号1＜i＜n时，根据水量平衡方程式(1)，如断面的流量和河段槽蓄增量用二维数组表示，则断面i 1与断面i之间、断面i与断面i+1之间的水量平衡关系式分别为Q(i 1，1)+Q(i 1，2) Q(i，1) Q(i，2)＝2ΔW(i 1，i)/Δt (4)Q(i，1)+Q(i，2) Q(i+1，1) Q(i+1，2)＝2ΔW(i，i+1)/Δt (5)联解方程式(4)、(5)Q(i，2)＝(ΔW(i，i+1) ΔW(i 1，i))/Δt+(Q(i 1，2)+Q(i+1，2))/2+(Q(i 1，1)+Q(i+1，1) 2Q(i，1))/2 (6)也可近似写成式(7)等形式Q(i，2)＝(ΔW(i，i+1) ΔW(i 1，i))/2/Δt+(Q(i 1，2)+Q(i+1，2))/2 (7)式中二维数组流量Q括号内的第一项为断面编号，第二项为时段Δt的初、末代号(初为1、末为2)；二维数组槽蓄增量ΔW括号内的第一项为河段上游断面编号，第二项为河段下游断面编号。因此，按水量平衡方程式计算的全河段累计槽蓄增量∑ΔW(水量平衡)为∑ΔW(水量平衡)＝Δt(Q(1，1)+Q(1，2) Q(n，1) Q(n，2))/2 (8)按动力方程式计算的全河段累计槽蓄增量∑ΔW(动力)为FSA00000086843000011.tif,FSA00000086843000014.tif,FSA00000086843000015.tif,FSA00000086843000021.tif
全文摘要
用槽蓄增量关系求解圣维南方程组模拟洪水演进的方法，是将圣维南基本微分方程组在Δt时段简化为水量平衡方程式和动力方程式后，根据水量平衡方程式和动力方程式分别计算河段槽增量一是根据水量平衡方程式可计算槽蓄增量，二是根据时段初、末的两个动力方程式也可计算槽蓄增量，由于两种方法计算的河段槽蓄增量应相等，因此本发明就为水力学方法求解圣维南方程组模拟洪水演进提供了一种新的、可靠的方法，在河道、渠道、分蓄洪区和水库合理确定防洪工程规模、洪水预报、科学调度洪水等方面发挥着重要作用。
文档编号E02B3/00GK101956382SQ20101016137
公开日2011年1月26日 申请日期2010年4月8日 优先权日2010年4月8日
发明者杜佐道, 杜曲 申请人:潜江市知识产权协会