胡瑞華，林 君，孫彩堂, 劉長勝, 周逢道

(吉林大學(xué) 儀器科學(xué)與電氣工程學(xué)院, 長春 130061)

0 引言

當前電磁法勘探已成為解決地質(zhì)問題的重要手段[1-2]，學(xué)者們通過正反演理論和技術(shù)來研究地質(zhì)模型的電磁響應(yīng)和地下結(jié)構(gòu)的復(fù)現(xiàn)，考慮到探測系統(tǒng)的發(fā)射功率和分辯能力，野外作業(yè)時常采用人工源激發(fā)達到獲取高信噪比的電磁數(shù)據(jù)的要求。磁性源和電性源是常被采用的兩種人工源，其中電性源使用一定長度的導(dǎo)線向地下供一定強度的交變電流來實現(xiàn)。在以電性源為場源的電磁探測方法中，其正演常需要研究基于電偶極子激發(fā)的場響應(yīng)，并根據(jù)響應(yīng)規(guī)律來解釋野外實測數(shù)據(jù)(即反演)[3～12]。這種正反演通常都是假定發(fā)射端與接收點距離遠大于供電兩端的長度，即要求導(dǎo)線長度滿足偶極子的條件。文獻[4]指出，當收發(fā)距大于等于10倍的導(dǎo)線長度時，導(dǎo)線可被看作是電偶極子。然而在實際應(yīng)用中由于受地形地貌、發(fā)射功率等多種因素的制約，導(dǎo)致長導(dǎo)線并不滿足電偶極子的條件，而使基于電偶極子的場計算公式不再適用，在這種情況下若仍然按照電偶極子激勵理論開展數(shù)據(jù)解釋，將造成較大偏差甚至是錯誤的地質(zhì)結(jié)構(gòu)的復(fù)現(xiàn)。為了解決這個問題，在正反演中可以先將長導(dǎo)線源切分成若干個電偶極子，然后對切分的電偶極子的響應(yīng)進行求和，以達到計算長導(dǎo)線源場響應(yīng)的目的。

針對如何將長導(dǎo)線源切分成多個電偶極子的方法，本次研究并實現(xiàn)了兩種切分算法，①借鑒二分法思想設(shè)計的遞歸切分算法；②是借鑒窮舉法思想設(shè)計的窮舉切分算法。通過對兩種切分算法的研究表明：遞歸切分算法切分出的電偶極子數(shù)目相對較多；窮舉切分算法在窮舉步長較小(如0.1 m)時能較準確地切分出真實電偶極子的位置。將兩種算法應(yīng)用有于CSAMT電場的計算中，計算結(jié)果表明：兩種切分算法都可應(yīng)用于實際場的計算。

1 切分算法原理

1.1 觀測點與長導(dǎo)線源的拓撲結(jié)構(gòu)

人工源野外測量時，觀測點離長導(dǎo)線源有一定的距離，接收和發(fā)射之間的是三維拓撲結(jié)構(gòu)[9,11]，常見的拓撲結(jié)構(gòu)有平地形的地面發(fā)射地面接收、起伏地形的地面發(fā)射地面接收、地面發(fā)射井中接收等(圖1)。

圖1 觀測點與長導(dǎo)線源的拓撲結(jié)構(gòu)Fig．1 Topology structure of survey point and long lead source

圖1中長方體代表大地，AB表示長度為AB的長導(dǎo)線發(fā)射源，A和B分別代表發(fā)射源的兩極，P1、P2、P3分別表示與AB共面、共線和地下的觀測點。P1、P2、P3與AB的三種位置關(guān)系與野外測量時的施工布局一致。偶極子的位置常用發(fā)射源兩極(導(dǎo)線兩端)的中點代表。若空間坐標系中A極、B極、測點P三點坐標分別取為A(xA，yA，zA)、B(xB，yB，zB)、P(xP，yP，zP)，其中x、y、z表示點的坐標分量，若用C表示AB的中點，則接收點與發(fā)射源的距離(發(fā)收距)表示為P點到C點的距離PC(常用r字符表示)。根據(jù)A、B坐標及代數(shù)關(guān)系很容易計算出C點的坐標、發(fā)收距PC和發(fā)射極距AB。另外，偶極子的條件是收發(fā)距大于等于N(通常最小值可取為10)倍的導(dǎo)線長度，即PC≥N*AB。實際應(yīng)用中A、B、P三點坐標已知，而N參數(shù)的確定可在綜合考慮工區(qū)、收發(fā)系統(tǒng)參數(shù)以及電磁數(shù)據(jù)質(zhì)量的情況下選定。

1.2 二分法遞歸切分算法

二分法算法是一種快速搜尋算法，算法的核心是從有序數(shù)據(jù)列中按照折半查找的方式尋找滿足條件的結(jié)果[12-14]。借鑒二分法思想研究切分長導(dǎo)線源為多個電偶極子的算法如下：

1) 根據(jù)長導(dǎo)線兩端A點、B點的坐標計算其中點C的坐標和A點到B點的距離AB；根據(jù)P點、C點的坐標計算P點到C點的距離PC，轉(zhuǎn)到步驟2)。

2) 如果PC≥N*AB則記錄下C點的位置，算法結(jié)束；否則轉(zhuǎn)到步驟3)。

3)用C點坐標替換B點的坐標，轉(zhuǎn)到步驟1)。

4)用C點坐標替換A點的坐標，轉(zhuǎn)到步驟1)。

此算法是遞歸算法，共分四步。第一步計算發(fā)射距AB和收發(fā)距PC，為算法的下一步作準備；第二步判斷AB是否為電偶極子，若是電偶極子，則AB的中點就是電偶極子的位置，算法結(jié)束，它是遞歸算法的出口；第三步當AB不是電偶極子時，搜尋CA段導(dǎo)線電偶極子的位置；第四步當搜尋完CA段導(dǎo)線的電偶極子時，再搜尋CB段導(dǎo)線的電偶極子。每遞歸一次就能找到一個電偶極子，當整個遞歸算法結(jié)束時，便找到了所有電偶極子的位置點(圖2)。

圖2 二分法遞歸切分算法流程圖Fig．2 Flow of binary segmentation method

1.3 窮舉切分算法

窮舉算法是一種試算遍歷算法，算法的核心是從取值區(qū)間對每個取到的值(或組合值)進行試算，通過判斷試算結(jié)果是否滿足要求，從而找到問題的解[12,15]，這種算法常被用在方程少而未知數(shù)多的問題中。將長導(dǎo)線源切分為多個電偶極子的問題，本質(zhì)上是查找AB上哪些點是偶極子的位置，但我們僅知道偶極子的條件PCnew≥N*AnewBnew(式中AnewBnew表示AB上的一段導(dǎo)線，Cnew為其中點)，這樣必須通過1個方程(PCnew=N*AnewBnew)求解AnewBnew的中點Cnew的位置(即要確定Cnew的3個坐標分量)，顯然方程少而未知數(shù)多，無法直接求解出方程的解，即無法定位偶極子的位置，因此這個問題屬于窮舉算法能解決的范疇。借鑒窮舉算法的思想，設(shè)計將長導(dǎo)線源切分為多個電偶極子的算法如下：

1)根據(jù)切分步長以A或B為參考點(算法描述中以A為參考點)對AB進行切分，得到切分點列表C=(C1,C2,…，Cn);設(shè)置兩個計數(shù)器變量i、j，并使i= 0，j= 0，轉(zhuǎn)到步驟2)。

2)改變i使i=i+1；計算點P到點Ci的距離PCi，點Ci到點B的距離CiB，點Ci到點A的距離CiA，轉(zhuǎn)到步驟3)。

3)如果PCi

4)改變i使i=i-1，如果CiB≥CiA，則記錄下Ci；用Ci+j的坐標替換A點的坐標；i=i+j，j=0，轉(zhuǎn)到步驟6)；否則轉(zhuǎn)到步驟5)。

5)記錄下AB的中點，結(jié)束算法。

6)改變j使j=j+1，轉(zhuǎn)到步驟2)。

算法共由六步組成。第一步是參考窮舉算法的原理根據(jù)切分步長對AB長導(dǎo)線源進行切分，切分的小段中除最后一小段外，其他的切分小段長度是等長的(若AB能被步長整除，則所有切分出的小段長度都是等長的)，并記下所有的切分點Ci，偶極子的位置就存在于這些切分點中和最后一個偶極子的中點，為了能找到偶極子的位置，用i變量跟蹤切分點的位置，用j變量跟蹤上一個偶極子B端到Ci的切分點個數(shù)；第二步計算PCi、CiB、CiA；第三步和第四步判斷PCi是否同時滿足偶極子的條件和CiB是否為原AB的最后一段，在兩個條件都滿足情況下，則記下Ci，然后計算出下一個待求偶極的A端位置，然后執(zhí)行第六步進入查找下一次偶極子的循環(huán)中；另外，在第四步中若CiB是原AB切分段的最后一段，則偶極子的位置直接取CiB的中點，算法結(jié)束(圖3)。

圖3 窮舉切分算法流程圖Fig．3 Flow of exhaustion segmentation method

2 算法應(yīng)用

使用兩種算法對地表長導(dǎo)線源AB，其中A極坐標為(0，0，0)，B極坐標為(500，0，0)，觀測點P坐標為(50，500，0)進行偶極子的切分，偶極子的條件取為收發(fā)距是10倍的發(fā)射距。切分的偶極子位置對比如圖4所示。

圖4 兩種算法切分結(jié)果對比Fig．4 Segmentation results comparison to two algorithms

窮舉切分算法切分步長分別選擇0.1 m、1 m和5 m，用這三種步長切分的偶極子數(shù)目都是9個。從圖4可以看出，用步長0.1 m和步長1 m切分出的偶極子位置差別很小，與理論上的準確偶極子位置基本一致；而用步長5m切分出的偶極子位置與用步長0.1 m切分出的偶極子位置差別明顯。如果將0.1 m步長切分出的偶極子看作是真實(理論上)的偶極子，那么切分步長越大切分出的偶極子的位置與真實偶極子的位置誤差就越大。采用遞歸算法切分出的偶極子數(shù)目為15個，多于0.1 m步長的窮舉法切分的偶極子數(shù)目，表明遞歸算法切分的偶極子長度要小于用窮舉法切分的偶極子的長度。

為了進一步了解場計算中需要將長導(dǎo)線源切分為多個電偶源的必要性及算法的有效性，本次研究以可控源音頻大在電磁法(CSAMT)一維電場計算為例對其說明。與切分算法相關(guān)的參數(shù)：發(fā)射源極距AB長 2 000 m，A點坐標為(0，0，0)，B點坐標為(2 000，0，0)，觀測點P位于離AB垂直距離 7 000 m 處，坐標為(0，7 000，0)，偶極子的條件取收發(fā)距是10 倍的發(fā)射距。使用文獻[1]中場的計算公式求解一個三層H型地電模型(電阻率：300 Ω·m，50 Ω·m，1 000 Ω·m；厚度：300 m，200 m)x方向的電場，并分別使用文中的兩種切分算法和不切分三種方式計算出的電場結(jié)果如圖5所示。

圖5是頻率-電場曲線圖，從圖5中可以看出，用兩種算法切分的偶極子數(shù)目都為 4 個，且用遞歸算法切分的偶極子計算出的電場曲線與用窮舉算法在步長分別取 0.1 m、1 m 切分的偶極子計算出的電場曲線基本重合，這表明了遞歸切分算法可用于實際應(yīng)用中；而步長取 50 m 切分的偶極子計算出的電場曲線與步長取 0.1 m 切分的偶極子計算出的電場曲線明顯分離，表明用較大步長切分偶極子計算出的電場會帶來較大的誤差；不做任何切分計算得到的場曲線與用0.1 m步長切分后計算得到的場曲線不重合，這表明在不滿足偶子條件的情況下，若仍然按照偶極子激發(fā)公式計算場將會得到不正確的結(jié)果。

3 結(jié)論

采用長導(dǎo)線源激勵時，為了避免收發(fā)距較小時源長度對場計算的影響，需要考慮到將長導(dǎo)線源切分為多個偶極子源，這樣長導(dǎo)線源場的計算就轉(zhuǎn)化為多個偶極子場的疊加。遞歸切分算法和窮舉切分算法都可以實現(xiàn)對長導(dǎo)線的切分。在實際物探應(yīng)用中，遞歸切分算法可直接使用，窮舉算法的切分步可取 0.1 m或 1 m。

圖5 兩種切分算法和不切分電場計算結(jié)果對比Fig．5 Electric field calculating comparison to two segmentation algorithms and no segmentation

參考文獻：

[1] 湯井田, 何繼善.可控源音頻大地電磁法及其應(yīng)用[M]．長沙: 中南大學(xué)出版社,2005．

[2] 何繼善．可控源音頻大地電磁法[M]．長沙:中南大學(xué)出版社,1990.

[3] 王剛,王書民, 雷達，等.CSAMT場源效應(yīng)試驗研究[J].物探化探計算技術(shù),2011,33(5):527-530.

[4] 趙廣茂,李志華,朱旭東，等.長導(dǎo)線源CSAMT一維正演研究[J].鐵道工程學(xué)報,2010,143(8):21-24.

[5] 殷長春. 可控源音頻磁大地電流法一維正演及精度評價[J].長春地質(zhì)學(xué)院學(xué)報,1994( 4):438- 453.

[6] 孫婭,何展翔,柳建新，等. 長導(dǎo)線源頻率域電磁測深場源靜態(tài)位移的模擬研究[J].石油地球物理勘探,2011,46( 1):149-154.

[7] 李勇,林品榮,徐寶利.電離層影響層狀介質(zhì)長導(dǎo)線源的電磁場[J].物探化探計算技術(shù),2009,31( 3):183-188.

[8] 劉昱,李云哲,黃子瑩.水平長導(dǎo)線源瞬變電磁三維積分方程模擬[J].現(xiàn)代地質(zhì),2010,24(2):397-402.

[9] 屈有恒,張貴賓,晉風明.傾斜線源的三維電場數(shù)值模擬研究[J].物探化探計算技術(shù),2007,29(5):431-435.

[10] 許振奎,張勝業(yè).層狀介質(zhì)中磁性源頻率測深的一種視電阻率算法[J].勘探地球物理進展,2004,27(5):333-336.

[11] 徐建華,朱德懷,胡文寶，等.多層介質(zhì)中偶極子場的系數(shù)遞推關(guān)系[J].石油地球物理勘探,1994,29(1):69-74.

[12] ALSUWAIYEL M H. Algorithms Design Techniques and Analysis[M].Beijing: Publishing House of Electronics Industry,2003.

[13] 王海濤,朱洪.改進的二分法查找[J].計算機工程,2006,32(10):60-62.

[14] 陳寶平.遞歸算法的設(shè)計模式與調(diào)試[J].電子科技,2011,24(9):28-33.

[15] 孫義欣,馮娜.窮舉法在程序設(shè)計中的應(yīng)用[J].計算機時代,2012,8:50-52.