页岩微观孔隙模型构建——以滇黔北地区志留系龙马溪组页岩为例

陈科洛，张廷山，陈晓慧，何映颉，梁兴

（1. 西南石油大学地球科学与技术学院，成都 610500；2. 四川省地质矿产勘查开发局测绘队，成都 610017；3. 中国石油浙江油田公司，杭州 310023）

0 引言

页岩气作为非常规油气的重要组成部分，以其资源量巨大、清洁、高效等特点而备受重视[1-2]。上扬子地区作为中国南方海相页岩气的有利勘探区，前人已对其岩相与沉积环境、微观孔隙结构和发育控制因素以及地球化学特征等方面进行了论述[3-6]，为页岩储集层评价奠定了基础。但页岩中微米—纳米级孔隙的复杂结构系统，给储集层评价与页岩气开发带来了困难。前人对页岩孔隙结构的评价方法多是在镜下研究的基础上结合压汞实验，以孔隙度、渗透率、孔喉中值等为参数进行分析[7]，但传统实验受仪器设备和操作因素影响大，很难达到地层温度和压力范围，且不能从分子级别的微观角度反映出页岩气吸附的本质。Tahmasebi等[8]基于交互式模拟（CCSIM）与信息熵的迭代算法对页岩进行了模型三维表征但缺少实验佐证。Yan等[9]提出基于双孔隙模型，并以TOC值作为主要建模参数对孔隙模型进行搭建的方法，对有机质孔隙进行了较好的表征，但该方法对复杂孔隙系统的验证仍需改进。以 CT扫描技术为基准建立的数字岩心表征对页岩孔隙结构在一定程度上具有较好的效果[10-12]，但是在模型验证上仍存在困难。本文基于分形理论，并根据扫描电镜和低温氮气吸附实验建立微观孔隙模型。基于分形Frenkel-Halsey-Hill公式（以下简称FHH公式）与数理统计理论，建立 Modified Frenkel- Halsey-Hill数学公式（以下简称 MFHH公式），对所建模型进行验证。模型验证中数学方法的提出，为页岩储集层模型表征提供了有利的依据，是模型定量化分析与验证的方向之一。

1 页岩微观孔隙特征

页岩微观孔隙具有不同的大小、形态与空间分布，其中微米—纳米级孔隙对页岩孔渗条件、烃类富集具有重要控制作用[13]。本文以滇黔北地区龙马溪组页岩为例（见图 1），对页岩微观孔隙特征进行分析。

图1 研究区位置与Y9井地层柱状图

1.1 微观孔隙类型

孔隙类型可以指示页岩所经历的沉积、成岩演化阶段，对孔隙模型的搭建起到重要的指导意义[14]。利用扫描电镜对研究区龙马溪组页岩样品进行系统观察分析后发现：目的层段孔隙类型多样，包括矿物铸模孔、粒间溶孔、粒内溶孔、黄铁矿晶间孔、黏土矿物层间孔、残余原生粒间孔与有机质孔等（见图2），分别代表了不同的孔径与成因类型，其中黏土矿物层间孔和有机质孔为主要孔隙类型。由于奥陶系五峰组—志留系龙马溪组页岩具有高演化、过成熟的特点，黏土矿物中高岭石、蒙脱石等不稳定矿物组分已大部分向伊利石转化，因此黏土矿物组成以伊利石与伊蒙混层为主。由于油气已经历了生烃高峰期，因而有机质孔隙较为发育。

图2 龙马溪组页岩孔隙与裂缝扫描电镜照片

页岩具有复杂的孔隙系统，孔径分布范围广，因划分依据不同对孔径分类方法存在差异。Loucks等[15]依据孔径大小将页岩孔隙划分为纳米孔（小于 0.75 μm）与微孔（大于等于0.75 μm）。目前对页岩孔隙分类被较为广泛接受的是 IUPAC（国际纯粹和应用化学联合会）孔隙分类方法[16]。本文参考IUPAC的分类方案将孔隙分为3类：微孔（孔径小于2 nm）、介孔（孔径为2～50 nm）和大孔（孔径大于50 nm）[17]。

1.2 微观孔隙结构特征

孔径大小与孔隙结构往往存在一定的对应关系。受限于扫描电镜的放大倍数，电镜下只能识别大孔和介孔。研究区不同孔隙类型的孔径差异较大，同一类型的孔径大小在不同深度也存在一定差异，但相同孔隙类型在不同井间孔径范围差异相对较小。研究区矿物铸模孔、粒内溶孔、粒间溶孔、黄铁矿晶间孔、黏土矿物层间孔与有机质孔孔径分别为 20.00～150.00 μm、0.08～15.00 μm、0.12～5.50 μm、0.02～0.16 μm、0.01～0.60 μm 与 0.01～0.50 μm（见表1）。

基于以上页岩孔隙特征，采用氮气吸附法对其微孔与微观孔隙结构进行识别分析。运用美国Quantachrome公司生产的Quadrasorb SI型比表面积与孔隙度分析仪进行测试。该仪器比表面积可测下限为0.000 5 m2/g，孔体积最小可测值为 0.000 1 cm3/g，孔径测量范围为0.35～400.00 nm，吸附-脱附相对压力范围为 0.004～0.995。不同孔径与氮气吸附机理具有不同对应关系：氮气在微孔材料上的主要吸附机理为单分子层吸附和微孔填充；介孔介质在低压区的吸附机理为单分子层吸附，在中等压力区为多分子层吸附，在较高压力区则发生毛细孔凝聚现象；大孔介质在低压区的吸附机理与介孔相同，但在相对压力较高时不发生毛细孔凝聚现象[16]。

表1 电镜下孔缝直径范围

孔隙比表面积、孔体积、平均孔径是孔隙结构描述中的重要表征参数。在液氮温度下，流经页岩的氮气被吸附饱和，在温度逐渐升至室温的过程中，被页岩吸附的氮气出现脱附现象[5]。依据国标 GB/T 19587—2004《气体吸附BET法测定固态物质比表面积》[17]，页岩样品比表面积计算采用BET方程，在相对压力为0.05～0.35时通过作BET图求得单分子层饱和吸附量，进而计算样品的比表面积[16]。由于氮气吸附量与页岩孔径大小有关，因此可采用BJH法测定并计算出孔体积和孔径分布，平均孔径则由相对压力约为0.993时的氮气吸附量计算获得。通过测试、计算发现：研究区页岩比表面积为4.365～23.015 m2/g，平均值为14.317 m2/g；孔体积为 0.005～0.049 cm3/g，平均值为 0.017 cm3/g；平均孔径为1.475～4.046 nm，平均值为2.200 nm（见表2）。较大的比表面积、孔体积与平均孔径将有利于页岩气的吸附聚集。

氮气吸附实验得出的吸附与脱附曲线组合形态可反映出不同的孔隙类型与孔隙结构特征，不重合的吸附与脱附曲线组成吸附回环（滞后环），吸附回环的形态与出现时间对孔隙结构的判定尤为重要。

IUPAC将吸附回环分为4类：H1型，吸附回环较为狭窄且出现时相对压力较大，吸附曲线与脱附曲线近于平行排列且在吸附回环出现的较窄压力范围内与压力轴近于垂直，指示样品孔径分布范围较小，孔隙类型为两端开放的圆筒或柱形毛细孔；H2型，吸附回环较宽，脱附曲线较吸附曲线斜率更大并出现快速下降的拐点，指示样品孔径分布范围较大，孔隙类型为细颈和广体孔或墨水瓶孔；H3型，吸附回环较窄，随着相对压力的增大，吸附与脱附曲线的斜率逐渐变大，上升速度较快，在平衡压力接近饱和蒸汽压时仍未出现吸附饱和现象，通常指示两端开口的楔形孔或似片状颗粒组成的槽状孔。H4型，吸附回环较窄，吸附与脱附曲线随相对压力的增大而缓慢上升，在平衡压力接近饱和蒸气压时出现吸附饱和现象，通常指示狭缝孔（微孔）（见图3）[18-19]。

表2 孔隙结构表征参数

研究区龙马溪组样品多在平衡压力接近饱和压力时出现吸附饱和现象，吸附曲线总体表现为低压处向上翘起、高压处陡增、中等压力缓慢上升的特点（见图 4a）。脱附曲线在低压与高压处特征与吸附曲线相同，但在中等压力处存在拐点（即吸附量随相对压力降低快速减少后呈缓慢减少的趋势）。脱附曲线随着相对压力降低而缓慢下降的过程证实了微孔的发育，而脱附曲线随着相对压力降低而快速下降的过程佐证了槽状孔占有较大比例。从孔径与孔体积关系图分析得出：对孔体积贡献最大的为3～5 nm的介孔，其次为微孔，大孔贡献的孔体积相对较小（见图4b）。综上分析，研究区样品以槽状孔为主并大量发育微孔。槽状孔对应黏土矿物层间孔、微孔对应有机质孔隙，电镜分析与实验测试结论具有较好的一致性。

图3 吸附回环分类及孔隙类型（据参考文献[18]修改）

图4 典型样品孔隙结构表征

2 页岩微观孔隙模型构建

微观孔隙模型的构建与验证基于扫描电镜（SEM）观察与氮气吸附实验测试结果。孔隙类型对孔隙模型构建有重要指示意义，实验测定的比表面积、孔体积与平均孔径对分形维数大小起到控制作用。

2.1 模型构建

储集层孔隙结构表征研究初期是通过欧几里得的几何学进行，但是经典几何学难以对不规则、具有强非均质性的孔隙系统进行恰当的描述，因此常常采用均值、分选、歪度等统计学参数对孔隙结构进行表征[20]。

自法国数学家曼德勃罗特提出分形几何之后，随着人们对世界认知的进一步加深，分形思想在研究中越来越受重视[21]。分形研究的是自然界中极其零碎复杂，但又具有自相似性、自反演性的体系，可以用分形维数进行表示。公式如下：

不少学者已经提出页岩具有分形特征[22-25]，并用分形维数的大小对页岩非均质性强弱进行判断。通常将绝对光滑表面的分形维数定为 2，将具有能提供大量吸附点位的极端粗糙表面的分形维数定为 3[26]。分形维数的确定方法有多种，其中气体吸附法中的分形 Langmuir公式、分形 BET公式和分形 FHH公式等已被广泛运用[27]，特别是分形 FHH公式已应用于多孔介质孔隙结构特征的分形维数计算。分形FHH公式为：

页岩具有规则分形特征，意味着无论怎样缩小（或放大）观测尺度去观察，其组成部分与原来没有区别，也就是说其具有无限的膨胀和收缩对称性。

分子大小与纳米尺度相当，因而分子模型在一定程度上能够从微观角度揭示纳米孔隙结构。基于分形理论，并通过等效模拟、空间上无限延展便能进一步表征微米—毫米级的孔隙系统。

有机质（干酪根）的组成和结构受有机物来源与沉积环境等诸多因素影响，往往没有严格固定的组成与结构[28]。国内外学者普遍认为石墨烯层状结构能在一定程度上对有机质（干酪根）孔隙进行表征[29-32]。鉴于研究区龙马溪组页岩黏土矿物以伊利石、伊蒙混层为主，孔隙类型以黏土矿物层间孔和有机质微孔为主，因此运用分子模拟软件 Materials Studio，以伊利石分子模型的空间排列表征黏土矿物层间孔、石墨烯分子模型的空间排列表征有机质孔，通过搭建分子层指令来建立页岩孔隙分子级模型（见图5）。建模中分形理论无限延展特性可将分子模型在空间上扩展到微米甚至毫米尺度。

2.2 模型验证

前人对模型的搭建多在定性分析或简易模型的基础上，基于一定的数学算法生成页岩模型[8-9]，虽然在一定尺度上具有较好的表征效果，但在模型验证方面往往缺乏足够的说服力。本文研究在氮气吸附实验的基础上，通过改进分形FHH公式并结合甲烷与二氧化碳竞争性吸附实验研究，提出MFHH验证公式：

图5 模型表征示意图

本文研究中p、V、V0、p0可通过氮气吸附实验获得，需求取的关键性参数为D、η与C1、C2。模型验证基本思路如下。

分形维数虽能定量表征孔隙非均质性的强弱，但是不同地区、不同模型对孔隙的表征仍然存在差异性与局限性，因此引入分形修正系数C1，C2。由于C1，C2用于表征孔隙系统，且研究区龙马溪组主要孔隙类型与孔隙结构已经确定，不同井间C1，C2数值应在较小范围内波动，即可假定C1，C2服从正态分布。选取一口标准井的C1，C2值分别作为 2个正态总体的数学期望，如果所建立的模型能够从整体角度在一定程度上较好地表征研究区页岩孔隙结构，则研究区其他井间C1、C2的数学期望应与标准井相同。

2.2.1 分形维数

本次研究中分形维数的计算基于分形FHH公式，前人指出（2）式中 lnV与 ln[ln(p0/p)]关系式的斜率A与D存在明确关系，但仍存争议[33-34]。主要计算方法有以下2种：

通过对（4）式、（5）式的计算发现，（4）式算出的分形维数存在偏差，而（5）式与前人的研究结果一致。比表面积、孔体积与分形维数具有较好的线性相关性，平均孔径与分形维数相关性不明显。比表面积、孔体积越大，分形维数越高，孔隙非均质性越强（见表2、见图 6）。

图6 比表面积、孔体积与分形维数关系

2.2.2 置换率

前人对甲烷与二氧化碳的竞争性吸附已有研究[35-36]，但是对二氧化碳在地层温压条件下对甲烷的置换定量表征研究较少[37-39]。以Materials Studio软件作为平台，通过模拟地层温压条件，在模型中先饱和甲烷并使结构优化稳定，之后注入二氧化碳与甲烷发生竞争性吸附。置换率是指被二氧化碳所置换出的甲烷体积与未注入二氧化碳之前甲烷体积的百分比。

结合滇黔北地区实际地质情况，设定地表温度为283 K、压力为0.1 MPa、地温梯度为30 K/km以及压力梯度为15 MPa/km。通过模拟后分别得出0，2，3 km埋深的地层条件下，孔隙壁面法线方向单吸附时甲烷密度曲线（见图7a—图7c）与混合吸附时甲烷密度曲线（见图7d—图7f），通过计算曲线与横轴所围成的面积便可得到置换率，如（6）式所示：

为便于计算，采用等面积割补法对复杂曲线进行等效替代（见图 8），依（6）式计算得出：埋深为 0时η为2%、2 km埋深η为26.1%、3 km埋深η为21.7%。

根据模拟与计算结果得出埋深与置换率的关系图。从埋深0到2 km，由于温度、压力的增加，甲烷与二氧化碳分子布朗运动加剧，置换率快速增大并达到峰值，随着埋深达到3 km，压力逐渐占据主导，虽然甲烷分子运动加剧但受到压力的制约反而置换率降低（见图9）。

图7 甲烷在孔隙壁面法线上的密度分布

图8 甲烷密度的等面积割补

2.2.3 分形修正系数

将所获得参数全部代入（3）式即可算得不同井的分形修正系数C1，C2。由于 Y9井电镜与实验测试资料最为详实，且分形维数为 2.6～2.8，与研究区总体范围相当，具有较好的代表性，因此以 Y9井C1，C2数值作为标准（数学期望），基于数理统计理论与模型验证思路，将模型验证问题转变为数学问题。

图9 埋深与置换率关系

已知两个结构表征参数C1，C2，其数学期望分别为-0.56，0.43，假定两个实际参数数值均服从正态分布。现已获得C1数值分别为-0.53，-0.46，-0.62，-0.48，-0.55，C2数值分别为 0.41，0.37，0.47，0.39，0.46，试判断C1，C2的数学期望E（C1）、E（C2）是否分别为-0.56，0.43。

假定E（C1）值为-0.56，E（C2）值为0.43。已知C1，C2均服从正态分布，即且未知。根据数理统计分布理论，在此假设条件下，有式（7）成立：

（7）式服从t（4）分布。将与C1相关的数值代入（7）式运算可得：判别值t1等于 1.14；将与C2相关的数值带入（7）式可得：判别值t2等于-0.51。本次研究假定显著性水平α为0.05，通过查表可得右侧边界值tα/2（4）为2.78，由于t分布函数关于y轴对称，因此左侧边界值为-2.78，那么拒绝域为X大于2.78或X小于-2.78（见图10）。计算所得数值1.14与-0.51均没有落在拒绝域上，接受假设，即E（C1）等于-0.56、E（C2）等于0.43，因此笔者认为所建立的模型在显著性水平α等于 0.05（正确率 95%）的条件下是可以接受的（见表3）。

图10 模型验证的t分布假设检验

表3 模型判别参数与判别结果

3 结论

研究区龙马溪组页岩孔隙类型多样，以黏土矿物层间孔和有机质微孔为主。氮气吸附实验表明孔隙内部结构以似片状颗粒组成的槽状孔和狭缝孔（微孔）为主。扫描电镜分析结果与实验结果具有一致性。

基于分形理论，运用Materials Studio软件搭建了以石墨烯和伊利石为分子单元的分子模型用以表征研究区页岩孔隙结构。

基于分形FHH公式，结合二氧化碳与甲烷竞争性吸附实验，提出MFHH公式用于模型验证的方法，并证实了分子模型基本可靠。

通过对页岩微观储集空间的直接观察、间接定量测定和分子模型构建，实现对页岩非均质微观孔隙系统的定性与定量表征。

符号注释：

A——lnV与 ln[ln(p0/p)]关系式的斜率，无因次；C1——截距修正系数，无因次；C2——斜率修正系数，无因次；D——分形维数，无因次；f(C)——检验函数的因变量，无因次；K——未知参数，无因次；n——单个参数对应数值的数目，无因次；N(λ) ——观测数目，无因次；p——平衡压力，MPa；p0——气体吸附的饱和蒸汽压，MPa；S——修正样本方差，无因次；S1——混合吸附时甲烷密度曲线与横轴围成的面积，cm2；S2——单吸附甲烷时甲烷密度曲线与横轴围成的面积，cm2；T(C)——检验函数的自变量，无因次；t——数学验证的判别值，无因次；V——平衡压力下所吸附的气体量，cm3/g；V0——单分子层吸附的气体量，cm3/g；X——随机变量，无因次；——参数的算术平均值，无因次；α——显著性水平，无因次；η——地层条件下二氧化碳对甲烷的置换率，%；λ——观测长度，m；μ——参数的数学期望，无因次；σ——参数的标准差，无因次。