直接确定地层有机碳含量测井方法的数值模拟

袁 超，周灿灿，薄亚杰，张 锋，李 霞，程相志（．中国石油勘探开发研究院测井与遥感技术研究所，北京 0008；．中国石油测井有限公司技术中心，陕西西安 70077；．中国石油大学（华东）地球科学与技术学院，山东青岛 66580）

直接确定地层有机碳含量
测井方法的数值模拟

袁 超1，周灿灿1，薄亚杰2，张 锋3，李 霞1，程相志1
（1．中国石油勘探开发研究院测井与遥感技术研究所，北京 100083；2．中国石油测井有限公司技术中心，陕西西安 710077；
3．中国石油大学（华东）地球科学与技术学院，山东青岛 266580）

摘要：非常规资源逐渐成为我国油气勘探的重要领域，其“源储一体”的特性使得对源岩特别是有机碳含量的连续定量评价成为关键。传统测井评价方法均是基于测井资料的间接评价，受多种因素影响，计算结果有一定的不确定性。本文利用快中子与地层元素原子核发生反应产生的非弹性散射伽马能谱和俘获伽马能谱，通过谱分析技术获取地层元素含量，根据与无机碳相关的钙、镁等元素含量确定无机碳含量，从总碳含量中扣除无机碳含量可直接提供高精度的有机碳含量信息。借助于蒙特卡罗数值模拟方法，建立了6种四川盆地龙马溪组和筇竹寺组页岩地层模型，模拟不同地层条件下的混合伽马能谱，并利用模拟元素的标准伽马能谱对模拟能谱进行谱解析。模拟结果表明，本文的方法可直接提供高精度的有机碳含量，有机碳含量大于2%时，其计算值相对误差小于6%。

关键词：有机碳含量直接确定；脉冲中子测井；伽马能谱；数值模拟

地层有机碳含量的确定有利于提高对储层产能和地下油气的评价，减小预测储层油气产能的不确定性［1］。虽然钻井取心分析可提供高精度的地层有机碳含量，但具有取心成本高、岩心样品有限、分析费用高、耗费时间长等不足，而且不能连续表征地层的有机碳含量信息［2］。从20世纪80年代开始，一些学者就尝试利用测井资料评价地层有机碳含量，Schmoker等［3－5］假设富有机质储层中地层密度变化是由于低密度的有机质存在引起的，利用密度测井评价地层有机碳含量；Fertl等［6］利用自然伽马能谱测井中获得的铀含量，根据区域测井值与岩心分析结果建立铀含量与有机碳含量的经验关系，推算未知地区有机碳含量；Passey等［7－8］将经过适当刻度的孔隙度曲线重叠到电阻率曲线上，由于孔隙度曲线响应于低密度、低声波速度的干酪根，电阻率曲线响应于地层流体，因此在富有机质层段两条曲线分开，根据其分开的距离可确定地层有机碳含量。随着技术的发展，测井新技术也融入到有机碳含量评价中，Jacobi等［9］结合核磁共振测井和密度测井评价有机碳含量，由于核磁共振测井响应于地层孔隙中的流体，而密度测井响应于整个岩石，因此两种测井方法确定岩石骨架密度的差值可反映地层有机碳含量；John等［10］对Jacobi等的方法做了改进，利用元素测井获得的元素含量精确估算骨架密度，使地层有机碳含量评价精度提高。但上述测井方法均是利用测井资料的间接评价，需要岩心资料刻度而建立经验公式，计算结果精度有待进一步提高。本文拟采用脉冲中子源和伽马探测器，同时获取快中子与地层元素原子核发生非弹性散射反应产生的非弹性散射伽马能谱以及快中子慢化到热中子与元素原子核发生辐射俘获反应产生的伽马能谱，通过能谱解析技术获得地层元素含量，根据与无机碳相关的钙、镁等元素含量确定无机碳含量，从总碳含量中扣除无机碳含量，以直接确定地层有机碳含量。

1　测井方法原理

脉冲中子源向地层发射的14MeV高能中子与地层元素原子核发生非弹性散射反应产生非弹性散射伽马射线，同时高能中子慢化到热中子与地层元素原子核发生辐射俘获反应产生俘获伽马射线。中子与地层元素原子核发生反应产生的伽马能谱可认为是地层中各单元素产生的伽马能谱（称为元素标准伽马能谱）的线性叠加［11］，因此测量伽马能谱第i道的计数xi为：

其中：aij为各元素标准伽马能谱的归一化谱组成的矩阵元；yj为第j种元素的相对产额；εi为误差；m为地层中元素种类数。

将式（1）写为矩阵形式，即为：

其中：X为记录的伽马能谱计数矩阵；A为元素标准伽马能谱矩阵；Y为元素相对产额矩阵。

其中：Si为元素i的探测灵敏度因子［12］；F为随深度变化的归一化因子。

俘获伽马能谱解析中，利用氧化物闭合模型［13］，即地层中所有元素的氧化物或碳酸盐矿物的质量分数之和为1，可以确定元素相对产额向元素含量转换的深度归一化因子Fc：

其中，Rj为元素j的氧化物或碳酸盐与元素j的质量比，称为元素j的氧化物指数。

由于碳元素俘获截面很小，从俘获伽马能谱中不能获得碳元素含量，因此需要利用非弹性散射伽马能谱。但由于非弹性散射伽马能谱解谱获得的元素相对产额种类有限，不能利用氧化物闭合模型获取深度归一化因子，因此需对俘获伽马能谱和非弹性散射伽马能谱进行联合谱解析获取元素含量。

俘获伽马能谱和非弹性散射伽马能谱联合谱解析中，采用参考元素法（即选取从俘获伽马能谱和非弹性散射伽马能谱中都能获得相对产额且地层中含量较多的元素作为参考元素），砂岩和碳酸盐岩地层中一般分别选用硅和钙元素。根据从俘获伽马能谱和非弹性散射伽马能谱中获取的参考元素含量在理论上是相等的，结合式（3）可得：

其中：WIR为由非弹性散射伽马能谱获取的参考元素含量；WCR为由俘获伽马能谱获取的参考元素含量；YIR为由非弹性散射伽马能谱获取的参考元素相对产额；SIR为由非弹性散射伽马能谱确定参考元素含量的探测灵敏度因子（由仪器刻度可得出的已知量）；FI为非弹性散射伽马能谱解析中元素相对产额向元素含量转换的深度归一化因子。

由式（5）可得到非弹性散射伽马能谱解析中的深度归一化因子FI：

因此，由非弹性散射伽马能谱解谱获取的待求解元素相对产额转换为元素含量可表示为：

其中：W

IUK

为由非弹性散射伽马能谱获取的待求解元素含量；Y

IUK

为由非弹性散射伽马能谱解谱获取的待求解元素相对产额；S

IUK

为由非弹性散射伽马能谱确定待求解元素含量的探测灵敏度因子（由仪器刻度可得出的已知量）。

通过对俘获伽马能谱和非弹性散射伽马能谱进行联合谱解析，可获取地层的总碳含量TC。地层中无机碳主要存在于方解石、白云石、菱铁矿、铁白云石和菱锰矿等矿物中［14］，利用与无机碳相关的钙、镁、铁、锰等元素含量可获得地层中无机碳含量TIC，从谱解析中获得的总碳含量TC中扣除无机碳含量TIC即可获得有机碳含量TOC：

TOC＝TC－TIC（8）

2　岩石物理模型及数值计算建模

2．1 岩石物理模型

根据四川盆地龙马溪组和筇竹寺组页岩的地层特征［15］，建立了岩石物理模型，如图1所示。岩石由骨架、干酪根、黏土矿物和孔隙组成，其中骨架由石英、正长石、斜长石、方解石、白云石和黄铁矿组成，黏土矿物为绿泥石和伊利石。

图1　 岩石物理模型Fig．1 Petrophysical model of rock

2．2 数值计算建模

利用蒙特卡罗方法［16］建立三维地层模型，采用MCNP5C软件包，首先模拟快中子与地层元素原子核发生非弹性散射反应以及快中子慢化到热中子后与元素原子核发生辐射俘获反应产生的伽马能谱的计数，进而用F8计数模拟记录伽马探测器的响应能谱。

数值计算建模条件如下：井眼直径为20cm，井眼充满淡水；地层为柱状模型，地层径向厚度为60cm、高度为150cm，将整个地层划分为1 470个栅元，地层充填物质按照图1所示岩石物理模型中的物质设定；伽马探测器采用高分辨率LaBr3∶Ce（5%）晶体探测器，能量分辨率设为5%，探测器直径和长度分别为5cm和10cm，探测器外包0．5cm的镉套；中子源采用D－T脉冲中子发生器，脉冲宽度为40μs；测井仪器直径为11．4cm，伽马探测器源距（中子源与伽马探测器之间的距离）为40cm，脉冲中子源与伽马探测器之间加屏蔽体，计算模型如图2所示。

3　数值模拟结果及分析

利用图2所示的蒙特卡罗计算模型，模拟计算设定地层的非弹性散射伽马能谱和俘获伽马能谱，并模拟用于能谱解析的元素标准伽马能谱，通过分析直接得出地层有机碳含量。

图2　 数值计算模型示意图Fig．2 Schematic drawing of numerical calculation model

3．1 混合能谱及元素标准能谱模拟

根据图1所示岩石物理模型设定6种地层模型，各模型组分的质量分数如表1所列。各地层模型中，正长石、斜长石、方解石、白云石、伊利石、绿泥石和黄铁矿的质量分数均分别为5%、15%、5%、2%、20%、5%、3%，干酪根、孔隙流体和石英的质量分数不同。

表1　 地层模型组分的质量分数Table 1 Mass fraction of components in different formation models

利用图2所示的蒙特卡罗计算模型，根据表1所列的6种地层模型的组分含量，利用中子截断方法模拟记录0～40μs时间窗内的非弹性散射伽马能谱计数，并利用非中子截断方法模拟记录50～100μs时间窗内的俘获伽马能谱计数，进而由得到的伽马能谱计数模拟伽马探测器响应能谱，如图3、4所示。

由式（1）可知，从测量的混合伽马能谱中获取元素相对产额需要已知元素标准伽马能谱。利用同样的计算模型，地层物质设定为及真空，利用中子截断方法模拟记录0～40μs时间窗内伽马能谱，获得Si、O、Ca、Al、C、Mg、S、Fe及仪器的非弹性散射标准伽马能谱计数；地层物质设定为

图3　 模拟地层非弹性散射伽马能谱Fig．3 Inelastic gamma energy spectrumof simulated formation

图4　 模拟地层俘获伽马能谱Fig．4 Capture gamma energy spectrumof simulated formation

3．2 结果分析及讨论

根据式（1），利用图5、6所示的非弹性散射标准伽马能谱和俘获标准伽马能谱，对图3和图4中模拟得到的混合伽马能谱，采用加权最小二乘法解谱可获得元素相对产额。进而，借助于式（4）氧化物闭合模型和式（7）联合谱解析方法，从俘获伽马能谱和非弹性散射伽马能谱中获取地层元素含量。利用与无机碳相关的钙、镁元素含量确定方解石、白云石等碳酸盐矿物中无机碳含量，从总碳含量中扣除无机碳含量得到地层有机碳含量：

图5　 非弹性散射标准伽马能谱Fig．5 Inelastic standard gammaenergy spectrum

图6　 俘获标准伽马能谱Fig．6 Capture standard gamma energy spectrum

模拟计算的6种地层模型中，部分元素解谱结果、有机碳含量计算值和理论值及计算值的绝对误差和相对误差列于表2。有机碳含量计算值与理论值的线性拟合方程为y＝1．065 5x＋0．001 28，R2＝0．993。

表2　 部分元素解谱结果及有机碳元素含量计算值与理论值对比Table 2 Unfolded results of some elements and comparison of theoretical and calculated value of TOC

由表2及线性拟合方程可看出，利用非弹性散射伽马能谱和俘获伽马能谱，通过联合谱解析方法获取元素含量，根据与无机碳相关的钙、镁等元素确定地层无机碳含量，从总碳含量中扣除无机碳含量可直接提供高精度的有机碳含量信息。6种地层模型模拟计算结果中，有机碳含量计算值绝对误差小于0．3%；有机碳含量小于2%时，由于计算相对误差基值较小，其相对误差大于10%；有机碳含量大于2%时，计算值相对误差均小于6%。

4　结论

1）传统有机碳含量确定方法是利用其他测井资料的间接评价，计算结果精度有待进一步提高。本文利用快中子与地层元素原子核发生核反应产生的非弹性散射伽马能谱和俘获伽马能谱，通过联合谱解析技术获取地层元素含量，从总碳含量中扣除无机碳含量直接得出地层有机碳含量。

2）借助于蒙特卡罗数值模拟方法，建立了6种四川盆地龙马溪组和筇竹寺组页岩地层模型，由有机碳含量计算结果知，本文的方法可直接提供高精度的有机碳含量信息，有机碳含量大于2%时，计算结果相对误差小于6%。

本文结果是利用数值模拟理论研究有机碳含量直接确定方法，但实际测井过程中井眼及地层条件会影响伽马能谱数据的测量质量，有机碳含量确定结果精度可能会受影响。

参考文献：

［1］ PEMPER R，HAN X G，MENDEZ F，et al． The direct measurement of carbon in wells containing oil and natural gas using apulsed neutron mineralogy tool［C］∥SPE Annual Technical Conference and Exhibition．New Orleans：［s．n．］，2009．

［2］ 霍秋立，曾华森，付丽，等．Δlg R测井源岩评价方法的改进及其在松辽盆地的应［J］．吉林大学学报：地球科学版，2011，41（2）：586－591．HUO Qiuli，ZENG Huasen，FU Li，et al．The advance ofΔlg R method and its application in Songliao basin［J］．Journal of Jilin University：Earth Science Edition，2011，41（2）：586－591（in Chinese）．

［3］ SCHMOKER J．Determination of organic content of appalachian devonian shales from formation－density logs［J］．AAPG Bulletin，1979，63：1 504－1 537．

［4］ SCHMOKER J．Organic content of devonian shale in western appalachian basin［J］．AAPG Bulletin，1980，64：2 156－2 165．

［5］ SCHMOKER J，HESTER T．Organic carbon in Bakken Formation，United States Portion of Williston basin［J］．AAPG Bulletin，1983，67：2 165－2 174．

［6］ FERTL W H，CHILANGAR G V．Total organic content determined from well logs［J］．SPE Formation Evaluation，1988，3（2）：407－419．

［7］ PASSEY Q R，CREANEY S，KULLA J B，et al．A practical model for organic richness from porosity and resistivity logs［J］．AAPG Bulletin，1990，74：1 777－1 794．

［8］ 朱振宇，刘洪，李幼铭．Δlog R技术在烃源岩识别中的应用与分析［J］．地球物理学进展，2003，18（4）：647－649．ZHU Zhenyu，LIU Hong，LI Youming．The

［9］ JACOBI D，GLADKIKH M，LECOMPTE B，et03－07．al．Integrated petrophysical evaluation of shale

［10］JOHN Q，ERIC M，GREG P，et al．A compari－estimation［C］∥AAPG Annual Convention and son of core and well log data to evaluate porosity，Exhibition．Pennsylvania：［s．n．］，2013．TOC，and hydrocarbon volume in the eagle ford

［11］张锋，刘军涛，冀秀文，等．地层元素测井技术最ZHANG Jingping，TANG Shuheng，GUO新进展及其应用［J］．同位素，2011，24（增刊）：Dongxin．Shale gas favorable area prediction of 21－28．the Qiongzhusi formation and Longmaxi forma－ZHANG Feng，LIU Juntao，JI Xiuwen，et al．tion of lower palaeozoic in Sichuan basin［J］．Development of formation element logging tech－Geological Bulletin of China，2011，31（2）：357－nique and its application［J］．Journal of Isotopes，363（in Chinese）．2011，24（Suppl．）：21－28（in Chinese）．

［12］HERRON S．Method and apparatus for determi－Carlo N－particle transport code［R］．US：Los ning elemental concentrations forγray spectros－Alamos National Laboratory，2000．copy tools：US，5471057［P］．1995－11－28．

［13］WRAIGHT P D，HERTZOG R C，DARWIN V analysis and application ofΔlog R method in theE，et al．Logging apparatus and method for source rock’s identification［J］．Progress in Geo－determining absolute elemental concentrations of physics，2003，18（4）：647－649（in Chinese）．subsurface formations：US，4810876［P］．1989－

［14］CHARSKY A，HERRON S．Accurate，direct gas reservoirs［C］∥CIPC／SPE Gas Technologytotal organic carbon（TOC）log from a new Symposium 2008Joint Conference．Calgary：［s．advanced geochemical spectroscopy tool：Comn．］，2008．parison with conventional approaches for TOC

［15］张静平，唐书恒，郭东鑫．四川盆地下古生界筇shale［C］∥SPE Annual Technical Conference and竹寺组与龙马溪组页岩气勘探优选区预测［J］．Exhibition．San Antonio：［s．n．］，2012．地质通报，2011，31（2）：357－363．

［16］BRIESMEISTER J F．MCNP—A general Monte

Numerical Simulation on Logging Method of Direct Determination of Total Organic Carbon Content in Formation

YUAN Chao1，ZHOU Can－can1，BO Ya－jie2，ZHANG Feng3，LI Xia1，CHENG Xiang－zhi1
（1．Research Department of Well－logging and Remote Sensing Technology，Research Institute of Petroleum Exploration and Development，PetroChina，Beijing100083，China；2．Technical Center，CNPC Logging Company，Xi’an710077，China；3．School of Geosciences，China University of Petroleum（East China），Qingdao 266580，China）

Abstract：Unconventional resources gradually become the important domain of oil and gas exploration in China，in which the evaluation of source rock，especially the total organic carbon content（TOC），is critical due to the“source and reservoir in one body”．Conventional evaluation methods of TOCbased on logging data are all indirect ones，and calculated results involve some uncertainties as they are affected by various influential究

factors．In this paper，inelastic gamma spectra and capture gamma spectra were recorded with the gamma ray detector，and elemental contents in formation were derived with spectral analysis technique．After that，the total inorganic carbon content（TIC），which is calculated with the TIC－related elemental contents，was deducted from total carbon content（TC）to directly determine TOC．According to the characteristics of shale in Longmaxi and Qiongzhusi formation of Sichuan basin，Monte Carlo method was employed to build 6formation models，gamma spectra under different formation conditions were simulated，and then the simulated gamma spectra were analyzed with elemental standard gamma spectra．The simulation results demonstrate that the method mentioned in this paper can directly provide high－accuracy results of TOC，and the relative error of calculated value is less than 6%when TOCis more than 2%．

Key words：direct determination of total organic carbon content；pulsed neutron logging；gamma spectrum；numerical simulation

作者简介：袁 超（1987—），男，山东临朐人，博士研究生，从事核测井方法、核测井数据处理及蒙特卡罗模拟研

基金项目：国家重大油气专项资助项目（2011ZX0520－008）；国家自然科学基金资助项目（41374125）

收稿日期：2013－12－25；修回日期：2014－11－21

doi：10．7538／yzk．2015．49．08．1497

文章编号：1000－6931（2015）08－1497－07

文献标志码：A

中图分类号：P631．817