车禹恒

（1.中煤科工集团重庆研究院有限公司，重庆 400037；2.瓦斯灾害监控与应急技术国家重点实验室，重庆 400037）

煤作为多孔介质材料，其孔隙结构特征决定了其物理性质，进而使煤层气在不同煤层构造中呈现不同的赋存状态和流动特性[1]。为了探索不同地质条件下煤层气在各煤阶储层中的赋存和流动性能，诸多学者从宏观角度出发对煤层孔隙分布、渗透性、吸附性及其相关影响因素进行研究[2-3]，鲜有学者利用煤微观孔隙结构揭示其宏观性质。随着高分辨率无损在线检测技术的发展，不仅能够实现煤微观孔隙结构可视化[4]，还能对其结构参数进行定量分析，对不同地质条件下煤层气的赋存、产气和控气机理作深入探索，成为目前的研究热点[5]。X-ray CT技术基于被检测样品的断层扫描图像，利用内置成像算法清晰重构出被测样品的内部结构，在实现无损检测的同时具备超高的分辨率以及三维数字化等优点，被广泛应用于煤微观结构探测领域。孙英峰[6-7]通过使用该技术将各种煤级的孔隙和喉道的空间形态进行了重构，随后将其处理和定量研究后的结果和液氮吸/脱附的实验结果进行了进一步的分析探讨。谢淑云等[8]利用X-rayμCT扫描技术对碳酸盐岩的微观孔隙结构的空间模型进行了重新构造，通过这个模型引入了一些多重分形指标来对孔隙系统的复杂程度和分布规律作进一步的阐述。陈昱林[9]在对泥页岩的微观孔隙的一些空间结构特征的研究与探索中也运用了这项技术来做定量分析。李伟等[10]通过利用CT扫描重构的技术把构造不同的各种煤样之间的微观孔隙进行了详细的辨析。姚艳斌等[11-15]利用了NMR和X-rayμCT技术相结合的方式将煤微观的孔隙空间结构做了更详细的阐述。这些研究对煤孔隙微观空间结构进行了丰富的表征及分析，然而大量微观孔隙中只有连通孔隙是煤体渗透性能的决定性因素，这些随机分布的孔隙相互连通形成结构迥异的孔隙团，其特殊的拓扑结构非均质性是煤体渗透各向异性的根本原因，目前关于这一特性的研究较少。因此，针对鄂尔多斯盆地低阶煤进行X-rayμCT扫描，重建煤孔隙微观空间结构模型。基于此，采用AVIZO软件提取最大连通孔隙团，并构建其等效孔隙网络拓扑结构，分析孔隙形状因子与孔隙等小直径和表面积的关系，探讨连通特性和分形维数在空间中的变化规律，从而揭示煤微观连通孔隙团拓扑结构的非均质特征。

1 X-rayμCT扫描及图像处理

实验煤样取自鄂尔多斯盆地西缘羊场湾矿褐煤（YCW）和鄂尔多斯盆地东缘斜沟矿气煤（XG）。原煤样品被打磨成5 mm×5 mm×10 mm，打磨过程中保留了煤体完整结构，表面光滑无划痕和裂隙，从而尽可能减少CT成像过程中尖锐的边角造成的边缘增强效应。煤样显微组分及工业分析见表1。

表1 煤样显微组分及工业分析Table 1 Microscopic composition and industrial analysis

实验采用Nano Voxel-3000系列高分辨率X射线计算机断层扫描仪。实验过程中测试电压120 kV，测试电流50μA，曝光时间1 000 ms，扫描模式为局部扫描，物理分辨率0.5μm，扫描时间52 min。

采用AVIZO软件对CT扫描图像进行处理，并重构煤微观孔隙空间结构模型。CT图像噪声的存在会是重构模型严重失真，因此使用中值滤波算法对图像进行降噪处理，CT图像滤波前后灰度值分布图情况如图1。

图1 CT图像滤波前后灰度值分布情况Fig.1 Distribution of gray value before and after CT image filtering

为了精准确定煤孔隙、矿物和基质的灰度阈值，依据孔隙度反演构建数学模型，确定不同介质的阈值区间。模型函数Ф（Tm）符合BiDoseResp增长型函数：

式中：Ф1、Ф2为函数极小值和极大值，分别取0和1；ζ1、ζ2为增长模型2个阶段的增长指数，表示2个阶段的线性增长速率；ω¯为增长模型发生转变前后的权重；γ1、γ2为增长模型2个阶段的状态参数。

基于建立的阈值模型，确定了最佳孔隙度和矿物含量阈值Tmp和Tmm，煤微观介质阈值结果见表2。

表2 煤微观介质阈值结果Table 2 Coal micro-medium threshold results

2 最大连通孔隙团拓扑形态构建

2.1 孔隙REV结构重建

代表性体积单元REV研究方法以研究内容在微观和宏观尺度的相似性为桥梁，实现微观领域物理性质表征宏观特性的跨越[16]。以该物理性质为基础选取能够代表煤宏观尺度物理性质的REV单元。为了确定每种煤样REV单元的尺寸，在Avizo重建的煤微观孔隙空间结构中选取4个目标点，以该点为中心选取不同尺寸的REV单元，并考察其孔隙率大小随单元尺寸的变化规律，从而选取与宏观孔隙率相符合的REV表征单元。经过计算，当煤微观孔隙REV边长大于200体素时，其孔隙率变化幅度小，与整体煤孔隙率接近。因此，将所有煤样微观孔隙的REV尺寸设置为200×200×200，物理尺寸为20 μm×20μm×20μm。REV灰度结构如图2，黑色表示孔隙，白色表示矿物质，灰色表示煤基质。可以看出，在该尺度下存在很强的非均质特征。

2.2 最大连通孔隙团提取

采用Axis Connectivity模块对REV单元中的最大连通孔隙团进行识别，然后采用Separate Objects模块中的分水岭算法对最大连通孔隙团进行分割，识别连通孔隙团中的单个孔隙并进行标记，REV孔隙结构如图3。

图3 REV孔隙结构Fig.3 REV pore structure

2.3 拓扑形态构建

等效孔隙网络模型PNM在提取的孔隙中轴线基础上，从拓扑学的角度表达了真实孔隙空间的对应关系，中轴线的节点即是孔隙中心，与此同时，采用最大球法对孔隙空间进行分割，从而简化了与真实煤微观孔隙系统具有等价关系的拓扑网络结构，PNM空间拓扑结构如图4。

图4 PNM空间拓扑结构Fig.4 PNM space topology

3 孔隙空间拓扑结构分析

3.1 形状因子

微观孔隙形状极其不规则，因此引入形状因子来考察孔隙形状的这种不规则程度。形状因子是三维结构球度的1个度量指标，定义为：

式中：η为孔隙形状因子；ν为孔隙体积，μm3；s为孔隙表面积，μm2。

为了进一步分析孔隙形状因子与等效直径之间的关系，对所有煤样各个REV区域孔隙形状因子进行统计，孔隙形状因子与等效直径的关系如图5。

图5 孔隙形状因子与等效直径的关系Fig.5 Relationship between pore shape factor and equivalent diameter

通过拟合分析发现，煤微观孔隙形状因子与等效直径具有如下指数关系：

式中：η0为常数；d为孔隙等效直径，μm；ζ为指数因子；γ为等效直径松弛度，%。

指数因子ζ和等效直径松弛度γ可以反映出孔隙形状因子对等效直径的敏感度。煤体微观尺度中孔隙形状因子不仅与其等效直径存在函数关系，与其表面积也存在较强的线性关系。对各煤样REV单元中孔隙的形状因子与表面积进行拟合并绘制曲线，孔隙形状因子与表面积的关系如图6。

图6 孔隙形状因子与表面积的关系Fig.6 Relationship between pore shape factor and surface area

形状因子与表面积之间存在如下线性关系：

式中：ηd为直线截距；k为直线斜率，同时表征了形状因子对表面积的敏感度。

3.2 连通特性

欧拉数分布如图7，图7对比了2个煤样REV最大孔隙团在xy切片上的欧拉数。

图7 欧拉数分布Fig.7 Distribution of Euler number

羊场湾褐煤平均欧拉数为37 271，A1、A2和A3单元欧拉数稳定在1 000～1 200范围内，但A4单元在切片距离10～16μm出现波谷，位于13μm的最低值仅350，原因是在A4中部边缘存在体积较大的矿物结核，此处孔隙团出现“缺口”。斜沟气煤平均欧拉数为7 344，B4单元在空间各处的欧拉数稳定在300附近。B1单元中部欧拉数小幅度升高，最高为500，B2单元欧拉数在5.8、14.5μm存在2个峰值，峰值欧拉数分别为550、700，表明该单元连通性在中部均一性较差，B3单元2～14μm范围内欧拉数稳定在650左右，当距离大于14μm时迅速降低至300附近，可见该单元孔隙团自下而上呈三角锥分布。

3.3 分形维数

分形维数分布如图8，图8具体给出了不同变质程度煤样最大连通孔隙团xy切片二维分形维数分布情况。

图8 分形维数分布Fig.8 Distribution of fractal dimension

二维分形维数通常介于1～2，能够表征孔隙团的破碎性，也是衡量孔隙团内部孔隙离散程度的有效指标。从图8可以看出，羊场湾褐煤A1～A3单元曲线相对稳定，说明这3个区域的最大连通孔隙团破碎程度相似，而A4单元在10～14μm出现波谷，最低值为1.3，该范围内孔隙团“缺口”的存在降低了整体的粗糙度。斜沟气煤B4单元内各xy切片的二维分形维数最为稳定，B1和B3单元的曲线变化有着一致性，距离越大孔隙团破碎程度越低，B2单元5、14.5μm处存在峰值，孔隙团伸出的侧链是分形维数急剧增加的主要原因。

4 结 语

针对鄂尔多斯盆地低阶煤进行X-rayμCT扫描，重建煤孔隙微观空间结构模型并提取最大连通孔隙团，构建其等效孔隙网络拓扑结构，揭示煤微观连通孔隙团拓扑结构的非均质特征：孔隙形状因子与等效直径存在指数关系，与孔隙表面积存在线性关系。表明连通孔隙团内部渗流孔隙空间越大，其形状越不规则，非均质越强；低阶煤连通孔隙团各断面的欧拉数变化较大，连通性具有更强的非均质性，孔隙离散程度更高，但粗糙度较低。