[发明专利]一种砂岩压实作用数值模拟方法

权利要求书

1.一种砂岩压实作用数值模拟方法，其特征在于含有以下步骤：基于全部地层的顶面构造图、地层厚度、地层年代、剥蚀厚度和岩性比率资料，重现目标层位砂岩的埋藏史和热史，作为压实作用数值模型；基于沉积微相平面分布和单井岩相分布，采用确定性建模方法和随机性建模方法编制出3D岩相模型，作为压实作用数值模型的输入参数；采用铸体薄片和阴极发光薄片定量统计数据评价出粒间体积-深度关系曲线和压溶石英胶结物含量，作为机械压实作用和化学压实作用数值模拟参数；基于此，采用压实作用数值模拟方法开展压实作用数值模拟研究，预测了3D空间上孔隙度高值区；

还包括以下步骤：

步骤1、收集资料：收集地质资料、测井资料、地震资料和分析测试资料；测井资料包括中子、密度、声波时差、自然伽马、电阻率、测井解释孔隙度和测井解释渗透率曲线；分析测试资料包括铸体薄片定量统计数据、阴极发光薄片定量统计数据、扫描电镜图像资料和常规物性分析数据；

步骤2、重现目标层位砂岩的3D埋藏史和热史：基于研究区内目标层位砂岩底界构造图到地表的全部地层的顶面构造图，结合地层厚度、地层年代、岩性比率、剥蚀厚度、地层温度、地层压力和现今热流资料，重现研究区目标层位砂岩的3D埋藏史和热史，作为3D砂岩压实作用数值模型；埋藏史要能反映出研究区目标层位砂岩经历的构造运动和构造演化阶段；并且，采用实测温度数据约束校正模拟的热史，达到地层温度的实测值与模拟值的一致性，确保热史的准确性；

步骤3、建立3D沉积微相模型：基于目标层位的顶面和底面构造图，建立研究区目标层位砂岩的构造模型；基于目标层位地震资料的振幅属性和砂岩时间厚度图编制出目标层位沉积微相图，采用确定性建模方法建立目标层位砂岩的3D沉积微相模型，评价出目标层位的沉积微相空间分布；

步骤4、建立3D岩相模型：在步骤3的3D沉积微相模型基础上，结合单井岩相的纵向分布，采用序贯高斯指示随机建模方法或者机器学习建立3D岩相模型，评价出研究区目标层位砂岩的3D岩相的空间分布，作为压实作用数值模型的输入参数；

步骤5、确定压实作用数值模拟参数：压实作用数值模拟参数包括机械压实作用参数和化学压实作用参数；采用铸体薄片图像资料定量统计数据评价出每种粒径岩性的粒间体积-深度关系曲线，粒间体积-深度关系曲线是机械压实作用的一个定量表征参数，作为机械压实作用数值模拟参数；采用阴极发光薄片图像资料定量统计数据评价出压溶石英胶结物含量，压溶石英胶结物含量是化学压实作用的一个定量表征参数，作为化学压实作用数值模拟参数；

步骤6、砂岩压实作用数值模拟：联合步骤2中3D埋藏史和热史、步骤4中3D岩相模型和步骤5中压实作用数值模拟参数，作为压实作用数值模型的输入参数，进行3D砂岩压实作用数值模拟研究，获得研究区目标层位砂岩的压力史和孔隙度演化史；

步骤7、采用地层压力数据的实测值与模拟值的一致性确保压实作用数值模拟结果的准确性：基于步骤6获得的目标层位砂岩的模拟压力史，采用实测的地层压力数据约束校正目标层位砂岩的压力史，直至达到地层压力的实测值与模拟值的一致性，确保目标层位砂岩的压力史和压实作用数值模拟结果的准确性；

步骤8、压实作用数值模拟结果及孔隙度高值区：采用压实作用数值模拟方法开展3D砂岩压实作用数值模拟研究，预测了3D空间上孔隙度高值区；

步骤2中重现目标层位砂岩的3D埋藏史和热史，埋藏史能够反映出目标层位砂岩的构造运动和构造演化阶段，采用地层温度的实测值与模拟值的一致性确保热史的准确性；具体如下；

基于研究区内从目标层位砂岩到地表的全部地层的顶面构造图，结合地层厚度、地层年代、岩性比率、剥蚀厚度、地层温度、地层压力和现今热流资料，重现研究区目标层位砂岩的埋藏史和热史，作为3D砂岩压实作用数值模型；

埋藏史重建采用回剥技术和超压技术；回剥技术属于反演法，适用于正常压实带，而超压技术属于正演法，适用于欠压实带；回剥技术的过程如下，根据质量守恒法，随着目标层位砂岩埋藏深度的增加，地层厚度逐渐减小，但是地层的骨架厚度始终保持不变；

埋藏史是基于回剥技术重建的，从现今研究区分层数据出发，按地质年代逐层剥去，直到全部剥完为止；基于回剥技术的数学模型评价出各套地层在不同地质时期的地层厚度h_i，然后，根据各套地层在不同地质时期的地层厚度，求取各套地层在不同地质时期的底界埋藏深度Z_ib，在此基础上，按照地层年代由老到新逐层叠加能够得到不同地质时期各套地层的埋藏深度，进而重建目标层位地层的3D埋藏史；

热史重建包括EASY％Ro和Mackenzie模型；建立镜质体成烃动力学标定模型；并且，构造各层沉积期间的古地温梯度DG_ij和古地表温度T_ij的目标函数；求解目标函数，目标函数为极小值时的古地温梯度DG_ij和古地表温度T_ij即为所求的各套地层沉积期间的古地温梯度DG_i和古地表温度T_0j，对目标函数求一阶偏导数后，采用优化算法求解出古地温梯度DG_i、古地表温度T_0j和剩余可反应潜量X_ij0；

基于上述古地表温度和古地温梯度的获取，评价出单井的温度场；采用EASY％Ro经验评价公式R_ij0＝exp(-1.6+3.7×X_ij0)评价出研究区内单井的镜质体反射率时空场R_ij0；

根据上述评价出来的单井的温度场和镜质体反射率时空场，对比分析单井镜质体反射率时空场、温度场中模拟值与实测值，判断是否满足精度要求，如果满足则继续往下运行；否则，重新优化调整古地温梯度，重新求解目标函数，评价单井的温度场和镜质体反射率时空场，直至评价出来的单井镜质体反射率和地层温度的模拟值与实测值满足精度；

在此基础上，重建研究区目标层位单井上热流时空场：

式中，Z_ijmax为第i层在第j层沉积时期的最大埋深；DG_ij为第i层在第j层沉积时期的古地温梯度；k(n)为不同粒径岩性的热导率；n＝1、2、3、4、5和6，分别对应着粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩；式中k(n)代表不同粒径岩性的岩石热导率；根据铸体薄片图像资料定量统计数据建立粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩的岩石热导率模型，包括孔隙空间和岩石骨架；孔隙空间分为天然气和地层水；岩石骨架分为石英、长石和粘土矿物；k(n)_地层水为地层水的热导率；k(n)_天然气为天然气的热导率；Φ为孔隙度；Sw为含水饱和度；Sg为含气饱和度；w(石英)为石英百分含量；w(长石)为长石百分含量；w(粘土矿物)为粘土矿物百分含量；n＝1、2、3、4、5和6，分别对应着粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩；

根据研究区目标层位单井热流时空场，重建研究区目标层位的三维热流时空场；采用双线性插值评价公式重建出盆地内的三维热流时空场内点(x，y)上热流值hf(x，y)；

式中，(x_m，y_m)、(x_m，y_n)，(x_n，y_m)，(x_n，y_n)分别为采样点1、2、3、4个点的横坐标和纵坐标值；hf_mm、hf_mn、hf_nm和hf_nn分别是采样点1、2、3和4的单井热流；

采用实测温度数据约束模拟的热史结果，直至达到地层温度的实测值与模拟值的一致性，满足模拟温度数据与实测温度数据的精度要求ε0.15和相关系数r≥0.925，见公式(3)；

式中，Temp(i)_meas是实测的温度数值点；是实测的温度的平均值；Temp(i)_cal是模拟的温度数据点；是模拟的温度的平均值；i代表不同井点位置；n是实测温度的数量；r是相关系数；

联合步骤3的3D沉积微相模型和步骤4的3D岩相模型，评价出3D岩相空间分布，具体如下；

基于收集到的地层顶面构造图，建立研究区目标层位的构造模型；基于研究区地震资料的目标层位振幅属性和砂岩时间厚度图，编制出沉积微相图；采用确定性建模方法建立基于地质、测井和地震多资料融合的3D沉积微相模型，评价出3D沉积微相空间分布；

在3D沉积微相模型基础上，结合单井岩相的纵向分布，采用序贯高斯指示随机建模方法或者机器学习建立3D岩相模型，确定出3D岩相的空间分布；

建立了低渗透-致密砂岩储层的3D岩相模型，这为研究区目标层位砂岩的压实作用数值模拟研究提供关键参数；

步骤5中确定压实作用数值模拟参数，压实作用数值模拟参数包括机械压实作用参数和化学压实作用参数，具体内容如下；

压实作用数值模拟参数包括机械压实作用参数和化学压实作用参数；采用铸体薄片图像资料定量统计数据评价出不同粒径岩性的粒间体积-深度关系曲线，粒间体积-深度关系曲线是机械压实作用的一个定量表征参数，作为机械压实作用数值模拟参数；采用阴极发光薄片图像资料定量统计数据评价出压溶石英胶结物含量，压溶石英胶结物含量是化学压实作用的一个定量表征参数，是化学压实作用数值模拟参数；

基于不同粒径岩性的粒间体积-深度关系曲线和压溶石英胶结物含量，构建出机械压实作用模拟参数和化学压实作用模拟参数的函数模型，将压实作用模拟参数的函数模型作为砂岩压实作用数值模拟研究的输入参数；

根据铸体薄片图像资料定量统计数据，构建出不同粒径岩性的粒间体积-深度关系曲线，见公式(4)；

Lithintervolcurve(i)＝f(Vol(i)，dep(i))              公式(4)

式中，Lithintervolcurve(i)代表不同粒径岩性的粒间体积-深度关系曲线；Vol(i)代表不同粒径岩性的粒间体积，％；dep(i)代表不同粒径岩性的埋藏深度，m；i＝1、2、3、4、5和6，分别对应着粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩；

基于不同粒径岩性的粒间体积-深度关系曲线，构建出不同粒径岩性机械压实作用参数的指数函数模型；根据粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩的粒间体积-深度关系曲线，采用Athy模型的公式(5)分别拟合出粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩的Athy指数模型的参数k(i)，作为粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩机械压实作用的函数模型参数；

φ(z)(i)＝φ(i)₁+(φ(i)₀-φ(i)₁)e^-k(i)z(i)    公式(5)

式中：φ(i)₀是初始孔隙度，％；k(i)是Athy指数模型的参数；φ(i)₁是最小孔隙度，为1％；z(i)是深度，m；i＝1、2、3、4、5和6，分别对应着粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩；

采用阴极发光薄片图像资料定量统计数据评价出压溶石英胶结物含量，压溶石英胶结物含量作为化学压实作用数值模拟参数；压溶石英胶结物含量作为石英胶结作用模型的约束参数；石英胶结作用被模拟为受温度和能够用于沉淀的石英表面积控制的沉淀作用限制的反应；石英胶结作用孔隙度损失率用公式(6)表示：

式中，C(i)是石英颗粒膜因子，覆盖黏土膜的、不可能沉淀的石英颗粒表面分数0-1，数值1阻止石英胶结作用；f(i)是石英颗粒体积分数，分布范围为0-1；D(i)是石英颗粒尺寸，代表石英颗粒的平均尺寸；φ(i)₀是初始孔隙度，％；φ(i)是地质时期内不同时期的孔隙度值，％；k(i)是Athy指数参数；Ac(i)和Eca(i)分别是石英胶结作用的频率因子和活化能；M(i)和ρ(i)分别是石英摩尔质量和密度；i＝1、2、3、4、5和6，分别代表粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩；这里，采用粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩的压溶石英胶结物含量作为约束参数，分别拟合出粗砂岩、砂质砾岩、中砂岩、细砂岩、粉砂岩和泥岩的频率因子和活化能。