基于扫描电镜的多孔介质迂曲度识别方法研究

尹婷婷 王超 张中标 孙剑楠 闫柯

摘要：迂曲度是研究多孔介质渗流特性的重要参数，反映了流体流经多孔介质的迂回曲折程度。目前主要通过试验以及仪器测量间接获取多孔介质的迂曲度参数，存在实验人员主观性影响较大、测量计算方法不明确的问题。为此，基于多孔介质材料的扫描电镜图像，利用水平集原理分割了多孔介质材料的孔隙通道，并计算了流体在孔隙通道里的流线，通过数字化识别方法提取了多孔介质的迂曲度参数。最后通过压汞仪测量对应多孔介质材料的迂曲度参数与识别结果作对比，验证了本文计算迂曲度方法的正确性，为多孔介质迂曲度的计算提供一种新的求解思路。

关键词：多孔介质；水平集原理；迂曲度；扫描电镜

中图分类号：TU441+.33                   文献标志码：A            doi：10.3969/j.issn.1006-0316.2024.05.002

文章编号：1006-0316 （2024） 05-0011-06

Tortuosity Identification Method of Porous Media Based on

Scanning Electron Microscope

YIN Tingting1，WANG Chao1，ZHANG Zhongbiao1，SUN Jiannan2，YAN Ke2

（ 1. Nuclear Power Institute of China， Chengdu 610213， China;

2. School of Mechanical Engineering， Xian Jiaotong University， Xian 710049， China ）

Abstract：Tortuosity is the key parameter to study the seepage characteristics of porous media， which reflects the tortuosity degree of the fluid flowing through the porous media. At present， the tortuosity parameters of porous media are mainly obtained indirectly by experiments and instrumental measurements， which have the problems of the subjective influence of the investigators on the experimenters and the unclear measurement and calculation methods. For this reason， on the basis of the scanning electron microscope images of porous media materials， this paper uses the level set principle to segment the pore channels of porous media materials and calculates the streamlines of the fluid in the pore channels， and extracts the tortuosity parameters of porous media through digital recognition methods. Finally， the tortuosity parameters of the corresponding porous media materials are measured by the mercury porosimeter are compared with the identification results to verify the correctness of the tortuosity calculation method in this paper， which provides a new solution idea for the calculation of the tortuosity of porous media.

Key words：porous media；level set principle；tortuosity；scanning electron microscope

作为影响多孔介质渗流特性的重要参数，迂曲度是弯曲孔道的实际路程与孔道首尾间可视距离即直线距离的比值，表示流体在多孔介质中流动路径的曲折程度[1]。目前，国内外针对迂曲度的研究方法主要包括试验法[2-3]、解析法[4-5]和模拟法[6-7]。在试验法方面，部分学者通过试验间接测量迂曲度，如吕道平[8]根据水力迂曲度与多孔介质中的孔隙半径、孔隙率、渗透率的关系，利用气测岩心渗透率求取了岩芯水力迂曲度。在解析法方面，将迂曲度的计算转化为简单的数学几何关系，如Yu等[9-10]基于简单的几何模型推导，间接得到了只与孔隙率有关的迂曲度表达式。颜瀚等[11]基于多孔介质颗粒间的位置关系和孔道的几何关系，提出迂曲度的上限和下限排列方式，并给出多孔模型的迂曲度取值范围。考虑到多孔介质中不同颗粒间随机堆叠的多样性，以某一确定多孔孔道的迂曲度来表征整体的迂曲度的精确度较低。在此基础之上，考虑几何模型的分形特征，提出分形毛细管力模型[12-13]。在模拟法方面，杨谦洪等[14]基于光滑流体粒子动力学对多孔介质的流动进行了模拟，研究了流体参数对多孔介质渗透率和迂曲度计算的影响。

综上所述，通过试验法和解析法间接计算迂曲度值时，存在实验人员主观性影响较大、测量计算方法不明确以及计算精确度较低等问题，难以直观精确地计算出迂曲度。而模拟法通过真实模拟孔道内的流体轨迹来计算多孔模型的迂曲度值，与迂曲度的定义完全符合，较直观地表达出迂曲度的真实影响因素，考虑了孔道形状、结构等真实因素，具备直观与精确的表征意义。为此，本文通过识别多孔聚酰亚胺（Polyimide，PI）材料扫描电镜图像内的孔隙通道，计算流体在孔隙通道里的流线，直观地提取了多孔PI材料的迂曲度参数。最后通过压汞仪测量对应PI材料的迂曲度参数与提取结果作对比，验证了本文计算迂曲度方法的正确性，为多孔介质迂曲度的计算提供一种新的求解思路。

1 孔隙识别

1.1 扫描电镜图像

为了能定性分析出多孔PI材料的微尺度孔隙参数，本文采用JEOL7800F场发射扫描电子显微镜（Scanning Electron Microscope，SEM），开展多孔材料截面孔隙分析。由于PI材料的非导电性，电子束在样品表面积累电子，样品温度急剧升高，破坏样品形貌，且因电流不能导通，致使显微镜不能成像。所以在测试前对样品进行喷金处理，在样品表面镀层形成导电通道，形成稳定的成像结果。

基于以上分析，选取体积为8 mm3的多孔PI材料样品开展扫描电镜分析。设置扫描电镜的工作电压为3 kV，放大倍数为1000倍，得到试样的SEM图像如图1所示。可以看出，多孔PI材料的孔隙由原料粉末颗粒经烧结融合后剩下的微量孔隙形成，排布复杂随机，孔隙之间通过微尺度的喉连通。

1.2 图像分割

基于扫描电镜的分析结果，采用水平集原理对获取的多孔PI材料的SEM图像进行孔隙分割，获取精确的多孔流道模型为后期的迂曲度计算做准备。

1.2.1 水平集原理

PI材料的SEM图像由像素点聚集而成。其能量函数，也称局部二元拟合能量函数式为：

（1）

式中：为局部拟合函数；x为像素点，所有x组成图像域I；C为图像轮廓；和为x附近的两个图像强度[15]；K为随着|x－y|减小并且接近零的定位属性核函数；和为正常数；in（C）和out（C）分别为轮廓C内部和外部的图像区域；y为接近目标轮廓的像素点。

对于每个x，当轮廓C与物体边界重合时，选择最佳的拟合值和，则被最小化。其原理如图2所示，黑色线条即为最小化后找寻的整个对象边界[16]，通过对在图像域Ω中对所有x求积分来实现，为：

（2）

式中：为局部拟合函数。

在此基础之上，将其转为等效水平集为：

（3）

式中：为水平集函数；为带宽＞0的高斯核函数；H为单位阶跃函数，由平滑函数u定义。

为了使稳定演化，减小与有符号距离函数之间的偏差，添加：

（4）

式中：P为距离正则项；为梯度算子。

同时，为了规范的零水平轮廓，定义：

（5）

式中：L为的零水平轮廓长度；为局部长度函数。

综上，定义：

（6）

式中：F为整个能量函数；和ν为非负常数，分别为距离正则化系数与零水平线规范化系数，可根据图像特征调整大小以达到较好的分割结果。

式（6）即为水平集函数的具体表达形式，通过最小化水平集能量函数的大小，找到u为零的点，即可得到分割结果。

1.2.2 孔隙分割结果分析

基于水平集原理，编写MATLAB程序，设置距离正则化系数和零水平线规范化系数，开展多孔保持架材料SEM图像中的孔隙分割。

运行MATLAB程序，设置迭代次数为500次，得到试样分割孔隙的图像识别结果如图3所示。图中红色线条为零水平曲线，红色曲线包围的部分为分割出来的孔隙。可以看出，SEM

图像的孔隙分割结果将孔隙大部分完整地识别出来，仅有极少部分灰度分布不明显的区域未识别。

为了更清晰地看出多孔PI材料内部的孔隙结构分布，在水平集函数图像识别结果的基础之上，设置图3中红色线条内孔隙部分的灰度值为255，红色线条外固体材料部分的灰度值为0，得到只包括0、255两种像素点的二值图如图4所示，其中白色部分为孔隙，黑色部分为PI材料。

由图4可以看出，多孔介质的大部分孔未连通，这是因为二维图只展现了多孔介质一个截面上的孔分布，缺少了在三维层面的连通孔分布。为了减少误差，将位图转化为由线连接点的矢量图如图5（a）所示，导入CAD软件中进行修改，将矢量图中较密集的孤立孔连接起来，模拟三维图中孔的连通状态如图5（b）所示，图中红色线条即为修改的连通孔。

2 迂曲度计算

2.1 流线数据分析

基于经过修改处理的多孔模型，利用COMSOL Multiphysics开展多孔模型的流动仿真，设置入口、出口以及0.05 m/s的初始速度为边界条件，得到多孔模型的流线数据如图6中的红色箭头所示。可以看出，多孔模型内具有多条流线，流线的走向与多孔模型的孔道形状有重要关系。

基于多孔模型的流动仿真结果，分段提取流线数据如图7所示，其中不同颜色代表不同孔道的连线数据。从图7红色框中流线可以看出，提取的流线数据存在首尾点特别近的流线，基于迂曲度的定义（孔道实际长度与首尾点直线距离的比值）进行计算易出现迂曲度无限大的情况。为解决此问题，编写程序计算迂曲度大小，删除异常的迂曲度值，筛选在合理范围的迂曲度。

2.2 迂曲度提取原理

基于流线数据计算迂曲度大小的流程如图8所示。为相邻两点实际距离。

首先对每一条孔道进行标记，判断是否为孔道的起始点，如果是起始点则计算上一孔道的可视距离，如果不是则累积相邻数据点的距离，计算孔道的实际距离。最后根据迂曲度参数的定义，利用实际孔道长度sn除以相应孔道的可视距离s′n计算迂曲度大小。

3 实验验证

3.1 压汞仪原理

压汞仪是集操作、测量、数据处理功能为一体的微型计算机，其孔径测定范围为0.003～150 μm，是测量多孔PI材料孔隙参数的常规方法。其原理是在真空条件下，施加外压力将非浸润液体压入多孔材料试样，记录施压过程进行孔隙参数的计算。

孔径计算：外压力与进汞孔的净宽成反比。采用圆柱孔模型，Washburn方程给出了压力与孔径间的关系为[17]：

（7）

式中：为测量的孔径大小；p为外压力；r为汞的表面张力，与样品的材料属性和温度有关；为接触角。

根据式（7）将压力数据转为测试样品孔径大小，并依据进汞体积与孔径之间的关系得到孔体积分布曲线。假设孔隙为圆柱孔，样品不含墨水形孔和在外压力下不变形，则可以从压力－体积曲线计算出进汞的比表面积。最后，根据进汞体积、进汞比表面积以及孔径分布等参数近似计算得到样品的迂曲度参数。

3.2 迂曲度参数验证

为了验证利用水平集图像分割算法计算多孔介质材料迂曲度参数的正确性，利用压汞仪对多孔PI材料样品进行测量，得到迂曲度参数大小为19.05，与2.2节中图像识别计算的迂曲度19.61相比，大小仅偏差2.9%。可以看出，图像识别计算的迂曲度具有较高准确性，进一步验证多孔结构分析方法的正确性。

4 结束语

本文基于水平集图像分割原理，识别了多孔聚酰亚胺材料扫描电镜图像的孔隙通道，并基于流体仿真软件计算了流体在孔隙通道里的流线，直观且高精度地提取了多孔介质材料的迂曲度参数。最后利用压汞仪测量对应多孔聚酰亚胺材料的迂曲度参数与提取结果作对比，两者之间的误差仅为2.9%，验证了本文计算迂曲度方法的精确性，为多孔介质迂曲度的计算提供一种新的求解思路。

参考文献：

[1]薛东杰，周宏伟，邓淋升，等. 低渗煤岩气液两相流分形运动方程[J]. 工程科学与技术，2018，50（4）：30-40.

[2]盛建龙，韩云飞，叶祖洋，等. 粗糙裂隙水、气两相流相对渗透系数模型与数值分析[J]. 岩土力学，2020，41（3）：1048-1055.

[3]李强，姚仰平，韩黎明，等. 土体的“锅盖效应”[J]. 工业建筑，2014，44（2）：69-71.

[4]SHENG D，ZHANG S，NIU F，et al. A potential new frost heave mechanism in high-speed railway embankments[J]. Géotechnique，2014，64（2）：144-154.

[5]ZHANG S，TENG J，HE Z，et al. Canopy effect caused by vapour transfer in covered freezing soil[J]. Géotechnique，2016，66（11）：927-940.

[6]ZHANG S，SHENG D，ZHAO G，et al. Analysis of frost heave mechanisms in high-speed railway embankment[J]. Canadian Geotechnical Journal，2015，53（3）：520-529.

[7]TENG J，ZHANG X，ZHANG S，et al. An analytical model for evaporation from unsaturated soil[J]. Computers and Geotechnics，2019（108）：107-116.

[8]吕道平. 多孔介质中水力学迂曲度因子的求取及应用[J]. 新疆石油地质，2000（6）：515-517，540.

[9]YU B M，LI J H. A geometry model for tortuosity of flow path in porous media[J]. Chinese Physics Letters，2004，21（8）：1569-1571.

[10]李滔，李闽，张烈辉，等. 微多孔介质迂曲度与孔隙结构关系[J]. 天然气地球科学，2018，29（8）：1181-1189.

[11]颜瀚，滕继东，张升，等. 考虑颗粒排列的土体几何迂曲度模型[J]. 岩土工程学报，2020，42（1）：195-200.

[12]员美娟. 分形毛细管中Reiner-Philippoff非牛顿流体的有效渗透率研究[J]. 武汉科技大学学报，2013，36（2）：158-160.

[13]袁哲，刘鹏程. 一个适用于煤岩的新的分形毛管力模型[J]. 科学技术与工程，2015，15（9）：63-67.

[14]杨谦洪，姚军，黄朝琴. 一种基于PF-SPH的多孔介质迂曲度和渗透率研究[C]. 第十一届全国流体力学学术会议论文摘要集，2020.

[15]LI C，XU C，GUI C. Level set evolution without re-initialization：A New Variational Formulation[C]. San Diego，CA，USA：Conference on Computer Vision and Pattern Recognition，IEEE，2005：430-436.

[16]阎贝．医学图像分割算法研究与实现[D]. 成都：电子科技大学，2012.

[17]WASHBURN HL，DRAKE LC. Macropore size distributions in some typical porous substances[J]. Ind. Eng. Chen. Anal. Ed，1945（17）：787-791.