自定义响应方程的测井最优化处理方法

冯周, 李心童, 武宏亮, 王华锋, 冯庆付, 王克文

(1.中国石油勘探开发研究院, 北京 100083; 2.中国石油长城钻探工程有限公司, 北京 100101; 3.北京航空航天大学, 北京 100191)

0 引 言

测井数据最优化处理方法最早由斯伦贝谢公司在地层矿物组分定量计算中提出[1-2],它首次引入广义地球物理反演的思想,将全部测井信息、误差及地区地质经验综合成一个多维信息复合体,借鉴数学中的最优化技术求取满足地层条件的最优解释结果。该方法最大的优点是可以综合利用各种测井曲线有针对性地建立符合地层实际情况的多矿物组分模型进行求解,提高计算结果的准确性,是目前常规测井计算地层岩石矿物组分通用的技术手段之一[3-7]。

最优化处理的基础是依据地层模型建立的各类测井响应方程,它表征了测井响应与地层特性参数之间的定量关系,响应方程的精度决定了最优化处理结果与实际地层的符合程度。在现有的最优化处理方法中,测井响应方程普遍只是采用预先定义的几类固定的形式,如声波测井方程采用Wyllie公式、Hunt-Raymer公式,电阻率方程包括Archie公式、Simandoux公式、Waxman-Smits公式等[8-10]。在实际应用中仍存在以下问题:首先,已有的测井响应方程都是基于传统经典模型建立的,因而其适用范围有限,在很多复杂类型储层中应用效果不理想;其次,预定的模型公式无法在所有油田区块、层系都适用,不同研究地区有其特定的形式;最后,在测井解释中通常还包括了很多根据岩石物理分析资料建立的经验公式,现有的最优化处理方法中欠缺对这些公式的处理能力,因而限制了方法的应用范围和处理效果。

对此,本文提出将计算机表达式解析技术[11-12]与测井最优化思想相结合的处理方法,利用计算机对自定义响应方程自动解析后再进行优化求解,实现对响应方程的体外添加处理,有效提升方法的计算精度和适用范围,实际井处理与岩心分析结果对比验证了方法的可靠性。

1 测井响应方程的构成

测井响应方程表征了测井响应和地层特性参数之间关系,通常由理论模型推导或实际资料统计获得。由于各测井方法原理的区别以及方程推导所采用地层模型的差异性,不同响应方程表达式的形式变化较大,复杂程度也各不相同。以Archie公式为例,转化为电导率响应方程形式为

(1)

式中,φ=Vw+Vo+Vg+Vb,Sw=(Vw+Vb)/φ;Ct为地层电导率值;Cw为地层水电导率;φ、Sw分别为地层孔隙度和含水饱和度;Vw、Vo、Vg、Vb分别表示地层水、油、气以及束缚水的体积含量;a、b、m、n是方程参数。

本文提出构成响应方程表达式的基本元素类型为运算符、数字、常量符以及变量符4种,其中变量符又可以分为基本变量和复合变量,基本变量是表示地层组分的基本单元,如式(1)中的Vw、Vo等;复合变量是由其他运算符、数字、常量符或变量符构成的表达式,如式(1)中的φ、Sw。利用上述4种元素符号,可实现对任意形式的测井响应方程进行表述。以用于复杂缝洞储层评价的李宁公式[13-14]为例,用于计算机解析的响应方程字符串形式可记为

CONDT=CW×FAIM/[A×(P1/SWN1+P2/SWN2)]

(2)

式中,CW、A、M、P1、N1、P2、N2是常量符号;FAI、SW是复合变量,其计算表达式可记为

FAI=UWAT+UGAS+UOIL+BWAT

(3)

SW=(UWAT+BWAT)/FAI

(4)

式中,UWAT、UGAS、UOIL、BWAT表示地层水、气、油、束缚水组分的体积含量,为基本变量(也就是最优化处理的待解变量)。

2 表达式解析与计算方法

通过运算符、数字、常量符及变量符4种元素可对测井响应方程形式进行表示,当表达式的定义完成后,可认为表达式自身的语法结构就不会再进行变更了,在表达式值计算过程中变化的部分就是常量与变量的数值大小。因此,对表达式后续处理的基础和核心是对其进行解析的过程,其重点在于对表达式元素(数值、括号、运算符和变量组合)的分离和抽象,通过解析将表达式进行运算层次的分析和存储,具体过程包括3个阶段:①表达式元素分析;②表达式分离和存储;③表达式偏导式处理及其值的计算。

阶段①中,对解释人员定义形成的响应方程表达式字符串中包含的常量符、数字、运算符、变量符元素进行分类匹配,检查表达式定义是否完整,运算符、括号等是否完备、合法,同时将表达式中复合变量元素逐次替换,直至表达式中只包含运算符、数字、常量符以及基本变量。

阶段②中,根据表达式元素匹配识别结果,对响应方程表达式进行语法分析,具体包括表达式中括号的作用、操作符的优先级、各个操作符出现的次序以及排除多义性等,转义为适用于计算机运算的后缀表达式形式(后缀表达式是将操作符放置在操作数的最后,可以避免不同运算符、括号位置等引起的运算优先级的问题,简化计算顺序)。同时,为了方便后续求导计算,针对测井响应方程表达式元素特点设计了1种新的数据结构体取代传统的栈式结构,实现了表达式解析和操作的完全分离。数据结构体包括3个字段,第1个字段为字符串型(string),用于存放表达式元素形式;第2个字段为双精度浮点型(double),用于存放该元素的值,数字元素的值即为其本身,变量符、常量符的值在后续计算中代入;第3个字段为整型(int),用于标识运算符和操作数的类型,常量为0,变量为2,定义的运算符为1。转义后的后缀表达式中每个元素采用该数据结构体记录并放入动态数组顺序存储。具体存储结构见图1。

图1 解析后的表达式存储示意图

阶段③中,根据式中包含的运算操作符类型,设置了求导字典表,方便查找对应的导数形式,同时方便了后续对运算操作符的扩充和修改。部分求导字典表定义为

……

#define d_sin(u) d(u)+″*cos(″+u+″)″ ∥du*cos(u)

#define d_cos(u) d(u)+″*-sin(″+u+″)″∥du*(-sin(u))

#define d_sqrt(u,stack) d(u,stack)+″/(2*sqrt(″+u+″))″

∥du/(2*sqrt(u))

#define d_exp(u,stack) d(u,stack)+″*exp(″+u+″)″

∥du*exp(u)

……

在此基础上,根据后缀表达式存储结果,遍历每个存储单元,判断每个元素与基本变量的关系,同时根据运算符类型,即可按照预定求导规则字典表确定元素间组合关系,直至全部替换完毕后完成表达式对各基本变量的偏导式求取。

在求导部分完成后,将对各偏导式进行重新解析,采用阶段②中所述存储结构进行记录。后续最优化计算中求取各偏导式值时,只需要将各常量符、变量符在对应深度点的值代入,依次遍历各存储单元,依据后缀表达式计算规则直接计算得到相应结果。

3 测井最优化问题求解方法

根据各测井曲线响应方程,采用最小二乘法建立最优化目标函数求解[15-17]。目标函数可用式(5)与式(6)表示

v*=arg min{F(v)}

(5)

(6)

式中,tci为地层理论测井响应值,根据响应方程计算获得;tmi为实际测量响应值;fi(v)为根据地层解释模型建立的测井响应方程;wi是不同测井方法标准化处理因子。

同时,上述建立的响应关系中地层组分还应满足一定的岩石物理以及地质的约束条件,约束条件方程组一般形式可记为

hk(v)=Ckj·vj-bk≤0

(7)

式中,Ckj代表约束条件系数矩阵;bk代表约束条件边界;vj是各地层组分含量。

由于用户定义的响应方程式通常都是关于地层流体组分的复杂非线性函数,因此目标函数式(6)是一个典型的带约束条件的非线性最小二乘问题。针对该类问题利用惩罚函数结合Levenberg-Marquardt[1819]的算法计算量较小、求解效率较高。针对目标函数的形式,惩罚项可设为约束条件的平方和形式,转化为无约束问题目标函数可表示为

(8)

对式(8)无约束非线性最小二乘问题,利用Levenberg-Marquardt算法求解的迭代增量计算式可表示为

(JTJ+μ·I)·h=-JTR

(9)

由于式(8)中,tmi、bk均为与地层矿物、流体组分含量无关的常量(对每个处理深度点而言),则

(10)

4 处理软件实现

综合上述表达式计算机解析技术和非线性最优化求解算法,建立了1套完整的测井最优化处理方法,包括地层解释模型建立、测井响应方程编辑、表达式解析与计算、最优化求解等步骤,其具体实现流程见图2。该流程与传统测井最优化处理最大区别在于通过响应方程编辑和表达式自动解析步骤,实现了对用户自定义形式响应方程的添加和处理。在此基础上,基于测井软件CIFLog平台[20-21]编制实现了完整的处理软件。解释人员可根据实际需要,通过在界面上点选表达式元素来自定义更符合地层条件的测井响应方程,从而进一步提高岩石矿物组分计算精度。

图2 基于计算机表达式解析的测井最优化处理流程

5 应用效果

利用上述处理软件,对实际井资料进行了处理,并与传统最优化方法计算结果以及岩心定量分析结果进行了对比。

龙×井储层位于6 052.60～6 077.40 m段,岩性为灰岩、灰质云岩,取心岩心描述和物性分析表明,该段储层溶蚀孔洞发育。为了进一步确定储层的岩石物理特性,该段储层全直径岩心还进行了驱替实验,表1中给出了该井段全直径岩心Sw—I实验测量结果及分别利用Archie公式与李宁公式拟合结果(见图3)。2种方程拟合曲线对比如图3右所示,其中蓝色实线是Archie公式拟合结果,红色实线是利用李宁公式拟合结果。由图3和拟合结果可见,2种方程均能较好反映岩心Sw—I的变化特征(相关系数均达到0.98以上),但李宁公式拟合曲线更准确,特别是在低含水饱和度时,Archie公式拟合结果与实测数据相差较大,这将直接影响油气饱和度解释的精度(见表1)。

表1 龙×井×号岩心Sw—I实验Archie公式与李宁公式拟合结果

图3 龙×井×号岩心照片及Sw—I实验Archie公式与李宁公式拟合曲线对比

基于上述Sw—I实验拟合确定的饱和度方程参数,分别使用Archie公式和李宁公式作为电阻率测井响应方程该井进行了最优化处理,其中Archie公式在传统最优化处理方法中可直接选用,李宁公式作为自定义形式方程输入到本文方法软件中处理,处理时电阻率响应方程参数采用表1中拟合结果,其他地层响应参数均保持一致。

图4是利用2种方法处理得到地层矿物、流体含量与岩心分析结果对比图。图4中第5至6道为Archie公式优化处理得到地层白云石含量、含油气饱和度与岩心分析结果对比;第7至8道为自定义的李宁公式优化处理得到地层白云石含量、含油气饱和度与岩心分析结果对比;第9道为自定义的李宁公式优化处理得到地层岩性剖面。图4各道中,横线+圆点符号表示相应的岩心分析结果。从对比结果可以看到,即使2种方法在输入电阻率响应方程精度非常接近的情况下,采用自定义李宁公式优化处理得到的地层白云石含量、含气饱和度与岩心分析结果一致性更好,Archie公式处理得到的白云石含量、含气饱和度值偏高。图5给出了2种方法计算结果与岩心分析结果的误差分布情况。显然,采用自定义李宁公式优化处理精度更高,矿物含量计算误差在5%左右,饱和度计算误差一般小于10%。

图4 龙×井地层矿物及流体组分含量计算结果与实验分析结果对比图*非法定计量单位,1 ft=12 in=0.304 8 m,下同

图5 2种方法处理得到的白云石含量及含油气饱和度与岩心分析结果误差分布

6 结 论

(1) 针对测井响应方程组成特点,建立了1套高效准确的用户自定义响应方程表达式解析与求导、求值计算方法,为目标函数优化求解奠定了基础。

(2) 建立了1套体外自定义响应方程的测井最优化处理方法并实现了相应软件,实际应用中,解释人员可根据研究区地层特征,采用更具针对性的解释模型和响应方程进行最优化处理。

(3) 实际井资料处理结果与岩心定量分析结果对比表明,新方法能够有效提升地层组分含量计算精度,从而验证了方法的可靠性。

参考文献:

[1] MAYER C, SIBBT A. GLOBAL. A New Approach to Computer-Processed Log Interpretation [C]∥SPE 9341, 1980.

[2] QUIREIN J, KIMMINAU S, LAVIGNE J, et al. A Coherent Framework for Developing and Applying Multiple Formation Evaluation Models [C]∥SPWLA 27th Annual Logging Symposium, 1986.

[3] QUIREIN J, WITKOWSKY J, TRUAX J A, et al. Integrating Core Data and Wireline Geochemical Data for Formation Evaluation and Characterization of Shale-Gas Reservoirs [C]∥SPE-134559-MS, 2010.

[4] QUIREIN J, GALFOHD J, WITKOWSKY J. Formation Evaluation and Characterization of Shale-gas Reservoirs by Means of Core and Wireline Data Integration [J]. The Leading Edge, 2013, 32(12): 1486-1492.

[5] 冯国庆, 陈军, 张烈辉, 等. 最优化测井解释的遗传算法实现 [J]. 天然气工业, 2002, 22(6): 48-51.

[6] 梁艳, 王燕, 李素兰, 等. 川东北海相碳酸盐岩ElanPlus测井解释模型及应用 [J]. 天然气技术与经济, 2011, 5(4): 39-41.

[7] 侯颉, 邹长春, 杨玉卿. 页岩气储层矿物组分测井分析方法 [J]. 工程地球物理学报, 2012, 9(5): 607-613.

[8] 高楚桥, 张超谟, 钟兴水. ELAN解释(程序)方法简析 [J]. 测井技术, 1995, 19(2): 135-139.

[9] 金勇, 张世刚, 顾列. FORWARD测井解释平台中使用的先进技术 [J]. 测井技术, 2000, 24(1): 64-78.

[10] PEETERS M, VISSER R. A Comparison of Petrophysical Evaluation Packages: LOGIC, FLAME, ELAN, OPTIMA and ULTRA [J]. The Log Analyst, 1991, 32(04): 350-357.

[11] 何云东, 黄昶. 复杂表达式解析和计算的研究实现 [J]. 中国科技信息, 2009, 8.

[12] 李世华, 刘晓娟, 姜晨, 等. 关于表达式求值的算法研究与实现 [J]. 甘肃科技, 2011, 27(1): 11-15.

[13] 李宁. 电阻率—孔隙度、电阻率—含油(气)饱和度关系的一般形式及其最佳逼近函数类型的确定(Ⅰ) [J]. 地球物理学报, 1989, 32(5): 580-591.

[14] 孙建孟, 王克文, 李伟. 测井饱和度解释模型发展及分析 [J]. 石油勘探与开发, 2008, 35(1): 101-107.

[15] 肖立志, 钟兴水. GLOBAL测井解释方法Incoherence函数性质研究 [J]. 石油学报, 1990, 11(2): 49-57.

[16] 雍世和, 孙建孟. 测井数字处理中最优化方法的选择 [J]. 石油大学学报(自然科学版), 1988, 12(4/5): 11-26.

[17] 高楚桥. 复杂储层测井评价方法 [M]. 北京: 石油工业出版社, 2003, 1-16.

[18] PHILIP E GILL, WALTER MURRAY, MARGARET H WRIGHT. Practical Optimization [M]. ACADEMIC PRESS, Harcourt Brace and Company, Publishers, 1981.

[19] MADSEN K, NIELSEN H B, TINGLEFF O. Methods for Non-linear Least Squares Problems [M]. Kgs. Lyngby, Copenhagen: Technical University of Denmark, 2004, 24-29.

[20] 李宁, 王才志, 刘英明, 等. 一体化网络测井处理解释软件平台CIFLog [J]. 石油科技论坛, 2013, 3: 6-10.

[21] 李宁, 王才志, 刘英明, 等. 基于Java-NetBeans的第三代测井软件CIFLog [J]. 石油学报, 2013, 34(1): 192-200.