污染物在包气带中迁移的HYDRUS-1D预测模型—以某焦化项目为例

邓强伟

(中国辐射防护研究院环境工程技术研究所， 太原 030006)

1 前 言

在自然界的水循环过程中，包气带是大气降水、地表水与地下水连接的纽带。无论对渗流本身，还是作为溶质迁移的载体，包气带水在水资源、农田灌溉、排水、生态平衡、工程地质和环境问题的研究中都占有重要地位[1]。在地下水系统中，污染物对地下水的危害主要是以“污染源-包气带-地下水”的途径污染地下水，污染物经过包气带(非饱和带)进入饱水带，在地下水流的作用下运移。由于地下水更新速度慢，埋藏情况复杂，一旦发生污染，依靠其自净能力很难降解污染物，需要投入大量的时间及成本进行治理[2]。因此，研究包气带水的污染过程并作出相应的防范措施就显得极其重要。

入渗过程是非饱和土壤水分的运动过程，属于广义渗流理论的研究范畴，其基础为法国工程师Darcy提出的达西定律[3]。自从Gardener和Bresler在对土壤与溶质间互相作用的研究过程中，依据菲克定律研究建立了一维土壤溶质运移方程[4]后，各种对非饱和土壤水分运动过程的研究层出不穷，后人在此基础上推导出各种土壤水入渗的基本公式和各类土壤水入渗过程的模型。其中HYDRUS模型被应用于分析水流和溶质在非饱和多孔隙媒介中的运移过程，它是用土壤物理参数模拟水、热及溶质在非饱和带水中运动的有限元计算机模型[5]。由于HYDRUS模型所需要的输入数据相对较少，模拟结果可靠，能够较为全面地对水分及溶质的运移进行模拟，在非饱和土壤水分运动过程的研究中发挥了很大的作用。

国内学者是从80年代起开始对水盐运移的理论及过程进行的研究，最初开始进行土壤水盐运移的模拟的是张蔚臻[6]。张效先[7]根据质量守恒方程推出了一种计算模拟溶质运移中水动力弥散系数的方法，黄康乐[8]按照质量守恒方程和利用瞬时剖面法求解非饱和的的思路提出了测定及推导非饱和弥散系数的方法。各学者将HYDRUS模型运用到土壤水分运移及溶质迁移的各过程及其各参数的求解过程中，并将研究结果推广到水利、农业等行业的实际运用中。

随着地下水环境影响评价的推广与发展，越来越多的项目开始关注到包气带对污染物运移的阻滞作用，预测特征因子在包气带中的迁移过程。因此，可以利用HYDRUS-1D软件来预测包气带对特征因子的影响过程，计算出特征因子通过包气带到达地下水面的时间和初始浓度值，避免了之前经常将污染物作为源强直接进入地下水造成的影响，更加科学准确。

2 材料与方法

2.1 Hydrus-1D软件

Hydrus系列软件是由美国农业部、美国盐碱实验室等机构在SUMATRA、WORM及SWMI等模型的基础上创建发展而来，模拟宏观及微观尺度上饱和及非饱和介质中一维水流、溶质、热和二氧化碳运移和反应的软件。它可以通过网格剖分精确刻画介质结构，根据反应物质及边界条件的不同，具有多种水流、溶质反应模型供用户选择，在农业、水利、环境学、微生物学等领域都得到了广泛应用[9]。近年来该软件主要应用于农业领域或室内试验模拟，如计算田间氮素流失及氮转化[10]、农药迁移转化[11]，灌溉水量周期性变化[12]，重金属离子运移[13]等。

Hydrus系列软件分为一维、二维和三维三种，分别命名为HYDRUS-1D、HYDRUS-2D和HYDRUS-3D，由于非饱和土壤水主要是一维垂向运动形式，因此HYDRUS-1D的应用非常广泛。

2.2 项目概况

研究目标位于山前倾斜平原区，根据项目厂址的岩土工程勘察报告，结合厂区及其周边村庄水井的柱状资料，项目区出露地层为第四系上更新统，下伏第三系、二叠系、石炭系及奥陶系等地层。地层岩性由新到老主要为：浅部第四系粉土、粉质粘土、细砂层和第三系红色粘土层；中部二叠系与石炭系的砂岩与泥岩层；下部为奥陶系石灰岩。地下水按赋存条件可分为第四系松散层孔隙水含水层和奥陶系岩溶裂隙水含水层两类含水层。

2.3 研究方法

2.3.1 概念模型

2.3.1.1 预测目标及其划分

本次的研究目标为该项目酚氰废水处理站调节池，根据其所在处的岩土工程勘察钻孔资料和厂址北侧村庄饮用水井柱状图，调节池底部0～20m为粉土层，20～30m为粉质粘土层，30～40m为细砂层，地下水水位埋深40m。

因此，项目区包气带厚度40m，包气带岩性从上至下主要为粉土、粉质粘土和细砂层。垂向上按50cm一格剖分，将包气带剖分为80格，剖分结果见图1。

图1 包气带剖分结果图Fig.1 Sectional profile of aeration zone

2.3.1.2 边界条件概化

将包气带水流概化为垂向一维流。调节池在非正常工况下发生渗漏，污染物随污水不断地渗入包气带，污染物在调节池破裂处浓度最高。HYDRUS-1D只考虑污染物在非饱和带的一维垂直迁移，因此模型的边界只有上边界和下边界，上边界为调节池场地的底断面，下边界为包气带底部(即潜水面)。

模型中的水流模拟采用经典Richards 方程[14]来描述水分运移的过程，模型上边界设定为大气边界，下边界设定为自由排水边界。而溶质运移模型采用对流-弥散方程，模型上边界设定为定浓度边界，下边界设定为零浓度边界。本次模拟的概念模型见图2。

图2 HYDRUS-1D概念模型Fig.2 Concept model of HYDRUS-1D

2.3.2 数学模型建立

HYDRUS-1D模型分为非饱和带水分运移模型和非饱和带溶质运移模型。

2.3.2.1 非饱和带水分运移模型

假定水分运移过程中气相作用很小，忽略温度梯度的影响，非饱和带水分运移采用Richards 方程的修改形式来表示，公式如下：

式中：θ—土壤体积含水率，L3/L3；

t—水分运移时间，T；

h—非饱和带压力水头，L；

K—土壤水的非饱和水力传导率，L/T；

θs—饱和含水率；

θ0—初始含水率；

V—渗透通量，L/T；

H—非饱和带深度，L。

2.3.2.2 非饱和带溶质运移模型

忽略污染物在气相中的扩散，不考虑液相中通过对流和弥散作用进行溶质运移时的化学反应，在固液相间的吸附作用采用线性平衡方程，公式如下：

式中：θ—土壤体积含水率，L3/L3；

t—水分运移时间，T；

z—空间坐标，向上为正；

Dw—水动力弥散系数，L2/T；

V—渗透通量，L/T；

C0(z，t)—初始浓度；

H—非饱和带深度，L；

Kd—分配系数，L3/M；

M—源汇项。

2.3.3 参数设置

水分迁移模型需要确定的水文地质参数包括：残余含水率θr、饱和含水率θs、垂直饱和渗透系数Ks以及曲线形状参数α、n，均采用HYDRUS-1D软件提供的土壤经验参数库中的数值。模型中采用的水文地质参数见表1。

2.3.4 污染物源强设定

根据工程分析，炼焦过程产生的熄焦废水、装煤、出焦除尘洗涤废水经沉淀处理后回用。其他废水如生活污水、化验废水、循环水排污水、化学水站中和废水入酚氰废水站一并处理。调节池为地下式钢筋混凝土结构，内侧刷防腐防渗涂层，底层、侧壁外侧铺设防渗土工膜，正常情况下污染物泄漏可以忽略不计。本次研究对酚氰废水站调节池设定为非正常工况。

2.3.4.1 情景设定

由于调节池为地下式水池，假定调节池污水发生泄漏，管理人员不能及时发现，致使少量污水通过泄漏部位渗入包气带，最终进入含水层。本次研究非正常工况设定为调节池因老化或者腐蚀局部出现小面积渗漏现象，将污染物泄漏位置看做连续注入的点源。

表1 模型中采用的水文地质参数Tab.1 Hydrogeological parameters used in the simulation

2.3.4.2 渗漏源强设定

调节池为钢筋混凝土构筑物，调节池因老化或者腐蚀产生的渗漏量按正常工况下最大允许渗漏量[15]的10倍考虑，渗漏量为20L/d·m2。调节池池体底面为水平矩形，容积为500m3(12m×8m×5m)，有效贮水深度3m，按涉池壁和池底的浸湿面积计算渗漏量。

研究目标的污染物组成中，可能造成地下水水质影响的主要有挥发酚、氰化物、氨氮和石油类，本次研究选取挥发酚作为预测因子，浓度为500mg/L。

因此，在发生以上景象的渗漏情况下，地下水环境污染物浓度及渗漏源强列于表2。

表2 非正常工况下污染源主要污染因子浓度和源强Tab.2 Concentration and source strength of the main pollution factors under abnormal conditions

3 结果与讨论

从环境安全角度考虑，不考虑吸附作用、化学反应作用等对溶质运移的延迟，采用连续注入模型预测挥发酚进入包气带后的迁移行为。

从预测结果来看，挥发酚进入包气带后，地表以下2.5m处(N1观测点)在15d时开始出现挥发酚，在924d时达到恒定浓度0.5mg/cm3。地表以下17.5m处(N3观测点)在118d时开始出现挥发酚，在710d时达到恒定浓度0.5mg/cm3。地表以下40m处(N6观测点)在306d时开始出现挥发酚，在979d时达到恒定浓度0.5mg/cm3，挥发酚在6个观测点的浓度随时间变化曲线见图3。

图4为100d、200d、300d、400d、500d、600d、800d和1 000d等时刻挥发酚在整个包气带剖面的浓度分布情况，随着时间的增加，污染物浓度逐渐接近浓度0.5mg/cm3，最后整个包气带剖面污染物达到平衡浓度。

图3 挥发酚浓度-时间曲线Fig.3 The concentration-time curves of volatile phenolic compounds

图4 挥发酚浓度-深度曲线Fig.4 The concentration-depth curves of volatile phenolic compounds

4 结 论

本文基于HYDRUS-1D软件对建设项目在非正常工况下调节池因老化或者腐蚀发生的渗漏情况进行预测，预测污染物泄漏以后在包气带中的运移情况。预测结果显示，挥发酚会在306d时进入含水层，说明包气带的阻滞作用使得污染物进入含水层的过程具有一定的滞后效应，在减轻污染影响结果的同时也增加了及时发现污染物泄漏的难度，因此在包气带厚度比较大的地区对涉及废污水的水池进行定期检测也就显得很有必要。

在获得污染物通过包气带到达地下水面时间的同时，也得到了污染物进入含水层的初始浓度。

但是，在进行环境影响评价报告的预测过程有以下需要注意的地方：

(1)由于模型没有考虑吸附，降解等因素，模拟结果不能反应真实情况，预测结果偏保守。

(2)由于某些环评项目实测数据的不足，导致模型的初始数据有所欠缺，某些情况下可以用年降水量来迭代预测稳定的初始含水率，准确率有所不足。

(3)模拟区地层的数据来源为岩土工程勘察报告或者水井柱状资料，其水文地质参数均为经验参数，相邻土层之间的衔接也会影响模型的收敛与最终结果。

参考文献：

[1] 许兆义. 包气带水文地质专论[M]. 北京:地震出版社, 1993.

[2] 李 纯, 武 强. 地下水有机污染的研究进展[J]. 工程勘察, 2007,(1):27-30.

[3] 赵西宁, 吴发启. 土壤水分入渗的研究进展和评述[J]. 西北林学院学报, 2004, 19(1):42-45.

[4] Nielsen D R, Th. V G M, Biggar J W. Water Flow and Solute Transport Processes in the Unsaturated Zone[J]. Water Resources Research, 1986, 22(9):89S-108S.

[5] 汤 英, 徐利岗, 张红玲,等. HYDRUS-1D/2D在土壤水分入渗过程模拟中的应用[J]. 安徽农业科学, 2011, 39(36):22390-22393.

[6] 张蔚臻,张瑜芳.土壤水盐运移数值模拟的初步研究[C].大连:土壤物理学术讨论会论文,1982.

[7] 张效先, 贾大林. 一维非饱和条件下水动力弥散系数的计算方法[J]. 灌溉排水学报, 1987，(2):8-15.

[8] 黄康乐. 野外条件下非饱和弥散系数的确定[J]. 土壤学报, 1988,(2):125-131.

[9] 李 玮, 何江涛, 刘丽雅,等. Hydrus-1D软件在地下水污染风险评价中的应用[J]. 中国环境科学, 2013, 33(4):639-648.

[10] 郝芳华, 孙 雯, 曾阿妍,等. HYDRUS-1D模型对河套灌区不同灌施情景下氮素迁移的模拟[J]. 环境科学学报, 2008, 28(5):853-858.

[11] Flury M. Experimental evidence of transport of pesticides through field soils - a review[J]. Journal of Environmental Quality, 1996, 25(1):25-45.

[12] 马 欢, 杨大文, 雷慧闽,等. Hydrus-1D模型在田间水循环规律分析中的应用及改进[J]. 农业工程学报, 2011, 27(3):6-12.

[13] 傅臣家, 刘洪禄, 吴文勇,等. 六价铬在土壤中吸持和迁移的试验研究[J]. 灌溉排水学报, 2008, 27(2):9-13.

[14] 杨 洋, 李 娟, 李鸣晓,等. HYDRUS-1D软件在地下水污染源强定量评价中的应用[J]. 环境工程学报, 2014, 8(12):5293-5298.

[15] GB 50141-2008, 给水排水构筑物工程施工及验收规范[S].