数值模拟在济宁铁矿矿坑涌水量预测中的应用

王仕昌,邱恺毅

(山东省物化探勘查院,山东 济南 250013)

济宁铁矿位于山东省济宁市兖州区、任城区、汶上县交界地带，面积约46km2。济宁铁矿赋存于新太古代济宁群变质岩中，埋深多在1000m以下，经勘查查明铁矿石资源量约18.34亿t[1]。因埋藏深，开采难度大，为开采可行性研究提供水文地质数据，施工了3个水文地质孔，进行了抽水试验，求取了水文地质参数，采用数值模拟法[2-3]，对矿坑涌水量进行了预测[4-5]。

1 矿区水文地质条件

矿床顶板为寒武纪馒头组地层，矿床周边为济宁群变质岩地层，岩石坚硬、完整，富水性、导水性较差，矿床底部仍为济宁群变质岩地层[6-8]，随着深度增加富水性及导水性趋于微弱，矿床底部可视为相对隔水层。

2 数值模拟预测矿坑涌水量

2.1 水文地质概念模型

2.1.1 模型模拟范围

模型平面范围选择依据①济宁矿区属于湖东汶泗河冲洪积平原中部水文地质单元；②铁矿资源储量范围远小于矿区范围，可以保证模拟范围在抽水试验时影响不到模拟边界；③能将区域内切割深入矿区的F1断裂较完整地包含于模型模拟范围内，方便模拟断层对抽水试验的影响[9]。

模型垂直范围模拟含矿层位济宁群变质岩裂隙含水层，模型上顶界为YK1第3次抽水试验起点-1113m，下底界定为YK1第6次抽水试验底界下移100m，即-1613m处。

2.1.2 地下水动力场及流动特征

济宁群变质岩裂隙水补给、径流、排泄条件差，根据钻孔钻探情况，地下水赋存于千枚岩破碎带中。千枚岩结构密实，孔隙、裂隙不发育[10-11]。根据抽水试验，断层上盘灰岩、白云岩裂隙岩溶含水层的地下水通过F1导水断层流向济宁群变质岩裂隙含水层。地下水运动为三维流；地下水水位是随时间变化的，由于介质的非均匀性造成水文地质参数随空间变化，体现了系统的非均质性。初始流场[12]根据抽水试验的静止水位给出，概化模拟含水层天然状态下水量保持稳定。

2.1.3 水文地质参数及其灵敏度分析

该次研究对象为深层的济宁群变质岩裂隙承压含水层，含水层的水文地质参数主要为渗透系数[13-14]。根据YK1钻孔四次抽水试验拟合求得的济宁群变质岩含水层的渗透系数分别为0.0055m/d，0.0083m/d，0.000021m/d，0.0046m/d。

水文地质参数的不确定性[15]往往导致计算水头或流量的不确定性，影响矿区涌水量的可靠性。因此应用灵敏度分析，对识别验证后的水文地质参数进行不确定性分析[8]，为研究矿区涌水量提供依据。

灵敏度可度量一种因子的变化对另一种因子的影响程度，某一模型因变量对模型输入参数的灵敏度可用因变量对该输入参数的偏导数来表示,即:

Xi,k=∂yi/∂ak

式中:Xi,k—模型的因变量对第i观测点,第k参数的灵敏度系数。由于不同的参数单位不同，这样不同参数的灵敏度系数就无可比性。此时,需对上式进行标准化处理,即:

根据上式计算得出灵敏度系数(表1)，从表1可以看出，渗透系数的变化对模型的运行结果有一定影响。

表1 灵敏度系数

2.2 数值模型

2.2.1 数学模型

对于建立的矿区水文地质概念模型，可用如下微分方程的定解问题来描述：

式中：Ω—渗流区域；h—地下水系统的水位标高(m)；K—含水介质的水平渗透系数(m/d)；Kz—含水介质垂向渗透系数(m/d)；ε—含水层的源汇项(1/d)；Ss—含水层储水率(1/m)。

2.2.2 子程序包的选择

运用基于有限差分法的GMS程序求解以上微分方程的定解问题。根据矿区的水文地质结构等特征，选择GMS中相应的子程序包来实现地下水流的模拟。

(1)LPF-层状特征流量程序包。该功能相当于计算单元间渗流子程序包，它对各个模型参数作出了相应的假定，选定该程序包可以通过参数分区方式赋值参数。

(2)WEL1-井流子程序包(WEL)。模型约定，在每个应力期，以指定流量从含水层抽水或向含水层注水；井流量不受井所在计算单元的大小及水头影响；负的流量值表示抽水井，而正的流量值则表示注水井。

(3)GHB-通用水头子程序包。将流向通用水头边界的流量项加进有限差分方程组。

2.2.3 数值模型

(1)空间离散[16]

将模拟区剖分为104行，81列，共计4层，各层均采用100m×100m的剖分格式。4层有效单元格共33696个。

(2)模拟期确定

该次地下水数值模拟模型的模拟期为YK1第一次抽水试验开始至最后一次抽水试验结束水位恢复。在抽水期和水位恢复期，每小时作为一个应力期，分层抽水的间隔时期，时间节点以小时的倍数逐步过渡增长：2h,4h,8h……共计233个应力期。

(3)定解条件的处理

初始条件：由于缺乏地下水流场，加之模拟含水层为深层承压水，循环条件差，遂根据分层抽水试验的静止水位分别给出统一水位作为每层初始流场。

边界条件：全部人为边界处理为通用水头边界。

2.2.4 模型拟合情况

根据YK1钻孔的第3～6次分层抽水试验的情况，拟合结果如图1～图4所示。

图1 YK1第3次抽水试验拟合图

图2 YK1第4次抽水试验拟合图

图3 YK1第5次抽水试验拟合图

图4 YK1第6次抽水试验拟合图

根据拟合图看出，YK1钻孔的第3～5次抽水试验的拟合效果均较好，第6次抽水试验的拟合有所偏差，分析原因是将模拟的含水层底板放低100m所致。

2.3 矿井涌水量预测

2.3.1 根据等效渗透系数原理合并模型为单层

矿层含水层为济宁群变质岩裂隙含水层1层，根据等效渗透系数原理将原本垂向划分为4层的模型合并为1层。

渗流方向与层面水平时，设其平行于层面的渗透系数为KP，当单一分层的渗透系数Ki和厚度Mi已知时，可求出KP。当水流平行于层面时，通过层状含水层总的单宽流量q等于各分层的单宽流量之和，总厚度M等于各分层厚度之和。对于每一分层而言，水力坡度J均为ΔH/l，因此，每一分层的流量为：

如果用等效的均质含水层代替层状岩层，显然等效层的厚度等于层状岩层的总厚度，并且在同一水力坡度ΔH/l作用下应当有相同的流量q：

由此得：

因而求得平行于层面方向的等效渗透系数为：

济宁群变质岩裂隙含水层的最终水文地质参数分区如图5所示。

图5 渗透系数分区图

2.3.2 涌水量预测

(1)疏干排水量

一般情况下，对于给定的疏干降深，可以分别求出不同的疏干流量所对应的疏干时间，然后作出不同疏干流量的降深与时间关系曲线(图6)和指定疏干水平的流量与时间关系曲线(图7～图10)，据此选择最佳的疏水量与疏干时间[17-19]。

图6 不同疏干流量的降深与时间关系曲线图

图7 开采水平-1200m的Q-t曲线

图8 开采水平-1300m的Q-t曲线

图9 开采水平-1400m的Q-t曲线

图10 开采水平-1500m的Q-t曲线

由图7～图10可以看出，开采水平为-1200m，-1300m，-1400m，-1500m所对应的最佳疏干时间为20.00h,25.23h,31.70h,39.73h；疏干排水量为3250m3/d。

(2)正常涌水量

将疏干工作面设置为定水头，用稳定流场求出矿井涌水量(表2)。

表2 正常涌水量预测

将模型各参数分别增减10%后进行稳定流模型正常涌水量的预测，可以得到正常涌水量的1个区间值(表3)。

表3 正常涌水量区间值

(3)最大涌水量

最大涌水量的预测考虑到矿体埋藏深度大，结合该次勘查阶段的精度低于C级、高于D级，而且未做群孔抽水试验，根据《地下水资源分类分级标准》(GB15218-94)附录C[20]，C级计算的允许误差为±35%， D级计算的允许误差为±50%，因此为提高安全程度，按D级允许误差±50%计算，最大涌水量为正常涌水量的1.5倍计算，结果见表4。

表4 矿井最大涌水量

3 结语

济宁铁矿埋藏于济宁群变质岩中，埋藏深度在1000m以下，该次通过数值模拟预测了-1200m,-1300m,-1400m,-1500m开采水平的正常矿坑涌水量分别为998.56m3/d,1079.27m3/d，1159.98m3/d,1240.68m3/d，但由于济宁铁矿赋存条件复杂，矿坑涌水量的预测较难，其预测数值与实际数值的一致性有待于进一步探讨研究。