1. 前言
矿山开采已发展成为我国主要的地下水污染源之一。由于矿山的建设和开采,矿坑涌水的抽排导致地下水位下降,地下水流动系统受到影响,严重情况下含水层会被输干从而导致地面沉降,严重破坏了居民生活质量及大自然的生态环境 [1] [2] ;普遍情况下,多个含水层被矿井贯穿,含有杂质的矿井水相互贯通并混合导致“串层”污染,水质产生严重恶化;另外,矿山建设中产生的废石堆含有大量的有害物质,矿产开采会使长期埋置于地下的有害物沉积物质暴露于地表,然后侵入含水层,造成地下水环境的污染,因此对矿山开采所引发的地下水环境问题进行研究是迫切需要的。目前国内外学者对矿山开采对地下水环境影响程度的研究日益增多,主要研究方法包括解析法、数值法、趋势外推法、类比预测法、回归分析法等。其中,将研究区地下水环境进行理想模型概化并进行数值模拟分析是行之有效的方法之一 [3] [4] [5] [6] [7] ,可以为矿山的地下水环境监控和管理提供依据。
目前的研究主要集中在将矿山开采过程中的水量问题与开采引起的水质问题进行单独研究,或者是仅研究单源污染物的迁移规律,但对于大型矿区污染风险评价及长观监测方案拟定在地下水水量和水质问题是同时发生和多污染源耦合污染方面研究还涉及不多 [8] 。
本文以某大型铜钼矿区为研究对象,详细调绘矿区水文地质条件并全面分析其水文地质环境,理论概化为与其实际相吻合的水文地质环境模型 [9] ,运用溶质运移模型和地下水流数学模型,以数值模拟软件FEFLOW为依托,建立其地下水流、溶质数值模拟模型,以在矿区有无抽排水条件下,多源污染物的迁移规律预测分析为目标,研究矿区多场耦合下污染物的迁移规律 [10] ,并根据预测模拟结果,并据此布设地下水监测井,为矿区地下水的监测和管理提供有力依据。
2. 矿区概况
2.1. 矿区简述
该铜钼矿床资源储量大、矿体规模大且埋藏浅,适于大规模开采。该铜钼矿体主要分布在二长花岗斑岩体内、外接触带中,并以二长花岗斑岩体周围砂岩地层矿化为主,斑岩体内接触带有部分矿化。经勘探,本铜钼矿矿石资源储量为:工业铜矿石4639万t,工业钼矿石14,305.8万t,低品位铜矿石19,532万t,低品位钼矿石18,600万t。
矿山基建期2.5年,矿山服务年限总共28年,其中首期12年。本工程主要由铜钼矿露天开采区、工业场地、排土场、尾矿库、矿区道路和矿区供水供电等部分组成。矿区的疏干降水方案考虑了矿山分期开采的特点,共设计两期分段疏干,第一期12年,疏干至35 m,第二期16年、疏干至−40 m,整个露采作业时间计28年。根据矿山开采的可行性研究设计,排土场下渗水、尾矿库下渗水及环保水库渗漏为该项目地下水污染的主要污染源,具体分布见图1。
![](//html.hanspub.org/file/12-2890287x9_hanspub.png)
Figure 1. Project distribution in the mining area
图1. 矿区平面示意图
2.2. 自然地理条件
矿区建设区属构造侵蚀中等山区,海拔以150~400 m为主,最高点592.0 m,位于圆珠顶山头,最低点93 m,位于建设区东北侧沙冲河一带,区内相对高差一般200 m,地形切割中等。地形总体北低南高,山势走向以北北西为主,其间发育南北向小沟壑。矿体分布范围内,大体上呈四周高中间低的小盆地,山脉呈环形,小沟壑呈树枝状。整体地形地貌见下图2。
矿区属亚热带季风气候,温暖潮湿,雨量充沛。据当地气象局1998~2007年多年观测资料统计,年平均气温21.6℃,最低气温−3.4℃,最高气温39.4℃;日最大降雨量为149.00 mm,年平均降雨量为1404.1 mm,区内50年内日最大降雨量出现在1966年6月11日~12日,为198.90 mm。降雨量多集中在4~9月,为当地丰水期,枯水期为12月至翌年2月,其余时间为平水期。
![](//html.hanspub.org/file/12-2890287x10_hanspub.png)
Figure 2. Topographic and geomorphological area
图2. 地形地貌示意图
2.3. 地质条件
区域内断裂发育,矿区距离沙头–夏郢大逆断层20 km,距离金装–都平大逆断层15 km,矿区范围内基本无构造影响,具体分布见图3。
矿区分布和出露最广的地层是寒武系水口群(∈sh)地层,边缘是不同时期的侵入岩体,以及沿河道分布的小范围第四系(Q)冲积物和局部的残坡积地层。
寒武系水口群(∈sh):分布广,由浅海相复理式砂页岩组成。层理较明显,多呈薄-厚层状,节理较发育。其下部夹硅质岩,中部夹细砾岩,上部夹粗砂岩或含细砾石粗砂岩。
第四系(Q):主要分布于东北部和西南部,由砾石、砂和粘土组成。
3. 矿区水文地质条件及概念模型
将矿区地下水系统实际的边界条件、渗透性质、水力特征、内部构造及补给排泄等条件合理概化为水文地质概念模型,从而进行数值模拟研究。科学并合理建立水文地质概念模型是矿区水文地质环境数值模拟研究的关键。
3.1. 模拟区范围的确定
根据研究区水文地质条件特点及区内环境保护目标和敏感区域分布,结合现有各类资料确定研究区模拟范围如下(见图4):模型北部以地形分水岭为界,全长约35 km,模型东部以贺江为界,全长约38 km,模型南部以地头—都平圩断层为界,全长约14 km,模型西部以东安江为界,全长约52 km,模拟区总面积约528 km2,基本构成一个较完整的水文地质单元。综上所述,研究区地下水系统概念模型为非均质、各向异性、三维非稳定流的地下水系统。
3.2. 矿区水文地质条件及结构的概化
研究区内含水层按埋藏条件及含水介质类型可分为松散岩类孔隙潜水含水层、基岩风化裂隙潜水含水层、基岩构造裂隙含水层(带)三类。松散岩类孔隙潜水分布范围小,厚度薄(0.9 m~2.3 m)且分布不均匀,没有形成连续稳定含水层,水量小,因此在进行的预测模拟研究中不予考虑。
研究区风化裂隙潜水主要通过大气降水补给,部分接受地表水的下渗补给,风化裂隙形成径流通道,补给的地下水沿其径流通道下渗与运移,从而形成地下径流,径流方向与地形有密切关系,总体由地形高处往低处径流,即在研究区内主体自西北向东南径流,局部由于地形变化,形成次级地下水汇流区。
![](//html.hanspub.org/file/12-2890287x12_hanspub.png)
Figure 4. Range of the simulation area
图4. 模拟区范围示意图
深部基岩裂隙水主要接受上部风化裂隙潜水的补给,径流方向与裂隙径流通道密切相关,总体自西北向东南流,往往沿构造破碎带或碳酸盐岩脉发育带运移,形成深部循环地下水。
基岩风化裂隙潜水和基岩构造裂隙水是本次模拟的主要对象,各自的基本特征如下:
1) 基岩风化裂隙潜水含水层
该含水层主要分布于研究区内寒武系水口群中的风化裂隙带中,岩石风化裂隙发育程度控制着含水层的厚度,含水层的底部边界与岩石风化带分界线保持一致,且与地形走势基本一致。据矿区钻孔揭露资料统计,风化裂隙潜水层的平均厚度55.06 m。
2) 隔水层
风化裂隙含水层以下的寒武系水口群上亚群中组(Єsh32)和燕山期二长花岗斑岩(ηγл)微风化层,岩体完整性较好,岩层厚度相对较薄,但后期因燕山期二长花岗斑岩的挤压和侵入,裂隙多为充填闭合状,透水性较差,仅在矿区及周边展布,该层在模拟研究中可视为相对隔水层。
3) 基岩构造裂隙含水层
研究区内微风化带以下变质石英粉砂岩等岩石的构造裂隙带为富水部位,为构造裂隙含水层。构造裂隙带和断裂破碎带含水性在垂向上变化较明显,总体上岩层含水性自上而下逐渐减弱。
3.3. 边界条件的概化
根据水文地质概念模型并结合已有各类水文地质资料,确定本次模拟研究区边界条件如下(见图5):
1) 四周边界
• 基岩风化裂隙潜水含水层
北部边界:模拟区北部边界为地形分水岭,将其定义为零流量边界,分水岭的提取主要利用模拟区DEM数据,在ArcGIS水文分析(Hydrology)模块下完成。
东侧边界:模拟区的东部为贺江,整体上贺江水自北向南流动,最后汇入长江。模拟区内河水平均水位标高50 m,区域上贺江离矿区较远,对项目建设区影响较小,故将其概化为第三类边界(交换边界)。
南侧边界:模拟区南侧发育和地头—都平圩断层,断层贯穿贺江和东安江,根据区域地下水的整体流向(自西北方向向东南方向流动)将该断层概化为流量边界。
西侧边界:模拟区的西部为东安江,呈南北走向,河水自北向南流动,最后汇集贺江。区内河水平均水位标高45 m,和东侧贺江类似,该河流概化为第三类边界(交换边界)。
• 构造裂隙含水层
北部为零流量边界(隔水边界),东部和西部因河水切割较浅,设为隔水边界,南部承压含水层被断层贯穿,同裂隙潜水含水层一样,设为流量边界。
2) 上边界为降水补给、蒸发和井排泄边界。
3) 下边界为构造裂隙含水层底部的寒武系水口群上亚群中组(Єsh32)和燕山期二长花岗斑岩(ηγл)新鲜基岩层,岩体完整性好且分布连续,可视为研究区内稳定的相对隔水层。
4. 地下水模型的建立
4.1. 数学模型
深入分析研究区地下水系统概念模型,基于达西定律和渗流连续性方程,建立与研究区水文地质概念模型相吻合的三维非稳定流数学模型:
![](//html.hanspub.org/file/12-2890287x13_hanspub.png)
Figure 5. Boundary types of the simulation area
图5. 模拟区边界条件示意图
![](//html.hanspub.org/file/12-2890287x14_hanspub.png)
, ![](//html.hanspub.org/file/12-2890287x16_hanspub.png)
, ![](//html.hanspub.org/file/12-2890287x18_hanspub.png)
, ![](//html.hanspub.org/file/12-2890287x20_hanspub.png)
![](//html.hanspub.org/file/12-2890287x21_hanspub.png)
![](//html.hanspub.org/file/12-2890287x22_hanspub.png)
式中,Ω:地下水渗流区域,量纲:L2;
μs:单位储水系数,量纲:L−1;
Kxx,Kyy,Kzz:分别为x、y、z主方向的渗透系数,量纲:LT−1;
H0:初始地下水位,量纲:L;
H1:指定水位,量纲:L;
S1:第一类边界;
S2:第二类边界;
w:源汇项,包括蒸发,降雨入渗补给,井的抽水量,量纲:T−1;
:表示在边界不同位置上不同时间的流量,量纲:L3T−1;
:表示水力梯度在边界法线上的分量。
qn:为正常的达西通量边界;
t:为时间(L);
:为转换系数(T-1)
溶质运移的三维水动力弥散方程的数学模型如下:
![](//html.hanspub.org/file/12-2890287x26_hanspub.png)
![](//html.hanspub.org/file/12-2890287x27_hanspub.png)
式中,等式右侧前三式是弥散项,后三式是对流项,末尾式为溶质的增量,溶质增量通过化学反应或吸附解析作用而产生;Dxx、Dyy、Dzz分别为x、y、z方向的弥散系数;μx,μy、μz为x、y、z方向的实际水流速度;c为溶质浓度,量纲:ML−3;Ω为溶质渗流的区域,量纲:L2;c0为初始浓度,量纲:ML−3。
4.2. 初始网格剖分及地质模型 [11] [12] [13]
本次研究区的模拟面积为528.71 km2,采用TMesh剖分方法将研究区地下水系统水文地质概念模型离散为不规则三角剖分网格,网格剖分严格遵循Delaunay法则 [11] ,同时根据矿区展布特点及两期开采范围,为提高模拟精度,对重点研究区(约92 km2)进行一级加密,对矿体一公里范围内进行二级加密,为防止因疏排水导致过高的水头差出现数值收敛性问题,对两期矿体边界进行多级加密,最终得到模拟区二维剖分结果,其中结点数30,402个,有限单元数60,137个。
本次地质概化模型共分为三层(layer)四片(slice)。
三层:第一层为寒武系水口群上亚群中组基岩风化裂隙潜水含水层,第二层为裂隙充填闭合且仅在矿区及周边展布的相对隔水层,第三层为寒武系水口群上亚群中组基岩构造裂隙含水层。四片:地面高程、上部隔水层顶板、构造裂隙含水层顶板和构造裂隙含水层底板。
研究区地面高程数据利用ASTER GDEM数据,通过ArcGIS软件处理得到。根据研究区钻探资料及区域地质资料,采用克里金插值法得到研究区其他各片高程,将高程数据导入FEFLOW后建立如图6所示的研究区三维地质模型,其中结点数121,608个,有限单元数180,411个。
![](//html.hanspub.org/file/12-2890287x28_hanspub.png)
Figure 6. Grid partition and model structure of the study area
图6. 研究区网格剖分与模型结构图
4.3. 水文地质参数的确定
矿区已完成的各类勘察及勘探工作中包含大量的水文地质试验,获得了较丰富的试验成果。本次模拟工作所用到的初始水文地质参数主要依据这些已有资料,同时根据研究区水文地质概念模型,将研究区渗透系数进行合理的统计分区,水文地质初始参数取值如表1所示。
![](Images/Table_Tmp.jpg)
Table 1. Initial hydrogeological parameters
表1. 评价区水文地质初始参数取值表
本次预测评价中初始流场采取的技术方法是将研究区水文地质参数分区并输入初始参数进行模拟,通过稳定流数值计算模型得到研究区稳定流条件下的天然流场,然后根据实际观测水位对天然流场下的水文地质参数进行参数反演,从而验证该模拟流场能否准确地体现研究区实际的水文地质环境 [11] 。
本次模型反演工作采用评价区工程勘察钻孔及水文地质钻孔实测数据,采用PEST (Parameter ESTimation)的参数估算模块,对各层的Kxx、Kyy、Kzz和入渗系数等共13组水文地质参数进行反演计算,得到最优水文地质参数(见表2),在此基础上确定模拟区地下水初始流场的水头(位)分布。
![](Images/Table_Tmp.jpg)
Table 2. Parameters after the PEST inversion
表2. PEST反演后参数取值表
由于水动力弥散尺度效应的原因,难以通过野外或室内弥散试验获得真实的弥散度。因此,根据Geihar等(1992)在世界范围59个大区域弥散数据资料统计分析的成果上 [14] ,本次评价区范围对应的弥散度应介于10~100之间,按照风险最大的评价原则,本次模拟的溶质运移模型参数如表3所示。
![](Images/Table_Tmp.jpg)
Table 3. Parameters for the solute transport mode
表3. 移模型参数表
5. 模型的验证
图7展示了校正后的地下水初始流场,将模拟水位值与矿区内具有的22个实测水位值进行拟合分析,结果见图8。通过图8的散点分布图可以看到,22个模拟值基本均匀分布在标准线两侧,说明数值模拟结果与实际值分布趋势基本相吻合,其中标准线上方水位点分布相对较多,体现了大部分模拟值比实际测量值略高,但所有数值模拟值基本都分布于90%置信区间内,反应了数值模拟值与实际测量值拟合情况较好,初始流场水位拟合折线图也完美体现了数值模拟值与实际测量值总体变化趋势的一致性。
通过上述试验数据对比,反应了本次建立的地下水系统概念模型基本吻合研究区实际水文地质环境,故利用该模型得到的水头值作为非稳定流的初始流场并以此为基础对研究区地下水环境影响进行预测评价是合理的。
![](//html.hanspub.org/file/12-2890287x30_hanspub.png)
(a) 第一层流场(基岩风化裂隙潜水含水层) (b) 第三层流场(基岩构造裂隙含水层)
Figure 7. Steady flow field in the simulation area
图7. 模拟区稳定流场示意图
![](//html.hanspub.org/file/12-2890287x31_hanspub.png)
Figure 8. Fitting scatter diagram of the initial flow field
图8. 初始流场水位拟合散点分布图
6. 情景设定
为了研究矿区在有无抽排水条件下多源污染物的迁移规律,设置以下两种情景:
情景一:天然防渗且在矿区无抽水条件下,沙冲沟、执茶沟、小段沟尾矿库,太金宝、思安口水库,排土场,1、2号浊水库同时发生渗漏,共渗漏28年,渗漏污染物均为Ni,污染物初始浓度为0.008 mg/L。
情景二:在天然防渗条件下,矿区分两期抽水:第一期12年,4380天,疏干面积为第一期露采区面积1.06 km2,第一期初始疏干强度6300 m3/d;第二期16年,5840天,第二期露采面积2.31 km2,第二期初始疏干强度9000 m3/d。在矿区抽水排水的同时,沙冲沟、执茶沟、小段沟尾矿库,太金宝、思安口水库,排土场,1、2号浊水库同时发生渗漏,共渗漏28年,渗漏污染物均为Ni,污染物初始浓度为0.008 mg/L。
7. 预测结果
研究区风化裂隙潜水主要接受大气降水的入渗补给,部分接受地表水的下渗补给,风化裂隙形成径流通道,补给的地下水沿其径流通道下渗与运移,从而形成地下径流,径流方向与地形有密切关系,总体由地形高处往低处径流,即在评价区内主体自西北向东南径流,局部由于地形变化,形成次级地下水汇流区。
(a)
(b)
Figure 9. Pollutant migration results without drainage
图9. 矿区无抽水条件下污染物迁移结果图
由图9可知,在矿区无抽水条件下,多源污染物的迁移基本随着天然地下水的径流方向迁移与扩散,污染物的迁移扩散主要受地形控制,由地形高的地方向低的地方迁移,且有效污染晕范围有限,对研究区地下水环境影响有限
由图10和表4的结果中可以看出,在矿区抽水条件下,与矿区临近的污染源不仅仅受到地下水流向的控制,主要还受控于矿区的抽排水,由于矿区的抽排水,使得矿区的地下水水位比临近的地方低,改变了矿区局部的地下水流场,使得临近的污染源的污染物向矿区迁移与扩散,从剖面图看出污染晕没有进入下层构造裂隙水含水层,仅仅污染上层风化裂隙含水层。由于合理地设计矿区的抽排水量,矿区开采对地下水水位和水质影响范围有限,但仍然需要针对模拟结果,对矿区做合理且有效的监控和管理。
8. 长效监测方案设计
长效监测方案布置的主要依据为:对矿区建成后可能造成水质污染的工程区,如尾矿库、浊水库、环保水库等,于其下游设置地下水水质监测孔;矿山建设过程中,由于采矿场疏干降水过程中对其周围水位带来影响,应于其主要影响范围内设置水位监测孔。
根据对前述对矿区水文地质条件的理解以及对污染物来源与迁移特征的认识,结合地下水模拟结果的分析,并根据在有抽水条件下,污染物迁移风险最大的位置下游设置监测井;并按照偏保守的原则,根据污染物在平面以及剖面上的迁移距离,定量设置监测井的位置以及井深,具体位置见图11,监测井的结构见表5。用于监测矿区内及影响范围内地下水,所有长期监测孔的监测项目都包括水位与水质。
(a)
(b)
(c)
Figure 10. Pollutant migration results with drainage
图10. 矿区有抽水条件下污染物迁移结果图
![](Images/Table_Tmp.jpg)
Table 4. Pollutant migration distance metre
表4. 污染物迁移距离表
![](//html.hanspub.org/file/12-2890287x37_hanspub.png)
Figure 11. Location of groundwater monitoring wells
图11. 地下水监测孔分布位置图
![](Images/Table_Tmp.jpg)
Table 5. Position and depth of monitoring wells
表5. 监测井位置及深度
据此,本监测方案分别于尾矿库、浊水库、环保水库下游设置SWZK01、SWZK02、SWZK09、SWZK04、SWZK10、SWZK11、SWZK12等7个钻孔开展水质监测工作,考虑到研究完整性等因素,对该批地下水监测孔均进行水质、水位监测,监测对象均为基岩裂隙水。
9. 结论
本文选取某铜钼大型矿区为研究对象,基于矿区水文地质条件,合理概化水文地质概念模型,采用FEFLOW软件进行数值建模,以在矿区有无抽排水条件下,多源污染物的迁移规律预测分析为目标,研究矿区多场耦合下污染物的迁移规律,得出以下结论:
1) 在矿区无抽水条件下,多源污染物的迁移基本随着天然地下水的径流方向迁移与扩散,污染物的迁移扩散主要受地形控制,由地形高的地方向低的地方迁移,且有效污染晕范围有限,对研究区地下水环境影响有限。
2) 在矿区抽水条件下,与矿区临近的污染源不仅仅受到地下水流向的控制,主要还受控于矿区的抽排水 [15] ,由于矿区的抽排水,使得矿区的地下水水位比临近的地方低,改变了矿区局部的地下水流场,使得临近的污染源的污染物向矿区迁移与扩散,且污染晕没有进入下层构造裂隙水含水层,仅仅污染上层风化裂隙含水层。说明由于合理地设计矿区的抽排水量,矿区开采对地下水水位和水质影响范围有限。
3) 根据建立的数值模型,对矿区在抽排水条件下多源污染物的迁移进行模拟,根据模拟结果,对地下水监测井进行定量布设,为矿山的监测和管理提供了有力的科学依据。
4) 对于类似矿山的多源污染风险评价的,建议根据研究区的具体水文地质条件,概化出合理的水文地质概念模型,利用数值模拟手段进行风险评价,并根据模拟结果,给出长效综合监控方案的建议。