自定義響應方程的測井最優化處理方法

馮周, 李心童, 武宏亮, 王華鋒, 馮慶付, 王克文

(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.