大型水库三维水质模型研究

刘中峰,李 然,陈明千,黄 翔,曹玲玲

(四川大学水力学与山区河流开发保护国家重点实验室,四川成都 610065)

近年来,二滩、溪洛渡、白鹤滩等一大批高坝工程已建、在建或待建。这些工程在带来巨大经济效益和社会效益的同时,不可避免地对河流水质产生了一定的影响。因此,对大型水库的水质预测模拟具有重大的理论和工程意义。

大型水库的水质要素分布有着与天然河道不同的特点[1]。与天然河道相比,水库水深显著增大,流速明显减小,水体水环境条件发生改变,同时由于地形、流场等的三维效应,温度、水质要素可能出现垂向分层现象。目前,一、二维水质预测模型已经比较成熟,三维模型的研究也正逐渐开展。聂晶[2]根据质量守恒定律和Fick定律,推导出水库水体总磷三维迁移转化生态数学模型,并将其应用到密云水库,取得了较好的模拟结果;MAO等[3]建立了一种物理过程和生化过程耦合的三维模型,用来模拟太湖中营养物质的年变化情况;George等[4-5]针对 Lake Washington提出了一个复杂的三维富营养化模型,成功分析了浮游植物的动力学特征。

以上模型能够较好地模拟水体中水质要素的基本分布规律,但它们大多针对特定的研究对象,且都没有考虑温度分层的影响,而在大型水库中温度分层对水质要素的分布有较大的影响[6]。另外,上述模型都是非稳态模型,这在提高模型预测精度的同时,却很大程度上提高了计算代价,尤其是大型水库的计算区域非常大,采用非稳态模型计算所花费的时间非常长。在实际工程中,大型水库在一个水文周期(丰、平、枯)内,水文条件和污染负荷往往变化不大,采用稳态模型就能满足工程要求。因此,笔者建立了一种考虑温度分层影响的大型水库稳态三维水质预测模型,并用二滩鱼感鱼河支库的实测资料对模型进行了验证。将模型应用到某大型水库,成功模拟出水质要素的三维分布特性。

1 数学模型

1.1 数学模型方程

数学模型方程包括连续方程、动量方程、k方程、ε方程及污染物输移方程。

连续方程

动量方程

k方程

ε方程)

污染物输移方程

其中

式中:i为坐标方向,i=1,2,3;ρ为水的密度,kg/m3;Ui为i方向的水流速度值,m/s;P为压强,Pa;μ为动力黏性系数,m2/s;μt为紊动涡黏系数,m2/s;k为紊动动能,m2/s2;ε为紊动动能的耗散率,m2/s3;σk,σε和 σC分别为紊动动能、紊动动能的耗散率和水质要素的Prandtl数,对应值分别为1.0,1.3和1.0;φ为DO,CODCr,TP,TN等水质要素的质量浓度,mg/L;Cμ,C1ε和 C2ε为模型常数 ,取值由基本实验确定[7],分别为0.09,1.44和1.92;Sφ为污染物输移方程的源项,具体形式见第1.2节。

1.2 污染物输移方程源项表达式

污染物输移方程源项表达式分为2种形式。

a.DO方程源项。DO方程源项分表层水体和下层水体2部分区别考虑。对于表层水体,DO方程源项包括表面复氧和降解耗氧2部分;对于下层水体,DO方程源项只有降解耗氧,而溶解氧的补充依靠输移扩散过程。表层水体和下层水体的DO方程源项 Sφ表及Sφ下的表达式分别为

其中

式中:k1(T)为温度T时CODCr的降解系数,s-1;k2(T)为温度T时的复氧系数,s-1;CS为水体中饱和溶解氧质量浓度(计算公式来源于文献[8]),mg/L;CDO为水体DO质量浓度,mg/L;LC为水体CODCr质量浓度,mg/L;θ为温度修正系数;k2(20)为温度在20℃时的复氧系数(计算公式来源于文献[8]),s-1;U10为水面上10m处的风速,m/s;h为表层水体厚度,计算中取0.5m。

b.一般污染物输移方程源项。假定污染物随时间变化为一阶动力学反应,源项表达式为

其中

式中:k(T)为温度T时水质要素的降解(衰减)系数,s-1;k(20)为温度在20℃时相应水质要素的降解(衰减)系数,s-1。

1.3 模型求解

控制方程组的离散采用有限体积法。离散后的代数方程组采用SIMPLEC算法求解。

2 参数率定与模型验证

采用二滩鱼感鱼河支库原型观测资料对模型进行参数率定与验证。

2.1 二滩水库概况

二滩水库位于金沙江支流雅砻江上。最大坝高240m,回水长度约140km。鱼感鱼河支库是二滩库区右岸的主要支库,回水长度约30km。2008年7月对鱼感鱼河支库水温和DO分布进行了原型观测,同时收集整理了有关部门的鱼感鱼河支库水文资料、TP与TN定点取样监测资料及相关污染源资料,用于对三维水质模型进行验证。

2.2 计算区域及网格划分

二滩鱼感鱼河支库上游由永兴河、惠民河和新坪河交汇而成。计算区域上游至永兴河汇口,下游至鱼感鱼河入二滩主库的汇口,全长约30km。计算区域上表面水位与水库正常蓄水位一致,为1200m。计算区域垂向最大水深为150m,横向最大水面宽度为1165m。计算区域平面如图1所示。

图1 计算区域平面示意图

计算区域共划分163740个网格,其中横向最大网格间距为55m,纵向最大网格间距为55m,垂向最大网格间距为10m。网格划分情况如图2所示。

图2 鱼感鱼河支库网格划分示意图

2.3 污染源

鱼感鱼河支库的点源污染主要来自渔门镇的生活污水排放。该点源的污水排放量为24.97万t/a,污水中各污染物的质量浓度分别为:ρ(CODCr)=363mg/L,ρ(TP)=4.05mg/L,ρ(TN)=61.6mg/L 。

经调查,库区所在位置为山区,耕地稀少,工农业很不发达,面源污染主要来自库区内的网箱养鱼。收集了库区渔业产量、饲料用量等渔业养殖资料,并根据文献调研成果[9-10],计算得到鱼感鱼河支库面源污染物排放量。在计算中,面源污染物随沿程汇流(末断面流量与起始断面流量之差,认为从两岸沿程均匀汇入)从两岸进入水体。面源污染排放特征值见表1,其中质量浓度值为污染物排放量与沿程汇流量之比。

2.4 边界条件

入流边界包括永兴河、惠民河和新坪河3条支流入口、渔门镇生活污水排放口及两岸沿程汇流,设置为速度入口。在入流边界给定速度分布和水质要素浓度分布,紊动特征量(紊动动能 k和紊动能耗散率ε)根据入口流速 u0和入口水深 H0按下述经验公式确定[11]:

表1 鱼感鱼河支库面源污染排放特征值

出流边界为鱼感鱼河支库进入二滩库区的汇口断面,设置为充分发展的紊流边界,即所有标量(k,ε,φ)法向梯度为零,平面的法向速度分量梯度为零。在水面上给定对称面边界条件,即没有垂向的速度分量。水底为固壁边界,水动力学条件按标准壁函数法处理,并认为通过壁面的污染物通量为零。计算中,水温采用鱼感鱼河支库2008年7月的实测数据。

2.5 验证结果与分析

模型中各水质要素的 k(20)及θ初值参照已有的研究成果[12-13]取值,见表2。

表2 各水质要素k(20)及 θ取值

采用表2中模型参数初值,模拟得到二滩水库各水质要素的空间分布。表3给出了Ⅰ断面(新坪河汇口)、Ⅱ断面(新坪河汇口下游3000m)、Ⅲ断面(新坪河汇口下游6000m)自由水面中间点上TP和TN模拟结果与实测结果的比较。由表中数据可以看出,数值模拟结果与实测结果吻合较好。

表3 自由水面中间点TP和TN模拟结果与实测结果对比

图3给出了惠民河汇口断面和新坪河汇口断面中垂线水面下15m范围内DO模拟结果与实测结果的比较。从图3可以看出,模拟结果与实测结果总体趋势接近,呈现上大下小的趋势。计算可得2个断面中垂线DO质量浓度的计算值和测量值的均方根误差分别为0.78mg/L和0.31mg/L,表明两者吻合较好。

图3 断面中垂线DO质量浓度模拟结果与实测结果对比

综合分析认为,基于上述模型参数的预测结果与监测结果误差较小,能够满足实际要求,可以用于对TP,TN,CODCr等污染物的数值模拟。

3 模型误差与适用性分析

在二滩鱼感鱼河支库水质模拟中,基于对网箱养鱼规模的调查以及国内关于网箱养鱼污染的研究成果,进行了面源计算。由于目前关于网箱养鱼污染物排放特征的研究资料和水平有限,因此面源污染负荷水平的计算误差可能直接影响到水质模型预测结果的精度。

分析DO模拟结果与观测结果的差别,认为除本模型中考虑的污染物耗氧和表面复氧过程外,实际水库中还存在底泥耗氧、光合作用产氧等过程。这些过程的综合作用可能造成DO沿水深的非线性变化,而笔者建立的模型中忽略了这些作用,因此存在一定误差。另外,本模型没有考虑风力作用对水库垂向混合的影响,因此造成表面DO浓度较观测值偏低。模型仅考虑了温度分层对模型中降解系数、复氧系数的影响,没有考虑秋、冬季节由温度差引起的垂向混合作用,因此模型对秋、冬季节DO模拟误差可能较大。

对于水库中TP和TN等污染物的模拟,由于它们在水库中的垂向分层效果不明显,同时相对于DO来说,它们受表层复氧过程的影响较小,因此忽略风力混合作用和温度差产生的误差相对于DO来说要小得多。另外,水质取样中的随机性以及测点代表性也可能影响到对模型精度的判断。

4 模型应用

将上述模型应用于国内某在建大型水库的水质预测研究,成功模拟出该水库BOD,DO,CODCr等水质要素的三维分布情况。受篇幅所限,在此仅给出该水库某库湾段丰水期CODCr的模拟结果。该库湾长约11km,共划分190×31×27(纵向×横向×垂向)个网格。

图4给出了自由水面CODCr质量浓度分布的模拟结果。可以看出,受污水排放影响,石灰窑沟汇口出现污染带,污染带横向宽度最大达到70m,纵向延伸至库湾出口。污染带横向影响范围较小,这主要是由于沟口下游河道束窄而形成回流区,抑制了CODCr在横向的扩散。

图4 丰水期自由水面ρ(CODCr)分布示意图(单位:mg/L)

图5 给出了黑水河汇口上游断面CODCr质量浓度分布的模拟结果。从图5中可以看出,CODCr质量浓度分布在垂向没有出现明显的分层现象。但是由于石灰窑沟污水汇入的影响,右岸CODCr质量浓度比左岸高。

图5 丰水期黑水河汇口上游断面ρ(CODCr)分布示意图(单位:mg/L)

5 结论与展望

建立了大型水库三维水质预测模型,并在污染物输移方程中采用了考虑温度分层影响的源项表达式。采用二滩鱼感鱼河支库原型观测资料对模型进行验证,模拟结果和监测数据吻合较好,表明该模型是可靠的。将模型应用于某在建大型水库的水质预测,成功模拟出水质三维分布特性,模拟结果对于开展水环境影响的定量预测评价、制定水库水污染防治规划、保护库区水环境具有一定的理论意义和应用价值。

需要指出的是,本文工作是对三维水质预测模型研究的初步探索。模型在污染源特别是面源估算、水温分布引起的密度流影响、水库底泥及光合作用影响等方面尚存不足。此外,没有模拟水温,而是将水温数据作为边界条件输入模型。因此,在深入污染源研究基础上,综合考虑各种影响因素,建立水温、水质耦合求解的三维水质预测模型值得进一步研究和探索。

[1]LI Ran,FU Xiao-ying,LI Jia,et al.Three-dimensional numerical simulation of BOD-DO for large reservoir[J].Journal of Hydrodynamics:Ser B,2002,3:50-53.

[2]聂晶.水库水体总磷三维数学模型及其应用[D].长春:吉林大学,2004.

[3]MAO Jing-qiao,CHEN Qiu-wen,CHEN Yong-can.Threedimensional eutrophication model and application to taihu lake,China[J].Journal of Environmental Sciences,2008,20(3):278-284.

[4]ARHONDITSIS G B,BRETT M T.Eutrophication model for Lake Washington(USA)PartⅠ:model description and sensitivity analysis[J].Ecological Modelling,2005,187(2/3):140-178.

[5]AR HONDITSIS G B,BRETT M T.Eutrophication model for lake washington(USA)PartⅡ:model calibration and system dynamic analysis[J].Ecological Modelling,2005,187(2/3):179-200.

[6]王冠,韩龙喜,常文婷.基于立面二维水动力-水温耦合模型的水库水温分布[J].水资源保护,2009,25(2):59-63.

[7]金忠青.N-S方程的数值解和紊流模型[M].南京:河海大学出版社,1989:54.

[8]JEAN-MARC T.Simulation of a mesotrophic reservoir(Lake Pareloup)over a long period(1983-1998)using ASTER2000 biological model[J].Water Research,2004,38(2):393-403.

[9]陈丁,郑爱榕.网箱养殖的氮、磷和有机物的污染及估算[J].福建农业学报,2005,20(B12):57-62.

[10]陈德春.投饵式网箱养鱼对水质的影响及保护措施[J].水利渔业,1993(2):23-26.

[11]邓云.大型深水库的水温预测研究[D].成都:四川大学,2003.

[12]饶群.大型水体富营养化数学模拟的研究[D].南京:河海大学,2002.

[13]黄真理,李玉梁,李锦秀,等.三峡水库水质预测和环境容量计算[M].北京:中国水利水电出版社,2006:211-221.