馮建新，杜 良

（1.湖南工程職業(yè)技術(shù)學(xué)院，長沙 410151；2.江蘇海洋大學(xué) 土木與港海工程學(xué)院，江蘇 連云港 222005）

大地電磁測深法[1-2]具有探測深度大、工作成本低等特點(diǎn)，最初主要應(yīng)用于地殼上地幔結(jié)構(gòu)探測。從1960—1965年開始，我國將大地電磁法應(yīng)用于大構(gòu)造區(qū)的深部結(jié)構(gòu)研究中，在探查深部礦產(chǎn)資源中發(fā)揮重要作用。隱伏地質(zhì)構(gòu)造是大地電磁法探測的難點(diǎn)，本文采用非均勻網(wǎng)格剖分的有限單元法和高效的數(shù)值算法對常見隱伏地質(zhì)構(gòu)造模型進(jìn)行了模擬分析，為解釋大地電磁法觀測資料提供重要參考依據(jù)。

1 大地電磁場特征與模擬策略

1.1 邊值問題

根據(jù)Maxwell方程組，角頻率為ω（時間因子為e-iωt）的定態(tài)電磁場方程為

在二維電性結(jié)構(gòu)中，假定構(gòu)造走向?yàn)閦軸，y軸垂直向上，x軸和y軸相垂直并保持水平。由式（1），可得出在橫電極化TE模式下的分量EZ和橫磁極化TM模式下的分量HZ，分別滿足以下偏微分方程

公式（2）和公式（3）可統(tǒng)一表示成

式中：?為二維哈密頓算子。在TE模式下，u=Ex，τ＝在TM模式下

TE和TM極化模式的研究區(qū)域如圖1所示。內(nèi)外邊界條件歸納為：①在上邊界AB上，取u=1；②設(shè)左右邊界AD、BC離電性不均勻區(qū)域足夠遠(yuǎn)，電磁場在AD，BC上左右對稱，取?u/?n=0；③設(shè)下邊界CD離電性不均勻區(qū)域足夠遠(yuǎn)，?。虎茉趦?nèi)部2種介質(zhì)的分界面Γ上，TE模式下取u1=u2；TM模式下取τ1?u1/?n=τ2?u2/?n。

圖1 TE和TM極化模式的研究區(qū)域

如果假定模型中存在一個電性不均勻區(qū)域，通過構(gòu)造泛函，可得到與上述邊值問題等價的變分問題

1.2 MT二維有限單元模擬策略

對于二維介質(zhì)的大地電磁正演模擬，假設(shè)研究區(qū)域中存在一個電性不均勻體，如圖1所示。可構(gòu)造泛函如下

對泛函式（6）進(jìn)行變分計(jì)算為

根據(jù)2種介質(zhì)分界面上電場切向分量具有連續(xù)性，即介質(zhì)分界面兩側(cè)是連續(xù)的，則上式中右側(cè)第二項(xiàng)和第四項(xiàng)為零，因此可推導(dǎo)得出式（8），Γ為研究區(qū)域的外邊界，Γ1為異常體內(nèi)邊界，n1和n2分別是研究區(qū)域Ω1和異常體Ω2的外法向，方向相反。式8如下

根據(jù)2種模式下的內(nèi)邊界條件，可得

因此，在F（u）的變分中不會出現(xiàn)內(nèi)邊界條件。在泛函去極值的過程中，內(nèi)邊界條件將自動滿足，屬于自然邊界條件，故當(dāng)研究區(qū)域存在多個異常體時內(nèi)邊界條件亦自動滿足。將外邊界條件帶入式（8）中，可得

對上式進(jìn)行移項(xiàng)化簡后，可得

因此，方程式▽·（η▽u）+λu=0與邊界條件所構(gòu)成的邊值問題等價于下面的變分問題

采用規(guī)則的矩形單元對研究區(qū)域進(jìn)行剖分，對電性均一的區(qū)域采用大網(wǎng)格剖分[3]，對電性變化的區(qū)域采用小網(wǎng)格剖分，給每個單元賦予電阻率值。以每個網(wǎng)格單元的4個角點(diǎn)為節(jié)點(diǎn)，在單元內(nèi)進(jìn)行雙線性插值計(jì)算生成單元內(nèi)的剛度矩陣，最后將每個單元內(nèi)的剛度矩陣合成總體剛度矩陣[4]。矩形網(wǎng)格剖分的單元編號和節(jié)點(diǎn)編號的順序原則如圖2所示。

圖2 矩形網(wǎng)格剖分單元圖

將公式（5）中的區(qū)域積分分解成各單元積分之和

求解積分方程，最后可得到一個線性代數(shù)方程組

式中：K為大型稀疏矩陣。

代入邊界條件，可形成方程右端項(xiàng)，上式變?yōu)?/p>

解線性方程組[5]，得到各節(jié)點(diǎn)的u值，它代表TE模式下的Ez或TM模式下的Hz。

2 模型正演模擬分析

圖3為一典型的正斷層構(gòu)造模型示意圖，該模型由上至下分為3個電性層。圖4、圖5為在ρ1＜ρ2＜ρ3與ρ1＞ρ2＜ρ3地層電性關(guān)系下的正演模擬結(jié)果。

圖3 隱伏正斷層構(gòu)造模型示意圖

圖4為隱伏正斷層在ρ1＜ρ2＜ρ3情況下正演模擬結(jié)果。對于TE模式，視電阻率-頻率擬斷面圖中可大致判斷斷層兩盤的相對位置關(guān)系，且在高頻部分對斷層與覆蓋層之間的上分界面表現(xiàn)比較明顯，而等值線在低頻部分對斷層與基地之間的下分界面表現(xiàn)相對平滑，兩盤之間的相對關(guān)系已比較模糊；相對于視電阻率-頻率擬斷面圖，阻抗相位-頻率擬斷面圖對巖層之間電性變化的反映較為明顯，可以更清楚地反映出斷層構(gòu)造的形態(tài)特征。對于TM模式，地下介質(zhì)的電性變化對視電阻率-頻率擬斷面圖造成的影響比較大，主要表現(xiàn)為在斷層面處形成了一個比較明顯的振蕩干擾，使頻率的穿透深度未達(dá)到基底，斷層高阻一側(cè)，低頻部分的視電阻率比真實(shí)電阻率高，低阻一側(cè)低頻部分的視電阻率比真實(shí)電阻率低；但在TM模式下的阻抗相位-頻率擬斷面圖的低頻部分并未受到干擾，因此可以彌補(bǔ)視電阻率-擬斷面圖在低頻部分對斷層形態(tài)的反映失真。

圖4 隱伏正斷層在ρ1＜ρ2＜ρ3情況下正演模擬結(jié)果（ρ1=60 Ω·m；ρ2=300 Ω·m；ρ3=1 000 Ω·m）

圖5為ρ1＞ρ2＜ρ3時的正演結(jié)果。2種極化模式下的視電阻率-頻率擬斷面圖和阻抗相位-頻率擬斷面圖對正斷層的響應(yīng)特征與高阻基底的基本一致?？蓪鄬訋Э醋饕粋€低阻板，比較2種模式下的視電阻率與阻抗相位對構(gòu)造體響應(yīng)特征，可知阻抗相位曲線的變化幅度較大，且在高頻部分的覆蓋層與斷層的分界面處的阻抗相位均大于45°。

圖5 隱伏正斷層在ρ1＞ρ2＜ρ3情況下正演模擬結(jié)果（ρ1=300 Ω·m；ρ2=60 Ω·m；ρ3=1 000 Ω·m）

3 應(yīng)用實(shí)例

MT1測線位于湖南省衡陽某勘探區(qū)的北部，總剖面長2 000 m。其中大地電磁法測深點(diǎn)19個，點(diǎn)距100 m；CSAMT測深點(diǎn)40個，點(diǎn)距50 m。根據(jù)已有的地質(zhì)資料，探測區(qū)域發(fā)育有正斷層構(gòu)造，地層地球物理參數(shù)特征見表1。

表1 地層地球物理參數(shù)特征表

據(jù)大地電磁法勘探成果可知，圖6中左側(cè)的虛線為地下隱伏正斷層的位置。推測在100號測點(diǎn)附近存在正斷層，斷層往東傾斜，傾角較大，傾角約70°。在下盤的斷裂邊沿的等值線出現(xiàn)了低阻振蕩，然后又恢復(fù)值相對高阻，這與TM極化模式的正演結(jié)果中在斷層面上出現(xiàn)的視電阻率振蕩干擾基本一致。兩者對正斷層的上下盤的空間位置關(guān)系的等值線響應(yīng)也大致相同。由此可得出正演模擬與實(shí)際觀測資料之間的互相驗(yàn)證關(guān)系。推測1 500號測點(diǎn)附近存在次級逆斷層，斷層往西傾斜，傾角約30°。

圖6 MT1測線大地電磁法勘探成果圖

根據(jù)綜合勘探成果（圖7）可知，MT1線開展了CSAMT勘探，該測線人文電磁干擾較少，中淺部的勘探成果可以參考CSAMT勘探成果，中深部勘探成果可以參考大地電磁法勘探成果。推測200號測點(diǎn)附近存在一正斷層，斷層往東傾斜，傾角較大，傾角約70°。推測該正斷層為樟樹灣斷裂。推測1 600號測點(diǎn)附近存在一次級逆斷層，斷層往西傾斜，傾角約40°。在現(xiàn)測中，阻抗相位對地電信息的響應(yīng)比較敏感，很難獲取有效的阻抗相位剖面圖。

圖7 MT1線CSAMT三維成果圖

4 結(jié)論

本文通過對隱伏正斷層地質(zhì)構(gòu)造體的正演模擬，總結(jié)了隱伏正斷層構(gòu)造體對大地電磁場的響應(yīng)特征，采用有效的信號處理方法提高了模擬精度。應(yīng)用工程實(shí)例[6]與正演模擬進(jìn)行了相互驗(yàn)證。比較隱伏正斷層模型2種極化模型下的正演結(jié)果，TE極化模式可以很好地反映出地質(zhì)構(gòu)造的形態(tài)特征，且橫向分辨率較高，可以用來判斷構(gòu)造類型；TM極化模式下受地下電性變化，畸變干擾比較嚴(yán)重，但在縱向上的分辨率較高，能放大不同巖性之間的電性差異，能更好地反映地下構(gòu)造體在縱向上的相對位置，可以用來判斷構(gòu)造的規(guī)模。通過正演模擬可以幫助認(rèn)識隱伏構(gòu)造體對大地電磁場的響應(yīng)規(guī)律，這些規(guī)律對觀測資料的解釋工作具有重要的指導(dǎo)意義。