基于Visual Modflow对某尾矿PRB设计的参数优化

颜秉源, 郑丹丹, 粟登峰, 李东阳

(重庆大学 煤矿灾害动力学与控制国家重点实验室,重庆 400044)

基于Visual Modflow对某尾矿PRB设计的参数优化

颜秉源, 郑丹丹, 粟登峰, 李东阳

(重庆大学 煤矿灾害动力学与控制国家重点实验室,重庆 400044)

渗透反应墙(PRB)的设计是在水文地质数据充分的情况下建立水动力模型进行试验来确定的。由于水文条件随季节变化给PRB设计参数造成了误差。以某尾矿为例,在获得模拟区域整体水文地质数据之后,通过Visual Modflow建立了合理的数值模型。在验证模型合理的情况下,结合不同补给排泄条件下的模拟结果,合理预测了渗透反应墙安装位置、走向以及上部埋深。得出在尾矿库南部南偏东9°方向548 m以内安装渗透反应墙为适宜,且走向垂直于粒子追踪路径。渗透反应墙上部埋深范围取91.11 m>h>87.60 m。通过与传统方案参数在成本及修复作用作对比研究,得知基于Visual Modflow选取的渗透反应墙上部埋深及位置优于传统方案所取,两种方案的走向偏差较小,对地下水的修复及成本影响相差不大。所以基于Visual Modflow选取的PRB参数比传统方案参数更优,可为PRB设计的参数选择提供依据。

Visual Modflow;渗透反应墙;地下水

渗透反应墙(PRB)技术是一种原位被动修复污染地下水技术,是一种新型的前景可观的地下水修复技术[1-3]。PRB的设计过程需要综合考虑多种因素。在通过前期的可行性调研和大量的水文地质勘查获得现场系统的水文地质参数后,进行实验室批量试验和圆柱试验,确定活性反应介质并测试其修复效果和反应动力学参数,再建立水动力学模型,得到系统的渗透反应墙的设计资料后才能进行设计。所以PRB的设计是一个相对复杂的过程[4-6]。且由于不同补给和排泄条件下地下水水流方位及水位会产生变化,所以地下水的运移和赋存通常存在季节性的变化。为了获得更充足的数据,在设计过程需要经过很长时间的观测。在观测时间不够充分的情况下,传统方案所获数据有限,对PRB设计定会存在一定偏差。

Visual Modflow是国际上最为流行且被各国同行一致认可的三维地下水流和溶质运移模拟评价的标准可视化专业软件系统[7]。Visual Modflow在国内得到了比较广泛的应用。马丛安等人用Visual Modflow对某露天矿地下水的水质点示踪和三维污染羽进行了模拟,为露天矿开采期间地下水位及水质变化的预测以及为监测和预保护措施的采用提供决策条件[8];马腾、王焰新运用Visual Modflow对中国南部某大型铀尾矿库进行了铀迁移的数值仿真模拟研究,对比研究了治理和不治理两种情况,得到治理后污染物可控结论[9];赵龙,张延军等人利用Modflow进行某露天煤矿开挖降水设计预测,为降水开挖方案提供依据。

本文基于Visual Modflow,在得到某尾矿系统水文地质资料的情况下,对该区域的丰水期和枯水期分别进行稳定流和瞬时流模拟。通过对比研究不同补给排泄条件下的模拟结果,合理预测了PRB设计的位置、走向及上部埋深,并通过与传统方案参数作对比研究,得出基于Visual Modflow选取的PRB参数比传统方案参数更优,可为PRB设计的参数选择提供依据。

1 Visual Modflow数学模型

Visual Modflow是一款用有限差分的方法来刻画地下水流运动规律的计算机程序。它通过把研究区在空间和时间上进行离散,建立研究区每个网格点的水均衡方程式,所有网格方程联立成为一组大型的线性方程组,通过迭代求解得到每个网格的水头值[10-11]。

1.1 MODFLOW数学模型

在不考虑水密度变化的条件下,孔隙介质中地下水在三维空间的流动可以用下面的偏微分方程来表示:

(1)

式(1)中:Kxx,Kyy,Kzz为渗透系数在x,y,z方向上的分量;h表示水头;W表示在非平衡状态下通过均质、各向同性土壤介质单位体积的流量,表示地下水的源和汇;Ss表示多孔介质的贮水率;t表示时间。

1.2 MODPATH数学模型

对于三维稳定流,MODPATH模型用关于有效孔隙率和渗流流速微分方程表示:

(2)

式(2)中:Vx,Vy,Vz分别表示线性流动流速矢量在各坐标轴方向的分量;n为含水层有效孔隙率(%);W表示由含水层内部单位体积源和汇产生的水量。

2 水文地质条件及应用模型

2.1 水文地质条件

该研究区域为平顶山近郊区,海拔高度80～100 m(图1)。据平顶山1966—1998年降水和蒸发资料,多年平均降水量为734.59 mm,年度降水量变化较大,丰水年最大降水量1 114.00 mm,枯水年最小降水量461.50 mm,7—9月降水量占全年的61.8%,多年蒸发量一般为1 133～1 488 mm,且该区域地层条件较好。地下水的补给主要方式是:大气降水、灌溉回渗。排泄包括蒸发和人工开采,其次为向区外径流排泄。该地形平坦,植被发育较好,包气带岩性结构松散,自地表向下地层分为:亚砂土、亚粘土层,厚度为7～9 m;中砂层和砂砾石层,厚度为13～15 m;淤泥、粉质粘土层,厚度为3～4 m;砂卵石层,厚度为17～19 m。

图1 模拟区域3D表面图Fig.1 3D surface map of simulation area

2.2 模型概化

由于水文地质条件所致,在垂向将模型概化为4层,由上到下依次为弱透水层、潜水含水层、相对隔水层、承压含水层。依据地下水水位数据,将模型北部边界概化为定水头边界,南部边界为常年有水河流,概化为河流边界。模型西部和东部概化为透水边界。区域水流走向整体趋于由西北向东南。模型上部边界为潜水水头,随时间的变化而变化,下部为承压含水层边界。在模型的西北角和南部共有三口抽水井。为了获取相应数据,在模型中共设置了7口观测井,其中1#～6#观测井是为了观测模型区域内最大水位值,7#观测井用来对模型的识别验证。

2.3 模型建立与模型参数

建模区域平面取为矩形,面积为5 400 m×3 900 m,将模型划分为60 m×60 m的单元格,垂向上划分为4个不同层,出于模拟运算的精确性,模型的潜水含水层和承压含水层网格分别细化为2层和4层,模型中共产生46 800个单元格(图2)。网格地面高程采用SRTM-DEM数据,经过软件处理后导入模型。

根据抽水试验和地层岩石性质以及结构特征,结合工程资料和经验,得出区域渗透系数和弹性释水系数(表1)。

图2 模型离散Fig.2 Discretization of model

表1 区域渗透系数与弹性释水系数范围表

Table 1 Regional permeability coefficient and elastic storativity range

P#岩性Kx/(m·s-1)Ky/(m·s-1)Kz/(m·s-1)S/(L·m-1)1粉土8．68e⁃68．68e⁃68．68e⁃70．000012中细砂1．16e⁃41．16e⁃41．16e⁃50．1～0．000013粘土1．04e⁃61．04e⁃61．04e⁃70．000014砂卵石1．67e⁃41．67e⁃41．67e⁃50．1～0．00001

2.4 初始水头的确定

在稳流状态下,将以实际观测井实测水头作为初始水头值,分别导入稳定模型和瞬变流,并在稳流模型中进行计算,得到水位分布状况和粒子追踪路径(图3、图4)。

图3 枯水期初始水位分布Fig.3 The distribution of the initial water level at dry season

2.5 模型验证

模型利用了研究区域钻孔在试验间的地下水观测资料进行调试验证。经参数调整,使研究区模拟地下水变化趋势与观测值基本一致。

在模型中,主要考虑含水层的渗透系数(k)、弹性释水系数(u*),以及边界条件,据经验参数和实测值对比校正,取中细砂含水层渗透系数为1.16e- 4 m/s,砂卵石含水层渗透系数为1.67e- 4 m/s。根据7#观测井分别在枯水期和丰水期的已有实际潜水层水位观测结果,结合模拟中的计算结果对模型进行验证。模拟中采用WHS计算,由图5、图6 7#井观测结果与模拟结果对比可以看出,误差较小,模拟结果较好。

图4 丰水期初始水位分布Fig.4 The distribution of the initial water level at wet season

图5 枯水期7#观测井实测值与模拟值对比图Fig.5 The comparison chart of simulated values and measured values for 7# observation well at dry season

图6 丰水期7#观测井实测值与模拟值对比图Fig.6 The comparison chart of simulated values and measured values for 7# observation well at wet season

2.6 结果分析及参数预测

考虑到优化参数为渗透反应墙的上部埋深(h)及安装方位及位置,所有水位参数取潜水层水位。在模拟过程中,采用PEST程序对模型参数进行校正,同时枯水期模型采用了稳定流和瞬时流两种模拟状态,在稳定流结果中主要获取粒子追踪路径,如图7用来分析渗透反应墙的安装方位以及走向;瞬时流模拟结果中主要获取相关观测井的观测水位,同时获得不同时刻下水位分布状态(图8-图11),以确定水位升降趋势,为渗透反应墙上部埋深选择提供依据。丰水期模型同样采用两种模拟状态,来获取相应的模拟结果。

图7 枯水期粒子追踪图Fig.7 The particle tracking figure at dry season

图8 枯水期第10天水位图Fig.8 The hydrograph of day 10 at dry season

从图7和图12可以看出,尾矿库所标记的示踪粒子从尾矿由北向南偏东向南部河流方向运移。以6#观测井为粒子路径上下游区分点,在追踪路径的上游,图7和图12表示出相对的一致性,重合度比较高,而在下游,丰水期路径向东的偏转角度较大,重合度较差,这是受补给和排泄作用影响的。可以得出,在补给排泄条件不同的情况下,地下水的运移和分布是有一定区别的,渗透反应墙的安装方位和走向不能根据单一的补给排泄条件下的地下水运移分布而确定,这也加大了渗透反应墙设计的复杂性。

图9 枯水期第120天水位图Fig.9 The hydrograph of day 120 at dry season

图10 丰水期第10天水位图Fig.10 The hydrograph of day 10 at wet season

图11 丰水期第120天水位图Fig.11 The hydrograph of day 120 at wet season

图12 丰水期粒子追踪图Fig.12 The particle tracking figure at wet season

选取观测井1#、3#、7#为参考,观察水位升降趋势。将图8和图9作对比,可以看出枯水期等水位线88 m、87 m、86 m均有不同程度向西北方向偏移,结合观测井水位随时间变化图13-图15可以得出,在枯水期受补给和排泄的影响下,模拟区域水位整体呈现下降趋势。同理,结合图10、图16-图18分析可知,丰水期模拟区域整体水位也呈现下降趋势。

图13 枯水期1#观测井水位时间图Fig.13 Water level change with time graph for 1#observation well at dry season

图14 枯水期3#观测井计算值Fig.14 Water level change with time graph for 3# observation well at dry season

图15 枯水期7#观测井计算值Fig.15 Water level change with time graph for 7# observation well at dry season

图16 丰水期1#观测井计算值Fig.16 Water level change with time graph for 1# observation well at wet season

图17 丰水期3#观测井计算值Fig.17 Water level change with time graph for 3# observation well at wet season

图18 丰水期7#观测井计算值Fig.18 Water level change with time graph for 7# observation well at wet season

渗透反应墙的安装一般布置在地下水流(污染羽)的下方,走向垂直于地下水流向[12-13]。由于补给及排泄的不同导致地下水运移和赋存也存在差异。结合图7、图12、图19可以看出丰水期稳定流模拟状态下的粒子追踪路径在下游出现了相对于枯水期较大角度的向东偏转,从节约成本和渗透墙的安装操作上考虑不适合于下游安装渗透反应墙;而在上游,丰水期和枯水期的粒子追踪路径基本吻合,所以适合于安装渗透反应墙。经计算丰水期粒子路径上游南偏东10°,下游南偏东46°;枯水期粒子路径上游南偏东8°,下游南偏东26°。因此设计时取粒子路径上游南偏东平均值(8°+10°)/2=9°。经测定,观测井6#和3#的垂向南北走向距离为822 m,为了缩小渗透墙布置宽度,取修正值822 m/1.5=548 m,设计时取位置为尾矿堆南部548 m以内。

图19 地下水矢量对比图Fig.19 Groundwater vector comparison chart

为了防止地下水溢出反应墙,出于对地下水位季节性波动情况的考虑,PRB的顶端需高于地下水最高水位[14-15]。为了能够监测到渗透反应墙安装区域的最高水位,在模型建立时设置了7口水位观测井。以渗透反应墙安装方位和走向的确定为依据,最优选择3#水位观测井作为最高水位监测井。从图14和图17得知3#观测井在丰水期最高水位低于87.60 m,枯水期最高水位低于87.0 m。因此可以判定渗透反应墙安装的上部埋深高于87.60 m。而上部埋深的最大值以该区域海拔高度为限,从SRTM-DEM数据可以得知为91.11 m。所以渗透反应墙的上部埋深(h)取91.11 m>h>87.60 m。

3 参数优化分析

不考虑其他因素影响,对比传统方案与基于Modflow方案的PRB设计参数,分析参数质量。

表2 参数对比表

从表2可以看出基于Visual Modflow方案对PRB的位置选取比传统方案更靠近尾矿堆(600-548=52 m)。而污染物的扩散有溶质一般迁移数学模型:

(3)

式中:C为溶质浓度(M·L-3);θ为有效孔隙度;Dij为弥散系数张量(L2·T-1);qi为达西速度(L·T-1);Cs为源汇项浓度(M·L-3);∑Rk为化学反应速率(T-1)。

可知污染物沿水流方向扩散快,以污染源为中心向外半径越大,扩散范围越广。所以基于VisualModflow选取PRB位置更靠近尾矿堆,这样有利于对污染物的控制修复,同时也可以减少修复成本。

基于两种方案的上部埋深都高于尾矿堆南端最高监测水位87.10m的情况下,传统方案的上部埋深下限高于VisualModflow模拟方案上部埋深下限1.6m。这样在能达到地下水同样修复效果条件下,在设计过程中传统方案会造成一定的资源浪费,附加成本增多。

在两种方案中PRB的走向选择,偏差只有2°,这对PRB地下水的修复效果及成本影响偏小,可以忽略不计。

VisualModflow对PRB的上部埋深及位置的预测结果比传统方案的预测结果更优。而在PRB走向选取上误差较小,优化效果不太明显。总结以上分析结果,得到基于VisualModflow选取的PRB参数比传统方案参数更优。

4 结论

基于MODFLOW对地下水三维流场进行模拟得出的水位图与粒子追踪图对渗透反应墙安装具有重要意义。

(1)MODFLOW对地下水的模拟结果能够很好地反应地下水的运移和赋存情况,通过与传统方案在修复效果和成本消耗作对比,得到基于VisualModflow选取的PRB参数比传统方案参数更优,可为PRB设计的位置、走向以及上部埋深提供依据。

(2) 为了预防补给和排泄不同导致最高水位高于上部埋深的情况,本文模拟采用两种典型的模型,丰水期和枯水期模型。这样能保证预测的上部埋深高于地下水各个时期的最高水位。由于在不同补给和排泄条件下地下水的运移和赋存有一定差异,所以在处理渗透反应墙的安装方位及走向能进行优化,但还不是最优值,有待在更多不同补给和排泄条件下更进一步的优化。

[1]McMurtyDC,EltonRO.Newapproachtoin-situtreatmentofcontaminatedgroundwater[J].EnvironProg,1985,4(3):68-70.

[2]GuerinaTF,HornerbS.AnapplicationofPermeablereactivebarriertechnologytoPetroleumhydrocarboncontaminatedgroundwater[J].WaterResearch,2002,36(1):15-24.

[3]BeitingerE.Permeabletreatmentwalls-design,constructionandcost[J].NATO/CCMSpilotstudy,1998,229(1):6-16.

[4] 赵勇胜.地下水污染场地污染的控制与修复[J].吉林大学学报:地球科学版,2007,37(2):303-310.

[5] 陈仲如,张澄博,李洪艺,等.可渗透反应墙的结构与设计研究[J].安全与环境学报,2012,12(4):56-61.

[6] 徐乐昌,周星火,詹旺生,等.某铀矿尾矿堆场受污染地下水的渗透反应墙修复初探[J].铀矿冶,2006,25(3):153-157.

[7] 丁继,周德亮,马生忠.国外地下水模拟软件的发展现状与趋势[J].勘察科学技术,2002,1(1):37-42.

[8] 马腾,王焰新.U(Ⅵ)在浅层地下水系统中迁移的反应-输运耦合模拟:以我国南方核工业某尾矿库为例[J].地球科学——中国地质大学学报,2000,25(5):456-461.

[9] 赵龙,张延军,巩立亮,等.利用Modflow进行某露天煤矿开挖降水设计预测[C]//第二届全国岩土与工程学术大会论文集编辑委员会.北京:科学出版社,2006:12-15.

[10] 马从安,李克民,石磊,等.VisualModflow在露天矿地下水模拟中的应用[J].环境工程,2011,29(1):98-101.

[11] 朱玉水,段存俊.MT3D—通用的三维地下水污染物运移数值模型[J].东华理工学院学报,2005,28(1):26-29.

[12]ARUNR,GAVASKAR.Designandconstructiontechniquesforpermeablereactivebarriers[J].JournalofHazardousMaterials,1999,68(1):41-71.

[13]CUNDYABA,HOPKINSONAL.Useofiron-basedtechnologiesincontaminatedlandandgroundwaterremediation:areview[J].ScienceoftheTotalEnvironment,2008,400(1):42-51.

[14] 胡黎明.地下水污染修复的活性渗滤墙技术[J].水利水电技术,2003,34(7):11-13.

[15]USEPA.Permeablereactivebarriertechnologiesforcontaminantremediation,EPA/600/R-98/125[R].Washington,DC:EPA,1998.

(责任编辑：陈文宝)

ParameterOptimizationofPRBDesignforaTailingBasedonVisualModflow

YAN Bingyuan, ZHENG Dandan, SU Dengfeng, LI Dongyang

(StateKeyLaboratoryforCoalMineDisasterDynamicsandControl,Chongqing400044)

Under the condition of the hydrogeology information is sufficient,the hydrodynamic model experiment was established and tested to design permeable reactive barrier(PRB).Due to the hydrological conditions change with the seasons caused PRB parameters have deviation.Taking a tailing as example,after the overall hydrogeology data of simulation district was obtained，the reasonable numerical model was established by Visual Modflow.After verifying the model is reasonable,combined with the simulation results in different supply and discharge,the location，strike，and the upper buried depth of permeable reactive barrier were predicted.The paper concluded that permeable reactive barrier was installed in the south of the tailings within 548m south by east 9° as appropriate.Permeable reactive barrier’s direction is perpendicular to the particle tracking path.The upper buried depth of permeable reactive barrier range from 91.11 m>h>87.60 m.By compare with the cost and repair function in the traditional scheme parameters,the location，and the upper buried depth of permeable reactive barrier based on the Visual Modflow are superior to the traditional scheme.The strike of two plans only have small deviation，which has small impact on the cost repair effect.So parameters selected by Visual Modflow than traditional scheme better，can provide the basis for parameters selection of PRB design.

Visual Modflow; permeable reactive barrier; groundwater

2014-12-25;改回日期:2015-03-19

颜秉源(1988-)，男，在读硕士研究生,采矿专业。E-mail:779190917@qq.com

P641.5

A

1671-1211(2015)03-0336-07

10.16536/j.cnki.issn.1671-1211.201503022

数字出版网址:http://www.cnki.net/kcms/detail/42.1736.X.20150422.1108.023.html 数字出版日期:2015-04-22 11:08