洱海主要污染物允许排放总量的控制分配*

卫志宏，唐雄飞，杨振祥，吕兴菊，孟 良，朱 江，窦嘉顺，杨四坤

(中国大理洱海湖泊研究中心，大理 671000)

洱海位于云南省大理州中部，属澜沧江流域黑惠江支流天然水域，是云南省第二大高原淡水湖泊，为国家重点保护的水域之一[1].湖泊形似耳状，略呈狭长形，南北长42.5 km，东西宽5.9 km，呈北北西-南南东向展开(图1).西洱河出口断面以上控制流域面积为2565 km2，湖面正常水位(1974 m)对应的湖面面积为249.3 km2，其中岛屿面积0.748 km2，相应的湖容量28.80×108m3，平均水深10.6 m，最大水深21.3 m，湖岸线长127.85 km.

洱海流域属澜沧江水系，境内有弥苴河、永安江、罗时江、波罗江、西洱河及苍山十八溪等大小河溪117条.洱海北部三江(弥苴河、永安江、罗时江)占洱海入湖总径流量的70%，其中弥苴河是洱海最大的入湖河流，汇水面积1233.58 km2，多年平均来水量为4.6×108m3，约占洱海入湖总径流量的57.1%;西部苍山十八溪约占洱海入湖总径流量的18.6%，南部波罗江约占洱海入湖总径流量的11.4%;东部的海潮河、凤尾阱、玉龙河等小溪流全年大部分时间处于断流状态.洱海多年平均入湖总水量约为8×108m3，2010年入湖水量为7.65×108m3.西洱河是洱海的主要出湖河道，其次是“引洱入宾”隧洞.

近二十多年来，随着洱海流域人口增加和经济快速发展，人类对自然资源的开发不断加剧，流域生态环境逐渐恶化，入湖污染负荷不断增加，洱海水质日益下降，水体呈现富营养化趋势.洱海水质已由1990s 的Ⅱ～Ⅲ类下降为近年来的整体Ⅲ类，其中2003年7-9月为Ⅳ类，洱海及湖湾水质整体严重下降[2].洱海面临丧失饮用水水源地功能的潜在危险，由富营养化导致的蓝藻水华也使洱海水生生态系统受到破坏，影响和制约了区域的社会经济发展.

从湖湾水质来看，洱海北部湖湾是北部三江主要入湖河流外源污染物进入洱海北部的第一驿站，是第一时间的纳污水体.北部的沙坪湾、沙村湾、海潮河湾、红山湾、长育湾、东南部向阳湾等的主要污染源物如总磷、总氮、COD 和溶解氧指标均各不相同，但水质类别一年中多数时段为Ⅳ类，有时甚至为Ⅴ类，红山湾、长育湾等湖湾为Ⅲ～Ⅳ类.洱海湖湾水环境-水生态特征存在明显的差异性和典型性，这一差异是外源污染通过湖湾向全湖扩散，又通过水动力向湖湾运动堆积的结果.因此，洱海北部湖湾已成为洱海水质好坏最直接的晴雨表，北部湖湾水环境质量和生态系统的好坏直接关系到洱海北部乃至整个洱海的水环境质量.

近年来，流域污染物允许排放总量控制的理念在国内逐渐成为研究的热点[3]，因为其综合考虑了环境保护目标、污染源负荷及其分布特征、环境背景值等条件，对流域的排放污染源从整体上进行优化从而有效地削减污染负荷，使得环境质量得以改善.污染物排放总量的分配则是总量控制管理的核心和关键[4]，其基本思想是以环境保护目标为基础，综合分析研究区域内各个污染源的特征，通过建立水质数学模型以及结合最优化技术理念，计算出每个污染源相应的允许排放量.

国外对于污染物排放总量分配的研究开展较早，最著名的如美国环保署为履行清洁水法的相关条例而制定的TMDL(Total Maximum Daily Loads)水污染控制计划方案，用于改善全国地表水的水质状况[5].Zou等[6]运用EFDC 结合WASP 建立Wissahickon 河的水质模型，探讨了流域污染物分配以及削减方案的可行性.Wool 等[7]通过建立流域水文水质模型，分析研究Neuse 河的TMDL 实施计划.

近年来，国内对于污染物总量分配的研究也取得了一系列成果.李适宇等[8]采用分区达标控制法，通过线性规划求解汕头海域的环境容量.林高松等[9]采用满意度和公平偏离度作为衡量指标，根据最大化最小值求取最佳协调解，在广州-佛山跨界河网地区开展了排污量公平分配的研究.吴悦颖等[10]通过借鉴评价收入公平性的指标——基尼系数的基本概念，构建了流域间水污染物总量分配方案合理性的评估方法.李如忠等[11]建立了一种定性与定量相结合描述判断矩阵的多指标决策的排污总量分配层次结构模型，用于合肥市水污染物排放总量的分配研究.

现阶段国内对于洱海的污染物允许排放总量的控制分配鲜有应用研究，因此本文利用环境流体动力学模型EFDC(Environmental Fluid Dynamic Code)建立洱海湖泊及湖湾的三维水动力水质模型，并在此基础上结合最优化技术计算洱海主要污染物的允许排放总量控制分配方案.

1 水动力和水质模型

1.1 模型概述

环境流体动力学模型EFDC(Environmental Fluid Dynamics Code)是美国威廉玛丽大学的佛吉尼亚海洋科学研究所(VIMS，Virginia Institute of Marine Science at the College of William and Mary)Hamrick 教授开发的三维地表水模型，是美国国家环保署推荐的三维地表水水动力水质模型，可实现河流、湖泊、水库、湿地系统、河口和海洋等水体的水动力学和水质模拟，是一个多参数的有限差分模型[12].EFDC 是一个公开的模型代码，已经应用于超过100 个实例，是一个非常成熟的模式，在高校、科研单位、政府部门和咨询公司得到了广泛的应用.

1.2 模型设置

1.2.1 模型计算网格 EFDC 平面采用贴体坐标网格，可以精确拟合洱海不规则的岸线.模型水平方向剖分为1013 个计算网格，平面空间分辨率最小为277 m，平均大小为409 m.网格在垂向上等距分5 个 σ 层，以便较好地拟合洱海的地形.计算区域的原始地形数据采用 1 ∶100000 的云南省大理洱海数字化地形图，采用克里金插值法得到各离散计算网格点的水深数据.洱海三维水动力水质模型的计算网格和监测站点如图1 所示.

洱海流域的主要入湖河流包括弥苴河、罗时江、永安江、波罗江，以及洱海西边入湖的阳溪、龙溪、茫涌溪、锦溪、万花溪等十八条小溪.洱海三维水动力水质模型包括了主要的23 条入湖河流以及西洱河出流和海东区“引洱入宾”取水口.

1.2.2 模型主要参数取值 水动力模型的底部摩擦系数取决于湖底海床以及底栖植物的分布情况，经过多次调试和测算后，在底栖植物分布较多的地方底部摩擦系数取为0.020 m，其他部分则取为0.015 m.

水质模型主要参考 Tetra Tech[13]和 Park 等[14]的研究成果，主要模型参数取值见表1.

1.2.3 风场处理 环湖布置了6 个风场监测站点(图2)，监测频率为30 min/次.EFDC 模型可同时考虑时间和空间变化的风场，因此采用泰森多边形法确定各个计算网格所采用的风场数据.将所有相邻的风场监测站连成三角形，作这些三角形各边的垂直平分线，于是每个监测站周围的若干垂直平分线便围成一个多边形，并且多边形内只存在唯一一个监测站，这个监测站的风场数据即可赋予多边形内的计算网格.

图1 洱海水动力水质模型计算网格及监测站点Fig.1 Computational grid of three-dimensional hydrodynamic and water quality model and the monitoring stations in Lake Erhai

表1 水质模型主要参数取值Tab.1 Main parameters of water quality model

1.2.4 初始条件 各水质因子初始条件参考实测资料进行给定.初始场赋值采取空间一致的原则给定一个较低的浓度值，水质模型在此基础上运行一段时间后可达平衡.为了最大程度上消除初始条件给模型带来的影响，先进行为期一年的计算热启动，计算结果作为正式计算的初始条件输入.

1.2.5 入流污染负荷 根据2002-2010年洱海入湖河流的水量和水质数据确定入湖污染负荷.主要的入湖河流如罗时江、弥苴河、永安江、波罗江、锦溪和茫涌溪均有逐日入流量数据，其他较小的河流如万花溪、龙溪等只有月均入流量，因此需转换成日均流量.水质监测频率为每月一次，利用每日入湖流量×当月的入湖水质浓度，即可得到该河流当日的入湖污染负荷.

由于缺乏洱海流域实测的干湿沉降资料，采用文献[15]的数据，干湿沉降入湖污染负荷总氮、总磷月平均值分别为28.3 和1.3 t，年营养盐入湖污染负荷总氮和总磷分别为339.6 和15.6 t.采用文献[1]的数据，洱海底泥总磷年释放量为9.7 t，总氮年释放量为333.2 t.

1.2.6 其他设置 与一般河流和海洋水体相比，湖泊和水库的流动性较差，流域的点源和非点源污染进入湖泊后由于水体滞留时间较长，需要长时间的模拟才能反映出水体水质的物理生化过程.因此洱海水动力水质数学模型的计算时间为10年(2001年1月1日至2011年4月30日).

2 模型验证和结果分析

2.1 水动力模型验证

2.1.1 水位和库容验证 湖泊的运行水位和库容是水环境容量的重要影响因素[16].通过对洱海模拟的水位、库容与相应实测值之间的比较表明，所建的水动力模型可以准确反映洱海水位和库容的季节性变化(图2).

图2 洱海三维水动力水质模型模拟水位(a)和库容(b)验证Fig.2 Simulated water level(a)and lake capacity(b)for verification of three-dimensional hydrodynamic and water quality model of Lake Erhai

2.1.2 流速验证 洱海部分监测站点的流速验证结果(图3)显示，水动力模型模拟值和观测值在同一流速范围，水动力模型能较好地反映洱海水动力特征.

图3 洱海三维水动力模型部分站点表层流速验证Fig.3 Flow velocity verification of three-dimensional hydrodynamic model of Lake Erhai with measured velocity data at selected monitoring stations

2.1.3 水温验证 洱海部分水质监测站点的水温时间序列和垂向水温分布验证结果显示，模拟和实测水温吻合度较好，洱海三维水动力模型可以很好地反映洱海湖体水温年际变化以及垂向水温分布情况.由实测和模拟的垂向水温分布情况可以看出，洱海垂向水体分层情况不明显，表层和底层水体温度差较小(图4).

图4 洱海三维水动力模型部分站点水温验证Fig.4 Water temperature verification of three-dimensional hydrodynamic model of Lake Erhai

2.2 水动力模拟结果分析

风是影响洱海水体流速的主导因素，表层水体流速可以达到风速的2%～4%左右.洱海表层水体受风拖曳力的直接影响，流速相对较大，流向基本上与风拖曳力方向一致;底层水体受底部摩擦力及地形等因素的影响，流速相对较小;中层水体为过渡层，流速介于表层和底层水体之间.洱海全湖各水域流域差异较大，近岸的流速大于远岸水域.湖心流速较小，一般均小于10 cm/s;而在近岸区域最大流速可以高于20 cm/s.这个差异主要是洱海本身地形地势造成的，洱海深水区位于湖心，而近岸边界水浅且较为平滑，易产生较强的沿岸流.

在主导风西南风和东南风控制时，北部湖湾形成顺时针流场;其他风向控制时，顺时针流场不明显.个别时刻北部湖湾会出现几处局地环流，都在靠近弥苴河、永安江和罗时江等区域，这是风场、地形和入湖水流共同作用的结果.洱海北部湖湾表层水体表现出明显的顺时针环流系统，中层和底层的水体因为表层水体的流动而表现为补偿上升流形态.三维流场也清楚地显示北部湖湾湖心处水深大而西岸水深浅，水深梯度导致西岸流速较大、湖心处流速小.顺时针环流系统可以将沿岸污染物快速带走，但湖心处则较难扩散.在沙坪湾、海潮湾、红山湾等多个局地海湾，局地流速和湖心流速形成相反的流态格局，这主要是因为湖心环流系统引起的水位压强梯度差.

洱海风生流形成过程的数值试验表明，早期阶段洱海水流由风拖曳力直接影响，过渡阶段由水位压强梯度力、地球自转柯氏力以及湖底地形共同作用，稳定阶段则在北部湖湾形成顺时针环流;中部靠西形成顺时针环流;中部靠东形成逆时针环流;南部形成顺时针环流.

2.3 水质模型验证

水质模型的验证效果直接关乎到主要水质因子环境容量计算结果的合理性.目前洱海全湖布设有多个水质连续监测站点(图1)，记录了从2001年至今的叶绿素a(藻类)、五日生化需氧量(BOD5)、氨氮、总氮和总磷等指标浓度.从部分监测站点水质参数验证结果可以看出，总磷、总氮、氨氮、叶绿素a(藻类)和BOD5的计算结果和实测结果吻合较好(图5)，较好地体现了污染物的时空分布规律.总体而言，水质模型的计算结果合理，与实测的水质情况及变化趋势一致，基本反映了洱海主要的水质变化过程，可用于实际的水质预测分析和主要污染物允许排放总量控制分配的计算.

2.4 水质模拟结果分析

氨氮、总氮、总磷和BOD5表现出了较为明显的时空差异性.洱海污染负荷基本上夏季比冬季高，这主要是因为夏季期间，降雨径流携带的入湖污染负荷比冬季高.在北部湖湾和南部湖区，因为入湖河流带来的含氮营养物质引起这些水域的氨氮浓度较高.在湖心区，污染物经过风生流引起的扩散降解作用，浓度较低.南部波罗江位于大理市市区，承载了众多的生活污染负荷，这些高氮、磷的生活污染源导致南部水域总氮、总磷浓度较高.北部湖湾集结了弥苴河、罗时江、永安江等较大的入湖河流，因此相应北部湖湾的总氮、总磷浓度会较湖心区高.

图5 洱海三维水质模型部分站点水质因子验证Fig.5 Water quality parameter verification of three-dimensional hydrodynamic model of Lake Erhai

洱海水体浮游植物的生长受到水温、营养盐及水动力因素的共同影响.夏季水温较高，有利于浮游植物的生长，同时北部湖湾的环流系统不利于污染物的扩散，在湾内累积了大量的氮、磷营养盐，这为浮游植物的生长提供了丰富的养料.而在冬季，一方面低水温抑制了浮游植物的生长;另一方面，冬季较强的水体交换能力使湾内的营养盐浓度下降，也不利于浮游植物的生长.这两方面的原因使得冬季洱海的叶绿素a 浓度明显低于夏季.洱海的垂向水温层结效应不明显导致藻类表、底层差别不大.叶绿素a 的浓度变化也表现出年际变化，当营养物质入湖负荷量较大时，叶绿素a 的峰值也会变大.

3 污染物允许排放总量控制分配

3.1 计算原理

污染物允许排放总量控制分配的思想可以定义为:在选定的各个水质控制点的主要污染物浓度不超过其对应的环境保护目标，以及入湖排污口不超过其相应的环境保护目标前提下，各个入湖排污口的污染负荷排放量之和最大，计算公式如下:

式中，L 为所有入湖口的污染物负荷排放总量(g/s);Qj为各个入湖口的流量(m3/s);Cj为各个入湖口的浓度(mg/L);Aij为第j 个入湖口的单位污染负荷量对第i 个水质控制点的污染贡献率;Cbi为第i 个水质控制点处的背景浓度(mg/L);Coi为第i 个水质控制点处的环境标准值(mg/L);Crj为第j 个入湖口的浓度上限值(mg/L).

此优化问题可以应用水质模型结合单纯形法得以解决，其本质是求解线性规划问题各个入湖排污口Cj的最优解.洱海主要污染物允许排放总量分配计算的流程见图6.

图6 洱海主要污染物允许排放总量分配计算的流程Fig.6 Flow chart of waste load allocation of Lake Erhai

图7 水域环境达标的分区控制Fig.7 Zone control of water environment

3.2 分区控制

本论文引进污染混合区的概念，当污染物由河流携带入湖时可以允许在污染物混合区出现水质超标现象，但在污染混合区的边界则要求水质达到相应的标准.例如入湖排污口附近的污染混合区(mixing zone)可以为Ⅳ类水质，但出了污染混合区后则要求水质变为Ⅲ类，经过一定的降解距离，水质可以达到Ⅱ类水标准.本文设置了两种水质控制点，在二级水质控制点处水质要求达到Ⅲ类水质标准;在一级水质控制点处水质要求达到Ⅱ类水质标准(图7).

根据洱海相关的水污染控制规划文献[2]和[15]得知，洱海入湖湖湾部分水域由于河流携带较高浓度污染物，可以允许出现水质超标现象，但经过一定距离后就要严格执行Ⅱ类水质标准.洱海主要的入湖排污口共19 个(图8)，一级水质控制点共228 个，二级水质控制点共36 个.水质控制点涵盖了洱海绝大部分区域，其中包括了北部湖湾、沙坪湾、海潮湾、红山湾、向阳湾等重要的湖湾，确保各区域水质达标.

表2 洱海部分入湖排污口对水质控制点的贡献率(×10 -2)Tab.2 Contribution rate of pollutant load inflows to lake water quality control points

3.3 贡献率和背景浓度

贡献率实际上会随着入湖排污口的浓度变化而变化.为了简化允许排放总量控制分配的计算，同时也确保计算成果的可操作性，本文采用的贡献率为一个模拟周期内的平均值.在研究中以弥苴河为例，设定弥苴河的污染物为1 个单位，其他入湖排污口污染负荷为0，进行连续3年的水质模拟，统计后两年各个水质控制点的平均值，进而求得平均贡献率.洱海主要入湖排污口共19 个，那么针对每个入湖排污口都需要单独进行一次水质模拟，求得的平均贡献率为一个19×264 的大型矩阵.因贡献率矩阵过于庞大，本文只列出部分结果(表2).

图8 洱海入湖排污口以及各水质控制点分布Fig.8 Distribution of the water quality control points and sewage outfalls in Lake Erhai

表3 洱海主要污染物允许排放总量及其分配Tab.3 Allowable waste loads allocated for different inflow sources with respect to major pollutants for Lake Erhai

各个水质控制点处的背景浓度理论上应该采用年均实测值，但由于洱海水质监测站点较少，而水质控制点较多，因此采用经过率定和验证的洱海三维水质模型计算结果提供.

3.4 排放总量分配结果及分析

根据上述的设定求解主要污染物总氮、总磷和CODMn的允许排放量及其分配方案.由于入湖排污口较多，因此本文根据洱海的北区、西区、南区和东区4 个分区进行污染物排放总量计算及分配方案统计(表3).由计算结果可知，北区的三条河流——弥苴河、永安江、罗时江因为本身流量相对较大，所以排放总量也较大，总氮、总磷和CODMn的允许排放总量分别占到整个洱海的47%、53%和49%.西区18 条小溪的排污总量相对北区要小，总氮、总磷和CODMn的允许排放总量分别占到整个洱海的34%、28%和30%.南区的允许污染物排放量主要由波罗江贡献，总氮、总磷和CODMn的允许排放总量分别占到整个洱海的13%、15%和16%.东区由于河流流量都非常小，导致污染物排放量也最小.

本研究所计算的主要污染物允排量及其分配方案是多年平均值，而环境标准值本身也是含有时间平均意义的，因而按照平均状态求环境容量及其分配是合理的.另外，如果想要求解使得各个水质控制点在任何时刻都不超标，只需要在计算公式上多加一组限制条件，即水质控制点处的水质在极端状况下浓度不可超过环境标准即可.

4 结论

1)为了研究洱海湖泊及其湖湾的主要污染物允许排放总量及其分配方案，本文利用环境流体动力学模型EFDC 建立洱海湖泊及湖湾的三维水动力水质模型.模型验证结果较好，表明可以用于模拟分析洱海的水流、水位和水质因子的变化趋势.

2)运用经过率定和验证的水质模型计算出各个入湖排污口对水质控制点的污染贡献率.通过设置一系列水质控制点并使其达到环境标准值，同时对入湖排污口的浓度也进行约束限制，运用单纯形法求解此最优化问题，得到了各个入湖排污口的允许排放量.

3)洱海北区的三条河流——弥苴河、永安江、罗时江因为本身流量相对较大，所以排放总量也较大，总氮、总磷和CODMn的允许排放总量分别占到整个洱海的47%、53%和49%.西区和南区的允许排放量相对北区要小，东区最小.

4)本文所建立的洱海三维水动力和水质模型，以及允许排放量分配方案求解的方法，为洱海流域的水环境综合整治提供了科学依据，具有重要的意义.

[1]彭文启.洱海水质评价.中国水利水电科学研究院学报，2005，3(3):192-198.

[2]洱海水污染防治规划(1996-2010).大理州城乡建设环境保护局，2001.

[3]梁 博，王晓燕.我国水环境污染物总量控制研究的现状与展望.首都师范大学学报:自然科学版，2005，26(1):93-98.

[4]张天柱.水污染物排放总量控制管理的经济原则.环境科学，1991，11(6):2-6.

[5]USEPA.Overview of current total maximum daily load-TMDL-program and regulations，2005.http://www.epa.gov/owow/tmdl/intro.html.

[6]Zou R，Carter S，Shoemaker L.Integrated hydrodynamic and water quality modeling system to support nutrient total maximum daily load development for Wissahickon Creek，Pennsylvania.Journal of Environmental Engineer，2006，132(4):555-566.

[7]Wool TA，Davie SR，Rodriguez HN.Development of three-dimensional hydrodynamic and water quality models to support total maximum daily load decision process for the Neuse River Estuary，North Carolina.Journal of Water Resources Planning and Management，2003，129(4):295-306.

[8]李适宇，李耀初，陈炳禄.分区达标控制法求解海域环境容量.环境科学，1999，20(4):96-99.

[9]林高松，李适宇，江 峰.河流允许排污量公平分配的多准则决策方法.环境科学学报，2007，27(3):494-500.

[10]吴悦颖，李云生，刘伟江.基于公平性的水污染物总量分配评估方法研究.环境科学研究，2006，19(2):66-70.

[11]李如忠，钱家忠，汪家权.水污染物允许排放总量分配方法研究.水利学报，2003，121(5):112-115.

[12]Hamrick JM.User's manual for the environmental fluid dynamics computer code，Special Rep.331 in Applied Marine Science and Ocean Engineering，Virginia of Marine Science，College of William and Mary，Gloucester，1996.

[13]Tetra Tech，Inc.The environmental fluid dynamics code theory and computation volume 3:Water quality module.Tetra Tech，Inc.10306 Eaton Place，Suite 340 Fairfax，VA 22030，2006.

[14]Park K，Kuo AY，Shen J.A three-dimensional hydrodynamic-eutrophication model HEM-3D:ion of water quality and sediment process submodels.VIMS VA:Sramsoe，1995.

[15]大理白族自治州人民政府.云南洱海流域水污染综合防治“十二五”规划，2012.

[16]Ji ZG.Hydrodynamics and water quality-modeling rivers，lakes，and estuaries.Hoboken:John Wiley ＆ Sons，Inc.，2009:510-511.