环境科学  2014, Vol.35 Issue (2): 555-561   PDF    
引用本文
贾海峰, 郭羽. 基于双向算法的湖库允许纳污负荷量计算及案例[J]. 环境科学, 2014, (2): 555-561.  
JIA Hai-feng, GUO Yu. Calculation of Allowable Pollution Loads for Lake and Reservoir Based on Bi-direction Algorithm and Its Case Study[J]. Environmental Science, 2014, (2): 555-561.  
基于双向算法的湖库允许纳污负荷量计算及案例
贾海峰, 郭羽    
清华大学环境学院, 北京 100084
收稿日期: 2013-6-6; 修订日期: 2013-8-27.
基金项目： 国家水体污染控制与治理科技重大专项（2009ZX07526-005）
作者简介：贾海峰(1967~),男,博士,博士生导师,主要研究方向为城市环境系统分析, E-mail:jhf@tsinghua.edu.cn
摘要：为了加强湖库流域水环境保护与污染负荷排放管理，支持湖库控制单元允许纳污负荷量的计算和污染负荷分配决策的制定，在总结国内外环境容量相关研究的基础上，针对现状水环境总量控制技术中允许纳污负荷计算与分配环节存在的问题，提出基于双向算法的湖库允许纳污负荷量计算方法. 并以EFDC（environmental fluid dynamic code）水环境模型为核心，通过模型概化、参数率定等建模技术与方法的研究，利用估算与精算相结合的方法制定排污情景，并应用模型完成情景分析. 应用多次情景试算、分析的方法，制定优化的排污情景，从而计算出湖库允许纳污负荷量. 最后以辽宁省柴河水库为实证区进行方法的应用，证明了方法具备科学性与计算精度.
关键词： 湖库控制单元     允许纳污量     双向算法     EDFC模型     柴河水库    
Calculation of Allowable Pollution Loads for Lake and Reservoir Based on Bi-direction Algorithm and Its Case Study
JIA Hai-feng, GUO Yu    
School of Environment, Tsinghua University, Beijing 100084, China
Abstract: Based on the reviews and summaries of water environment carrying capacity researches and practices, the main problems of allowable pollutants load estimation and its allocation in China were analyzed. Then a bi-direction algorithm for allowable pollutants loads calculation was proposed to support the pollutants loads management in the lake and reservoir control units. It was the combination of forward algorithm and backward algorithm. The two major steps were modeling and scenario analysis. Firstly, the basic scenario was proposed using the estimation model. Then the basic scenario was analyzed using the water quality simulation model to assess the compliance of water quality objectives. The allowable pollutant loads were calculated after several loops of scenario simulation, result analysis and scenario optimization. Finally, the Chaihe Reservoir in Liaoning Province, China was used as a case study using Environmental Fluid Dynamic Code (EFDC) model as the kernel model. The results demonstrated that the algorithm proposed provided an efficient and appropriate methodology for allowable pollutant load calculation.
Key words: lake and reservoir control unit     allowable pollution loads     bi-direction algorithm     environmental fluid dynamic code (EFDC)     Chaihe Reservoir    

在经济发展与环境保护并重的今天，进行湖库流域管理时，既要允许湖库接纳一定量的污染物负荷，又要保证湖库水体的环境质量，维持湖、 库的环境功能. 在此目标下，各类湖库水环境模型以及基于水环境容量的环境管理方法得到广泛的研究与应用[1, 2, 3]. 随着对流域水环境管理要求的逐步提高，人们已经不满足于仅使用零维或一维水环境模型进行环境容量的估算，希望利用机理模型，根据不同污染源排放的时空特征以及水体水文、 水质时空变化特征，刻画出以保护湖库水体环境保护敏感区环境目标为前提的允许纳污量，并分配到相应的污染源[4, 5, 6].

EFDC模型是美国EPA推荐使用的湖库三维水动力-水质模型[7]. 它是一种在设定了上游来水水文、 水质条件以及排污口位置和排污特征状况下，求解水体控制断面水质指标的正向算法模型. 而目前常用的水环境容量计算模型，则是基于受纳水体的结构、 水文特征以及水质保护目标，反推水体的环境容量和允许排放量[3]. 常用的水环境容量模型属于静态模型，模拟因子少，优点是简单实用，适合用于数据情况较差，水环境容量计算精度要求不高的湖库水环境管理中. 而利用EFDC模型构建的湖库水环境容量模型是一个考虑多水质因子的动态水环境模型，它对水质因子间的交互关系以及水质的时空变化刻画的较为精确，但对数据的要求高，需要通过情景分析法等方法，寻找达标情景. 本研究将集成这两类方法，以EFDC模型为核心，建立基于双向算法的湖库允许纳污负荷量计算方法，并以柴河水库为例进行实证分析.

1 基于双向算法的湖库允许纳污负荷量计算方法
1.1 技术路线

针对现状水环境总量控制技术中允许纳污负荷计算与分配环节存在的问题，借鉴美国TMDL技术[8， 9]，提出基于双向算法的湖库允许纳污负荷量计算方法，其技术路线可分为两条主线：①允许纳污负荷计算研究； ②污染负荷分配研究. 在传统的容量总量控制中，允许纳污量计算与污染负荷分配是相互独立的，而在本研究中，将两条主线进行整合，让容量计算、 负荷分配、 模型模拟工作有机结合，互相参考、 互为支持. 技术路线如图 1所示.

图 1基于双向算法的水体允许纳污负荷量计算技术路线Fig.1Roadmap of allowable pollution loads calculation based on bi-direction algorithm

图 1中反算法是根据控制断面或控制区域的水质目标，反推水体纳污总量的估算方法. 通过反算法可以得到允许纳污负荷总量的初算结果，以此为基础，结合污染负荷的分配方法，制定不同的负荷分配情景，然后通过正算法进行分配方案的达标校核与其对水质影响的分析. 在负荷分配方案优化中，利用模型对污染源贡献度的分析，计算污染源权重，并结合分配方案的环境影响程度、 经济因素、 技术可行性因素等，对方案进行筛选和优化.

根据此技术路线图，具体的技术环节可以分为：①水环境问题分析与水质目标制定； ②允许纳污负荷总量估算； ③负荷分配方法与技术； ④水动力-水质模型对情景方案的模拟； ⑤模型不确定性分析与安全余量计算； ⑥污染源-水质响应关系分析与分配方法的优化.

1.2 基于EFDC模型的湖库动态水质模型构建

包括模型建立和情景制定分析两个环节，如图2所示.

1.2.1 模型概化与数据输入

(1)空间概化 利用EFDC构建湖库三维水环境模型，首先要进行空间概化. 基本过程是将湖库水体划分为模型计算单元网格的组合. 网格划分后，按照模型的数据需求，将网格坐标、 网格间连接方式、 网格内的初始水动力条件、 初始水质指标等进行整理，作为模型输入.

(2)时间概化 时间步长的划分关系到模型求解的时间与运算的稳定性. EFDC需要输入多种时间序列的边界条件(例如：入流、 出流、 蒸发、 降水、 气温、 负荷输入等)，由于数据条件限制，一般需要按照时间步长要求，将已有数据进行平滑、 插值等处理，完成时间序列相关的数据与参数的输入.

1.2.2 参数率定和模型验证

EFDC模型常规水质指标模拟中，相关反应主要包括4类：藻类生长动力学、 氮反应动力学、 磷反应动力学、 溶解氧-有机碳相关反应. 基于4类动力学反应过程，EFDC模型可模拟22种水质指标在水体中的相互转换与影响. EFDC模型反应机理较为复杂，控制各种反应的水质参数众多，因此，模型的参数率定是保证模拟结果精度的重要技术环节.

图 2基于EFDC的湖库动态水质模型建立Fig.2Setup process of the EFDC dynamic water quality model

如果率定结果满意，再利用另外一套与用于参数率定的数据独立实际监测数据进行模型验证. 如果验证不能满足要求，则重新进行参数率定，直到验证结果满足要求.

1.3 情景分析法计算允许纳污负荷量

EFDC模型是在已知排污口及其排放特征的情况，求解控制断面水质的正向算法模型. 可以使用情景分析法，进行允许纳污负荷计算. 情景分析法通过设定排污情景、 模拟排污情景下控制断面水质、 分析模拟结果并计算排污口对控制断面水质的贡献度、 调整排污情景等步骤，多次循环寻找优化排污方式以及此种排污方式对应的湖库允许纳污负荷量.

1.3.1 估算模型与初始情景设置

合理的设置初始排污情景，有利于缩短利用情景分析法寻找优化排放的时间. 可以利用传统水环境容量模型作为估算模型，利用其由控制断面水质推算排污口允许排放量的特点，粗略的估算允许排放量，并将此计算结果作为EFDC模型情景分析中的初始情景.

1.3.2 EFDC模型与情景分析

情景设定为模型确定了输入条件，通过对水质模拟结果进行分析，判断设定的排污情景是否达到水环境管理的标准，并通过多次模拟寻找排污口与控制断面之间的水质响应关系，通过线性拟合、 多元回归等方法，建立关系函数，求解单个排污口对控制断面水质的贡献度，有利于指导更合理的情景方案的设定.

2 柴河水库实证研究

选择柴河水库作为研究的实证区，利用EFDC模型对柴河水库进行三维建模，并设置多种水文条件，计算不同条件下柴河水库最主要的入库河流——柴河在不同月份允许带入水库的最大污染负荷量，通过计算结果，对柴河水库流域污染源排放进行管理.

2.1 研究区概况

柴河水库是辽河一级支流柴河上建设的一座大型水库，位于辽宁省铁岭市东南，是铁岭市主要饮用水源地之一. 水库总库容为6.36亿m3，控制流域面积1355 km2. 柴河水库流域以非点源污染为主，重点污染期出现在丰水期，超标的主要水质指标为氨氮和总磷[11].

根据管理需求，模型模拟的主要水质指标为氨氮和总磷. 主要控制断面为柴河水库坝前断面，以柴河水库的水环境功能分区要求(《地表水环境质量标准(GB 3838-2002)》Ⅱ类标准)为控制断面水质目标，并对污染物混合区(污染物入库后允许超标的区域)面积、 藻类浓度等进行考虑. 具体控制目标如下.

(1)控制断面水质目标. 控制断面在全年范围内，月均值达到Ⅱ类水标，氨氮为0.5mg ·L-1，总磷为0.025mg ·L-1.

(2)水库在控制断面达到峰值浓度时，要求全库混合区面积小于30%，即至少有70%网格的模拟结果在全年90%的时间内的污染物峰值浓度低于0.5 mg ·L-1.

2.2 模型概化
2.2.1 模型空间概化

将柴河水库多年平均蓄水水位下区域概化为模型模拟区域，对此模拟区域进行三维空间网格划分，在垂直方向上将水库分为上下两层，每层由480个网格组成，每个网格为200 m×200 m. 模拟区域总面积为19.2 km2.

2.2.2 模型时间概化

为保证模型的运算的稳定性，并考虑模型的计算效率，确定本次水动力模拟的基本时间单元为1 d，计算时间步长为1 min. 水质模拟的时间步长为水动力模拟的2倍.

2.3 模型的参数研究
2.3.1 参数灵敏度分析

根据模拟对象的设置，需要确定EFDC的34个水质反应参数. 首先对参数进行灵敏度分析，选取需要精确率定的敏感参数[12]. 下面以氮反应动力学过程中的主要参数为例，进行参数灵敏度分析，表1为分析计算的结果.

表1 参数灵敏度分析结果 Table 1 Parameters sensitive analysis results

表1中显示，在柴河水库EFDC模型对氨氮的模拟中，最大硝化速率系数、 水温大于最适硝化作用温度时的温度系数、 有机氮矿化速率系数、 有机氮相间分配系数是相对敏感的参数. 对于其他反应参数也可以使用类似的灵敏度分析方法.

2.3.2 参数率定与模型验证

参考湖库水质模型的参数研究成果[13, 14, 15, 16]，使用2005年柴河水库入库口、 库中、 坝前的实际监测数据，分别完成了对EFDC模型的水动力、 水质模拟的参数率定； 使用2009年同样上述监测点的实际监测数据进行模型验证. 以坝前监测点为例，参数率定与模型验证结果如图 3和图 4所示.

图 3EFDC水力学模拟参数率定和模型验证结果Fig.3Calibration and verification result of hydrodynamic simulation

图 4水质参数率定与模型验证水质指标图

Fig.4Calibration and verification result of water quality simulation

采用中值误差对模型参数率定与模型验证的模拟结果进行误差分析，如表2所示. 可见模型可以从宏观和趋势上反映柴河的精度，可以支持柴河水库的水环境管理和允许纳污量计算.

表2 模拟结果中值误差分析表 /% Table 2 Median error of simulation/%

2.4 允许纳污量计算
2.4.1 初始情景设定与分析

使用式(1)所示的简单营养盐容量模型进行了氨氮与总磷的允许纳污负荷量的估算. 

式中，W为水库典型污染因子的环境容量(t ·a-1)； P为水质目标(mg ·L-1)； σ为典型污染物的沉降速率系数(m ·a-1)； Z为湖库水面面积(km2)； H为湖库平均深度(m)； ξ为湖库的冲刷系数(1 ·a-1)； V


为湖库的蓄水量(万m3)； Qout为水库的出流量(m3 ·s-1)； 各系数为单位转换系数.

在初始情景设置中，参考估算模型计算的年柴河入库口允许纳污负荷量，按照入库水量等比例将其分配到各月.

在模型水动力模拟上，选取2005年柴河水库水动力条件作为模拟对象，将2005年8月作为柴河水库典型丰水月，在此条件上进行污染排放控制，可以达到控制最不利条件，保证其他水文条件下水质达标得目的. 将柴河水库2005年出、 入流量数据，蒸发、 降水、 渗漏、 取水工程取水量等数据，估算模型得出的允许排污量数据等，输入EFDC模型进行模拟. 模拟结果如图 5所示.

图 5初始情景模拟结果

Fig.5 Simulation results of basic scenario

可以看出在初始情景的计算结果中，污染峰值浓度较为接近水质标准，具备一定的估算价值. 不过高估了氨氮的降解速度，低估了总磷的沉降速率. 因此，需要在初始情景基础上进行调整.

2.4.2 排污口-控制断面水质响应关系

分别调整氨氮与总磷的入库浓度进行调整，设置情景后，输入EFDC模型进行模拟，模拟结果如表3所示.

表3 排污口-控制断面水质响应关系 /mg ·L-1Table 3 Response relationship of pollution load and water quality at the control point/mg ·L-1

通过多次模拟，分析污染物排放口与控制断面之间的水质响应关系，可以看出排放口污染物浓度与控制断面水质之间具有较好的线性关系，可以用线性模型进行拟合. 基于此线性关系，可以求解在典型丰水月条件下的柴河入库口的污染物允许入库负荷量.

2.4.3 混合区面积校核

根据控制目标的要求，水库在控制断面达到峰值浓度时，要求全库混合区面积小于30%，即至少有70%网格的模拟结果在全年90%的时间内的污染物峰值浓度低于0.5mg ·L-1. 对于两个主要控制指标氨氮与总磷，氨氮在柴河水库的控制难度较大，因此在进行面积校核时，选取氨氮模拟结果进行.

根据入库污染源-控制断面水质相应关系分析，以及其超标区域面积分析，在氨氮入流浓度为0.5 mg ·L-1，库中区氨氮超标面积仍较大，超标区内氨氮峰值浓度接近0.8 mg ·L-1. 其主要原因是入流总氮浓度较高，总氮沉降后底泥释放出氨氮，以及溶解态有机氮转化为氨氮. 因此，需要通过限制总氮入流浓度，以满足氨氮超标区域面积限制的控制目标.

通过模型多次试算，计算总氮入流浓度限制为1.37 mg ·L-1. 在此情景下，氨氮超标区域大于30%的天数被控制在36 d之内.

2.4.4 模型不确定性与最终允许入库负荷量

安全余量(margin of safety，MOS)的概念在美国TMDL计划中被提出，是指在水环境管理中由于计算污染负荷与受纳水体之间关系时存在的不确定性，而需要在水环境管理过程中，在模型计算结果的基础上，留出余量，保证水质达标率[17, 18, 19].

在柴河水库EFDC模型进行允许纳污负荷计算过程中，使用FOEA(First-Order Error Analysis)法对模型进行不确定性分析[20]，并计算安全余量.

FOEA法中，将模型的参数不确定性简化假设为单个参数不确定性的线性加和，通过泰勒多项式对模型在参数取值点进行展开，通过求解单个参数灵敏度以及参数经验取值范围，而求解模型的整体参数不确定性.
FOEA法求解的安全余量结果如图 6、 表4所示，这样就可计算出最终的允许入库负荷量.

图 6氨氮的安全余量与允许入库负荷量

Fig.6MOS and allowable pollution loads of NH3

图 4典型丰水月水文条件下允许入库负荷量及安全余量计算结果 /t ·month-1 Table 4 Calculation results of allowable pollution loads and MOS in typical wet month/t ·month-1
3 结论

针对我国在水环境总量控制中允许纳污负荷计 算与分配环节存在的问题，借鉴以美国TMDL技术为主的国外先进水环境管理技术，提出基于双向算法的湖库允许纳污负荷量计算方法，让容量计算、 负荷分配、 模型模拟工作有机结合，互相参考、 互为支持，支持水环境的水质目标管理. 并基于EFDC机制模型的特征，建立已EFDC模型为核心的湖库允许纳污负荷量计算模型，可通过“初步情景设置-计算模拟-情景分析与修改-再模拟分析”过程，加快达标情景与优化情景的求解速度. 最后以柴河水库为研究实证区，通过模型概化、 参数灵敏度分析和率定与模型验证，完成了柴河水库允许入库负荷量的计算. 该方法可用于其它水体的允许纳污量计算和污染负荷管理.

参考文献
[1] 宋国君. 论中国污染物排放总量控制和浓度控制[J]. 环境保护, 2000, 6 (1): 11-13.
[2] 刘洁, 冯银厂, 朱坦. 总量控制在环境管理中应用[J]. 城市环境与城市生态, 2003, 16 (1): 59-61.
[3] 程声通. 河流环境容量与允许排放量[J]. 水资源保护, 2003, 19 (2): 8-10.
[4] 孟伟, 张楠, 张远, 等. 流域水质目标管理技术研究(Ⅰ)——控制单元的总量控制技术[J]. 环境科学研究, 2007, 20 (4): 1-8.
[5] Armstrong N E. Water quality modeling and total maximum daily loads[J]. Water Environment Research, 2001, 73 (2): 131.
[6] DePinto J V, Paul L, Freedman P E, et al. Models quantify the total maximum daily load process[J]. Journal of Environmental Engineering, 2004, 130 (6): 703-713.
[7] Hamrick J M. A Three-Dimensional Environmental Fluid Dynamics Computer Code: Theoretical and Computational Aspects[R]. Gloucester Point: The College of William and Mary, Virginia Institute of Marine Science. Special Report 317, 1992.
[8] 叶兴平, 张玉超. TMDL计划在污染物总量控制中的应用初探[J]. 环境科学与管理, 2008, 33 (8): 13-16.
[9] Elshorbagy A, Teegaravatapu R S V, Ormsbee L. Total maximum daily load (TMDL) approach to surface water quality management: concepts, issues, and applications[J]. Canadian Journal of Civil Engineering, 2005, 32 (2): 442-448.
[10] 王建平, 程声通, 贾海峰, 等. 环境模型参数识别方法研究综述[J]. 水科学进展, 2006, 17 (4): 574-580.
[11] 元红. 柴河水库饮用水源地水质现状调查及评价[J]. 内蒙古环境科学, 2008, 20 (3): 48-50.
[12] 于顺东, 尤学一. WASP水质模型检验及参数敏感度分析[J]. 水资源与水工程学报, 2007, 18 (6): 41-47.
[13] 贾海峰, 程声通, 杜文涛. GIS与地表水水质模型WASP5 的集成[J]. 清华大学学报(自然科学版), 2001, 41 (8): 125-128.
[14] 王建平, 苏保林, 贾海峰, 等. 密云水库及其流域营养物集成模拟的模型体系研究[J]. 环境科学, 2006, 27 (7): 1286-1291.
[15] 贾海峰, 张岩松, 何苗. 北京水系多藻类生态动力学模型[J]. 清华大学学报(自然科学版), 2009, 49 (12): 1992-1996.
[16] Zou R, Lung W S, Wu J. Multiple-pattern parameter identification and uncertainty analysis approach for water quality modeling[J]. Ecological Modeling, 2009, 220 (5): 621-629.
[17] Dilks D W, Freedman P L. Improved consideration of the margin of safety in total maximum daily load development[J]. Journal of Environmental Engineering, 2004, 130 (6): 690-694.
[18] Benaman J, Shoemaker C A. Methodology for analyzing ranges of uncertain model parameters and their impact on total maximum daily load process[J]. Journal of Environmental Engineering, 2004, 130 (6): 648-656.
[19] Melching C S, Yoon C G. Key source of uncertainty in QUAL2E model of Passaic River[J]. Journal of Water Resource Planning and Management, 1996, 112 (2): 105-113.
[20] Zhang H X, Yu S L. Applying the first-order error analysis in determining the margin of safety for total maximum daily load computations[J]. Journal of Environmental Engineering, 2004, 130 (6): 664-673.