不同水压加载下煤体裂隙渗流数值模拟

李治豪, 陈世江*, 冯宇迪, 郭文彬, 姬长兴

(1.内蒙古科技大学矿业研究院, 包头 014010； 2.呼伦贝尔学院， 内蒙古自治区矿山压力重点实验室, 呼伦贝尔 021008； 3.内蒙古福城矿业有限公司, 鄂尔多斯 016299)

煤矿地下水丰富且复杂,随着煤矿开采的深入,水的渗流对煤层开采的影响愈发明显,关系着煤田的勘探和煤矿的安全生产。地下水在煤岩的渗流中,水压会影响煤基质本身的收缩效应和应力效应,对渗流特性产生重要的影响[1]。王刚等[2]应用电子计算机断层扫描(CT)三维重构技术和 ANSYS 软件对煤体的微观孔隙结构进行了数值模拟,并得出温度载荷下孔隙结构的变形差异较大。郑司建等[3]通过高速离心和低场核磁共振实验,得出了低阶煤渗流实验束缚水的最佳离心力和测得低阶煤的平均孔隙半径。王开德等[4]以广义达西定律为基础的数学模型,利用COMSOL软件对煤层注水流固耦合模型进行计算,得出了孔隙压力的分布情况及渗流速度的最大值。白若男[5]采用X射线三维显微镜和应用三维重构技术和FLUENT有限元软件模拟不同压力的水渗情况,得到了煤体微孔隙内的渗流速度和通道压力分布规律。张永将等[6]提出了高压水射流环切割缝煤层自卸压增透技术,并且建立普通钻孔及环割钻孔瓦斯流动微分方程,应用FLAC3D软件分析了高压水射流割缝后的煤体的应力演化规律。李波波等[7]利用热-流-固耦合三轴伺服渗流装置,开展了不同含水条件下孔隙压力升高过程中煤岩渗透特性的试验研究,得到了含水率恒定,孔隙压力升高时,煤壁介质的吸附能力先是缓慢增大,增大到一定程度后趋于平缓。王浩等[8]通过不同气体压力下的煤岩瓦斯渗流实验,得出外部压力和孔隙压力是引起煤体渗透率发生变化的主要原因。

目前众多学者对于煤体渗透特性规律开展了很多研究,并且取得了丰硕的成果[9-10]。在前人的研究基础上,进一步在细观方面研究了煤体的渗透特性,借助纽迈核磁共振和COMSOL软件研究在不同水压加载下煤体内部渗流变化以及煤体内部渗流运动物理场的耦合变化情况。

1 渗流运动方程

在煤体介质内,由于渗透系数较小,在数值模拟中认为渗流过程水是仅在裂隙内运动的。在地下水渗流力学中,煤体的渗流运动符合达西定律,水头的损失和流速呈线性关系。

煤体裂隙水流层流的计算公式为[1]

(1)

相应的线性定律上线雷诺数为

Rekp=600[1-0.96(Δ/b)0.4]1.5

(2)

式中:q为渗流量，cm3；g为重力加速度，m/s2；J为与裂隙平面平行的水力梯度；b为裂隙宽度，cm；Δ为裂隙壁面绝对粗糙度，cm；μ为水的运动黏滞系数，cm2/s。

在水渗入煤体裂隙的过程中,由于竖直向下的水压作用,会引起煤岩内部裂隙的体积变化。渗透过程中,水压的变化会引起煤体裂隙介质体积的变形,煤体在发生破坏的情况下,宏观变形并不明显,微观变形则是裂隙孔隙的大小和基质体积的变化[4]。假设煤体为理想的线弹性体,忽略地温条件下,考虑水应力、煤体内壁吸附解吸膨胀引起的体积应变。煤体从结构可以划分为煤骨架、宏观裂隙、微观孔隙和节理,数值模拟主要研究宏观贯通裂隙的渗流运动规律。

2 数值模型建立

COMSOL Multiphysics软件是一款基于有限元的大型高级数值仿真软件,将采用该软件进行煤体渗流模拟分析。该软件广泛应用于各个领域的科学研究以及工程计算,模拟工程和科学领域的多数物理过程。COMSOL Multiphysics软件集前处理器、网格剖分、求解器和后处理于一体,具有用途广泛、灵活易用的特点,较其他有限元分析软件的优势在于,利用附加的功能模块,能够很容易的进行扩展。被誉为“第一款真正的任意多物理场直接耦合分析软件”。为了使所建模型更加逼近煤体真实细观结构特征,应用纽迈核磁共振扫描与图像重构技术建立煤体微孔隙渗流模型。

2.1 煤体裂隙细观几何模型

为便于计算和分析,对煤体渗流实验的试样进行扫描,截取较为典型的一部分裂隙通道,并且将纽迈核磁共振扫描裂隙图像导入CAD中进行图像处理,将调整好的裂隙图像导入COMSOL数值软件的几何模块中,导入的图形为高207 μm,宽 190 μm 的矩形,完成几何模型的建立。如图1、图2所示。由于裂隙渗流是在垂直方向表示,故可以用横截面的渗流二维模型来演示。将截取的裂隙外部边框设定为固定边界条件,渗流场入口处和出口处的边界设定为自由边界条件。

图1 核磁共振扫描裂隙图像Fig.1 Fissure image scanned by MRI

图2 裂隙图像CAD截图Fig.2 CAD screenshot of crack image

2.2 模拟条件

在软件模拟中选择COMSOL中的稳态物理场进行研究。在成功建立模型之后对贯通裂隙采用自由三角形网格进行划分,如图3所示,然后进行整体几何模型的网格划分,如图4所示,模型的尺寸为207 μm×190 μm的矩形。

图3 裂隙通道网格划分Fig.3 Meshing of fracture channels

图4 整体几何模型自由三角形网格划分Fig.4 Free triangle meshing of the overall geometric model

根据地应力的研究成果,其变化的一般规律为:地应力随深度正比例增长,1 MPa地应力约为100 m水深的压力。模拟地下煤层上部分别存在不同承压地下水的渗流情况。地下水以煤体的裂隙为通道发生渗流,裂隙上部施加压力为水压,且施加水压以变量为1 MPa逐级递增,裂隙模型两侧为固定边界,将两侧边界设定为固定载荷边界代表固定地层边界,不设定围压,模拟产生的围压效果仅为煤岩介质内渗流运动对固定边界产生的反作用力,如图5所示。煤体渗流数值模拟的各项参数如表1所示。

表1 数值模拟计算相关参数Table1 Numerical simulation calculation of relevant parameter

图5 裂隙渗流模型外部压力分布Fig.5 External pressure distribution of fracture seepage model

3 水压加载下煤体渗流数值结果及分析

根据扫描截取的渗流几何模型,设置上部仅有的一个裂隙口为渗流入口,下部的3个裂隙口为渗流出口,如图2所示。分别计算水压从1 MPa逐级加压到10 MPa时每级水压加载下的渗流结果。图6 为1 MPa水压下的渗流速度场结果。

图6 水压为1 MPa时的渗流速度场结果Fig.6 Seepage velocity field results when the water pressure is 1 MPa

3.1 不同水压下煤体渗流导水通道的分布

通过对渗流模拟结果的分析,裂隙通道存在贯通通道,同时也存在很多的非贯通通道。裂隙呈现出不规则且互相连通的状态,形成了较为良好的导水通道,从渗流速度场的速度大小和水流模拟箭头来看,非贯通通道的速度场在水压变化的过程中没有发生变化且渗流速度基本为零,所以认为非贯通通道在充满水形成饱和空间之后,不参与渗流运动,对煤体内部裂隙整体的渗流速度场没有影响。

3.2 不同水压下渗流场出口平均速度和内部最大流速分析

渗流完成后依次监测3个出口处的渗流速度的平均值,如表2所示。从绘制曲线(图7)可以看出,在水压逐渐增大的同时出口处的渗流平均流速也是在逐渐增大的,但是平均流速的增长率却是在逐渐减小,说明在水压逐渐增大的同时,渗流过程存在某种阻力影响了水渗流的速度,阻力的产生原因可以通过渗流场流速变化和裂隙介质通道的变形情况来进行分析。

图7 水压与渗流出口平均速度的关系Fig.7 Relationship between water pressure and average velocity of seepage outlet

表2 各级水压加载下渗流出口的平均速度Table 2 Average velocity of seepage outlet under pressure at all levels

对渗流过程中不同水压下的渗流最大流速分析,渗流通道的最大流速是随着初始水压的增长而增大的(表3),但是从绘制的折线图(图8)来看,水压和裂隙通道的最大流速的变化关系存在不同趋势的3个阶段：第一阶段,在水压达到4 MPa之前,水压和最大流速的增长曲线是递增的直线且呈现出正比例关系；第二阶段,在水压为4～5 MPa这个短暂的水压区间内,最大流速随着水压的增大曲线呈现出近水平的直线,这一阶段可以作为最大流速的增长率发生突变的过渡阶段；第三阶段,在水压达到5 MPa之后,水压和最大流速的增长曲线依然是成直线递增的,但是增长曲线的斜率明显小于第一阶段,也就是最大流速的增长速度是明显小于第一阶段。

表3 各级水压加载下渗流通道的最大流速Table 3 Maximum velocity of seepage channel under pressure loading at all levels

图8 煤体渗流场最大流速和水压的关系Fig.8 Relationship between maximum flow rate and water pressure in coal seepage field

在渗流过程中,考虑重力场的影响在内,在水压的作用下煤体上下存在水头差的影响,主导裂隙和孔隙中的水发生了渗流运动。渗流运动的渗流场产生了渗透体积力,渗透体积力产生后随之在煤体内部长时间存在,相当于内部载荷长时间作用在内部孔隙介质上,使得裂隙内部的应力场发生了改变,也就是流固耦合的作用效果。在第一阶段水压未达到4 MPa之前,渗流场所产生应力场不足够大,故而应力场对渗流运动的反作用力不够明显。而在第二阶段,在水压逐级增大至4 MPa之后,由于渗流场所产生的应力场对渗流运动的反作用力蓄势已达临界值,表现出突然显现的状态,但在4～5 MPa,由于此时为裂隙内部应力场产生反作用力效果显现较为迅速,表现为明显的临界反应,此时的渗流场和应力场需要短暂的时间来达到一个应力平衡,故此区间内的最大流速则会基本处于一个稳定值,表现为一个短暂的过渡期。而在第三阶段,水压超过5 MPa之后,水力耦合作用愈加明显且趋于稳定,应力场的阻力作用也逐渐趋于稳定,最大流速随着水压的增大仍然在逐渐增大,但是由于应力场产生的阻力也在随着渗流场稳定增长,故此阶段内最大流速的增长率是小于第一阶段的。

3.3 不同水压下煤体裂隙介质最大变形与渗流最大流速比较分析

煤体在渗流过程中的空间网格位移代表煤体介质的变形特征,通过模拟计算可以分别得到煤体裂隙介质在水压区间为1～10 MPa下变形的最大位移。表4为水压1 MPa时裂隙介质在X、Y方向的最大变形位移。从X方向最大位移来看,在水压达到3 MPa时的X方向最大位移会出现一个峰值,达到90.89×10-18μm,在此之后,水压逐渐增大,X方向的变形突然减小然后经过短暂的稳定之后又继续缓慢增大直至趋于基本稳定的过程。而要使裂隙介质的变形在达到峰值之后再下降,煤体内必然存在阻碍煤体变形的作用力。在水压逐级增大的过程中,Y方向的变形位移总体上呈现出逐渐增长的趋势,故推理在水压增长的过程中裂隙介质在纵向的变形是在不断加剧的,直至煤体发生破坏。

表4 不同水压下裂隙介质内X和Y方向最大位移变形Table 4 Maximum displacement deformation in X and Y directions in the fractured medium under different water pressure

通过裂隙介质最大形变与裂隙内部最大渗流速度的曲线(图9)对比分析,在水压载荷为1～3 MPa,水压逐级增大,内部裂隙介质在X方向的最大位移变形也在增大,同时最大渗流速度也在逐渐增大。在水压到达3 MPa之后,裂隙介质的最大变形突然出现了大幅度减小,并且在水压为4～5 MPa,裂隙介质的最大位移变形则变得较为平缓,同时裂隙通道内最大渗流速度的变化也变得较为平缓。在水压超过5 MPa后,裂隙介质在X方向的最大变形又开始缓慢增大,但是增长的速率却明显变小,且通过观察曲线图,在水压达到8 MPa后,裂隙介质的变形基本趋于稳定。

图9 裂隙的最大应变和渗流最大速度之间的关系Fig.9 Relationship between the water pressure and the maximum deformation displacement of the fractured medium

通过对水压产生的煤体介质的横向最大位移变形和渗流场最大流速分析,煤体变形与最大流速的变化趋势基本一致,但是两者之间存在一个短暂的时间差。裂隙内部的应力场发生改变之后,煤体介质内部结构也会随之发生改变,介质颗粒的空间位置和裂隙孔隙的大小就会存在一定程度的改变,进而影响渗流场的渗流速度。

综合分析,在渗流运动过程中渗流场由于流固耦合作用产生的应力场就是影响渗流的主要因素。

3.4 不同水压下煤体介质内部应力分布变化分析

根据在1～10 MPa水压下煤体内部渗流运动内部介质的应力场分布云图,如图10所示,可以看出,在水压从1 MPa增大至4 MPa时,应力云图分布显示应力较大的区域面积在逐渐减小,说明在此过程中流固耦合作用产生的应力场的作用效果在逐渐降低。在水压达到5 MPa之后,应力云图分布应力较大区域面积又在逐渐扩大,说明在水压达到 5 MPa 之后,应力场所产生的作用效果越来越明显,在水压达到8 MPa之后,应力场分布面积较为稳定。由此可以看出,应力场的分布情况变化过程与裂隙变形过程相似,故说明渗流运动过程中水力耦合作用引起的应力场的应力变化是裂隙介质发生变形的主导因素。

4 结论

(1)煤体渗流过程中,只有贯通裂隙通道参与渗流运动,非贯通通道不参与渗流运动。

(2)在水压逐渐增大的过程中,出口平均速度和内部通道最大流速都在逐渐增大。但是平均速度的增长趋势类似与对数函数,而最大流速的增长则分为:增长速率较大的正比例阶段、近水平阶段、增长速率较小的正比例阶段。

(3)在水压增大的过程中,裂隙最大变形是先增大到一个峰值又突然反弹降低,经过短暂的平缓期后又缓慢增大直至基本趋于稳定的过程。

(4)煤体渗流过程中裂隙最大位移变形的变化情况与最大流速的变化情况基本一致,而最大位移的变化又与应力场的应力变化基本相似。在水压变化的同时,渗流场内部的水力耦合作用是影响渗流规律的重要内在因素。