基于CT 成像的煤岩孔裂隙结构重建及渗流模拟研究

刘维福，丁序海，周 光

（神华国能集团有限公司 煤炭管理部，北京100033）

随着井工煤炭生产开采的需要[1]，地下矿井高效挖掘开采过程截割齿截割破碎过程总是伴随着大量煤尘的产生[2]，长期暴露于工作面上高浓度粉尘条件下的煤矿工人及易罹患严重肺尘病（CWP）[3]，高浓度粉尘也是诱发粉尘爆炸灾害的主要原因，严重制约井下煤炭安全生产[4]。煤尘注水作为井下高效防尘防突、驱排瓦斯主要方法[5]，已被越来越多的应用于生产过程中[6]。在煤层内存在孔隙和裂隙2个系统[7]：孔隙系统发育在煤岩基块中，天然裂隙是煤中流体渗透通道。煤的孔隙结构和形态直接关系着孔隙率、渗透率。连通性是孔隙的重要特征，它直接关系到储层的渗透作用。目前绝大多数的学者主要聚焦于煤的内部孔隙裂隙分布、CT 扫描复原上，关于煤的细观表征建模、水在煤体中渗流流动规律上的研究很少。蔺丽君[8]等通过COMSOL 构建一种煤体二维孔隙尺度流动与变形数值模型，分析了在应力作用下岩石孔隙结构演化趋势及对渗流影响。郎颖娴[9]等通过RFPA 建立了可以反映玄武岩内部孔隙结构的巴西盘数值模型，研究了孔隙对试样的破坏机制和抗拉强度的影响。李相臣[10]运用FDK 算法重建煤岩的灰度图像，对煤岩中的孔隙尺度以及空间分布进行综合表征。白若男[11]等通过CT 扫描，重构了煤体模型，分析了流场的流线、压力、速度及渗流方向的演化规律。陈俊国[12]等通过图像分析技术，分析了渗透率随围压的变化规律。为了准确描述煤层注水时煤样渗流过程，通过对煤体进行CT 扫描，经过MATLAB 对扫描图像进行中值滤波处理，对孔隙变化准确描述，并运用数值模拟软件对煤体进行了重构，分析研究水在煤体流动的机理及运动规律，为煤层注水降尘工艺提供指导意义。

1 实验煤样及实验方法

研究煤样取自阜新市雷家区4 号煤层，煤样种类为长焰煤。按照GB/T 482—2008 国家标准进行采样。在实验室用标准煤样切割机将煤样切割成5 cm×5 cm×10 cm 的长方体，获得3 个标准煤样，选取1 个裂隙发育完好煤样进行探究。

工业CT 采用先进的高频恒压X 射线源、数字图像探测器以及高精度机械检测平台。工业CT 原理是利用X 射线穿过被检测物体的截面，提取衰减后的信息。被扫描物体的密度由X 射线衰减系数来体现。物体的密度越大，对X 射线吸收能力越强。通过衰减系数转换成CT 数，CT 数越大相当于物体密度越高。通过对比CT 图像中的黑白分布展现出被测物体的密度分布情况。通过CT 扫描后的物体，可以真实无损地将煤体内部孔隙结构表征出。

采用丹东奥龙射线仪器有限公司AL-CT-225工业X 射线CT 检测系统，采用锥束扫描和DR 实时成像多种检测方式，可扫描具有缺陷、孔隙分析和被检测工件区扫功能。层析识别系统的CT 空间最小分辨率为6 μm，密度分辨率为0.5%，最大成像素为1 024×1 024，实验共扫描出800 个切片图像。

2 CT 图像数据处理

2.1 图像中值滤波处理

采用MATLAB 数字图像技术对CT 扫描图像中的孔隙结构信息进行识别处理，将煤岩的孔隙裂隙和基质进行分割。利用X 射线CT 扫描得到的图像存在系统噪声，噪声定义为不可预测，只能用概率统计方法来认识的随机误差。即噪声在图像上是随机分布的，这些噪声点位在图像上会增加或减小图像像素点位的真实像素值，常表现为孤立像素点或像素块。所以降低了图像的质量，通过相关的滤波函数处理噪声，选取中值滤波进行处理。

中值滤波对脉冲噪声有良好的滤除作用，在滤除噪声的同时还可以保护信号的边缘，使之不被模糊。因此，为了将每个点的信息准确提取，将图片中噪声点使用中值滤波函数进行处理，并利用MATLAB 将图像上存在的每个像素点的值导出为矩阵的模式。中值滤波的基本原理为把数字图像或数字序列中的1 点值，用其邻域中各点值进行代替。对1 个数字信号序列进行滤波处理时，首先定义1 个长度L 为奇数的窗口，L＝2n＋1。窗口信号样本为（x1，x2，x3，…，xn），将n 个数据大小进行排列：x1＜x2＜x3＜…＜xn。则有：

其中y 为窗口序列（x1，x2，x3，…，xn）的中值，将窗口在数据上滑动，将窗口正中所对的像素值用窗口各个像素的中值代替，即：

y（i）＝Med［x（i－N），…，x（i），…，x（i＋N）］

MATLAB 中值滤波算法前后原始图像和中值滤波滤除噪声比较如图1。

图1 原始图像和中值滤波滤除噪声比较Fig.1 Comparison of noise removal between original images and median filtering

利用MATLAB 中值滤波对原图像进行滤除噪声处理，通过中值滤波法保护了图像的原有信息，滤除了原始图像的噪声点位，使得处理后的图像信息更加准确，并作为后续二值化处理的图像基础。

2.2 图像二值化处理

二值图像是1 个取值只有0 和1 的逻辑类型的二维数组。所以二值图像可以看成是一种特殊的亮度图像。1 种只有2 个亮度级（最暗和最亮）的亮度图像。CT 扫描后的图片为灰度图像，其中的微小裂隙较难提取。选取OStu 法选确定最佳阈值，通过处理的图像可以提取图像的孔隙结构信息特征。二值化处理的图像如图2。

图2 二值化处理的图像Fig.2 Binary images processing

原始灰度图像经MATLAB 的中值滤波处理后，去除了图像的噪声点位。清晰地分割了煤体中的孔隙和基质。有效的保留了煤体的孔隙特征和边界信息。而通过二值化处理后的原始图像清晰地分割了煤体孔隙和基质。为下一步的二维模型建模提供了基础，并以此作为有限元方法构建数值模型。

2.3 煤体微观精细模型建立

选取了煤体试样2 个孔隙发育完好煤体图片，煤体微观孔隙结构模型如图3。由图3 可以看出，主裂隙贯穿整个煤样，其中次生裂隙在煤体中分布。部分次生裂隙和主干裂隙相互连接，孔隙在试样内的分布均匀。

图3 煤体微观孔隙结构模型Fig.3 Microscopic pore structure models of coal

将中值滤波、二值化处理后的图像信息进行提取，运用COMSOL 中的插值函数将煤样切片建模，以此为数值模拟的模型基础。

2.4 CT 图像的孔隙率计算

孔隙率指散粒状材料表现体积中材料内部的孔隙占总体积的比例。CT 扫描图像中包括了煤体的骨架孔隙信息特征，其中包含了主裂隙和次裂隙可以将煤体内部的孔隙结构展示出来。

Taud H[13]在CT 扫描图像内，把CT 灰度级类比成为三维地形图的高低，灰度水平等效于三维地形图的地形高低。灰度值越大的点位地形图所在位置越低，结合数理统计方法获得孔隙率的分布函数，通过灰度值的大小和统计灰度值出现的次数可以计算出试样的孔隙率，孔隙率φ 公式为：

式中：ri为各像素点的灰度值，取值范围［rmin，rmax］；H（ri）为［rmin，rmax］范围内的灰度直方图；rmin为图像灰度最小值；rmax为图像灰度最大值。

由此可以得到图像中每1 个像素点的孔隙大小，煤样每层孔隙率如图4。

图4 煤样每层孔隙率Fig.4 Porosity of each layer of sample coal

图4 中试样煤样的孔隙率平均值为14.82%，试样煤样每个切片层数的孔隙率波动变化在0.02 范围内。说明试样内孔隙数据具有一定的可靠性。

3 煤体内部水渗流通道的数值模拟

3.1 孔隙度的分布

通过对比经过中值滤波、二值化、COMSOL 建模的模型，可以看出煤体内部的主裂隙分布明显，在中值滤波和二值化处理后的图像，煤的裂隙分布在黑色线条内。通过COMSOL 建模后得到的孔隙度分布图如图5。可以看出，煤体的主裂隙呈现红色区域，而蓝色区域代表的是煤基质。其他区域代表的是次生裂隙和煤基质夹杂部分。图5（a）中平均孔隙度为14.37%，图5（b）中平均孔隙度为15.26%。

图5 孔隙度分布图Fig.5 Porosity distribution map

3.2 水在煤体中流动的数值模拟

煤体内部的孔隙结构错综复杂，这也就导致了水在煤体中运动过程的复杂性。流体介质在多孔介质流动中流动呈向的是多孔流动流出，当水在煤体中运动过程中会通过主孔隙流动，分裂隙连同基质处水部分流动。考虑到水流动在煤体中的复杂性，现简化水在煤体的流动过程。煤体的上边界设置20 mm 的水头，在下边界设置为0 的水头。这样产生了20 mm 的水头差，驱动流体流动,左右边界不可渗透，流体密度为1 000 kg/m3,动力黏度0.001 Pa·s。渗流水模拟结果图如图6。

图6 渗流水模拟结果图Fig.6 Diagram of seepage water simulation results

从图6 可以看出，煤体含有主裂隙和次生裂隙，水在经过煤体的过程中，在主裂隙流动过程明显，水也会在次生裂隙中流动。水在流动过程中，速度会呈现逐渐减小后的趋势。水在主裂隙流动过程中，速度变化的趋势不大。水流在主裂隙水流流速如图7。

在图6（a）中水流速度从0.017 5 m/s 初始速度衰减到出口速度0.017 2 m/s。这主要是因为在此煤样主裂隙中，裂隙发育较为平缓，水在流动过程中受到的阻力较小。在图6（b）中由于主裂隙发育的较为完好且弯折不大，水在流动过程中受到阻力变大，在这一过程中水的流速会出现剧减的状态，之后速度趋于稳定减小。水在经过次生裂隙的过程中水流动速度会出现局部增大的趋势，这主要是因为水在经过次生裂隙的过程中，因为次生裂隙的孔口半径狭小且弯曲曲折方式急剧，水会出现流动速度激增的状态。虽然煤体中的次生裂隙很多，但水的流向不会通过所有的次生裂隙。这主要是因为水的流动具有选择性，水会选择从孔隙结构发育完好的地方渗透。在次生裂隙和基质夹杂部分流动不是很明显。

图7 主裂隙水流流速图Fig.7 Velocity diagram of main fractures

4 结 论

1）通过CT 扫描煤体图像，通过MATLAB 中值滤波对图像进行处理进而得到图像的每一点灰度值，并以二值化处理得到了煤体的孔隙结构。运用COMSOL 插值函数将煤体进行二维模型重构，将煤体上的孔裂隙分布特征进行了精细表征。

2）运用Taud H 的数理统计方法获得了孔隙率的分布函数，对比COMSOL 模拟做出了孔隙度的分布情况，两相对比反演了煤的孔隙率。运用数理统计方法得到煤体的平均孔隙率为14.82%，在COMSOL中模拟数值的结果为14.37%、15.26%。证明了这种建模的可靠性和精准性。

3）模拟结果表明，水在煤体流动过程中，裂隙的发育程度会影响水在煤体中的运动。整体的趋势为水流速呈现减小的趋势。沿着主裂隙流动的流速衰减趋势很小，而在次裂隙流动的水跟裂隙的发育完好程度以及孔裂隙的弯曲程度有关。在煤体基质的流动速度基本为0。