基于数字图像重构的含瓦斯煤失稳模拟研究

张雪强，刘利娜，李鑫鹏

(山东科技大学 矿业与安全工程学院，山东 青岛 266590)

近年来，国家经济发展对煤炭的需求量逐渐增大。煤与瓦斯突出是煤炭开采过程中一种异常动力现象，具有突发性，严重制约着煤矿的安全生产[1-4]。煤与瓦斯突出灾害主要受地应力、瓦斯压力和煤体力学特性等因素影响[5]，因此厘清上述因素对煤体力学行为的影响规律对防治煤与瓦斯突出灾害具有重要意义。

目前，国内外学者针对瓦斯压力影响煤体力学行为展开了大量的研究。WANG S等[6]分析了气体解吸、加载速率及应力等对含瓦斯煤变形的影响；张东明等[7]针对煤体力学特性与能量特征受瓦斯压力影响规律展开了研究；高保彬等[8]研究了不同瓦斯压力影响煤体力学特性和声发射特征规律；尹万蕾等[9]通过实验对比分析煤样弹性模量、峰值强度与瓦斯含量之间的关系，得到了瓦斯对煤力学特性的影响规律；何俊江等[10]通过含瓦斯煤吸附实验与全应力应变实验研究了瓦斯压力和煤体力学特性之间的关系；程远平等[11]研究了地应力演化特征对煤与瓦斯突出的控制作用规律；王维忠等[12]研究了煤体内瓦斯压力对煤与瓦斯突出强度的影响规律；王家臣等[13]通过含瓦斯煤加卸载实验，从宏观角度研究了瓦斯压力及围压对煤体力学特性的影响规律；刘延保[14]从细观角度研究了不同瓦斯压力和围压下煤样的动态破裂演化过程，并采用分形维数评价应力引起煤样损伤的破裂程度。

目前对含瓦斯煤样力学特性的研究，主要集中在不同应力路径下的煤样力学特性演变规律及影响因素方面，但对煤样结构变形及破坏过程的研究较少。且煤在自然界中是一种非均质材料，裂隙较为发育，目前数值模拟对煤的真实结构考虑较少，与现实煤样组成差别较大，结果会有一定偏差。笔者基于数字图像技术，对煤样图像进行数字处理，并采用FLAC3D有限差分软件进行模型重构，研究不同瓦斯压力下煤样的变形破坏过程及塑性区扩展规律。

1 煤样断面细观结构分析

1.1 煤样数字图像处理原理

通过数字图像技术可有效获取煤样表面细观结构，煤样的图像由大量的像素点构成。彩色图像中任何颜色都是由红色(R)、绿色(G)、蓝色(B)3种颜色按照不同的比例混合而成的，所以真彩色图像数据可以由3个离散的函数fk(i，j)表示，其中k=1，2，3。对于灰度图像，k=1。

在煤样图像处理时通常使用灰度图像，组成灰度图像的每个像素点都对应一个整数，用来表示该像素点的亮度，即灰度。灰度值范围为0～255，亮度由深至浅，对应灰度图像中的颜色由黑至白。灰度图像中每个像素点都为256个灰度值之间的一种，整个图像由不同灰度值的像素点阵(数字图像矩阵)组成，二元图像与对应像素点阵实例如图1所示。

图1 二元图像与对应像素点阵实例

上述数字图像矩阵的灰度构成了一个离散函数，该函数表达式为[15-16]：

(1)

式中：f为灰度值；i、j为像素点所对应的行与列；m、n为数字图像所包含像素点的行数与列数。

在该离散函数中，灰度值代表了组成数字图像像素点的信息，这些像素点坐标及灰度值大小是采用FLAC3D进行数值模型重构的基础。

煤样三维图形数据由多个二维图形数据叠加而成，笔者选取煤样中的一个切面做裂纹扫描，扫描图像如图2所示。

图2 煤样断面真实图像

煤样主要由煤、裂隙、裂隙填充物、颗粒胶结物等构成，为简化计算，现只考虑该煤样由煤和裂隙填充物构成。

采用Matlab软件读取图像，并采用rgb2gray函数对图像进行二值化处理，读取灰度值，其灰度直方图如图3所示。

图3 煤样扫描图像灰度直方图

1.2 煤样细观结构重构

采用FLAC3D有限差分软件将数字图像重构为可进行数值计算的模型，该模型由一系列网格构成，每个网格对应灰度图像中一个像素点。数字图像中灰度值分布较分散，经对比，提取0～20的灰度值能更真实表征煤样表面裂隙。煤样裂隙分布如图4 所示。

图4 煤样裂隙分布图

使用Matlab软件提取图4中像素点灰度值与坐标值，在FLAC3D中建立一个尺寸为50 mm×100 mm的模型，网格数为100×200，每个网格对应图4中 1个像素点坐标。

选择应变—软化模型，对模型整体定义煤样参数，然后Fish程序开始遍历每个像素点的坐标，若该坐标的灰度值处于0～20内，则对该坐标所对应的网格重新定义材料种类，并给予相应材料参数如体积模量、黏聚力、内摩擦角、抗拉强度等，相关模型力学参数如表1所示。

表1 模型力学参数

像素点处理完毕后，建立煤样的真实结构模型，真实结构网格如图5所示。

图5 煤样真实结构网格图

2 含瓦斯煤样失稳数值模拟研究

数值模拟采用的边界条件为：下端施加垂直约束，左右两端自由，上端加载应力。通过控制位移实现应力加载，加载速率为1×10-4mm/s，根据计算步数控制加载位移大小，直至煤样破坏。

2.1 不含瓦斯煤样失稳数值模拟结果

通过模拟煤样无侧限压缩过程，得到不含瓦斯煤样塑性区扩展过程如图6所示。

图6 不含瓦斯煤样塑性区扩展过程

由图6可知，在加载应力过程中，煤样裂隙部位先产生塑性区，随后塑性区由裂隙向外延伸；继续加载，塑性区随应力增大继续扩展，连接成主塑性区。主塑性区尖部两端存在一些次塑性区；继续加载，次塑性区会逐渐扩展为主塑性区，并在主塑性区两端再产生次塑性区，塑性区扩展随加载过程循环前进，直至煤样失稳破坏。

塑性区沿裂隙逐渐贯通如图6(a)～(c)所示；然后破坏由顶部和底部沿边缘向下发展，同时向煤样内部发展，如图6(d)、(e)所示；最后塑性区贯通，煤样发生破坏，如图6(f)所示。

2.2 含瓦斯煤样失稳数值模拟结果

为考察瓦斯压力对煤样失稳破坏过程的影响，数值模拟时在煤样裂隙中施加不同大小的气体压力，以代替孔隙压力及内部膨胀应力。数值模拟所用模型、边界条件及力学参数与不含瓦斯煤样失稳数值模型相同，气体压力分别为0.5、1.0、1.5、2.0 MPa。

由理想气体状态方程式(2)可知，在气体物质的量一定的情况下，气体压力和体积成反比：

pV=nRT

(2)

式中：p为气体压力，Pa；V为裂隙体积，m3；n为气体物质的量，mol；R为气体常数，J/(mol·K)；T为热力学温度，K。

采用FLAC3D中Fish语言实现数值模拟过程中气体压力的变化过程，如图7所示。为提高计算效率，并考虑计算精度，每50步进行1次循环，达到目的时间步后，计算结束。

图7 气体压力变化流程图

在气体压力为1.0 MPa和2.0 Ma条件下，得到含瓦斯煤样塑性区扩展过程如图8所示。

图8 不同气体压力下含瓦斯煤样塑性区扩展过程

由图8可知，含瓦斯煤样在剪切塑性区的基础上，增加了拉伸塑性区。当加载0.4 mm位移时，气体压力为2.0 MPa时煤样塑性区已完全贯通；气体压力为1.0 MPa时煤样仍在扩展阶段，塑性区扩展速度随气体压力的增大而增大。

不同气体压力下煤样应力—应变曲线如图9所示，含瓦斯煤样无侧限抗压峰值强度及峰值点应变与气体压力的关系如图10所示。

图9 不同气体压力下煤样应力—应变曲线

图10 不同气体压力下峰值强度和应变图

由图9、图10可以看出，随着气体压力的增大，峰值强度逐渐减小，峰值点应变也逐渐减小，呈线性关系，拟合公式分别为y=-1.460x+7.596和y=-0.113x+0.634。

煤样的失稳破坏是外部载荷及内部瓦斯压力二者共同作用的结果。当气体压力较小时，瓦斯压力可以促进煤样骨架变形，使其变形模量降低；当有效应力超过极限应力时，瓦斯压力增大了裂隙端部受力，促进裂隙扩展。瓦斯压力使煤的横向应力随着瓦斯压力的增大而增大，因此煤样的破坏形态在剪切破坏的基础上增加了拉伸破坏。

3 结论

1)基于数字图像处理方法，分析处理煤样断面图像，编写Fish语言模型重构程序，在FLAC3D中建立了煤样的真实结构数值模型。

2)采用Fish语言实现了瓦斯压力随塑性区扩展逐渐变化，模拟结果表明：无瓦斯煤样以剪切破坏为主，存在一个贯穿整个模型的主剪切破坏面。含瓦斯煤样存在2种破坏形式，在剪切破坏的基础上增加了拉伸破坏。

3)煤样中瓦斯压力可以促进煤样骨架变形，进而加快煤样裂隙扩展，煤样的失稳破坏是外部载荷及内部瓦斯压力二者共同作用的结果。含瓦斯煤样的压缩强度与峰值点应变均随着瓦斯压力的增大而减小，呈线性关系，分别符合公式y=-1.460x+7.596和y=-0.113x+0.634。