岩溶区某磷石膏堆场地下水渗流与污染物运移数值模拟

张科正，陈 舟，陈长生，赵贵清

(1.河海大学 地球科学与工程学院，南京 210000； 2.长江三峡勘测研究院有限公司，武汉 430074)

0 序言

我国西南岩溶地区的部分磷石膏堆放场存在选址不当、堆场设计不合理、防渗措施不彻底等造成的渗漏污染问题[1-2]，主要污染物有磷、氟、硫酸盐等。污染物对地下水的污染及对地表水体的排放对水环境和居民生活造成了较大影响。针对此类问题，已有众多学者对渗漏类型、特征和定量分析做出了相关研究，其中，赵伟丽等通过10个工业废渣堆场的岩溶含水层介质特征及污染质运移特征的分析,划出4种岩溶渗漏污染的水动力弥散类型，分别为管道系统、管道裂隙、孔隙裂隙、采空区渗漏污染类型[3]。彭展翔等对本研究堆场渗漏区、渗漏类型、渗漏途径、渗漏主方向等渗漏特征进行了分析，确定渗漏主方向为发财洞方向[4]。在定量评价方面，由于岩溶地区地质条件的复杂性，对于岩溶渗漏污染的数值模拟与预测研究还不够完善，现主要预测模型有水动力弥散叠加模型、稀释混合模型、基于MATLAB BP神经网络模型、一维解析解法等。赵良杰等采用MODFLOW中排水沟(Drain)和河流(River)模块概化模拟了岩溶管道特征，对比分析并指出Drain模块中岩溶管道仅起排水作用，不允许管道内水流流向周围含水层，而River模块允许岩溶管道与周围含水层进行水量交换，因而相对于Drain模块，应用River模块概化模拟岩溶管道更加精确[5]。总结可知，岩溶区地下水污染物迁移的数值模拟鲜有报道，本文拟以某磷石膏堆放场的污染渗漏为例，运用Groundwater Modeling System(GMS)软件搭建了该区域的三维地下水渗流场和污染物迁移的数值模型，以期查明污染物(磷)的分布迁移特征，为堆场的污染治理工作和岩溶水污染预测评价提供依据。

1 研究区基本地质条件

2 堆场渗漏的数值模拟

2.1 模型建立

2.1.1 水文地质概念模型

图1 研究区模型边界及类型Figure 1 Study area model boundary and type

图2 堆场—发财洞剖面图Figure 2 Dump-Facaidong section

2)含水层结构。研究区属裸露型岩溶山区，碳酸盐岩出露面积占70%以上，地下水类型为裂隙水和岩溶水为主，区域内的断层、天窗使得上下含水层相互串通，成为统一的混合含水层。

3)含水层的水力特征。结合研究区地质资料、地下水补、径、排特征及岩溶发育特征，将研究区概化为非均质各向异性三维地下水流场。研究区出露泉点水位多年动态变化较稳定，所以本次模拟将建立枯水期研究区地下水稳定流场。地下水溶质主要通过对流-弥散作用随地下水流一起运移。

2.1.2 地下水渗流与溶质运移耦合数学模型

1)地下水渗流数学模型。通过概化得到的非均质各向异性等效连续介质模型，地下水稳定运动数学模型为

(1)

式中：H—地下水水头，m；

Kx，Ky，Kz—各向异性主渗透系数，m/d；

Γ1—模拟区域第一类边界；

Γ2—模拟区域第二类边界；

H0(x,y,z)—含水层初始水头m；

Hг(x,y,z)—第一类边界条件边界水头，m；

q0(x,y,z)—第二类边界单位面积过水断面补给流量，m2/d；

ε—源汇项强度(包括开采强度等)，1/d；

Ω—渗流区域。

(2)地下水溶质运移数学模型。

(2)

式中：n—介质孔隙度，无量纲；

C—组分的浓度，mg/L；

t—时间，d；

Dij—水动力弥散系数张量，m2/d；

υi—地下水渗流速度张量，m/d；

qs—源和汇，1/d；

Cs—源或汇水流中组分的浓度，mg/L；

C0(x,y,z,t0)—已知浓度分布，mg/d；

Ω—模拟区域；

Γ1—通量边界；

fi(x,y,z,t)— 边界上已知的弥散通量函数。

2.2 参数选择

2.2.1 渗透系数

本模拟结合钻孔压水试验，并对比了相关的工程实例，进行了渗透系数的选取。考虑岩体各向异性，其渗透系数按照三维取值，各岩性对应渗透系数见表1。

2.2.2 降雨入渗系数

模型输入时的参数Recharge rate项等于降水强度与降雨入渗系数的积，降雨入渗系数在枯水年平均有效降雨入渗系数为0.35，参考其他文献和区域的地形和岩性特征，区域的降入入渗分区及参数取值可见表2。

表1 研究区岩体渗透系数赋值

表2 不同岩性对应降雨入渗系数

2.3 模型识别与校正

研究区地下水数值模型采用GMS10.2软件进行数值离散，采用等间距有限差分的离散方法对模拟区进行自动剖分，网格单元行×列×层为200×200×5，模拟区的三维尺度在X方向上长度为8 349.66 m，Y方向上长度为7 492.21 m，Z方向的长度为500m，网格剖分及岩性三维分区见图3。初始水头设置为网格顶部高程。基于GMS创建各岩溶通道并将其属性设置为River，通道的两个端点高程设置为岩溶管道实际高程。

图3 网格剖分及岩性三维分区Figure 3 Mesh generation and lithologic 3D partitioning

在数值模拟过程中，由于地质条件的概化、模型参数的选取等都会造成模拟地下水流场与实际地下水流场不符。因此，对模型进行识别与校正是数值模拟的一个重要环节。模型识别与校正不能片面追求拟合效果而任意调整模型参数，甚至改变含水层的边界条件及源(汇)项等原始数据，这种做法不仅容易使模型脱离实际的水文地质条件，而且还造成模型的多解性。研究区的面积很大，钻孔水位资料相对缺乏，但岩溶泉点有50个，根据泉水点位的出露高程插值可得研究区的初始流场图，见图4。此外，模型的识别利用了2016年踏勘所调查的16个岩溶泉出露水位作为监测井，模型的初始水头设置为网格顶板高程[8]。

图4 研究区地下水初始流场图Figure 4 Study area groundwater initial flow field

利用GMS软件中的参数自动识别模块(PEST)与手动调参相结合的方法进行模型校正与参数反演，识别后模型的地下水流场和监测井水位计算误差图可见图5和图6。

图5 研究区地下水场拟合图Figure 5 Study area groundwater field fitting diagram

图6 各观测井实测值与计算值的拟合效果图Figure 6 Fitting effect diagram of measured and computed values from observation wells

由结果可知，16口观测井的观测水位与计算水位的平均残差值为-0.253 7，拟合误差的平均绝对值为0.586 4，拟合误差小于0.8m的观测井占总的井数量的87%以上，可见实测水位与计算水位之间达到了较好的拟合效果，参数校正结果比较可靠，参数的校正结果满足精度要求。此外，模型River模块反应出了岩溶通道中水流集中排泄的特征，通道区域周围的地下水汇集到管道中以较快的流量进行排泄，所建立的模型切合实际。

2.4 堆场渗漏预测

2.4.1 溶质运移模型

本研究利用GMS软件中MT3DMS模块进行污染物迁移的模拟和预测。模拟参数中的孔隙度可见表3。

纵向弥散度由Xu和Eckstein(1995)[7]回归分析得到公式：

aL=0.83×(logLs)2.414

=0.83×(log 2 000)2.414

=14.8m

(3)

式中：aL—纵向弥散系数；Ls—观测尺度。横向及垂向弥散度为纵向弥散度的十分之一。

表3 不同岩性对应的初始和识别后孔隙度

污染渗漏工况设置如下：由文献[3]可知，摆纪渣场的渗漏量约为300 m3/d，渗漏液中磷的浓度约为3 000 mg/L，本次模拟按照渗漏液的10%会入渗地下水的工况进行堆场渗漏模拟。

污染物迁移模型的识别和校正的数据利用了2016年发财洞处污染物(磷)监测数据。观测值与计算值的误差图可见图7。孔隙度识别后的参数值可见表4，弥散度识别后的数值为12.5 m。

图7 磷浓度计算值与观测值比较Figure 7 Comparison of phosphorus measured and observed concentration values

如图7所示，发财洞处一年内磷浓度的计算值与观测值较为一致，数据间相关系数为0.997，可认为模型较为真实地反映出地下水中污染物的运移情况。

2.4.2 渗漏预测

为预测持续渗漏情况下污染物的浓度变化，利用识别后的模型预测了堆场连续渗漏后1 a、5 a、10 a和20 a后污染物在含水层的运移情况，具体可见图8。此外，为了对比含水层污染物迁移的垂直分布和三维分布情况，以污染物迁移10 a后为例展示了相关图件，具体可见图9。

由预测可知，污染物持续渗漏十年后发财洞处总磷浓度可达400 mg/L左右，磷石膏堆场～发财洞方向为堆场的主要渗漏方向，老落凼岩溶管道与摆郎发财洞岩溶管道为主要渗漏通道，渗漏污染物随地下水水流方向河流处排泄， 摆郎发财洞岩溶管道水系统为主导方向，形成堆场至发财洞污染物浓度梯次递减形分布的规律。

图8 不同时段污染物迁移模拟图Figure 8 Pollutant migration simulated diagram of different time intervals

图9 10 a后污染物运移预测三维及剖面图Figure 9 Predicted pollutant migration 3D map and section after 10 years

3 结语

为了定量研究岩溶区某磷石膏堆放场地下水渗流和污染物渗漏污染问题，本文利用地下水数值模拟软件GMS建立了岩溶区地下水流场及污染物(磷)在岩溶含水层的运移情况，在模型的识别和校正后预测了1～20 a后污染物的浓度分布。

由结果可知：River模块能较好的模拟岩溶区的岩溶管道渗流。堆场至发财洞方向为研究区的渗漏主方向，老落凼岩溶管道与摆郎发财洞岩溶管道为主要渗漏通道，形成堆场至发财洞污染物浓度梯次递减分布的规律。预测持续渗漏情况下，10 a后发财洞处总磷浓度可达400mg/L左右。

由于基础资料不够充足且岩溶区的部分特征未能在GMS模型中体现，数值模拟的结果还存在很大的完善空间；溶质运移模拟中，本文只考虑了污染物运移过程中的对流、弥散作用，并未考虑吸附作用及不同管道特征的影响，这也是相关工作需进一步深入之处。

参考文献:

[1]段先前,韦俊发,丁坚平.贵州某磷石膏堆场渗漏污染评价[J].资源环境与工程, 2008, 22(2):218-221.

DUAN X Q,WEI J F, et al.Evaluation on seepage pollution of a phosphogypsum yard in GuiZhou Province[J].Resources Environment And Engineering, 2008, 22(2):218-221.(in Chinese)

[2]郑文成,贾洪彪,郭明,等.平坝磷石膏堆场渗漏分析与污染防治对策[J].中国岩溶, 2010, 29(1):75-80.

ZHENG W C, JIA H B, GUO M, et al. Analysis and countermeasures of pollution prevention and control of Pingba phosphogypsum yard leakage[J].Carsologica Sinica,2010,29(1):75-80.(in Chinese)

[3]赵伟丽,褚学伟,董毓,等.岩溶含水层渗漏污染弥散类型分析——以贵州废渣堆场污染为例[J].地下水, 2011, 33(2):6-7.

ZHAO W L,CHU X W, el al. Analysis on leakage pollution dispersion type in karst aquifer[J].Groundwater, 2011,33(2):6-8.(in Chinese)

[4]彭展翔, 褚学伟.摆纪磷石膏堆场渗漏分析[J].地下水, 2012, 34(5).

PENG Z X, ZHU X W. Analysis on the leakage of the phosphorus gypsum storage yard in Baiji[J].Groundwater,2012,34(5).(in Chinese)

[5]赵良杰,夏日元,杨杨,等.基于MODFLOW的岩溶管道水流模拟方法探讨与应用[J].中国岩溶, 2017, 36(3):346-351.

ZHAO L J,XIA R Y, el al. Discussion and application of simulation methods for karst conduit flow based on MODFLOW[J].Carsolgica Sinica,2017,36(3):346-351.(in Chinese)

[6]张景辉, 孙曙光. 地下水非稳定运动的数学模型及应用[J].上海国土资源, 1981(4):43-52.

ZHANG J H,SUN S G. The mathematical model of unsteady groundwater movement and application[J].Shanghai geological,1981, 41-50.(in Chinese)

[7]M Xu,Y & Eckstein. Use of weighted least-squares method in evaluation of the relationship between dispersivity and field scale[J].Groundwater, 1995, 33(6):905-908.

[8]陈舟,赵贵清,王志光,等.岩溶区某磷石膏堆放场渗漏特征分析[J].水文地质工程地质,2017,44(2):144-150.

CHEN Z,ZHAO G Q, el al. Leakage Characteristics of a Phosphorus Gypsum Storage Site in Karst Area[J].Hydrogeology And Engineering Geology,2017,44(2):144-150.(in Chinese)

[9]祝晓彬.地下水模拟系统(GMS)软件[J].水文地质工程地质, 2003, 30(5):53-55.

ZHU X L. Groundwater modeling system(GMS) software[J].Hydrogeology and engineering geology,2003,5:53-55.(in Chinese)

[10]谭文清,孙春,胡婧敏,等.GMS在地下水污染质运移数值模拟预测中的应用[J].东北水利水电, 2008, 26(5):54-55.

TAN W Q,SUN C, el al. Application of GMS in simulation of pollutants migration for groundwater[J].Water Resources And Hydropower of Northeast China,2008,5:54-55.(in Chinese)

[11]骆祖江,于丹丹,刘艳,等.江苏长强钢铁厂地下水环境影响评价[J].中国煤炭地质, 2014, 26(5):34-37.

LUO Z J,YU D D, el al. Jiangsu Changqiang Steel Mill Groundwater Environment Impact Assessment[J].Coal Geology of China,2014,26(5):34-37.(in Chinese)