基于数值模拟模型研究滨海化工企业地下水环境影响预测评价

，，

(南京大学环境规划设计研究院股份公司，江苏 南京 210093)

随着我国经济的快速发展，沿海地区加大了滨海地区的开发强度，化工、医药、农药等行业在滨海地区迅猛发展，人工填海地区成为该行业重要的选择[1-3]。滨海地区的过度开发导致污染问题日益突出，并逐渐成为了建设生态文明和谐社会的一个主要矛盾。尤其是对于地下水环境而言，地下水污染具有隐蔽性、滞后性以及一旦污染便难以清除的特点，因此备受关注[4]。而针对建设项目的地下水环境影响预测评价，能预先评估和预防项目可能造成的地下水污染，起到源头防治和管理的作用，对保护地下水环境意义重大[5-7]。此外，滨海地区地下水、地表水以及海水之间的水力联系密切且复杂，滨海地区地下水环境评价成为项目环评中的难点。这类项目地下水环境评价研究较少，目前可参考的文献非常有限。本文以某化工企业为例，基于数值模拟模型研究滨海地区地下水环境影响评价。

1 研究思路与背景

1.1 研究思路

通过资料收集和分析，结合现场水文地质环境调查进行综合研究，确定合理的模拟范围；分析研究区地下水流动特征、含隔水层系统结构和边界条件等，建立研究区水文地质概念模型；通过对研究区进行空间离散、水文地质参数赋值、地下水流场识别与校正，构建研究区地下水水流模型；研究项目污染源特征，确定合理的污染因子，设计地下水泄露污染情景，在地下水水流模型的基础上耦合污染物运移模型，对研究区地下水进行数值模拟预测，选择合理的水质标准，针对模拟结果开展评价[8-11]。

1.2 研究背景

1.2.1 研究区概况

研究区为南通某沿海经济开发区，累计落户项目100多个，形成集医药、农药、高分子材料为主导的三大产业板块群。开发区东枕黄海，占地面积约9.7km2，区位条件得天独厚。研究区地貌特征为典型的滨海平原，场地较平缓开阔，地形坡度3°以内，自然地面标高3.79～5.69 m(黄海高程系)，总体呈西南高东北低状。

1.2.2 水文地质条件

本区大地构造处于扬子准地台东部，地层属于扬子地层区。地表第四系广泛覆盖，无基岩出露。区内地下水类型主要为松散岩类孔隙水，具有分布广、层次多、水量丰富，水质复杂等特征。区内潜水含水层地层岩性主要为粉土、粉砂，粉土层厚1.35～1.60 m，平均1.51 m，粉砂层厚19.65～25.60 m，平均22.61 m；潜水含水层下隔水层层厚7.00～17.00 m，平均12.06 m，岩性主要为粉质粘土，渗透性较弱。区内浅层地下水埋深较浅，一般在1.55～2.83 m左右，地下水流向主要由西南向东北，研究区地下水流场图见图1。

图1 研究区地下水流场图

2 数值模型

2.1 水文地质概念模型

研究区水文地质条件较单一，主要含水层为第四系松散岩类含水层，含水层分布较均匀，在潜水含水层和承压含水层之间有一层良好的粉质粘土隔水层，使得潜水含水层与承压含水层之间水力联系较弱。

结合研究区水资源分区、水系分布，考虑区域地质、水文地质、环境水文地质条件，研究区西部以振洋河为河流界，北部以黄海为水头界，东部以匡河为河流界，南部以一内陆河为河流界，整个研究区范围面积约16 km2。含水层上边界为地面，通过该边界含水层系统与大气降水、地表水等产生垂向上的水量交换；下边界为透水性差的粉质粘土弱透水层，该层阻断了潜水含水层与下伏承压含水层的水力联系，故定义为隔水边界。根据模拟区地层条件，污染进入地下主要污染潜水含水层。

2.2 数学模型

根据水文地质概念模型，评价范围内地下水流运动的数学模型可以表示为潜水含水层非均质、各向异性三维非稳定流数学模型，其控制方程及定解条件如下：

H(x，y，z，0)=H0(x，y，z)(x，y，z)∈Ω

H(x，y，z，t)|Γ1=H1(x，y，z，t)

式中：Kxx,Kyy,Kzz为主坐标轴方向多孔介质的渗透系数，[LT-1]；h为水头，[L]；w为单位面积垂向流量，[LT-1]，用以表示源汇项；μ为多孔介质的给水度(或饱和差)；z为潜水含水层的底板标高，[L]；t为时间，[T]；Ω表示渗流区域；Γ1表示第一类给定水头边界。

污染物在地下水中的运移包括对流、弥散以及溶质本身的物理、化学变化等过程，可表示为：

本文采用GMS软件求解，用MODFLOW计算模块求解地下水水流运动数学模型，用MT3DMS模块求解地下水污染物运移数学模型。

2.3 模型离散

本文数值模拟求解采用有限差分法，需将评价区剖分为四边形网格，在厂区污水处理站处进行加密处理，最小网格空间长度达到2 m，网格垂向上剖分依据场区建设特点以及评价区内潜水含水层特征划分为两层：第一层为粉土，厚度1.5 m左右；第二层为粉砂，含水层厚度22.6 m左右；在含水层下部设置一层相对隔水层，厚度为12.1 m左右；整个模型在垂向上一共三层。一共剖分148 512个网格，其中有效网格为89 733个。详细见图2。

图2 研究区三维网格剖分示意图

2.4 水文地质参数

水文地质参数的选择主要依据该项目水文地质调查所进行的各类野外及室内水文地质试验的结果，并结合以往在附近区域各类水文地质试验数据资料确定，其水文地质参数取值如表1所示。

表1 研究区水文地质参数取值

弥散度是研究污染物在地下水中迁移转化规律的重要参数之一，由于弥散度的不确定性，本文将通过对比不同弥散度取值条件下污染物迁移规律，帮助企业科学管理地下水。

2.5 模型校正与识别

将数值模型进行求解，对比模型计算结果与实际观测数据差异程度，从而校正检验模型的合理性、准确性和稳定性。

2.5.1 地下水水位拟合

模拟计算含水层地下水水位与实测地下水水位关系见图3。图中可以看出实际观测水位与计算水位误差均在0.2 m以内，模拟误差小，说明模型是合理的。

图3 研究区地下水水位拟合图

2.5.2 水均衡

模拟计算得到的模拟范围内水均衡结果如表2所示。

表2 研究区水均衡结果表 m3/a

根据水均衡结果，评价区每年地下水排泄进入地表水1 757 080.8 m3，地表水补给地下水的量为785 582.2 m3，表明地下水和地表水存在密切的水力联系。综上，根据对地下水水位及水均衡计算结果的分析，模型能较好反映该地区地下水流运动特征，可以用于地下水环境影响的预测评价。

3 地下水模拟预测与评价

3.1 地下水污染源

建设项目工程可能存在的地下水污染源主要是各生产装置产生的生产废水、各生产废水的输送管线、污水处理池等，根据对本次改扩建项目的生产工艺、设备设施分析后，项目装置区内的污水处理池防治措施较薄弱，也是最重要的地下水污染源。非正常状况下污水处理池池壁因系统老化、防渗措施失效等原因发生渗漏可能对地下水产生影响。本工程排放的废水主要产生于噻虫胺、乙螨唑、虫胺等原药生产过程中的废水。污水处理池废水源强进水情况见表3。

表3 污水处理池进水源强情况表

非正常状况下，污水处理池发生渗漏，废水经包气带进入潜水含水层。挑选污水处理池调节池底部及侧壁面积约为220 m2，渗漏面积按面积的5‰计算，根据《给水排水构筑物工程施工及验收规范》(GB50141-2008)，钢筋混凝土结构水池渗水量不得超过2 L/(m2·d)，非正常状况按照正常状况的100倍考虑，则污水处理池调解池渗水量为0.22 m3/d。

3.2 预测因子及时段

项目地下水的主要污染源来自污水处理池中的废水，根据对项目废水的主要成分分析，确定废水中主要种类和标准值属较大的COD和甲苯作为预测的对象，此类物质有较大危害性，如进入地下水易对环境造成危害，因此选取该指标进行地下水溶质迁移预测。

考虑项目建设、运营和退役期，将地下水环境影响预测时段拟定为10 000 d。结合工程特征与环境特征，预测污染发生100 d、1 000 d及10 000 d后污染物迁移情况。

3.3 评价标准

COD标准限值参照《地下水质量标准》(GB/T14848-93)Ⅲ类标准值[12]，甲苯超标范围参照《地下水质标准》(DZ/T0290-2015)Ⅲ类标准值[13]。因COD与高锰酸钾指数存在换算系数，根据南通地区地表水监测数据显示换算系数在2.5～3之间，为保守起见，本文COD浓度与高锰酸钾指数按4倍进行换算。各指标具体情况见表4。

表4 地下水污染物评价标准

图4 不同时刻下COD运移分布图(D=10m)

3.4 预测结果分析

3.4.1 弥散度为10 m时，污染物迁移规律

将污染物泄露量和水文地质参数带入校正模型进行模拟预测，弥散度为10 m时，100 d、1 000 d、10 000 d时刻下COD和甲苯运移规律见图4和图5。

从图中可以看出，随着时间不断增加，污染物在中心处的浓度逐渐升高，往东北方向不断迁移，距离逐渐增加，污染范围不断扩大；污染物详细迁移规律见表5。

表5 D=10 m时不同时刻污染物运移特征表

3.4.2 不同弥散度条件下污染物迁移规律

溶质运移模型的弥散度参数具有“尺度效应”(随研究尺度的增大而增大)，因此弥散度具有很大的不确定性。为了更好的了解污染物可能的运移规律，可以通过分析不同弥散度条件下污染物迁移规律，得到污染物迁移距离、污染范围、污染物浓度等参数，从而方便企业后期科学管理判断可能存在的污染现象。本文选取弥散度为2 m、10 m、50 m进行模拟预测，COD和甲苯模拟结果分别见图6和图7。

图5 不同时刻下甲苯运移分布图(D=10m)

图6 不同弥散度下COD运移分布图(T=10 000 d)

图7 不同弥散度下甲苯运移分布图(T=10 000 d)

图8 不同弥散度下污染物中心浓度变化趋势图

图9 不同弥散度下迁移距离与超标范围变化趋势图

模拟结果显示，模拟时间为10 000 d时，随着弥散度的不断增加，污水处理池中心浓度逐渐下降，迁移距离和超标范围不断扩大，见图8和图9。弥散度为2 m、10 m、50 m时，污染物COD和甲苯中心处浓度与初始浓度的比值从分别为0.97、0.94、0.90和0.96、0.91、0.87，迁移距离分别为32.39 m、52.12 m、88.24 m和36.89 m、58.38 m、102.49 m，超标范围分别为2 263.3 m2、4 007.4 m2、10 080.3 m2和2 585.5 m2、4 763.1 m2、12 851.1 m2。COD比甲苯的浓度下降的更慢，而甲苯比COD的迁移距离和超标范围更大。

发生这一现象的主要原因是，弥散度的增加有利于污染物的扩散，同时随着对流作用，污染物越容易向外扩散，从而导致迁移距离不断增加，污染物扩散范围及污染范围变大；与此同时，中心处的污染物浓度随着扩散作用的增强会逐渐下降。

由于弥散度的尺度效应，在目前无法有效解决尺度效应情况下，通过模拟不同弥散度下污染物迁移规律，确定弥散度对溶质运移的影响。与确定的弥散度模型相比，本文可以更好的弄清污染物可能的迁移距离和超标范围，有利于企业加强地下水管理和判断可能发生的污染现象。

4 结语

(1)本文在对滨海地区水文地质条件详细了解和调查的情况下，构建了地下水数值模型，选择合理的水文地质参数校正模型，所建立的模型能够真实地反应研究区实际的水文地质条件。

(2)针对本项目特点，设计了合理的污水处理池泄露情景，分析了非正常工况下污染物对地下水的影响，模拟预测了污染物的时空迁移与扩散规律，结果表明随时间的增加，中心处污染物浓度逐渐升高，迁移距离逐渐增加，污染范围不断扩大。

(3)由于弥散度的尺度效应，研究了不同弥散度条件下污染物迁移规律，结果表明弥散度的不断增加，污水处理池中心浓度逐渐下降，迁移距离和超标范围不断扩大。通过模拟不同弥散度下污染物迁移规律，有利于企业加强地下水管理和判断可能发生的污染现象。

[1]薛禹群，谢春红，吴吉春，等．山东龙口-莱州地区的海水入侵[J]．地质学报.1992.66(3):280-291.

[2]薛禹群，谢春红，吴吉春．海水入侵研究[J]．水文地质工程地质.1992.19(6):29-33.

[3]左锐，王金生，杨洁，等．滨海石化项目地下水环境影响评价的关键问题[J]．水文地质工程地质.2010.37(3)：97-101.

[4]余期冲，祝晓彬，吴吉春，等．死端孔隙对溶质运移影响的实验研究[J]．水文地质工程地质.2017.44(4)：160-164.

[5]卢晓华．基于数值模拟的企业地下水重金属污染的环境影响预测评价[J]．安全与环境工程.2014.21(1)：93-97.

[6]郑玉虎，吴明洲，徐爱兰，等．考虑土壤吸附作用的地下水污染物运移特征研究[J]．地下水.2017.39(3)：4-7.

[7]岳强，范亚民，耿磊，等．地下水环境影响评价导则执行过程中遇到的问题及建议[J]．环境科学与管理.2012.37(10):174-177.

[8]薛禹群，吴吉春．地下水动力学[M]．北京:水利水电出版社．2009．

[9]王焰新．地下水污染与防治[M]．北京：高等教育出版社.2007.

[10]朱学愚，钱孝星．地下水环境影响评价的工作要点[J]．水资源保护.1998(4):48-53.

[11]余期冲，祝晓彬，吴吉春，等．场地地下水源热泵系统抽灌实验及数值模拟研究[J]．地下水.2015.37(3)：24-27.

[12]GB/T14848-93 地下水质量标准 [S].

[13]DZ/T0290-2015 地下水质标准 [S].