一维水量水质模型

1、.第七章一维非恒定河流和河网水量水质模型对于中小型河流， 通常其宽度及水深相对于长度数量较小， 扩散质（污染物质、热量）很容易在垂向及横向上达到均匀混合， 即扩散质浓度在断面上基本达到均匀状态。 这种情况下， 我们只需要知道扩散质在断面内的平均分配状况， 就可以把握整个河道的扩散质空间分布特征， 这是我们可以采用一维圣维南方程描述河流水动力特征或水量特征（水位、流量、槽蓄量等） ；用一维纵向分散方程描述扩散质在时间及河流纵向上的变化状况。 特别地，对于稳态水流， 可以采用常规水动力学方法推算水位、 断面平均流速的沿程变化； 采用分段解析解法计算扩散质浓度沿纵向的变化特征。 但是，在非稳态情况下

2、 （水流随时间变化或扩散质源强随时间变化） 解析解法将无能为力 （水流非恒定） 或十分繁琐（水流稳态、源强非恒定），这时通常采用数值解法求解河道水量、水质的时间、空间分布。在模拟方法上， 无论是单一河道还是由众多单一河道构成的河网， 若采用空间一维手段求解， 描述水流、水质空间分布规律的控制方程是相同的， 只不过在具体求解方法上有所差异而已。7.1单一河道的控制方程7.1.1水量控制方程采用一维圣维南方程组描述水流的运动，基本控制方程为：QB WZq(1)xtQ2uQgAZu 2Agn2 u Q0(2)R4 / 3txxx式中 t 为时间坐标， x 为空间坐标， Q 为断面流量， Z 为断面平

3、均水位， u 为断面平均流速， n 为河段的糙率， A 为过流断面面积， BW 为水面宽度（包括主流宽度及仅起调蓄作用的附加宽度） ，R 为水力半径， q 为旁侧入流流量（单位河长上旁侧入流场）。此方程组属于二元一阶双曲型拟线性方程组，对于非恒定问题，现阶段尚无法直接求出其解析解， 通常用有限差分法或其它数学离散方法求其数值解。在水流稳态、 棱柱形河道条件下， 上述控制方程组退化为水力学的谢才公式，可采用相应的方法求解水流特征。7.1.2扩散质输运控制方程描述河道扩散物质运动及浓度变化规律的控制方程为： 带源的一维对流分散（弥散）方程，形式如下：( AC)(QC)xAExcKACA Sr S(

4、3)txxh式中， C 为污染物质的断面平均浓度，Q 为流量， E x 为纵向分散系数， S 为单位时间内、单位河长上的污染物质排放量， K 为污染物降解系数， Sr 为河床底泥释放污染物的速率。.此方程属于一元二阶偏微分方程， 对于非恒定水流问题， 微分方程位变系数的偏微分方程， 现阶段尚无法直接求出其解析解， 通常用有限差分法或其它数学离散方法求其数值解。 在水流稳态、 污染源源强恒定条件下， 可按水动力特征将河道分为若干子段，在每个分段上，上述控制方程简化为常系数的常微分方程，可采用解析方法秋初起理论解。7.2 单一河道一维水量水质模型7.2.1 单一河道一维水量模型（ 1）控制方程的离

5、散采用四点隐式差分格式离散方程组。如图1 所示，河道被 (n+ 1)个断面分为 n个子河段，在第 i 个子河段 M(i,i+ 1)上，对任一变量取：( M )jj1 ) / 2(4)( ii( M )j !j1(1)jji 1ii 1ixx(5)( M )j 1j 1jjii1ii1(6)t2t第i 个子河段2ii+1nn+1图 1 计算断面示意图式中，上角标表示时间坐标，下脚标表示空间坐标。为空间差商的权重系数 (01), =0 时，此格式为显式格式， 而当0 时，此格式具有隐式差分的特征。为使差分方程保持无条件稳定， 必须05. 。采用下式进行阻力项的线性化：n 2Q un2ujn2jQi

6、j 1uQij 1(7)ggR4 / 31 0.5R4 / 3R4 / 3i1i将式 (4)-(6)代入连续方程得第 i 个子河段的差分方程：C i Z iC iZ i 1Q iQ i 1D i(8)式中， C iBwi 1/ 2xi， D i(1) (QijQi 1j )C i ( Z ijZi j1 )qixi2t下角标 i+1/2 表示断面 i 与断面 i+1 河段的均值。按照同样的方法，可得动量方程的差分方程：Ei QiGi Qi 1 Fi Z i Fi Z i 1 H i(9).式中， EiFiGiH i对任一河段xigxin2 uj2u j1 / 22ti2R4 / 3igABu2

7、ji1/ 2jxijgxin 2 u2t2ui 1/ 22R4 / 31ixi Qi j 1/ 22 1uij1 / 2 (Qi jQi j 1 )1gA Bu 2jj 1 Z1j )i 1/ 2 ( Z iti (i1,2,n) ，可得到方程组：Ci Z iCi Z i 1QiQi 1D i(10)Ei QiGi Qi 1Fi Z iFi Zi 1H i对每一河段可列出两个线性代数方程， 再加上上下游边界条件， 构成完备的封闭方程组，采用追赶法可求得各个断面的水位流量。（ 2）边界条件根据上有下游边界条件类型的不同可以写成如下两种追赶形式：上游水位边界条件 Z1Z1* (t ) ；下游水位（

8、或流量）边界条件 Z n 1Z n 1* (t )（或Qn 1 Qn 1* (t ) ），追赶形式为：Z1P1R1Q1Q1L2M 2Q2Z 2P2R2 Q2QiLi 1M i 1 Qi 1Zi 1Pi 1(11)Ri 1Qi 1QnLn 1M n 1Qi 1Z n 1Pn 1Rn 1Qn 1Z n1Z n*1 (t ) (or , Qn 1 Qn*1 (t )式中， P, R, L, M 为已知系数，依据上述方程组，可逐步由下边界水位或者流量，推算得到上游各个断面水位流量值。上游流量边界条件 QQ *(t )；下游水位边界条件 Zn 1*(t ) ，追赶形式为：11Z n 1.Q1P1R1 Z

9、1Z1L2M 2 Z 2Q2P2R2 Z 2Z iLi 1M i 1 Zi 1QiPi 1(12)1Ri 1 Z i 1Z nLn 1M n 1Z i 1Qn 1Pn 1Rn 1 Z n 1Z n1Z n*1 (t)式中， P, R, L, M 为已知系数，依据上述方程组，可逐步由下边界水位，推算得到上游各个断面水位流量值。7.2.2单一河道一维水质模型（ 1）控制方程的离散与求解对方程 (3)进行离散，空间差分采用隐式迎风差分格式。顺流时 (从断面 i 流向 i+ 1) 有：( AC)( AC )i ( AC ) intt( AUC )(QC)i(QC) i 1xxi 1( AE xC1(

10、AE x)iCi 1CiCiCi 1x)xi( AE x )i 1xi 1xxi1KACSKi 1 Ai1/ 2 CiSi 1得到统一形式的差分方程：ai Cibi C i1ci Ci1di(13)式中， ai , bi ,ci , di 为系数，分别表示为：bi Qi 1/ 2 ,0( AEx)i1/ 2 txi 1aiViQi 1/ 2 ,0 t Qi 1/ 2( AEx)i 1 / 2( AEx )i1/ 2t K iVit,0 ttxixi 1ciQi 1/ 2 ,0( AEx )i 1/ 2 txi.d iVi n CinQi ,0 Sixit Qi ,0 Si 1xi 1 t方程（

11、 13）两边同时除以 ai 得到：CiEi Ci 1 Fi Ci 1 Gi(14)在顺流情况下，各河段差分方程可写成：C 2E2C1F2 C 3G 2C3E3C 2F3C 4G3CiEi Ci 1 Fi C i 1 Gi(15)C n1En 1C n 2Fn 1C n G n 1C nEn C n 1Fn C n 1 G n对首断面给定第一类边界条件， 对末断面给定第二类边界条件， 可得到如下封闭的方程组：C1C1* (t )C2E2C1F2 C3G2C3E3C2F3C4G3CiEi Ci 1 Fi Ci 1 Gi(16)Cn 1En 1Cn 2Fn 1Cn Gn 1CnEn Cn 1Fn C

12、n 1GnCn1En 1CnGn 1对方程组（16）采用追赶法可容易求得n, n1, n2,3,2 等断面的扩散质的浓度。（ 2）参数确定 纵向分散系数 EX 的确定EX与水流流速、水面宽度成正比，与水深成反比，常采用下面的经验公式：E xC02 q式中， C0c /g 是无尺度谢才系数， c 为谢才系数，=B/ h 为宽深比， q 为单宽流量，0.01 为经验常数。 降解系数 K 的确定.可采用监测资料对降解系数进行率定，或根据经验得到。7.2.2应用实例三峡大坝位于宜昌县三斗坪中堡岛 , 葛洲坝位于南津关下游的宜昌市境内 . 两坝间水域处于鄂西山区向平原的过渡地带 , 周围地形地貌呈西北高

13、东南低之势 . 三峡大坝至葛洲坝河段长 38Km,两坝间江面宽 210m至 1500m,大部分处于西陵峡谷中。该实例建立了两坝间水量水质模型 , 分别运用大坝一期围堰及二期围堰施工期间的同步水文水质实测资料对模型进行了率定和验证, 取得了较好效果 . 此模型可用以预测大坝施工期间及投入运行后两坝间水流及水质特性的变化.（）水量模型率定利用 1996 年三斗坪、白庙子及黄陵庙等断面水文观测资料率定糙率，得到各子河段的糙率， 率定结果显示糙率分布取值范围为 0.036 至 0.050 。1996 年白庙子及黄陵庙等断面实测水位过程线及计算水位过程线图略。（）水量模型验证采用太平溪断面 1998 年

14、实测流量作模型验证的上边界条件，采用葛洲坝坝前断面 1998 年实测水位作模型验证的下边界条件， 对 1998 年两坝间水流进行模拟，通过对各水文观测断面的水文要素的观测值及计算值进行比较可见，吻合程度较好，因此此水量模型可用来模拟两坝间的一维水流情况，1998 年白庙子及黄陵庙等断面实测水位过程线及计算水位过程线见图。H(m)三 斗 坪 断 面H(m)黄陵 庙 断 面74.0072.0072.0070.0070.0068.0068.0066.0066.0064.0064.00060120180240300360060120180240300360T(daY)T(daY)图 2 模型验证各断面

15、实测及计算水位过程线( )CODMn模型验证高锰酸盐指数的降解系数由实测资料，取经验值，不采用模型率定。采用1998 年 712 月的 CODMn 浓度监测资料对 CODMn 模型参数可靠性进行验证，模拟时 CODMn 污染源分点源和面源两种情况。 部份断面 CODMn 水质因子实测及计算浓度值见图 3。.三斗 坪 断 面COD(mg/L)白 庙 子 断 面COD(mg/L)3.00实 测 值3.00计 算 值实 测 值计 算 值2.002.001.001.0003060901201501800306090120150180T(day)T(day)COD(mg/L)黄 陵 庙 断面COD(mg

16、/L)南 津关断 面3.00实 测 值3.00实 测 值计 算 值计 算 值2.002.001.001.0003060901201501800306090120150180T(day)T(day)图 3 模型验证 CODMn 浓度实测及计算值7.3平原河网水量水质模型河网地区是中国社会、经济、文化的发达地区， 在国家经济发展中占有举足轻重的地位。随着地区经济的进一步发展、居民物质精神生活水准的进一步提高，水资源问题日益突出。如中国著名的长江三角洲、珠江三角洲地区，已成为水质型缺水区，水资源、水环境问题已成为制约经济社会发展、 事关区域可持续发展的重要因素。生产实践的迫切需要使得河网地区的水环境

17、保护研究呈不断深入、系统之势。 作为基础， 河网水流水质模拟方法是进行区域环境规划、 环境管理等的必备工具，在水环境问题研究中占有重要地位。因此，该领域的研究一直是环境科学研究人员、 环境管理决策部门十分关注的重要问题。在物理模型、 数学模型这两大手段中，因为数学模型具有经济、快捷、实用等优点，加之水网地区河道密布等客观条件的限制， 现阶段只能采用数值方法模拟水网地区的水流运动及污染物输运规律。但由于理论、技术及各种客观条件的制约， 目前水流特别是水质数值模拟精度不是十分令人满意。7.3.1 平原河网水量模型按河网水流的控制方程及对河网的概化处理方式不同，河网地区水流数值模拟方法可分为两大类：

18、 第一类为常用的一维圣维南方程组数值解法；第二类为所谓的 “组合单元解法” 。，其中， 一维圣维南方程组数值解法又可分为直接解法和间接解法两种。在直接解法中， 较有代表性的有文献 G.Noseda, Mathematical Model of Unsteady Flow in Open Channels Networks, Proceedings of the International Symposium on Unsteady Flow in Open Channels,1976 ，J.J.R. Williams, T.R.E. Chidley, Nonliear Analysis of

19、Unsteady Flow in Open提出的方法， 该类方Channel Networks, Proceedings of the International Symposium on Unsteady Flow in Open Channels,1976法将计算断面交替取为水位和流量断面，对河网在所有计算断面上统一对一维圣维南方程组差分离散并求解。 但该方法未知数数量较多， 在河网规模较大的情况下，因为河道的交叉衔接，形成的矩阵是一个不规则、不对称的大型稀疏矩阵。为减少存贮，J.J.Dronkers 于 1976年提出间接解法 5 的思想，以后又有许多学者对其作了进一步完善。为提高计算效

20、率，需缩小矩阵规模， 中山大学数力系 1976年提出了河网非恒定流隐式方程组稀疏矩阵解法6 李岳生.等，河网不恒定流隐式方程组稀疏矩阵解法，中山大学学报（自然科学版），1979 年 3 月 ，该方法从河网矩阵本身的特点出发，能够有效地节省存储并提高计算速度，但矩阵中需包含所有断面的未知数，方程规模仍然较大，其实际使用也受到限制。 而间接解法是将断面未知数往交汊点集中，待求出交汊点未知数后，再求解各单一河道未知变量，计算效率较高。间接解法的思想首先由荷兰水力学专家Dronkers 于 1976 年提出，以后又有许多学者对此方法进行进一步的完善，相继提出了河网非恒定流的二级解法徐正凡，明渠非恒定流

21、 M. 武汉水利电力学院， 1983 年. 、三级解法 8张二骏，河网非恒定流的三级联合解法，华东水利学院学报， 1982 年第 1 卷： 112和四级解法 9 吴寿红，河网非恒定流四级解法，水利学报，1985 年，第8 期： 4250 。此类间接解法中， 以三级解法最为常用 。其基本求解思路可概括为“单一河道交汊点单一河道” ，即（ 1）先将单一河道划分为若干子河段，在计算断面上对一维Saint-Venant 方程组进行有限差分运算，得到各单一河道差分方程组，进行消元计算，得到单一河道首、末断面间流量与水位的相互关系；（ 2）根据河道交汊点水量守恒方程，得到并求解交汊点水位方程组，得到所有交

22、汊点水位；（ 3）根据交汊点水位返回各单一河道，求得各计算断面水位、流量值。“组合单元解法” 10 由法国水力学专家Jean A. Cunge 于 1975 年首次提出， 国内也有研究者采用此方法进行了水网地区的水力模拟1112 。此方法的基本思想是：将河网地区水力特性相似、 水位变幅不大的水体概化成单元。取单元中心的水位为代表水位，采用谢才公式模拟单元间流量交换， 根据水量守恒建立每一单元的微分形式的水量守恒方程，离散并得到以单元水位为自变量的代数方程，辅以边界条件，可求得各单元水位、单元间流量。在上述两类方法中， 组合单元解法对河道进行了简单概化，以单元为计算单位，计算相对简单，但模拟精度

23、相对较低，仅实用于大尺度水域的水力模拟。而Saint-Venat 方程组数值解法可以精确计算每一条河道的水流状况，所以成为目前河网水力模拟的主流方法，其中又以三级解法最为常用 13141516 。 当采用三级解法模拟长系列大范围河网水流特性时，需对河道进行概化处理，即将等级较小的河道概化为一条“概化”河道，要求该概化河道与被概化河道的过流能力、调蓄量相当。 同时，为考虑降雨对河网水力特性的影响，还必须对河道包围的陆域面积进行产汇流计算， 将产流量以包围陆域的河道长度为权重分配到周围河道17 。由于三级解法以单一河道为模拟对象，同时考虑了降雨等的影响，计算精度较高，可以满足水利、航运、环境保护等

24、的需要。（ 1）控制方程的离散与求解采用四点隐式差分格式离散方程，得子河段差分方程：C i Z iC i Z i 1Q iQ i 1D i(17)Ei QiGi Qi 1Fi ZiFi Zi 1H i(18)式中 Ci , Di , Ei , Fi , Gi , H i 由时段初值及河道特征求得。（ 2）节点连接条件水流运动在河网各节点上应满足质量守恒及能量守恒，即满足以下两个连接条件：质量守恒条件进出某一节点的流量与该节点内水量蓄量的增减相平衡，定量表示为：mj 1jQi (19)KK ) / ti 1式中 K 为节点编号， m 为流入 (流出 ) 节点 K 的河道数，K 为节点蓄量， Q

25、为流量 (流入为.正；流出为负 )。若节点汇合区容积与子河段容积相比可忽略不计，则此节点称为无调节节m点，方程 (19)可简化为：Q i0i1能量守恒条件不计节点处能量损失有：Z iZ ji =1,2, m,j 1,2, m(20)(3) 界条件有三种 型的 界条件：水位 界条件即在 界河道上 定水位随 的 化 程：ZZ ( t ) ；流量 界条件即在 界河道上 定流量随 的 化 程：QQ ( t )；水位流量关系当 界河道上有水工建筑物(如水 、堰、堤 等 ) ，通常 定水位流量关系： QQ ( Z ) 。（ 4） 方程的求解利用消元法，方程(3) 、(4) 经递推运算，写成如下形式：Q i

26、ii Z ii Z n1(21)Q i1i 1i 1 Z i 1i 1 Z1(22)式中， i , i ,i , i ,i ,i 为追赶系数，可逐步递推求得。利用公式(7)依次由末断面向首断面递推，可将各断面流量表为该断面水位及末断面水位的函数。利用公式(22)依次由首断面向末断面递推，可将各断面流量表为该断面水位及首断面水位的函数。特殊的， 单一河道首末断面分别与节点相连，分别对应如下追赶方程：Q 111 Z 11 Z n1(23)Q n 1n1n1 Z n 1n 1 Z 1(24)将 点各支流相 的 (23)或 (24)代入公式 (19)，并将 点各相 断面水位 一表示成 点水位，得 点方

27、程：f i ( Z i , Z 邻 1, Z 邻 2 ,Z 邻 KI )0i=1,2,M式中， M 点 数， KI 第 i 个 点相 河道数。当某 一河道 界河道 ，分三种 型分 有： 于第一 界条件， 式(24)中 Z1 已知； 于第二 界条件，式(23)中 Q1 已知； 于第三 界条件，已知 性关系 Q1 =Q( Z1 )，代入 (23)消去 Q1 。无 何种 界条件，最 可增加一.水位方程，封闭节点方程组f i ( Z i , Z 邻 1 ,Z 邻 2 ,Z 邻 KI )0 ，可解得各节点水位，再返回单一河道方程，最终求得各断面水位及流量。7.3.2平原河网水质模型(1) 单一河道控制方

28、程的离散采用如下网格对方程 (3)进行离散，空间差分采用隐式迎风差分格式，由于河网地区水流流向可能谁随时间变化， 根据流向的不同， 分顺流及逆流两种状况给出差分形式：Ci 1CiCi 1123 i 1ii 1 n 1 nxi1xi图 4 水质模拟网格布置示意图顺流时 (从断面 i 流向 i+ 1)有：( AC )( AC )i( AC) intt( AUC)(QC)i( QC) i1xxi1( AE xC )1( AE x ) iCi 1Ci( AE x ) i 1CiCi 1 xxxi1xixi 1KACSKi1 Ai1/2 CiSi 1逆流时 (从断面 i+ 1 流向 i)有：( AC )

29、( AC )i( AC) intt( AUC)(QC)i1(QC) ixxi( AE xC1( AE x)iCi 1CiCiCi 1)xixi( AE x )i 1xi 1xxKACSKi Ai1/ 2CiSi为考虑流向顺逆变化的影响，引入流向调节因子rc 及 rd ，对于断面 i 表示为：.Qw(QiQi 1 ) / 2Qe(QiQi 1 ) / 2rcw(QwQw ) / 2Qwrce(QeQe ) / 2Qerd w(QwQw ) / 2Qwrd e(QeQe ) / 2Qercrd 0,（当 Qw , Qe0 ）得到任意流向下， 一形式的差分方程：i Ci 1i Cii Ci 1Zi式

30、中， i ,i , i, Zi 系数， 于一般断面 (i=2,3,n-1)分 表示 ：i rcw D wwrd w D pwFcw t / Vi rcw Dwwrd w D pwrce D p prd e De pFcp Fd p t / Vrcw Ki 1rd e Ki t10.i rce D p prd e De pFd e t / Vn11Z iCircw Si 1 Ai 1/ 2rd e SiAi 1/ 2 t 于首断面 (i =1)，逆流 有：1 0t1 rd e De pFd p rd e Ki trd exi Ai1/2 rd ed ept1Fd e xi Ai1/ 2Z1C1n

31、trd eS1Ai 1/ 2 于末断面 (i =n)， 流 有：tn rcw D wwFcw An 1/ 2rcw D wwtrcw Ki 1 trcwnFc p xi 1 Ai1/ 2n0ZnCnnrcw Sn 1tAi 1/ 2式中，Ai 1/ 2( Ai 1Ai ) / 2 ,Ai 1/ 2( Ai 1Ai ) / 2.(25).V1Ai 1/ 2 xi 1 ,V2Ai 1/ 2 xiVrcw V1rd eV2Dww( AE x ) i 1 / xi 1 ,D p p( AE x ) i / xiD pw( AE x ) i / xi 1 ,De p( AE x )i 1 / xiFcw

32、(Qi 1Qi 1 ) / 2 ,Fc p(QiQi ) / 2Fd p(QiQi 1 ) / 2 ,Fd e(Qi 1Qi 1 ) / 2将差分方程 (25)用向量形式表示为：A c C cB c(26)式中 Ac 为一 (n 1)n 阶的三对角的系数矩阵， C c 为此河道断面平均浓度的 n维列向量， Bc 为已知的 (n-1)维列向量。 (26)式中含 (n-1)个方程， n 个未知数，方程组不闭合，需引入节点方程及边界条件。（ 2）河网节点方程对于污染物充分混合的节点，从时刻 n t 到 (n 1) t ，根据节点有无调蓄作用，分别给出相应的节点方程。若节点有 m 条单一河道，其中流入

33、节点的河道m1 条，流出节点的河道 m2 条，m1 及 m2 随流场变化而变化 ( m1 + m2 =m)。若节点本身具有调蓄作用，有：C out , iC N， i= m11, m1 2 , m(27)m1mQin , i Cin , iQout , i Cout ,ii 1i m1 1(28)1N CNNj CNj e K tSNt式中， CN 为节点的浓度， Cout , i 为流出节点的第i 条河道与该节点相邻断面的污染物浓度， Ci n , i 为流入节点的第i 条河道与该节点相邻断面的污染物浓度，Q 分别为与之相应的流量，N为节点的蓄量， S为节点的污染源加入项，方程 (28)N右端括号中末项为时段初节点污染物降解至时段末的残留量。 (27)式称为充分混合假定， (28)式反映节点内污染物的质量守恒。当节点调蓄可忽略时，称为无调蓄节点，有：.N 0（ 3）模型的动态耦合数值计算在流场已知时，由 (26)-(28)式加上边界条件及初始条件，构成封闭的线性代数方程组，即可求解。但对于大型复杂河网来说，其规模巨大，较为经济的做法是：将方程 (26)作适当的运算，得单一河道首末断面浓度关系方程，并与充分混合假定 (27