基于Matlab的大连地区突发水污染事件预报模型研究及应用

孙大明

(辽宁省大连水文局，辽宁 大连 116023)

水质安全关系到区域社会经济发展，特别是北方山地丘陵地区，基本上面临资源型和污染型并重的缺水状况。大连市境内多山地丘陵，少平原低地，河流均独流入海，属于典型的北方季节性山区河流[1,2]。为了能够及时掌握突发水污染事件的水质污染程度、污染物持续时间以及污染物扩散范围等数据，为后续的应急处理提供详细的数据支撑，很多学者在突发水污染事件模拟和预警方面作了大量研究。钟名军[3]结合水质数学模型和GIS技术，对西江局部河段进行了水质计算和结果评价，并将结果在开发的系统中实现可视化。邵世保等[4]选择EFDC水动力-水质模型作为慈湖河水环境模型，通过CAD、GIS软件，建立二维模型，并将其应用于水质预警，模拟污染事故发生后慈湖河水质变化情况。夏雨晴[5]基于HYDROINFO水利信息系统建立了松花江流域一维水动力及水质模型，并对该流域洪水及水污染事件进行模拟预报，探究污染物在水体中的输移扩散、演变规律和各断面水力特征的变化过程。本文主要在目前研究工作的基础上，基于Matlab软件的强大数值计算、图形显示功能，建立突发水污染事故的预报模型，预测发生突发水污染事件时可能造成的水环境影响，为决策部门启动应急预案和污染事件处理提供科学依据。

1 水质模型的确定

1.1 模型基本方程

本文针对大连地区独流入海季节性河流特点，利用二维非恒定流水质扩散模型对突发性水污染事故的水质进行了数值模拟，模型按照有无岸边限制分列如下。

1.1.1 无岸边限制的河流瞬时点源模型

如果污染源不是稳定连续排放，则任意一点的浓度是时间的函数。如果在t=0的那一时刻，在x=0处向河流中投放了质量为M的污染物，则称为瞬时点源排放。当河面较宽时，水质扩散近似于无岸边限制，污染物在水体中浓度为C(x，y，t)，则无岸边限制瞬时点源水质模型为[6]

(1)

1.1.2 有岸边限制的河流瞬时点源模型

若研究区河流较小，宽度有限，向河流释放总量为M的污染物。污染物扩散到岸边时，会被岸边反射，如镜子般，会呈现另外一个像源，这样污染程度会增加，若发生全发射，则相当于再次污染，且污染强度叠加。

如果河宽为B(B>0)，点源离岸边距离为a(0≤a≤B)，考虑一次虚源反射可得有岸边限制瞬时点源二维扩散模型[7]：

(2)

上二式中：x为预测点距污染物排放点的距离，m；y为预测点距污染物排放口的横向距离，m；t为时间，s；c为预测点(x，y)处污染物的浓度，mg/L；M为单位时间内污染物排放量，g/s；a为污染物排放口离河岸距离(0≤a≤B)，m；h为平均水深，m；Dx为纵向扩散系数，m2/s；Dy为横向扩散系数，m2/s；ux为河流纵向平均流速，m/s；uy为河流横向平均流速，m/s；K为污染物衰减系数，1/d；C0为水体中污染物的本底浓度，mg/L；n为污染物反射次数，当a=0或a=B时，相当于点源在岸边排放，a=B/2，相当于河中心排放。在计算时，如果水体存在本底浓度，需在上述公式计算结束后加上C0exp(-Kt)。

1.2 模型参数的确定

大连山地丘陵地区季节性河流模型需要确定的参数主要有横向扩散系数(Dy)、纵向扩散系数(Dx)、衰减系数K以及流速(u)和水深(H)。

1.2.1 横向扩散系数

天然河流的横向扩散系数表示为Dy，通过经验公式估算法确定，常用的经验式公式为：

Dy=αHu*

(3)

(4)

上二式中：α为横向扩散系数；g为重力加速度，m/s2；I为水力坡降；H为水深,m；h为平均水深，m。

大连地区河流主要属于季节性河流，横向扩散系数会因汛期大流量和非汛期小流量状态下河流的水力坡降不同，而发生较大变化，详见表1。

表1 大连地区主要河流非汛期和汛期水力坡降

1.2.2 纵向扩散系数

河流纵向扩散系数Dx采用Elder经验公式计算:

Dx=5.93Hu*

(5)

式中：H为水深,m；u*为摩阻系数。

1.2.3 断面流速和水深

模型中的河流流量、流速和水面宽度等参数主要基于水文监测站确定，这些水文指标主要用于计算河流污染物在水体中的传输情况。断面流速的确定主要有以下几种方法：

a.对于有实测河流水文资料的断面，其水位、流量和流速的资料比较丰富，通过绘制水位流量关系曲线、水位面积关系曲线，根据公式V=Q/A计算断面的平均流速[8]。

b.上下游有实测资料的断面：

断面没有水文监测资料，但其上游或者下游存在水文站，可以采用曼宁公式利用河道比降、糙率、平均水深推算，最后利用公式V=Q/A求出断面流速。

对于上下游有实测河流水文资料的断面，还可以利用断面流速、水深与断面流量具有正相关关系，采用流量系数法求出流速：

u=aQb

(6)

H=αQβ

(7)

上二式中：Q为流量，m3/s；u为河流平均流速,m/s；H为水深,m。

对式(6)做线性转化：

y=B0+B1x

式中：y=lnU，x=lnQ，B0=lna，B1=b。

据此，断面平均流速与流量的关系就归结为y依x的线性回归方程。

本文利用大连市冰峪沟水文站、沙里涂水文站、茧场水文站、关家屯水文站以及登沙河水文站多年汛期和非汛期流量、流速和水深的实测资料进行相关性分析，取得了大连市主要河流断面流量流速以及水深流速经验关系式，见表2。

表2 大连市河流断面经验公式

c.无资料的河段。对没有实测流量流速资料的河段，可以采用现场实测或类比法推求断面流速和水深。

1.2.3 综合衰减系数

本文在构建水质污染物扩散模型时，综合衰减系数确定主要采用以下两种方法：

经验公式法，可以根据各种经验公式推求综合衰减系数。

实测法，对于河流，选取一个河道顺直、水流比较稳定、无支流和排污口的河段，在河段内设置多个监测断面监测污染物的浓度值，通过分析计算污染物的综合衰减系数。

2 实例应用

2.1 模型应用情境设定

本文研究采用假设模拟的方法对大连市突发水污染事件预报模型进行应用。假设有一瞬时点源向碧流河水库上游茧场水文站断面附近释放1t挥发酚污染物质，模型选择在有岸边限制条件下，考虑一次虚源反射。

表3 汛期、非汛期茧场水文站水文参数

苯酚是一种重要的化学原料，广泛用于农药、杀菌剂、化学药品、医药、合成纤维等行业。本研究确定挥发酚污染物质的降解速率常数K=0，释放点位于河中心，挥发酚污染物的水质标准值为0.005mg/L，本底值C0=0.0025mg/L。表3为汛期、非汛期茧场水文站水文参数。

2.2 基于Matlab软件模拟结果分析

2.2.1 污染物扩散情况模拟

根据上述模型对碧流河水库上游污染物浓度变化的实际情况进行模拟。同时借助Matlab软件强大的计算与图形界面功能得到模拟结果，并取得了良好的模拟效果。输入参数界面见图1。

图1 输入参数界面

图2为运用Matlab软件模拟汛期挥发酚类物质污染发生后，5min、10min、20min和30min时，水中的挥发酚浓度三维空间变化图。图2显示，污染团在降解的过程中随着水流向下游漂移，污染范围有所增加，但浓度也在不断降低。同时，用该软件模拟了非汛期挥发酚浓度变化趋势，得出的结论是：非汛期污染物扩散情况与汛期类似，区别在于汛期污染物浓度降低与扩散范围远大于非汛期。由此可知，污染物漂移速度与水流流速有关，水顺流方向的流速是导致污染团向下漂移的决定性因素。而污染范围的增加，主要与横向和纵向扩散系数有关。浓度的降低是污染物自身降解与水体的自净能力直接导致的。通过不断模拟分析可知，时间不断增加，水体中污染物浓度不断减小。

2.2.2 特定点污染物浓度随时间变化曲线

基于Matlab软件构建突发水污染事件预报模型，点击参数输入，选择污染物浓度随时间变化曲线。输入预报模型相关参数后，点击结果按钮，得到特定的水质控制点污染物扩散情况。某时刻该点污染物浓度超过该污染物的水质标准值即可认为该点受到了事故危害的影响，此时预警开始。本文模拟污染点下游1000m处应急入连取水口处的污染物变化情况。运用软件模拟可以得到应急入连取水口处污染物超标，即开始预警的时间，污染物浓度达到峰值的时间和最大浓度，以及污染预警解除的时间，并分汛期以及非汛期两种情况讨论，见图3。

由图3可知，汛期挥发酚超过水质标准值的起始时刻为污染事故发生后的498s，即开始预警；应急入连取水口处挥发酚浓度达到峰值的时间为854s，最大浓度为2.14万mg/L；应急入连取水口处挥发酚浓度重新回到标准值之下的时间为1469s，即预警解除。挥发酚污染持续时间971s。

非汛期应急入连取水口处挥发酚浓度预警开始时间为3954s；预警解除时间为6307s；污染持续时间2353s，挥发酚浓度峰值出现的时间为4993s，峰值为121万mg/L。相比汛期，非汛期应急入连取水口处挥发酚污染持续时间延长1382s，预警开始时间延后3456s，挥发酚浓度峰值出现时间延后4139s。

图3 应急入连取水口处挥发酚浓度随时间变化曲线

3 结 论

a.基于二维非恒定流水质扩散模型构建大连山地丘陵地区季节性河流水质扩散模型，并对模型主要参数横向扩散系数、纵向扩散系数、衰减系数以及流速和水深等的确定方法进行综述。

b.基于Matlab软件对大连地区水质扩散模型进行研究，并以挥发酚为污染物，模拟大连地区突发水污染事件污染物汛期和非汛期扩散情况。结果表明：汛期时，随着时间的变化，污染团在降解的过程中随着水流向下游漂移，污染范围有所增加，但浓度也在不断降低；非汛期污染物扩散情况与汛期类似，区别在于汛期污染物浓度降低与扩散范围远大于非汛期。

c.基于Matlab软件模拟污染源下游1000m处应急入连取水口处污染物变化情况，结果表明：该处污染物浓度变化为先增加后减少，并存在最大峰值，由此可得到此处污染物预警的开始与解除时间，并且相比汛期，非汛期应急入连取水口处挥发酚污染持续时间延长1380s，预警开始时间延后3456s，挥发酚浓度峰值出现时间延后4140s。

本文基于Matlab软件对构建的山地丘陵地区季节性河流发生的水污染事故水质扩散模型进行编程，进一步建立水污染动态时空可视化展示方法，通过研究，能够为水资源管理部门及时了解汛期以及非汛期的水中污染物扩散变化情况，提供数据基础和科学依据，显著提高相关部门的应急处理效率，为水质安全提供技术保障，最大限度避免或减少突发水污染事件造成的损失。