刘庆龙 金斌全 赵伟 李一卓

摘要：本文以广西龙江河镉污染事件的数据为依据，从重金属镉传播特性角度出发，结合系统动力学理论，率定了重金属镉污染扩散系数，建立了镉污染一维传递模型，探讨了在不同流速下镉在各个监测点之间的传播过程。利用模型预测了洛东电站镉污染物浓度，与实际污染数据对比其[R2]=0.867，并进行可视化展示。

关键词：水体污染;迁移过程;数据分析;一维水质模型;龙江河

中图分类号：TP311      文献标识码：A

文章编号：1009-3044（2021）32-0011-05

Design and Implementation of Cadmium Pollution System Model in Longjiang Based on MATLAB

LIU Qing-long， JIN Bin-quan， ZHAO Wei， LI Yi-zhuo

（College of Engineering and Design， Hunan Normal University， Changsha 410081， China）

Abstract： Based on the data of cadmium pollution event in Longjiang river of Guangxi， the diffusion coefficient of heavy metal cadmium pollution was calibrated from the perspective of transmission characteristics of heavy metal cadmium， combined with the theory of system dynamics. The one-dimensional transfer model of cadmium pollution was established， and the transmission process of cadmium between monitoring points at different flow rates was discussed. The model was used to predict the concentration of cadmium pollutants in Luodong power station， and compared with the actual pollution data， the R2 = 0.867， and the visualization was carried out.

Key words： water pollution; migration process; data analysis; one dimensional water quality model; Longjiang river

水资源是人类赖以生存的基础。然而随着经济的发展，近年我国突发水污染事故频率增加，其中流域突发水污染事故约占突发环境事故的94.7%[1，2]。根據中华人民共和国生态环境部在2019年5月份发布的《2018中国生态环境状况公报》显示，2018年，全国地表水检测的1935个水质断面（点位）中，Ⅰ～Ⅲ类比例为74.3%;劣Ⅴ类比例为6.9%[3]，如图1所示。重金属镉在自然界中多以化合物的形态微量存在，我国规定居民用水中镉的含量不得超过0.02mg/L[4]。由于镉元素难以被生物体内分解，镉元素浓度高于标准的生活用水会导致镉元素在生物体内聚集，长期使用被污染水体、食物等会严重影响居民的身体健康。水体中重金属污染物的污染程度取决于重金属的类型、含量和水体的形态。随着重金属废水的排放，进入水体的重金属离子会发生吸附、扩散、沉降等作用，但重金属元素本身的衰减极其微弱。重金属离子能够通过水体自身的迁移作用进行传播，并在传播过程中不断被水体中的悬浮物以及河床中的沉积物吸附，从而达到新的平衡状态[5]。相对于有机物污染，重金属镉的衰减性很弱，在研究中常称之为保守型污染物质或惰性污染物。根据重金属污染物在水体中迁移的方向不同，可以分为水平、垂直方向的扩散和沿水流方向的迁移[6]。

现在学界广泛采用的水体重金属污染传播模型，大部分根据质量平衡原理构建。这类模型的优势在于可以分析河床对重金属的吸附，能够有效提高模型的精度，但其弊端在于建模过程比简单模型更复杂，模型建立所需的时间较长。突发型水污染事件往往在应急方面有很高的要求，模拟的精度方面可以通过后期的应急监测来不断调整和完善，一般的水相污染模型可以满足突发型水污染事件的要求，故本设计讨论的污染模型不考虑底泥对重金属的吸附和降解作用，主要考虑重金属污染在水体中的水平输移过程。

本文通过分析广西龙江河基本水体特征，从重金属镉传播特性角度出发[7]，对龙江河2012年水质数据进行处理，利用系统动力学理论建立一维龙江水质模型，体现出镉在各个监测点之间的传播过程，并进行可视化展示。

1水污染事故水质模拟案例分析

1.1研究背景及概况

广西龙江河地处黔桂珠江水系典型的亚热带喀斯特岩溶区域，具有独特的水质特征。广西河池市作为中国著名的有色金属产地，各类有色金属矿产资源储备量高居全国前列。龙江河发源于拉浪水库，截止于红花电站。龙江水域地图如图2所示。

龙江河总水域面积约68km2，多年平均流量为67.2m3/s，河床宽100～200m，龙江最大流量为4240m3/s，最大流速2.7m/s，最大水面宽172m。正常流量147m3/s，正常流速1.03m/s，正常水面宽105m。枯水期流量1963年仅13.1m3/s[8]。

根据龙江水域的基本水文数据，采样点附近的水力学特征参数如表1所示。

1.2实验数据筛选

由于镉污染事发突然，受采样站点的布置等因素的影响，原始采集数据中存在一定的异常值（负值、异常突变值等），所以需要对龙江水污染原始数据中的异常值进行缺失值替换处理，预处理前后数据如图3所示。

对龙江水污染数据进行筛选和处理，经处理后的数据较处理前的数据相比，更能够清晰、直观地体现出龙江水污染各个点镉浓度峰值及其变化曲线，处理后的数据也更加准确，便于进行后续龙江水污染模型的建立与分析。

2 系统模型简介

2.1一维水质模型的建立

对于突发的水污染事故而言，进入水体的一般为守恒污染物，即在自然条件下不会进行主动降解的污染物（如重金属、高分子有机物等），故一般不考虑污染物的衰减过程。而污染物的分散作用过程中，一般情况下河流污染物的弥散过程影响要远大于污染物的分子扩散与湍流扩散，因此对于污染物的分散作用，一般只考虑弥散过程。河流污染物的迁移原理图如图4所示。

系统动力学以系统分析为基础，强调动态的、复杂的（非线性、反馈回路、延迟和随机性等）系统结构，认为系统的行为模式与特征主要根植于系统内部的反馈结构与机制，因此适合于研究复杂系统的结构、功能与动态行为之间的关系[10]。

根据龙江水污染事故的监测以及镉元素传播途径特点分析，可以建立一维水质模型如式（1）所示;

[∂c∂t=Ex∂2c∂x2−ux∂c∂x−kc]                 （1）

其中，假设在体积微元处[（x）t]时刻的污染物质量浓度为[c]，[ux]为流速，[Ex]为弥散系数，污染物的衰减速率系数为[k]，i和j分别表示当前位置和当前时刻，i-1和i+1分别表示上一位置和下一位置，j-1和j+1分别表示上一时刻和下一时刻，[∆x]表示中间站点到两边站点距离的平均值，因此式（1）变形为：

[cj+1i−cji∆t=Eicji+1−2cji+cji−1∆x2−uicji−cji−1∆x−kicji]   （2）

当污染物为镉这类难溶性的守恒污染物时，设定式（2）中的[ki=0]，则整理可得弥散系数如式（3）所示：

[Ei=（cj+1i−cji∆t+uicji−cji−1∆x）/（cji+1−2cji+cji−1∆x2）]  （3）

2.2弥散系数求解

将龙江水域水力学特征数据带入一维水质模型。本文取糯米滩、西门崖、露塘断面三个测量点的镉浓度差值以及流速[ui]为0.28m/s代入公式（3），来求取各个点的弥散系数[Ex]值，如表2所示。

然后通过多组数据求出弥散系数[Ex]关于镉浓度差[Δc]的函数关系，最终确定龙江一维水质模型参数[Ex]的非线性拟合，R2=0.5153，如图5所示。

从图5可以看出：

l当浓度[Δc]为负值时，求得弥散系数[Ex]也为负值，也就是点落在第三象限;同理，当浓度差[Δc]为正值时，求得弥散系数[Ex]也为正值，也就是点落在第一象限。彌散系数[Ex]的绝对值随浓度差[Δc]绝对值的增大而增大，即弥散程度随浓度差扩大而提高;

l当浓度差[Δc]为0mg/L时，弥散系数[Ex]理论上也为0m2/s;当浓度差Δc为0mg/L时，弥散系数[Ex]为216.4m2/s，是由于污染物浓度的数据测量存在一定误差。

l当浓度差[Δc]为正值时，有的点会取到负值，即有点落在第二象限;当浓度差[Δc]为负值时，有的点会取到正值，即有点落在第四象限。简而言之，若有点落在第二象限或者第四象限，原因是每一个位置到达峰值的时间不同，会导致公式中各个位置之间的浓度差不同，出现正负性不一致的情况。

采用二次拟合方法，得出的龙江水质模型的弥散系数：

[Ex=8313x2+8120x+216.4]               （4）

由于模型受到地理环境等外界不可抗因素的影响，会导致河流流速发生改变，从而影响弥散系数的大小。根据实际的理论分析，可以得出当流速[ui]增大时，会使得弥散系数[Ex']增大;当流速[ui]减小时，会使得弥散系数[Ex'']减小，如图6和图7所示。

采用二次拟合方法，得出流速增大到0.36m/s时龙江水质模型的弥散系数：

[Ex'=1.88×105x2+8.595×104x+2133]   （5）

采用二次拟合方法，得出流速减小到0.2m/s时龙江水质模型的弥散系数：

[Ex''=3850x2+6404x+15.93]             （6）

取浓度差[Δc]的值为0.0341、0.0083、-0.0191、-0.0285、-0.1581，通过其相对的[Ex]、[Ex']和[Ex'']可以看出弥散系数的变化，和[Ex]的绝对值相比，[Ex']的绝对值是增大的，而[Ex'']的绝对值是减小的，如表3所示。

表3   不同流速下的弥散系数变化

[浓度差[Δc]（mg/L） [Ex]（m2/s） [Ex']（m2/s） [Ex'']（m2/s） 0.0341 405.60 794.03 234.24 0.0083 1553.97 1735.30 1372.64 -0.0191 1263.78 1844.07 683.49 -0.0285 1493.81 2537.57 450.06 -0.1581 -4925.68 -7752.26 -2099.10 ]

将模型参数求解出后，为了验证模型是否能在镉污染的实际复杂情况下保持合理与准确性，结合数据，分段模拟污染物传播预测浓度，将预测实现可视化展示，并与实际监测数据进行比对。

当采用拟合[Ex]的数据组时，污染物的实际数值相对于预测值偏低，考虑可能是因为河流底泥对于镉的吸附作用。从整体预测效果来看，污染物浓度的变化趋势、变化速率均能较好匹配，数据比对如图8所示。但是，图8中某些实际数据会出现0值，这是因为实际测量中存在一定测量误差，在预测中模型能将这些误差值计算出来进行补足。

9数据比对（去除0值）

将实际数据中0值点进行去除处理，然后再次进行数据比对，如图9所示。计算实际数据及预测数据之间差距，其中均方差为：[7.53×10−4]，[R2]为：0.867。

3模型可视化

龙江水污染模型的预测系统整体设计满足了系统要求，可视化系统所实现主要功能模块包含了用户登录，读取已知数据，设置参数等模块。具体实现形式如流程图10所示。

系统主要包含了用户登录界面和水污染数值预测界面两大部分，其中用户登录界面如图11所示，水污染数值预测界面如图12所示。

水污染预测系统界面主要通过了GUI的文本框按钮，图形窗，编辑文本，以及静态文本来构成。静态文本主要用作展示各模塊名称，编辑文本用于导入数据，设定参数。点击函数回调按钮会将设定点的预测结果在图形axes窗口中画出图来，达到预测的目的。

4 结论

本文以龙江镉污染的实测数据来建立基于系统动力学的镉污染预测模型，求取了[Ex]，并对[ux]的变化进行定性分析，实现了对洛东电站下一时刻镉污染浓度变化的预测，[R2]=0.867，为突发性的镉污染提供理论预测模型。

参考文献：

[1] 艾恒雨，刘同威.2000—2011年国内重大突发性水污染事件统计分析[J].安全与环境学报，2013，13（4）：284-288.

[2] 刘丹，黄俊.流域突发水污染事件应急能力体系建设[J].人民长江，2015，46（19）：71-74.

[3] 2018年《中国生态环境状况公报》（摘录二）[J].环境保护，2019，47（12）：50-55.

[4] 林立丰，潘尚霞，王凤岩，等.基准剂量法在突发性镉污染饮用水应急阈值研究中的应用[J].环境与健康杂志，2018，35（4）：295-298.

[5] 宫凯悦.松花江哈尔滨段河流底泥重金属污染及内源释放规律研究[D].哈尔滨：哈尔滨工业大学，2014.

[6] 赵华林.城市半封闭河道水体重金属污染分配特征研究[J].四川环境，2019，38（6）：120-124.

[7] 刘潭秋，沈新平，王汉华.基于一种时间序列模型的河流重金属污染浓度预测研究[J].计算技术与自动化，2012，31（3）：29-33.

[8] 童娟.珠江流域概况及水文特性分析[J].水利科技与经济，2007，13（1）：31-33.

[9] 张永祥，蔡德所，吴攀高.广西龙江河沉积物重金属污染特征及生态风险评价[J].中国农村水利水电，2014（6）：50-53.

[10] 王其藩.系统动力学理论与方法的新进展[J].系统工程理论方法应用，1995（4）：6-12.

【通联编辑：唐一东】