王 超， 吴 丰， 陈义国， 李广涛， 吴 凤

(1.陕西延长石油(集团)有限责任公司研究院， 西安 710075； 2.西南石油大学地球科学与技术学院， 成都 610500)

黄土是一种广泛分布的区域性土状堆积物，在中国中部干旱和半干旱地区均有不同程度的发育。由于其具有整体排列松散，黏土内部孔隙发育、弱胶结、强度高，压缩性小，含水饱和较低等特征[1-4]，在遇水浸湿后黄土内部结构发生变化，强度降低，易发生下沉现象，从而引发黄土崩塌等宏观的地质灾害现象。而这些现象的发生与其黄土内部的微观结构、孔隙分布及渗流特征等密切相关[5-8]，因此精细评价黄土微观结构与渗流情况对工程地质与水文地质研究均有非常重要的意义。

目前，针对黄土中颗粒类型，微观孔隙分布以及渗透率各向异性，前人均开展了一定研究。其中，雷祥义[9-10]综合分析黄土中存在的多种孔隙类型，探讨了黄土的孔隙结构与工程性质之间的关系；李兰等[11]采用电镜扫描手段和计算机图像处理，定量分析了黄土中微结构孔隙性；李林翠等[12]采用“注胶法”制成的黄土试样分析了不同埋深的各类形状孔隙分布。王梅等[13]提出了黄土微观结构的制备方法，分析出图像中的孔隙、颗粒展布状态；胡海军等[14]采用Autocad软件获得颗粒和孔隙的边界，综合分析了不同区域的黄土颗粒以及孔隙大小、形状和排列的定量差异；而苏英等[15]分析了陕西关中地区晚更新世黄土在粒度、厚度、湿陷性等方面的宏观差异，但尚未对其微观特征开展研究。在渗透性研究方面，张宗祜[16]于1962年即明确了渗透性分析在黄土类研究中的重要性。之后李喜安等[17]基于黄土渗气及饱和渗透试验，对孔隙分布与渗透性关系进行了相关研究。此外，随着新技术不断发展完善，文献[18-21]分别基于核磁共振、扫描电镜、CT(computed tomography)扫描等技术分析了黄土内部孔径分布并对渗透性进行了一定探讨。但相比于新技术在砂岩孔隙结构中的广泛应用[22-23]。在黄土孔隙结构、孔隙分布以及渗透性能等各方面的研究仍然有待完善。综合而言，黄土岩中不同矿物颗粒形状、大小开展具体研究较少，而矿物形态对孔隙赋存成因、黄土岩中孔隙网络变化、渗透率各向异性的相关研究仍亟需不断推进。

借鉴前人研究思路与成果，基于A地区的黄土岩微米CT图像，模拟其微观孔隙网络，在黄土成分、矿物类型的研究基础上，综合分析不同矿物与孔隙赋存状态之间的关系，讨论了黄土颗粒中微观特征对渗透率各向异性的影响，旨在为黄土岩以及各类岩性的孔隙空间分布及渗流能力研究提供一个可供参考和借鉴的思路与方法。

1 实验和模拟原理

1.1 实验准备

黄土样品选自A地区，处理方法：从原始黄土样品上切出一个大约5 mm×20 mm的柱塞样，柱塞样的中轴线与井眼竖直方向一致，利用ZEISS VersaXRM-410CT设备，扫描完成得到二维灰度图像(实验电压及功率分别为140 kV和9.9 W)，重构图像像素比例为1.2 μm/px，叠加扫描后的灰度图像即可得到样品的三维图像。由于样品切割中存在误差，使得获取到的样品外围数据及图像难以反映真实情况，因此为了避免样品切割面影响原有数据体，对CT数据及图像进行处理，提取出一个直径1.2 mm，长1.2 mm的微米CT图像圆柱体(图1)。

图1 黄土样品及微米CT实验过程图Fig.1 Loess sample and original scan image of micro-CT

1.2 模拟原理

1.2.1 孔隙网络模拟方法

黄土样品孔隙结构模拟原理参考砂岩中数字重构及渗流模拟[24-26]。为了确定CT样品的孔隙度值，使用图像处理软件ImageJ及仿真模拟Avizo 9.0软件对提取出的微米CT图像进行分割，分割过程中的阈值采取氮气所测的孔隙度值，之后利用Avizo软件提取出样品中的孔隙结构。此外，利用最大球方法提取微米CT中的数据可定量分析孔隙结构和网络[27]，通过孔隙或者喉道的横截面面积A及横截面周长P可确定其形状因子G[28]为

(1)

1.2.2 渗透率模拟方法

分析发现，流体在孔隙网络中流动时具有非常低的雷诺数(表征流体流动情况，Re≪1)，流体流速较低，相对比较稳定，多相流加速度基本可以忽略，流体流动可被视为层流，即可利用Navier-Stokes描述样品中流体流动[26]。此外，考虑到单相流体流动形式，Sciffer等[29]对Navier-Stokes方程简化，提出非压缩性Stokes等式，即

(2)

式(2)中：∇为梯度运算符；∇·为散度运算符；∇2为拉普拉斯运算符；u为流体速度；μ为流体动态黏滞度，常数，Pa·s；p为流体压力，Pa。

该式中流体为不可压缩的牛顿流体，基本呈现一种流速不随时间变化的稳态流动。

利用软件对微米CT图像中的孔隙立方体中进行渗流模拟以确定绝对渗透率时，需求解方程组1[式(2)]，并综合达西定律进行估算。而模拟过程需分别设置输入压力及输出压力(两者压力差为ΔP)，通过计算即可得到流体的速度，而流速可进一步确定流量Q，便可依据达西渗流计算样品绝对渗透率值[30]，即

(3)

式(3)中：k为绝对渗透率，m2；Q是单位时间流量，m3/s；S为孔隙立方体的横截面面积，m2；L为孔隙立方体长度；ΔP为孔隙立方体模拟出入口的压力差。

注：式(1)和式(2)计算均可通过Avizo软件模拟的绝对渗透率实验完成。

2 结果

2.1 黄土成分分析

基于微米CT实验，获取岩样三维立体结构图像。通过分析，可在不损坏岩样外观特征及孔隙网络前提下定量表征黄土岩中不同类型矿物含量及孔隙大小。分析表明，A区黄土成分由孔隙、泥粒、石英、长石、云母以及重矿物组成(图1)。孔隙具有最低的相对密度，在图像上显示为黑色的斑点。重矿物具有较高的相对密度，在CT图像上显示为亮色的斑点。利用CT图像可自动识别这两类组分。同时，由于泥粒、石英、长石、云母其相对密度大小范围接近[图2(a)]，利用CT图像无法自动识别这些矿物。通过手动自动和手动图像分割配合，提取了该样品中57个体积较大的石英及长石颗粒[图2(b)]、58个云母颗粒[图2(c)]以及多个重矿物颗粒[图2(d)]。由于泥粒边界太模糊，手动也难以识别。

图2 相对密度分布及微米CT数据提取图件Fig.2 Relative intensity distribution and segmentation of micro-CT data

2.2 矿物形状及方位研究

一般颗粒形状可被分为球形、板状、长形、叶片状等形状[31]。基于微米CT图像矿物识别结果，分别测量了黄土样品中石英、长石、云母以及重矿物的长轴长度L、短轴长度I以及厚度S[图3(a)]。结果表明，云母颗粒厚度S与长轴长度L比值较小，形状接近于板状和叶片状[图3(b)、图3(e)]。而石英、长石、重矿物具有较大的S/L[图3(c)～图3(e)]，三种矿物形状更接近于圆球状，通过不同颗粒三维长度及横截面面积进一步确定矿物的磨圆度，即

图3 矿物形状定量描述图Fig.3 Quantitative description of the shape of grains in loessite

(4)

式(4)中：S为磨圆度；V为样品体积；A为表面积。

其中，黄土矿物成分中重矿物磨圆度最大，石英、长石颗粒磨圆度相对较小，而云母颗粒磨圆度最小。虽然云母颗粒长轴较大，但相比体积占比仅次于石英长石含量，而重矿物的长轴长度及体积均较小(图4)。结果表明，黄土矿物中不同矿物的长轴长度与其磨圆度之间有较好的对应关系，一般矿物长轴长度小，其磨圆度也较小，即表明矿物整体也较小，所占体积组分在黄土整体体积中占比较小。

图4 黄土中矿物磨圆度和累计体积大小分布Fig.4 Distribution of mineral roundness and cumulative volume in loess

此外，提取不同矿物排列方位，总结不同矿物在平面上的分布情况，石英和长石颗粒方位主要在250°～0°之间，云母颗粒主要以0°、180°南北朝向分布为主[图5(a)、图5(b)]，而重矿物颗粒在平面上各个方向均有分布[图5(c)]。整体来看，石英和长石、云母、重矿物等3种成分在平面分布上规律性均比较弱。

图5 3种矿物颗粒的长轴倾角与倾向分布Fig.5 Long axis dip and dip distribution of three mineral grains

2.3 孔隙网络分析

利用图2(a)中各矿物的阈值分析，将原始微米CT图像中[图6(a)]提取出来的孔隙分割成12个小块[图6(b)]，采用最大球方法提取各部分孔隙网络空间-球棍模型图[图6(c)]。对比原始图像与提取出的孔隙网络，并分析孔隙网络中孔隙、喉道、孔隙形状因子以及配位数分布。大孔隙的空间分布非均质性较强，黄土孔隙网络中孔隙半径主要分布在2～6 μm，喉道半径主要分布在1 μm左右(该区原始压汞数据显示喉道半径的主要分布在0.037～0.547 μm，即受微米CT本身分辨率的限制，从微米CT图像中提取出的喉道数据整体有一定程度偏大)，而孔隙形状因子主要分布在0.027左右[图7(b)]，对应的孔喉配位数z=3[图7(d)][32]。

图6 微米CT孔隙网络提取过程图Fig.6 Raw image, extracted pore subvolumes and pore network of micro-CT

图7 黄土岩孔隙、喉道半径等微观参数特征Fig.7 Characteristics of microcosmic parameters such as pore and throat radius of loess rock

2.4 渗透率各向异性特征

为有效确定黄土内部孔隙度及渗透率各向异性特征，将模拟的孔隙网络分割成12个方形孔隙体(单个孔隙体尺寸：360 μm×360 μm×360 μm)，并对其三个不同方位的渗透率进行模拟[图6(b)]，采用Navier-Stokes方程及达西渗流公式确定孔隙体渗透率值，模拟结果表明：孔隙度分布在14.35%～19.28%之间，各个孔隙体的渗透率基本分布在33.95～91.22 mD(图8)，对比黄土样品原始实验分析数据(气测孔隙度：13.35%～18.93%，渗透率：38.24～86.75 mD)发现，模拟数据基本接近实测数据。一般来说，孔隙体的渗透率值应随着孔隙度的不断增大而增大，但是，不同孔隙体模拟的部分结果并不符合这种模式，例如，在X、Z方向上，孔隙度14.35%对应的渗透率高于孔隙度16.96%对应的渗透率值[图8(b)]，这个现象表明：黄土中的孔隙分布的很分散，不同孔隙体中的不同方位孔隙连接方式变化也很大。因此，孔隙的空间分布和其相互之间连接方式的差异是造成渗透率各向异性的主要因素之一。例如，不同孔隙中不同方向的渗透率均不一致，对于大多数孔隙体，Z方向的渗透率相比X、Y方向的均比较低，有一小部分渗透率值大小介于X、Y方向之间(图8)。

图8 渗透率各向异性特征Fig.8 Permeability anisotropy

3 讨论

3.1 孔隙赋存形式

图5(c)球棍模型图中大孔隙分布规律不明显，与之连接的喉道数目(配位数z)也有较大不同，显示出较强的非均质性，而图6也暗示不同孔隙体中的孔隙的空间分布和其相互之间连接方式也有较大差异。微米CT图像中3个切片(切面为XZ平面)显示(图9)：大孔隙的聚集与云母矿物的赋存形式有较好的对应关系，图中所示，云母矿物在切面中显示为长方形[图9(a)～图9(c)]，其正如在图1(c)中显示的扁平状云母颗粒一样。而云母颗粒与大孔隙之间连接比较紧密，相比其他小孔隙，被云母覆盖或者被云母侧向顶起的孔隙相对较大。以微米CT照片为基础，提取孔隙与云母二者分布之间的形态并建立相关模式图(图10)：云母矿物在黄土后期沉积过程中可直接覆盖在原有泥粒、石英以及长石颗粒中间的空隙上面，可有效防止其他较小的颗粒以及黏土矿物充填空隙，在后期黄土逐渐沉积胶结后，很大一部分大孔隙便分布在片叶状的云母颗粒周围，由于云母的赋存并没展现出明显规律，由云母赋存导致的孔隙非均质性也相对较强。

图9 扁平状云母与孔隙赋存方式、空间分布之间的关系图Fig.9 Relationship between pore reservation and spatial distribution of tabular shaped mica grains

图10 云母对孔隙留存的影响模式图Fig.10 Positive effect of mica grains in preserving pores

3.2 微观特征与渗透率各向异性关系

分析表明，孔隙空间分布、孔隙相互之间的连接是造成渗透率各向异性的重要因素。考虑到云母颗粒在孔隙形成及保存的重要作用，云母颗粒间接对渗透率各向异性的产生起着至关重要的作用。

分析提取出的孔隙体尺寸大小，云母颗粒长轴的长度为80～206 μm，石英、长石颗粒长轴的长度为50～156 μm，对比分析孔隙、喉道尺寸大小，云母、石英、长石以及泥粒似乎均可以对不同方位渗透率均可产生足够大的影响(图8)。然而，考虑到这些颗粒的不同形态，扁平状的云母颗粒对渗透率各向异性影响更大。

对于泥粒，由于其塑性较强，在沉积期间受到上覆岩层的压力，其形状易于发生变化，泥粒周围的喉道半径均较小(喉道半径主要分布范围：0.037～0.547 μm)，且利用微米CT难以识别，因此，泥粒并不是渗透率各向异性产生的主要因素。对于石英和长石颗粒，尽管两者在平面分布上有较强的方向性[图5(a)]，但由于两者具有相对高的磨圆度[主要分布区间：0.53～0.85，图4(a)]，石英和长石更易被进一步磨圆，对渗透率各向异性产生的影响也较小。而云母颗粒除在平面分布上有较强的方向性[图5(b)]，其磨圆度也较低[主要分布区间：0.27～0.69，图4(a)]。当扁平状的云母颗粒在水平方向上大量展布，即可保留大量的大孔隙，一方面提高了X、Y方向上的孔隙连通性，另一方面由于云母颗粒覆盖在孔隙之上，对垂向的孔隙连通性有一定的阻碍作用，进而会影响其垂向渗透率的变化。

4 结论

针对黄土岩微米CT图像，识别出了不同的矿物组分，分析不同矿物形状，模拟其微观孔隙网络，综合分析了黄土岩孔隙赋存状态及申渗透率各向异性，共得到以下结论。

(1)A区黄土成分主要由孔隙、泥粒、石英、长石、云母以及重矿物组成，各类矿物分布规律性不强，孔隙喉道之间连接方式变化大，孔隙空间分布非均质性强。

(2)相比磨圆度高的石英、长石、重矿物等成分，扁平状的云母颗粒在水平方向上大量展布，由于其呈板状形态、拥有较低的磨圆度，对黄土岩中孔隙后期赋存有较大影响。

(3)孔隙空间分布、孔隙相互之间的连接是造成渗透率各向异性的重要因素。当岩石中孔隙被如云母之类的矿物在垂向或侧向遮挡之后，孔隙连通性严重降低，并会间接影响黄土岩中不同方位渗透率的变化。

(4)微米CT图像及孔隙网络模拟方法相结合，可有效表征黄土岩矿物类型、大小及内部孔隙网络、渗流等特征，实验过程及模拟方法两者可操作性强，且结果与原始实验分析结果相近，表明该方法在研究岩石矿物成分、孔隙网络、渗透特性方面具有较好的效果，具有一定的可推广性。