北方干旱地区水库突发污染事件应急模拟分析

孙凯迪，徐明德，安 静

(太原理工大学环境科学与工程学院，太原 030024)

我国北方干旱地区缺水问题严重，水库是很多城市的主要供水水源，汾河水库呈峡谷形、是一座具有调节径流、防汛、灌溉、城市供水等多项功能的北方水库，水库控制流域大，承担着工业、农业及城市生活用水的任务，是山西省内最大的集中式饮用水源地[1,2]，其水质好坏与太原市超过400 万人口的饮水安全息息相关。一旦突发污染事件将严重影响水质，对供水终端水处理技术提出了更高的要求，增加了处理成本，威胁供水区域居民饮水安全，因此对其进行突发污染事件研究具有重要意义。

国内外关于突发污染事件的研究起步较晚[3,4]，已有研究主要针对海域及河流[5-7]，水域突发污染事件数值模拟研究重点关注海上溢油事件，针对干旱地区中小型水库的研究相对较少[8]，研究所使用的模型从简单的回归分析、营养物平衡模块发展到生态动力学水质模型[9,10]。鉴于此，本文以汾河水库为研究载体，识别其移动风险源，结合水动力及物质输移扩散模型，设置不同情景，模拟水库突发环境事件时污染物的迁移扩散规律。

1 研究区域概况

1.1 水库概况

汾河水库南北长15 km，东西宽5 km，最大回水面积32.0 km2，水库库容7.22 亿m3，但经过多年的淤积，截止到2008年，汾河水库的库容已经减少为3.6 亿m3。近10 a，年均来水量2.14 亿m3，年均供水1.79 亿m3，入库河流主要为水库北部的汾河和西部的涧河。水库地理位置见图1。

图1 汾河水库地理位置Fig.1 The location map of Fenhe Reservoir

1.2 风险源

研究区域内固定风险源主要位于汾河水库上游区域，各工业企业水处理措施比较完善，故本研究主要考虑移动风险源引起的污染事件。如果汾河水库周边道路运输有毒有害危险物质的车辆发生交通事故造成泄露物进入汾河水库，将对水库水源地构成严重威胁。交通风险路段包括：

(1)省道252(岚古线)，处于水库库区东边，车流量为450～750 辆/d。

(2)省道217(岚马线)，处于水库库区西边，车流量为300～550 辆/d。

危险物质包括汽油、柴油、粗苯、液氨等污染物。

2 模型建立

2.1 水动力模型

汾河水库是典型的干旱地区中型水库，水深较浅，浓度、流速、水深等水力参数在垂直方向变化较小，故水库水动力模拟研究建立浅水流动的二维水动力-水质耦合模型[11]：

连续性方程：

(1)

x方向的动量方程：

(2)

y方向的动量方程：

(3)

式中：ζ为自由水面水位；H为水深，H=h+ζ；u、v为x、y方向上的速度分量；f为科氏力系数；fw为风阻力系数；g为重力加速度；μ为黏性系数；Wx、Wy为风速在x、y方向上的速度分量；t为时间；p为静压强。

2.2 水质模型

(4)

2.3 模拟范围确定

模拟区域为水库全部水域范围，考察水库实际地形，基于DEM模型、遥感影像资料及GPS实测信息，分析确定水库边缘点位坐标，利用GIS地理信息系统对其进行数字化处理，生成汾河水库水陆边界。采用三角网格系统，将汾河水库库区概化约3 160个三角形网格[12,13]。

2.4 参数及定解条件

(1)模拟时间步长。根据CFL条件取模拟时间步长Δt=1 800 s。

(2)涡黏系数。采用Smagorinsky公式[14]，涡黏系数取0.4。

(3)底床阻力。根据曼宁系数(Manning)，底床阻力取32 m1/3/s。

(4)摩擦力。根据Wu[15](1980)的经验公式，取风场摩擦力ca=1.255×10-3；底床摩擦力选取曼宁系数32 m1/3/s。

(5)降雨量与蒸发量。根据气象部门相关资料生成对应的降雨和蒸发时间序列文件。

(6)源与汇。将西部涧河入库口、北部上游入库口及南部坝址处的进出水口作为水库源汇项。各源汇项生成对应的时间序列文件。

(7)初始条件。假设模拟开始前水库处于静止状态，x、y方向上的速度分量均为0，水库水位初始条件定为边界水位平均值，取1 126.2 m。

(8)边界条件。闭边界采取岸壁法，取法线方向为不可入条件。开边界为北部、西部进水口和南部出水口，为避免过强的浅水效应[16]，采用了“干湿点判别法”，即增水时，水深大于0.1 m，视为水域处理；退水时，水深小于0.005 m，不作为水域处理。

(9)扩散系数。根据Smagorinsky公式，扩散系数与涡黏系数之间的比例系数取0.1，即扩散系数取0.04。

(10)降解系数。出于安全性，降解系数取0。

2.5 模型率定与验证

在水库西北部和东南部布设A、B2个监测点，见图2。采用2016年4月5日至2016年8月5日实测A、B点的溶解氧浓度值进行率定和验证，A、B点的溶解氧浓度值见表1。

图2 汾河水库水质监测点位Fig.2 The monitor position map of Fenhe reservoir

(1)参数率定。反复调整模型参数使模拟结果与A点实测数据相吻合，设置模型最终计算参数，率定结果见图3。

(2)模型验证。在已确定模型参数的基础上，用B点实测数据重新计算，模拟值接近实测值，两者拟合程度较好，验证结果见图4。根据偏差统计法，最大绝对值偏差为0.2 mg/l,相对误差2.1%，满足精度要求，因此本文所构建的模型是有效的。

图3 参数率定Fig.3 The parameter calibration

图4 模型验证Fig.4 The model validation

3 水质模拟

3.1 情景设置

本次研究以苯胺泄漏为例，考虑移动污染源下排放入口、风等因素对污染物的扩散情况及浓度分布规律产生的影响，主要设置了10种情景进行模拟。若汾河水库周边道路上运输苯胺的车辆突发交通事故，造成苯胺泄漏，考虑从水库上游北部进水口、西部进水口和东侧排放口不同入口流入水库，具体情景设置见表2。

表2 突发事件情景设置Tab.2 Classifications of scenarios for emergencies

3.2 水质模拟结果分析

3.2.1 模拟结果

水库突发污染事件48 h后，各情景下污染物的扩散面积及浓度分布情况见图5。

图5 苯胺的浓度分布Fig.5 Concentration distribution of Aniline

由情景1的模拟结果可以看出，苯胺从上游北部进水口进入水库，由于西南风的作用，苯胺从上游北部进水口顺着水库东岸自西北向东南不断扩散，48 h扩散了16.07 km2(苯胺浓度大于0的库区面积)，平均扩散速率约93 m2/s，超标面积为3.43 km2(苯胺浓度大于0.1 mg/L的库区面积)；48 h后距离入口1 300 m区域苯胺浓度最大，浓度值21.5 mg/L，污染物超标215倍。

由情景2的模拟结果可以看出，苯胺从西部涧河进水口进入水库，由于西南风的影响，苯胺顺着水库北岸自西南向东北扩散，48 h扩散了15.98 km2，平均扩散速率约92.48 m2/s，超标面积为0.60 km2；48 h后距离入口西北部300 m区域苯胺浓度最大，浓度值3.65 mg/L，污染物超标36.5倍。

由情景3的模拟结果可以看出，苯胺从东侧排放口进入水库，由于西南风的影响，苯胺从进水口处逐步向外扩散，48 h扩散了25.79 km2，平均扩散速率约149.25 m2/s，超标面积为0.55 km2；48 h后东侧排放口周围区域苯胺浓度最大，浓度值1.21 mg/L，污染物超标12.1倍。

由情景4的模拟结果可以看出，苯胺从上游北部进水口进入水库，由于西南风的影响，苯胺顺着水库南岸自西南向东北扩散，48 h扩散了为31.74 km2，平均扩散速率约183.68 m2/s，超标面积为7.86 m2；48 h后距离北部进水口东部1 800 m区域苯胺浓度最大，浓度值8.2 mg/L，污染物超标82倍。

由情景5的模拟结果可以看出，苯胺从上游北部进水口进入水库，由于西北风的影响，苯胺从北部进水口顺着库区西岸自西北向东南扩散，48 h扩散了17.36 km2，平均扩散速率约100.46 m2/s，超标面积为3.69 km2；48 h后在入口东南部1 500 m区域苯胺浓度最大，浓度值22.2 mg/L，污染物超标222倍。

由情景6的模拟结果可以看出，苯胺从上游北部进水口进入水库，由于西北风的影响，苯胺从进水口顺着库区西岸自西北向东南扩散，48 h扩散了29.89 km2，平均扩散速率约105.6 m2/s，超标面积为7.26 km2；48 h后距入口南部1 000 m区域苯胺浓度最大，浓度值9.53 mg/L，污染物超标95.3倍。

由情景7的模拟结果可以看出，苯胺从西部涧河进水口进入水库，由于西北风的影响，苯胺顺着库区南岸自西南向东北扩散，48 h扩散了16.44 km2，平均扩散速率约95.14 m2/s，超标面积为0.61 km2；48 h后距入口西南侧500 m区域苯胺浓度最大，浓度值1.55 mg/L，污染物超标15.5倍。

由情景8的模拟结果可以看出，苯胺从北部进水口进入水库，沿着水库中心区域由西北向东南逐步扩散，48 h扩散了12.79 km2，平均扩散速率约74 m2/s，超标面积为3.92 km2；48 h后入口东南部1 500 m区域苯胺浓度最大，浓度值25.2 mg/L，污染物超标252倍。

由情景9的模拟结果可以看出，苯胺从西部涧河进水口进入水库，从进水口自西南向东北扩散，48 h扩散了6.99 km2，平均扩散速率约40.39 m2/s，超标面积为0.46 km2；48 h后在排放入口泄漏点周边区域苯胺浓度最大，浓度值3.7 mg/L，污染物超标37倍。

由情景10的模拟结果可以看出，苯胺从东侧排放口进入水库，从东北向西南不断迁移，48 h扩散了13.81 km2，平均扩散速率约79.91 m2/s，超标面积为0.62 km2；48 h入口处邻近区域苯胺浓度最大，浓度值2.82 mg/L，污染物超标28.2倍。

对比情景1和5，情景4和6污染物的扩散情况可知，当排放入口、初始浓度及风速相同时，西北风条件下污染物扩散速率比西南风快，汾河水库呈峡谷形，进水情况为西北至东南走向，由此可得风向是影响污染物扩散的重要因素。

对比情景1、2、3，情景8、9、10污染物的扩散情况可知，当事件源、初始浓度、风力条件相同时，东侧排放口污染物的扩散范围最大，北部上游干流入口扩散面积大于西部涧河入口。汾河水库东侧入口的水面明显较北部和西侧涧河进水口宽广，水体流量也大，由此可得排放入口也是影响污染物扩散的重要因素。

对比情景3和4、情景5和6污染物的扩散情况可知，若初

始浓度、排放入口及风向相同时，风速对污染物扩散规律影响较大，风速增加，污染物的扩散及超标面积随之增加，当西北风风速由1.4 m/s增加到5.6 m/s时，苯胺扩散面积从17.36 km2变化至29.89 km2，苯胺超标面积从3.69 km2变化至7.26 km2，扩散面积增加了72%，超标面积增加了49%。

以此水质模型为基础，对不同风速苯胺的浓度分布进行模拟，(浓度取常量1 000 mg/L)，将风速与污染物超标及扩散面积的模拟结果进行拟合，得到扩散面积及超标面积随风速的变化曲线，给出西北风条件下不同风速西部涧河入口苯胺的浓度分布见图6，此方法可以应用在水质控制管理上。

图6 污染物扩散及超标情况Fig.6 The proliferation of pollutants and excessive situation

3.2.2 应急预估

汾河水库突发污染事件时，为了能及时准确地确定污染物的扩散情况及浓度分布规律，结合汾河水库实际地形，以DEM模型、遥感影像资料、GPS实测信息及GIS地理信息系统为基础，模拟输出了不同时间点汾河水库突发污染事件时的污染物浓度分布与迁移扩散规律，模拟结果与实际相符，建立应急预估表(见表3)，可准确直观的反映突发水污染事件的影响范围及发展趋势，应急预估具有一定的可靠性。

表3 汾河水库突发污染事件时污染物扩散预估表Tab.3 The diffusion law of pollution emergencies in Fenhe reservoir

续表3 汾河水库突发污染事件时污染物扩散预估表

注：移动源，初始浓度1 000 mg/L。

4 结 论

(1)当排放入口、事件源及初始浓度相同时，风向决定污染物的扩散路径，风速决定污染物的扩散速率。静风条件下污染物以泄露点为中心逐步向外扩散，有风时污染物扩散路径随风向改变，扩散速率增加，污染物扩散面积、超标面积与风速近似服从线性关系；风速相同时，污染物的扩散速率西北风>西南风；排放入口污染物扩散面积大小为东侧排放口最大，北部上游干流入口大于西部涧河入口。

(2)情景4即污染物排放入口为北部上游干流入口，西北风为5.6 m/s时，污染物扩散48 h的扩散面积和平均扩散速率最大，面积为31.74 km2，速率为183.68 m2/s；情景9即排放入口为西侧涧河，静风条件下48 h的扩散面积和平均扩散速率最小，面积为6.99 km2，速率为40.39 m2/s。

(3)以汾河水库为研究载体，识别其移动风险源，构建水动力及物质输移扩散耦合模型对不同情景汾河水库突发环境事件进行模拟分析，建立水库突发环境事件应急预估表，一旦污染事故真正发生时，可快速确定污染物的扩散路径、扩散范围以及浓度分布情况，为合理进行应急监测，制定应急处置方案提供数据支持，从而及时采取应急防范措施。

□