呂玉增,,彭蘇萍
(1.中國礦業大學 煤炭資源與安全開采國家重點實驗室,北京 100083;2.桂林理工大學 地球科學學院,桂林 541004)
地~井方位激電(IP)人機交互解釋軟件開發及應用
呂玉增1,2,,彭蘇萍1
(1.中國礦業大學 煤炭資源與安全開采國家重點實驗室,北京 100083;2.桂林理工大學 地球科學學院,桂林 541004)
基于MSR(Modified Sparse Row)非零元素存儲和SSOR-PCG方程求解技術,編制了三維IP的有限元法快速正演程序,在此基礎上,利用Delphi等開發了可視化地~井方位IP人機交互解釋軟件。在地~井方位IP的數據解釋中,遵循“先井口,后其他方位”的分析方法,先通過井口方位觀測曲線初步斷定異常體的電性特征、大致埋深和傾向;然后對比分析其它方位曲線的異常特征,推斷出異常的賦存方位。對于異常復雜情況,該軟件可利用測井數據、鉆孔柱狀圖等,形成分層背景,通過差值法和比值法數據預處理技術以突出井旁異常,提高異常分辨能力。實測數據的解釋結果表明,針對地~井方位IP觀測數據量少,充分利用測井、鉆孔等信息,用人機交互正演擬合反演方式解釋是可行的。
地~井方位觀測;有限單元法;IP;人機交互解釋;軟件
地~井方位激電觀測是井中激電的工作方式之一,指在地面供電,井中測量的井中激發極化法,主要用來發現井旁盲礦,確定其埋藏深度及方位,追蹤和圈定礦化帶等。由于地~井方位觀測是在井中觀測,能更好地接近深部的目標異常體,因此,在深部隱伏礦產的勘查以及老礦山的“探邊摸底”等方面有著廣泛的應用。
近年來,國內、外對井地、井井電法的研究非常重視,取得了豐富的研究成果[1~13],但對地~井電法的研究成果偏少[14、15]。一方面由于地~井方位觀測野外施工難度大;另一方面地~井方位觀測正反演所用的模型單元數多,而可利用的實測數據少,正演模擬和數據解釋難度相對較大[16]。針對地~井方位IP的實際情況,作者提出人機交互解釋思路,并開發相應軟件[17~19]。
地~井方位IP的人機交互解釋需建立在快速、準確的三維正演基礎之上。有限單元法是研究復雜三維地電模型正演的最有效數值模擬方法之一,但它需要在區域內進行網格剖分。對于地~井方位觀測而言,由于沿井方向的網格剖分數量大導致整個網格剖分數量巨大,數值計算上存在存儲、計算速度和精度等多個矛盾。
下面介紹解決上述矛盾,實現三維地-井方位IP有限元快速正演的技術。
有限元法求解地~井方位IP的正演問題,通過網格剖分、單元分析、矩陣合成等步驟,最終形成系線性方程組:

其中 K是系數矩陣;u是全部節點的電位ui組成的列向量;p是與點電源項有關。
有限元法正演形成的系數矩陣K是對稱、稀疏矩陣,其中絕大部份元素值為零,但每行零元素的分布和個數不同,與網格剖分和節點編號有關例如對于直角坐標系中網格節點(x×y×z)=(5 ×5×5)的三維四面體剖分,網格節點個數125(5 ×5×5),系數矩陣K的大小125×125=15 625,常規變帶寬存儲這個下三角陣需要3 225個元素,而非零元素個數僅有810個,可見非零元素存儲能明顯減少內存存儲。
MSR(Modified Sparse Row)是一種改進的行壓縮索引存儲技術,是以行為單位,順序存儲系數矩陣中的非零元素。考慮到在方程組的求解過程中,系數矩陣對角線元素均不為零,且對角線元素的訪問和運算最頻繁,把對角線元素提出來單獨存儲(見圖1)。該存儲方式僅需要兩個數組:一個實型數組AA和一個整型數據ICOL,具體數組說明見表1。
SSOR-PCG迭代法求解方程組(1)的過程為[19]:對稱系數矩陣K可寫為對角陣D和嚴格下三角陣E的和形式:K=D-E-ET,而M為K的 SSOR預條件矩陣,則:

其中 ω∈[0,2]為松弛因子,通過試算確定SSOR-PCG求解的整個迭代過程如下:

由于預條件矩陣M與系數矩陣K具有相同的稀疏性,同樣也采用MSR壓縮存儲,整個迭代過程只是系數矩陣非零元素的運算,計算量小,計算速度快。

圖1 四面體網格剖分和系數矩陣下三角非零元素分布Fig.1 Tetrahedral mesh and lower triangular distribution of nonzero elements in the coefficient matrix

表1 MSR存儲說明Tab.1 MSR storage instruction
該軟件主要有兩大功能:①地~井方位IP三維有限元正演模擬;②地~井方位IP數據的人機交互反演。軟件是基于三維有限元正演基礎上開發的,適合于地~井多方位觀測方式下,任意的三維地電體的直流激電正演模擬和正演擬合交互反演解釋。三維有限元正演軟件包含測井基本信息輸入、三維異常體模型構建、有限元正演、顯示結果曲線、保存結果數據等功能,軟件執行流程見圖2。野外數據的反演擬合步驟大致分為實測數據輸入、圖示數據、測井數據分層、建立層狀模型、三維模型構建、曲線對比擬合等步驟,詳見圖3流程圖。
地~井方位IP人機交互解釋軟件是在Delphi 7.0和Fortran 4.0開發平臺下開發完成的可視化軟件,軟件設置了文件、數據輸入、數據切換、模型構建等菜單項,各菜單項功能和對應的快捷方式見圖4和下頁表2。
建立測井工程文件是地~井方位IP正演和交互解釋的第一步,用戶根據地~井方位的模型和實際工作情況,通過“工作參數輸入”對話框輸入(編輯)方位測井的工作參數(見后面圖5),包括地面方位點的供電位置、井孔深度、測點距等信息。當完成工作參數輸入后,軟件系統會根據工作參數信息,繪出地~井方位測井示意圖。這時用戶可以進行三維建模等操作,還可以通過功能項保存當前的模型到工程模型文件(.GP)中,便于下一次的直接打開和編輯。


表2 菜單選項、工具條及其功能介紹Tab.2 Menu options,toolbars and features
測井工作參數輸入完成后,系統自動對三維空間進行水平分層。三維建模的思路是:將每個深度層看作是一定厚度的板狀體,通過對各不同的深度層建模,建立復雜的三維模型。具體編輯模型操作步驟如下:
(1)選擇要編輯的深度層。
(2)用不同的顏色區分不同的異常體,通過鼠的顏色示意在主界面右側深度層選擇區域上。

數據預處理的目的是為下一步的分析、解釋提供幫助。根據方位測井的實際經驗,軟件提供了比值法、差值法二種預處理方法,以突出旁側異常。比值法(差值法)預處理就是把所有方位觀測數據統一除以(減去)井孔方位觀測數據值,以壓制孔壁局部影響,突出旁側異常信息。
根據預處理后的曲線結果,結合原始數據曲線特征和鉆探等地質信息,遵循“先井口、后其他方位”的分析方法,推斷異常體的大致深度、方位和異常體電性特征。并初步建立預測初始模型,進行正演試算。將模型正演結果與實測數據對比,不斷修正模型并再次正演,直到結果滿足要求為止。
廣西某礦山鉆孔深度493m,孔徑約100mm,采用了地~井方位激電觀測,目的是了解已知低阻高極化礦帶在該鉆孔附近的位置和隱伏情況??睖y單位進行了鉆孔四周四個方位探測(A1,A2,A3和A4),各方位距離相等,井中觀測電極距MN=10m,相鄰測點距10m,井中實際觀測位置從深度為100m~480m。各方位點相對坐標、高程和方位距離見表3,視電阻率、視極化率等實測數據保存在數據文件中。為了便于對比分析,我們用鉆孔四周四個方位探測(A1、A2、A3和A4)數據的平均值當作是井口方位(A0)的探測結果。

表3 各方位點相對坐標、高程和方位距離Tab.3 Relative coordinates and azimuth distances of all sites
人機交互正演擬合反演解釋過程如下。
(1)工作參數和實測數據導入。進入“工作參數輸入”對話框,導入供電點坐標、輸入井深、直徑和觀測參數(見圖8(a))。當基本參數輸入完成后,導入視電阻率、視極化率等實測數據到主界面左側相應數據表格中,并繪制出四個方位的觀測曲線(見圖8(b))。

圖8 工作參數與實測數據輸入和觀測曲線顯示Fig.8 Working parameters(a),the measured data input an result curves show(b)
(2)根據實測曲線異常特征建立初始模型。從圖8方位測井視電阻率和視極化率曲線上,有二處較明顯的異常段:①深度300m-350m區段,四個方位都有明顯低阻高極化異常,且視電阻率和視極化率異常對應較好,推斷礦體大致埋深320m另外,四個方位激電異常不對稱,對比發現A1和A2方位視電阻率和視極化率異常相近,而A3和A4方位觀測的視電阻率和視極化率異常相近且幅值較大,初步斷定礦化體靠近A3、A4方位,結合地質資料,可初步推斷該異常體為近似水平板狀體②在鉆孔的底部450m附近,低阻高極化異常,初步推斷引起該激電異常的是位于鉆孔底部附近另一個礦化體。
通過上述分析,建立初始模型,背景電阻率和極化率近似取觀測數據平均值,分別為410Ω·m和8%;依據該地區的礦石標本物性參數,異常體電阻率和極化率取值分別為50Ω·m和20%。
(3)模型修改和交互擬合解釋。在初始模型基礎上,利用“共軛梯度法”進行模型快速正演計算并反復修改模型、正演計算,對比模型正演與實測曲線的擬合情況,正演擬合以曲線形態擬合為主圖9為修改后最終模型有限元數值解與實測數據曲線的對比情況。

圖9 模型正演(紅)與實測數據(藍)擬合對比Fig.9 Fitting comparison forward modeling(red)with the measured data(blue)
(4)解釋結果。最終的擬合結果顯示,在該鉆孔330m~360m深度上,礦帶比較靠近鉆孔,且在A3、A4方位方向有一定延伸,延伸長度約40m;另外,推斷在鉆孔底部下方深度450m以下,賦存另一高阻高極化礦化體。
根據勘測單位鉆探情況,該解釋結果基本證實了地質、物探人員的判斷,即該鉆孔已經穿過淺部的礦脈,但還沒有揭露深部的礦化帶。同時,礦體埋深、賦存方位等推斷解釋結果也為勘測單位提供了幫助。

圖10 正演擬合交互解釋結果Fig.10 The result of forward interactive interpretation
(1)利用MSR非零元素存儲和SSOR-PCG方程求解技術實現了三維IP的快速正演,為實現地-井方位激電(IP)人機交互解釋奠定了基礎。
(2)利用Delphi等軟件平臺開發了可視化地~井方位IP人機交互解釋軟件,在適合于地~井多方位觀測方式下,任意的三維地電體的直流激電正演模擬和正演擬合交互反演解釋。
(3)在地~井五方位IP實測數據的推斷解釋中,可遵循“先井口,后其他方位”的分析方法,先通過井口方位觀測曲線初步斷定異常體的電性特征、大致埋深和傾向;然后對比分析其它方位曲線異常特征,推斷異常的賦存方位。對于異常復雜情況,可利用測井數據、鉆孔柱狀圖等,形成分層背景,通過差值法和比值法對比分析不同方位觀測曲線特征,以突出井旁異常,提高異常分辨能力。
(4)實測數據的解釋結果表明,針對地~井方位IP觀測數據量少,充分利用測井、鉆孔等信息,用人機交互正演擬合反演方式解釋是可行的。
致謝:
感謝上海應用技術學院黃俊革教授給予的幫助。
[1] ZHDANOV M S,YOSHIOKA K.Cross-well electromagnetic imaging in three dimensions[J].Exploration Geophysics,2003(34):34.
[2] ZHOU B,GREENHALGH S A.Rapid 2-D/3-D crosshole resistivity imaging using the analytic sensitivity function.Geophysics[J],2002,67(2):755.
[3] WILKINSON P B,CHAMBERS J E.Extreme sensitivity of crosshole electrical resistivity tomography measurements to geometric errors[J].Geophys.J Int.,2008,173:49.
[4] PARDO D,CALO V M,TORRES-VERDIN C,e al.Fourier series expansion in a non-orthogona system of coordinates for the simulation of 3D-DC borehole resistivity measurements[J].Comput Methods Appl.Mech.Engrg.2008(197):1906.
[5]NAM M J,PARDO D,TOORES-VERDIN C Simulation of DC dual-laterolog measurements in complex formations:A Fourier-series approach with nonorthogonal coordinates and self-adapting finit elements[J].Geophysics,2009,74(1):E31.
[6] 安然,李桐林,徐凱軍.井地三維電阻率反演研究[J]地球物理學進展,2007,22(1):247.
[7] 王志剛,何展翔,劉昱.井地直流電法三維數值模擬及異常規律研究[J].工程地球物理學報,2006,3(2)87.
[8] 柯敢攀,黃清華.井地電法的三維正反演研究[J].北京大學學報:自然科學版,2009,45(2):294.
[10]底青云,王妙月.積分法三維電阻率成像[J].地球物理學報,2001,44(6):844.
[11]沈平,強建科,李永軍,等.井間視電阻率的幾何成像方法[J].中南大學學報:自然科學版,2010,41(3)1079.
[12]潘紀順,葛為中,折京平.地面/井地/井間超高密度電阻率成像技術[J].華北水利水電學院學報,2010,31(2):74.
[14]郭文波,宋建平,李貅,等.層狀介質井中電測數值計算及其應用研究[J].地球物理學報,2006,49(5)1561.
[15]蔡柏林,黃智輝,谷守民.井中激發極化法[M].北京地質出版社,1983.
[18]呂玉增.地-井、井-地IP三維快速正反演研究[D].長沙:中南大學,2008.
book=7,ebook=7
1001—1749(2012)03—0358—07
P 631.3+24
A
10.3969/j.issn.1001-1749.2012.03.22
呂玉增(1978-),男,河北邢臺人,博士后副教授,主要從事電法數值模擬與反演成像研究。
國家自然科學基金(40774057);全國危機礦山接替找礦項目新技術新方法示范項目(200799086);廣西自然科學基金(桂科自0832263)
2011-07-27改回日期:2011-10-31