新疆北部某煤矿持续采掘条件下矿井涌水量动态模拟预测

杨 晟，贾 超，袁 涵，张建配，陈 阳

(1.山东大学土建与水利学院，山东 济南 250061；2.山东大学海洋研究院，山东 青岛 266237)

矿井涌水量是衡量矿井开采方式合理性和安全性的重要指标，直接关系到矿井的安全生产。因此矿井涌水量的准确预测对于矿井开采和生产具有重要意义[1]。矿井开采是一个动态的过程，随着开采面的变化，对应的开采区疏水面积和疏水层位也发生变化，因此矿井的涌水量也是动态变化的。而以往预测矿井涌水量时往往较少考虑矿井的开采进度和开采面的变化，因此导致预测矿井涌水量与实测值差距较大，且无法反映矿井涌水量的动态规律，因此考虑矿井采掘进度的动态涌水量预测会更加真实准确地反映矿井开采过程中涌水量实际情况。

目前矿井涌水量的预测方法主要有水文地质比拟法[2]、时间序列法[3-4]、解析法[5]、大井法[6]、神经网络法[7-8]、数值法[11-15]等。其中数值法可以模拟矿区复杂的三维地质结构，并且可以概化仿真采掘进度，使涌水量的模拟预测更为精确。因此以新疆北部某矿区为例，在利用数值法的基础上考虑矿区的采掘进度，运用GMS数值模拟软件对矿区涌水量进行动态模拟预测[14]，预测结果与实际情况较为符合。

1 研究区概况

研究区位于新疆北部吐鲁番盆地一带，地形属于西北边缘低山丘陵地带，矿区内绝大部分为斜平地，坡度50‰，北高南低。区内干旱少雨，降雨对地下水补给微弱，侏罗系基岩裂隙含水层主要接受来自北部山前雪山融水的深层补给，并向南部以地下径流的方式排泄。地下水流动方向主要为东北向西南。区域广泛覆盖第四系松散堆积物，主要断层导水性较弱，地层构造较为简单。区内主要地层自下而上分别为石炭系中统、三叠系上统、侏罗系下统八道湾组、侏罗系下统三工河组、侏罗系中统西山窑组、侏罗系中统头屯河组、第三系及第四系。矿区的含水层组主要为第四系松散孔隙含水层、下第三系陆相半胶结孔隙裂隙含水层、侏罗系基岩裂隙含水层和石炭系含水层。主采煤层赋存于侏罗系西山窑组含水层，含水层岩性为灰白色粗砂岩-细砂岩与灰色粉砂岩、泥岩互层，偶夹砂砾岩，含5号煤组，平均总厚度12.1 m，含煤系数5.38%。矿区主要分北一采区、北二采区和北三采区三个采区，采区内开采顺序为后退式，即由井田边界向上山方向回采，区段及水平开采顺序为下行式，即先采上水平区段，后采下水平区段。各采区分布如图1所示。

图1 采区分布图Fig.1 Mining area distribution map

2 水文地质模型的建立

2.1 模拟范围及边界条件的概化

矿区北部边界为山前地下水径流补给边界，作为定流量边界；南部边界为地下水排泄边界作为定流量边界；东西边界为与流线垂直的边界，作为隔水边界处理；模型顶部覆盖大面积第四系松散堆积物，接受大气降雨入渗，概化为补给边界；煤层所在西山窑组含水层以下紧邻着三工河组黑色粉砂岩致密隔水层，岩层致密，隔水性能极好，将其作为模型的底部边界，概化为隔水边界。

2.2 模拟区源汇项的概化

1) 降雨蒸发。矿区降雨稀少，根据矿区附近气象站的监测数据，矿区年平均降水量为41.4 mm，地表第四系岩性主要为中粗砂，根据地下水流数值模拟技术要求，降雨入渗系数取经验值0.3，通过Recharge模块进行概化(表1)。

矿区蒸发强度较大，年均蒸发量为1 500 mm，地下水水位埋深普遍大于5 m，根据地下水流数值模拟技术要求，蒸发系数取经验值0.01，通过Evaporation模块进行概化(表2)。

表1 区域多年平均降雨量表Table 1 Regional average annual rainfall

表2 区域多年平均蒸发量表Table 2 Regional average annual evaporation

2) 边界补给和排泄。研究区地下水补给主要来自于北部山前雪山融水的地下侧渗补给，并由南部地下水流出边界排泄。根据断面法结合达西定律初步估计北部补给边界流量约为q补给=1 489.64 m3/d；南部排泄边界的流量约为q排泄=-1 242.15 m3/d，同时流量值在后期流场校核检验过程中再不断调整优化以使其符合实际情况。

2.3 地层结构概化

依据矿区水文地质资料及钻孔勘探数据，将矿区水文地质模型地层结构划分为6层含(隔)水层(图2)，自上而下主要为第四系松散孔隙含水层(H-I-1)、下第三系陆相半胶结孔隙裂隙含水层(H-I-2)、侏罗系头屯河组承压裂隙含水层(H-Ⅱ-1)、西山窑组上段5-3煤层上部至5-1煤层底板间承压裂隙含水层(H-Ⅱ-2)、西山窑组下段5-1煤层底板以下粗粒相承压裂隙含水层(H-Ⅱ-3)、三工河组粉砂岩致密隔水层(G-Ⅱ-1)。各地层在空间上由于岩性的差异，水文地质参数表现出差异性，因此将各地层概化为非均质各向异性地层，模型概化为6层非均质各向异性三维非稳定流模型。

图2 区域地层综合柱状图Fig.2 Regional stratigraphic comprehensive borehole histogram

3 数学模型

根据矿区水文地质概念模型，矿区地下水渗流符合达西定律，结合连续性方程可建立矿区地下水非稳定流数学模型见式(1)[16]。

式中：H为水头，m；Kxx、Kyy、Kzz为x、y、z方向上的渗透系数，m/d；n为法向向量；Ss为弹性释水率，1/d；W为降雨入渗补给强度和蒸散发强度，m3/d；H0为模拟区初始流场，m；q为第二类定流量边界流量，m3/d，流入为正，流出为负；Ω为模拟区。

4 数值模型建立

4.1 模拟区域网格剖分及时间步长的确定

在满足模型精度要求前提下，将模型在横向剖分300个网格，纵向剖分280个网格，垂向有6层含(隔)水层，因此按6层剖分，最终得到300×280×6个单元格(图3)。

模拟过程中选取2007年1月～2007年12月作为模型的识别期，识别调整模型水文地质参数；以2008年1月～2008年12月作为模型的检验期，检验计算流场和水位的准确性；矿区三个采区采掘期为2009年1月～2019年1月，将2009年1月～2020年1月作为涌水量动态模拟期，时间划分以月为步长。

4.2 水文地质参数分区

根据矿区水文地质资料，将各地层参数分区概化如图4所示，同时各分区的水文地质参数初始值由钻孔抽水试验等资料初步给定，后期根据各井点地下水位实测值进行校核调参确定最终校正值。

图3 模型网格剖分图Fig.3 Mesh of the numerical simulation model

图4 各含(隔)水层水文地质参数分区图Fig.4 Aquifers hydrogeological parameter zoning map

4.3 采掘空间的概化

为了较为准确地刻画出矿井采掘过程中的涌水，运用drain模块刻画矿井的排水[17-18]，三个采区分为10个工作面，每个工作面按计划在一年内完成开采。开采前先对工作面单元格进行短期疏水，将水位迅速下降至工作面的安全开采平面，然后进行开采。对于已开采和正在开采的工作面均设置排水(drain)边界，保持地下水位低于该工作面的安全开采平面，未开采的工作面则不进行排水，通过这样的方式可以较为合理准确地刻画矿区开采过程中地下水流场的变化。

模型排水(drain)边界的计算见式(2)。

(2)

式中：QD为排水量，m3/d，即由含水层流入drain的水量；CD为采取渗透性系数，m2/d；HD为排水标高，m。

4.4 模型的识别与验证

运用矿区水位观测孔2007年的水位监测数据对模型计算的流场进行识别，最终调整得到矿区水文地质模型各分区水文地质参数见表3。

将识别后的水文地质参数和概化的水文地质条件代入模型，最终计算得到天然状态下地下水流场计算结果如图5所示，计算流场和实测流场在水位、流向、流态上基本一致，计算流场能基本反映实际流场的宏观形态。

表3 模型水文地质参数一览表Table 3 Hydrogeological parameters of the hydrogeology model

图5 计算流场与实测流场对比图Fig.5 Comparison diagram of calculated and observed flow field

选取矿区6号观测孔为例，用2008年检验期地下水位动态监测数据对模型计算水位进行动态检验，水位拟合曲线如图6所示，结果显示模型计算水位与观测孔实测动态水位基本拟合，计算水位可以基本反映实际水位的动态变化规律，模型整体可靠度较高，可作下一步预测分析。

图6 地下水水位计算值与实测值拟合图Fig.6 Fitting curves of observed and calculated groundwater level

5 矿井涌水量的动态预测

5.1 考虑采掘进度的矿井涌水量动态模拟预测

运用识别检验后的数值模型对矿区三个采区的涌水量进行下一步的模拟预测。三个采区按照北一采区、北二采区、北三采区顺序依次开采，北一采区安全开采面底板高度+450 m，北二采区安全开采面底板高度+430 m，北三采区安全开采面底板高度+400 m，按照由高到低、自上而下的原则开采。

北一采区分三个工作面，按计划三年开采完毕；北二采区分为四个工作面，按计划四年开采完毕；北三采区分为三个工作面，按计划三年开采完毕。对于已开采及正在开采的工作面，地下水位疏干至安全开采平面以下，还未开采的工作面则不设置排水。按照采掘先后的顺序模拟预测，煤层所在西山窑组含水层的水位变化如图7～9所示。

图7 北一采区疏干后水位等值线图Fig.7 Contour map of groundwater after dredging in north 1 mining area

图8 北一采区、北二采区疏干后水位等值线图Fig.8 Contour map of groundwater after dredging in north 1 and north 2 mining area

图9 三采区全疏干后水位等值线图Fig.9 Contour map of groundwater after dredging in three mining areas

通过Flow budget模块计算出了三个采区的矿井涌水量随时间的动态变化曲线如图10所示。从计算结果来看，在一个采区的开采过程中，当一个工作面开采完成，开始开采下一个工作面时，由于开采新的工作面需要迅速将水位疏降至安全开采面以下，因此短时间大量排水，涌水量突增，随后将水位维持在安全开采面以下，涌水量会慢慢降至平稳。最终随着开采面的不断扩大，总体上涌水量在不断增大。此外由于北二采区、北三采区邻近北一采区，因此其疏水过程中将分担北一采区部分涌水量，从而使得北一采区涌水量进一步减小。而由于北三采区与北二采区间隔一定距离，因此北三采区疏水对北二采区的涌水量影响较小，北二采区涌水量下降幅度也较小。由此说明相邻采区在疏水过程中由于流场的相互干扰和影响，因此涌水量也会互相影响，邻近采区的疏水会分担一部分先前开采区的静储水量，使得先采区的涌水量减小。而随着流场的稳定，各采区的涌水量也逐渐趋于稳定。最终计算得到模拟期末北一采区涌水量为211 m3/d、北二采区涌水量为613 m3/d、北三采区涌水量为563 m3/d。

图10 采掘期内三采区计算涌水量动态变化曲线Fig.10 Curves of calculated dynamic mine inflow in three mining areas during the mining period

表4 采掘期逐年平均涌水量统计表
Table 4 Statistical table of average annual mine inflow during the mining period

年份20092010201120122013平均涌水量/(m3/d)309457556763899年份20142015201620172018平均涌水量/(m3/d)9761 0531 2791 4011 456

同时为了探究矿区开采过程中总体涌水量的变化规律，列出开采过程中逐年的平均涌水量见表4，并作变化曲线如图11所示，结果表明矿区涌水量随着开采的深入，总体呈现上升趋势，但中间各个采区过渡时涌水量有较大的增幅，呈现阶梯状上升，原因是由于北一采区至北三采区的安全开采平面高度呈梯度下降，因此由已采区向未采区过渡时，疏水降深也随之增大，涌水量会有短期较大幅度的增长，但随着新采区的逐步开采，涌水量的上升幅度逐渐变缓，涌水量慢慢趋于稳定。

5.2 未考虑采掘进度的矿井涌水量预测

不考虑采掘进度的情况下，则不考虑开采的实际顺序，将三个采区的地下水位同时疏降至各采区的安全开采平面以下，三个采区涌水量动态变化曲线如图12所示。三个采区初期均大量排水，后期随着疏水的稳定，涌水量慢慢减小并趋于稳定。最终计算得到模拟期末北一采区涌水量为236 m3/d、北二采区涌水量为671 m3/d、北三采区涌水量为611 m3/d。

图11 采掘期内逐年平均涌水量变化曲线Fig.11 Curve of annual average mine inflow during the mining period

图12 采掘期内三采区计算涌水量动态变化曲线Fig.12 Curves of calculated dynamic mine inflow in three mining areas during the mining period

5.3 大井法计算矿井涌水量

根据矿井的水文地质资料，采用大井法计算三个采区涌水量。将各采区巷道系统的排水视为理想大井，周界为统一的降落漏斗，由于侏罗系西山窑组含水层本为承压含水层，随着开采疏水作用，采区水位慢慢下降，含水层水位下降至顶板以下，开采区域由承压转无压，因此采用承压转无压完整井公式(式(3))。最终得到大井法计算的矿井涌水量,计算成果见表5。

(3)

式中：K为渗透系数，m/d；H为自然水位至含水层底板距离，m；M为含水层厚度，m；h0为疏干水位至含水层底板距离，m；R0为引用影响半径，m；r0为引用半径，m；S为水位降深，即自然水位至疏干水位的距离，m；P为预测区面积，km2。

5.4 分析结果

最终将考虑采掘进度预测矿井涌水量、未考虑采掘进度预测矿井涌水量、大井法计算矿井涌水量和实测涌水量进行对比，将各种方法计算涌水量及其与实测涌水量的相对误差列于表6。由对比结果分析，考虑采掘进度模拟矿井涌水量最为接近矿井的实测涌水量，相对误差最小，因此考虑采掘进度模拟矿井涌水量相对其他方法精度更高，并且可以反映矿区涌水量的动态变化规律，相比其他方法具有一定的优越性。

表5 大井法矿井涌水量计算成果表Table 5 Table of mine inflow calculated by large well method

表6 各方法预测矿井涌水量统计对比表Table 6 Statistical comparison table of the mine inflow predicted by various methods

6 结 论

1) 开采新工作面时由于需短期将水位疏排并维持在安全开采面以下，涌水量先期会突增，而后保持工作面水位不变进行排水，涌水量慢慢降至平缓。新采区的开采疏水会分担相邻采区的静储水量，导致邻近采区的涌水量减小，采区距离越近，互相影响越为明显。

2) 矿区涌水量随着开采的进行总体上是逐年上升的。在旧采区向新采区过渡时由于安全开采面的下降，涌水量会有短期较大幅度的增长，后期随着该新采区的开采，涌水量增长变缓，慢慢趋于稳定。

3) 最终将各方法计算的涌水量和实测涌水量进行对比分析，发现考虑采掘进度计算得到的矿井涌水量与实测值相对误差最小，相比于其他方法准确性要高，能更真实反映矿井涌水量情况，为矿井后期的水灾防治和安全生产提供一定的参考依据。