河流交汇区污染物扩散数值模拟研究

杨志山，米潭，姚建，张弢

(1.四川大学 建筑与环境学院，四川 成都 610065；2.成都市生态环境科研监测所，四川 成都 610072)

河流交汇现象普遍存在于自然的河流水系中，交汇河道的边界岸线、河床地形等具有不断变化的特点[1]，导致交汇区污染物的扩散规律难以明晰。不同的交汇流量是影响污染物扩散的主要因素[2]，研究不同交汇流量下污染物的扩散特征，对工程实践中涉及的污水排放方式、排污口位置的选择等具有重要指导意义。目前，水质数值模型逐渐成为研究水体中污染物扩散规律的主要工具[3]，其中，平面二维自由表面流MIKE 21 FM模型在河流、湖泊、水库的污染物扩散研究中已得到广泛应用。如：张守平等[4]通过MIKE 21 FM模型研究了入河排污口对水体水质的影响，结果表明该水质模型能准确反映污染物在水体中的扩散规律；贾瑞鹏等[5]建立了万宝湖水动力水质二维耦合模型，分析了万宝湖各类污染物的迁移转化机制；杨晨等[6]考虑风场对污染物扩散的影响，模拟了汾河水库11种污染物的迁移扩散规律。

关于河流交汇区不同汇流比情况下污染物扩散特征方面的研究，多为建立水槽物理模型进行探讨。如：魏娟等[7]选用水气两相流数学模型，通过物理实验，对河流交汇区不同流量情况进行实验和数值模拟研究，结果表明污染带的狭长程度受流量比的影响较大；袁航[8]通过建立45°支流斜接干流型的交汇水槽物理模型，探究了交汇河道不同汇流比对污染物混合特性的影响，结果表明水流交汇区会产生二次环流，对污染物混合有很大的影响；HAN S等[9]对室内交汇明渠的污染物混合过程进行了数值模拟，发现交汇口附近污染物的混合速度主要由上游来流流量决定，来流流量较小时混合速度较快。但是，上述研究中物理实验模型往往将河道简化为顺直河道，难以反映自然河流交汇区边界、地形等复杂多变的特点。

为此，本研究基于平面二维自由表面流MIKE 21 FM模型，充分考虑自然河流交汇区的特点，选取中国西南中小河流石亭江与其支流鸭子江的交汇区，研究交汇区COD、氨氮质量浓度在不同汇流比情况下的扩散规律，为交汇区污染物扩散机理研究提供借鉴。

1 模型的建立与验证

1.1 模型建立

平面二维自由表面流MIKE 21 FM模型的控制方程为基于Boussinesq假定和流体静压假定的二维不可压缩雷诺平均N-S方程。建立的模型包括二维水动力学模型MIKE 21 HD和水质模型MIKE 21 AD。

网格划分：本研究计算区域范围为东西向2.00 km、南北向1.95 km。采用无结构三角形网格对研究区域进行网格划分，最后共生成网格1 962个，计算节点1 300个。在搭建计算网格的基础上，利用实测天然地形，对网格地形进行插值，得到计算区域的地形模拟图(图1)。

图1 地形模拟图

边界设置：本研究计算区域共设置3个开边界：交汇区上游约1.0 km处石亭江省控断面L设为上边界1；交汇区上游约800 m处鸭子江国控断面S设为上边界2；交汇区下游约1.2 km处石亭江出境断面W设为下边界。上边界采用实测的流量数据作为控制条件，下边界采用实测的水位数据作为控制条件。

断面设置：为方便模型验证及污染物扩散特征分析，交汇区域共设置6个断面：Z1、Z2、Z3共3个横向断面，反映污染物横向扩散特征；C1、C2、C3共3个纵向断面，反映污染物纵向扩散特征。边界及断面设置如图2所示。

图2 边界及断面设置图

计算条件设定：①初始条件：二维水动力模型选取2015年1月1日的水动力学指标实测资料；水质模型选取该日各监测断面COD、氨氮的质量浓度实测资料。②边界条件：二维水动力模型选取2015年1月1—30日L、S断面的实测流量过程及W断面的实测水位过程；水质模型选取2015年1月1—30日L、S及W断面COD、氨氮的质量浓度实测资料。

1.2 模型验证

采用2015年1月1—30日Z1、Z2、Z3断面的流量及COD、氨氮质量浓度的模拟值，分别计算其平均绝对误差(Mean Absolute Error，MAE)、决定系数(R2)、均方根误差(Root Mean Square Error，RMSE)、平均相对误差(Mean Relative Error，MRE)4种拟合优度指标来评价模型的性能。通过查阅文献与模型调试，确定石亭江、鸭子江汇流前河床糙率为0.03，汇流后河床糙率为0.05；涡黏系数采用Smagorinsky公式估算，取值为0.28 m2/s；干水深为0.005 m，淹没水深为0.05 m，湿水深为0.1 m；风阻系数为0.001；污染物水平扩散系数为0.05 m2/s；COD、氨氮质量浓度的衰减系数分别为0.5、0.1 d-1。

Z1至Z3断面上的流量、COD和氨氮质量浓度的MAE、R2、RMSE、MRE值见表1。由表1可知，各监测指标在不同断面处的MAE、R2、RMSE值均比较接近，说明各监测指标测量值与模拟值的误差在同一数量级；各监测指标的R2值均在Z3断面处最大，MRE值在Z3断面处均为最小，表明Z3断面处测量值与模拟值拟合最好。MRE值越小，模型的精度越高，本文的MRE值均在30%以内，说明用实测值标定的MIKE 21 FM模型可以保证精度，满足模拟需求。

表1 各变量测量值与模拟值比较的指标计算

2 结果与分析

采用本文建立的MIKE 21 FM模型对研究区域不同汇流比情况下的COD、氨氮扩散特征进行数值模拟，模拟工况设置见表2，其中Q在3种工况下分别选取石亭江、鸭子江在枯水期、平水期及丰水期时的平均流量，汇流比R=Q鸭子江/Q石亭江，模拟初始污染物质量浓度C选取2015年1月1日L、S断面实测的COD、氨氮质量浓度。

表2 数值模拟工况参数

2.1 污染物浓度场分析

采用水质模型MIKE 21 AD对不同工况下的COD、氨氮质量浓度进行模拟计算，获得3种工况下COD、氨氮的质量浓度场如图3所示。根据《地表水环境质量标准》(GB 3838—2002)[10]中Ⅲ、Ⅳ类水质的COD、氨氮质量浓度范围，得到COD、氨氮质量浓度的分布面积见表3。

图3 不同工况下COD、氨氮质量浓度场

表3 COD、氨氮质量浓度分布面积统计

由图3可知，在石亭江和鸭子江交汇前，两条河流的污染物质量浓度均逐渐降低。石亭江COD的质量浓度在不同汇流比下均大于22 mg/L，汇流比越大，COD质量浓度分布在22～24 mg/L范围内的面积越小。鸭子江COD质量浓度变化受汇流比变化的影响较大，当R=0.33时，COD质量浓度自上游至下游依次分布在18～20 mg/L、16～18 mg/L的范围内，在交汇口附近COD质量浓度低于16 mg/L；当R=0.50时，COD质量浓度分布在18～20 mg/L的范围内，在交汇口附近分布在16～18 mg/L的范围内；当R=0.67时，COD质量浓度仅分布在18～20 mg/L范围内。由此可见，在两条河流交汇前，汇流比越大，污染物越不易降解。石亭江、鸭子江氨氮质量浓度在交汇之前随汇流比变化不明显，均分布在1.4～1.5 mg/L、低于1.0 mg/L的范围内。在两条河流交汇后，污染物混合情况各有不同，当R=0.33时，COD质量浓度自上游至下游依次分布在20～22 mg/L、18～20 mg/L、16～18 mg/L的范围内；当R=0.50时，COD质量浓度分布依次为20～22 mg/L、18～20 mg/L的范围内；当R=0.67时，COD质量浓度分布依次为22～24 mg/L、20～22 mg/L、18～20 mg/L的范围内。可见，在相同的混合河段距离内，汇流比越小，污染物混合质量浓度也越低。不同汇流比下，在交汇后氨氮质量浓度均分布在1.4～1.5 mg/L、1.2～1.4 mg/L、1.0～1.2 mg/L的范围内。由于初始浓度CM>CN，交汇后，受石亭江影响，右岸(鸭子江侧)污染物的质量浓度沿程呈上升趋势，河流左岸至右岸污染物的质量浓度逐渐降低。

研究区域所在水功能区水质管理目标为Ⅲ类。由表3可知，随着汇流比的增加，COD、氨氮质量浓度超过《地表水环境质量标准》Ⅲ类水(COD质量浓度>20 mg/L，氨氮质量浓度>1 mg/L)标准部分在研究区域内形成的污染带面积明显增加，当R=0.33时，COD形成的污染带面积为159 844 m2，氨氮形成的污染带面积为351 155 m2；R=0.50时，COD形成的污染带面积为207 650 m2，增加30%，氨氮形成的污染带面积为366 708 m2，增加4.4%；当R=0.67时，COD形成的污染带面积为309 728 m2，增加93.8%，氨氮形成的污染带面积为388 465 m2，增加10.6%。可见，汇流比越大，污染带分布越广，适当降低汇流比，对交汇区水质可起到一定的改善作用。

2.2 汇流区污染物扩散特征分析

图4为河流交汇后COD、氨氮质量浓度的分布情况，取不同横向、纵向断面分析其污染物扩散特征。

图4 汇流区COD、氨氮质量浓度场

2.2.1 交汇区污染物纵向扩散特征分析

图5为不同汇流比下COD、氨氮质量浓度在同一纵向断面上随不同横向断面的变化情况，选取C1、C3断面进行纵向扩散分析。

由图5(a)可知：由于石亭江初始污染物的质量浓度较高，COD、氨氮沿C1纵向断面的质量浓度逐渐下降；当R=0.33时，Z1至Z2断面、Z2至Z3断面COD的质量浓度分别降低14.2%、10.4%，氨氮质量浓度分别降低12.3%、6.7%；当R=0.50时，COD质量浓度分别降低11.7%、7.4%，氨氮质量浓度分别降低11.4%、5.0%；当R=0.67时，COD质量浓度分别降低11.3%、3.0%，氨氮质量浓度分别降低11.0%、3.4%。由图5(b)可知：由于鸭子江初始污染物的质量浓度较低，受石亭江污染物扩散的影响，COD、氨氮沿C3纵向断面的质量浓度逐渐升高；当R=0.33时，Z1至Z2断面、Z2至Z3断面COD质量浓度分别增加6.5%、2.4%，氨氮质量浓度分别增加20.2%、8.0%；当R=0.50时，COD质量浓度分别增加5.5%、1.7%，氨氮质量浓度分别增加18.0%、5.0%；当R=0.67时，COD质量浓度分别增加4.9%、1.5%，氨氮质量浓度分别增加15.4%、4.2%。根据C1、C3断面污染物纵向质量浓度变化情况，当汇流比分别为0.33、0.50、0.67时，COD纵向质量浓度的梯度变化分别为8.4%、6.6%、5.2%，氨氮纵向质量浓度的梯度变化分别为11.8%、9.9%、8.5%。模拟发现，在交汇后未出现COD明显的混合均匀位置；当R=0.33时，氨氮在1 820 m处均匀混合；当R=0.50时，氨氮在1 960 m处均匀混合；当R=0.67时，氨氮在2 050 m处均匀混合。可见，汇流比越大，污染物纵向质量浓度梯度越小，在交汇区出现的污染带越长，出现混合均匀的位置越远离交汇口。

图5 不同汇流比下交汇区污染物质量浓度纵向变化

2.2.2 交汇区污染物横向扩散特征分析

图6为不同汇流比下COD、氨氮质量浓度在同一横向断面上随不同纵向断面的变化情况，选取Z1、Z3断面进行横向扩散分析。

由图6(a)可知，C1至C3断面污染物质量浓度逐渐降低：当R=0.33时，C1至C2断面、C2至C3断面COD质量浓度分别降低13.6%、5.3%，氨氮质量浓度分别降低14.8%、7.8%；当R=0.50时，COD质量浓度分别降低14.1%、6.6%，氨氮质量浓度分别降低19.3%、10.7%；当R=0.67时，COD质量浓度分别降低15.2%、8.0%，氨氮质量浓度分别降低21.7%、11.1%。由图6(b)可知，C1至C3断面污染物质量浓度也逐渐降低：当R=0.33时，C1至C2断面、C2至C3断面COD质量浓度分别降低3.4%、1.8%，氨氮质量浓度分别降低2.3%、1.6%；当R=0.50时，COD质量浓度分别降低4.6%、3.2%，氨氮质量浓度分别降低3.1%、2.4%；当R=0.67时，COD质量浓度分别降低5.7%、4.1%，氨氮质量浓度分别降低3.7%、3.1%。根据Z1、Z3断面污染物横向质量浓度变化情况，当汇流比分别为0.33、0.50、0.67时，COD横向质量浓度的梯度变化分别为6.0%、7.1%、8.3%，氨氮横向质量浓度的梯度变化分别为6.6%、8.9%、9.9%。模拟发现，横向断面越靠近交汇口，污染物横向质量浓度梯度越大，汇流比的增加会增大污染物的横向质量浓度梯度，尤其对于交汇口附近断面的影响更为明显，故而会加剧污染物在交汇口的横向扩散趋势。

图6 不同汇流比下交汇区污染物质量浓度横向变化

3 结论

较大的汇流比会导致污染带分布变广、污染物纵向质量浓度梯度变小、纵向扩散距离变长，同时，汇流比的增加会增大污染物横向质量浓度的梯度变化，加剧污染物在交汇口的扩散趋势；MIKE 21 FM模型可较好地反映自然河流交汇区的特点，可应用于自然河流交汇区污染物扩散机理的研究；其他影响自然河流交汇区污染物扩散规律的因素如流速、污染物质量浓度等有待开展深入研究。