基于Lorenz規范條件下磁矢勢和標勢耦合方程的頻率域電磁法三維正演

陳輝,殷長春,鄧居智

1 吉林大學地球探測科學與技術學院，長春　130026 2 東華理工大學放射性地質與勘探技術國防重點學科實驗室，南昌　330013

?

基于Lorenz規范條件下磁矢勢和標勢耦合方程的頻率域電磁法三維正演

陳輝1,2,殷長春1*,鄧居智2

1 吉林大學地球探測科學與技術學院，長春130026 2 東華理工大學放射性地質與勘探技術國防重點學科實驗室，南昌330013

為了克服空氣層和地表耦合以及避免一次場計算，開發適合不同類型場源、不同應用范圍的頻率域三維正演模擬統一平臺，本文從麥克斯韋基本方程出發，推導基于Lorenz規范條件的磁矢勢和標勢耦合方程；通過將不同類型場源分解成一系列短導線(電性)源組合，采用交錯網格采樣和有限體積技術對方程進行離散得到對稱大型稀疏線性方程組，并采用Jacobi迭代預處理QMR(Quasi-Minimum-Residual，擬最小殘差)算法進行求解，我們成功實現不同類型場源、不同應用范圍的頻率域電磁法三維正演模擬.通過層狀模型下大地電磁法以及有限長接地導線和大回線磁性源激發下的電磁場響應模擬，并與一維解析解對比驗證算法的有效性.進而，我們利用該算法平臺的模擬結果對典型地電模型在不同場源激發下頻率域電磁法響應特征進行對比分析.本文算法研究及實現為建立頻率域電磁法三維正反演統一框架打下基礎.

頻率域電磁法；有限體積法；正演模擬；Lorenz規范；耦合勢方程

1　引言

頻率域電磁法(Frequency-Domain Electromagnetic,FDEM)通過觀測和分析不同頻率的人工場源或天然場源激發地下介質產生的二次場或總電磁場分布規律以探明地下電性分布特征.該方法種類繁多，按照場源性質(Nabighian，1988)可分為平面波場的大地電磁法(MT)、音頻大地電磁法(AMT)、甚低頻電磁法(VLF)；接地導線激發的可控源電磁法(包括地面可控源音頻大地電磁CSAMT和海洋可控源電磁法MCSEM)；不接地回線激法的電磁剖面法(包括地面電磁和航空電磁法AEM)等.這些方法在深部地球結構探測 (Selway,2014)、地熱勘探(Muoz,2014)、礦產勘查(Smith,2014)、油氣勘探(Strack,2014)以及壞境與工程勘探(Sheard et al.,2005)等各領域起到關鍵作用.

頻率域電磁法三維正、反演已成為目前電磁勘探研究的熱點和難點問題.其中，MT三維正反演技術已基本成熟，并開始走向實用化(Siripunvaraporn，2012；Miensopust et al.,2013).Avdeev(2005)，B?rner(2010)和Newman(2014)對電磁法三維正反演技術發展和挑戰進行綜述，指出在正演算法中積分方程法(Wannamaker,1991;Song et al.,2002;Farquharson et al.,2006;Zhdanov et al.,2006;Avdeev and Knizhnik,2009;Bakr and Mannseth,2009)通常利用Green函數求解關于二次電場的Fredholm積分式，方法已基本成熟，但對復雜模型的正演精度不高；有限單元法(Yoshimura and Oshiman,2002;Mukherjee and Everett,2011;Pardo et al.,2011;Yahya et al.,2012;Schwarzbach and Haber,2013；Cai et al.,2014)通常利用變分原理或加權余量法求解關于電場、磁場或耦合勢的微分方程，能夠模擬復雜地形和地下結構，但計算量巨大、求解速度慢；有限差分或有限體積法既能模擬復雜模型，又具有較快計算速度，已成為三維電磁正反演中主流算法.Mackie等(1993)使用有限差分法實現了三維大地電磁正演模擬并將計算結果與積分方程法進行了對比，驗證了方法的正確性；Smith(1996a,1996b)對交錯采樣有限差分法的原理、誤差分析以及加快計算速度的雙共軛梯度法進行了詳細闡述；譚捍東等(2003),Tan等(2006)系統論述了消去電場分量得到關于磁場的大地電磁三維交錯采樣有限差分數值模擬算法，并對實現過程中交錯網格剖分、積分公式離散化、邊界條件、方程組求解、三維張量阻抗的計算等內容做出研究，并實現了并行計算；陳輝等(2011)在此基礎上提出利用磁場散度實校正方法(RRCM)提高三維正演精度和加快正演速度；沈金松(2003)利用交錯網格有限差分法消去磁場分量得到關于電場的線性方程組，實現了三維頻率域電磁響應正演模擬，并且在文中提出在迭代過程中施加了散度校正的方法；Weiss和Constable(2006)利用交錯網格有限差分求解二次電場的雙旋度方程實現海洋電磁法三維數值模擬；鄧居智等(2011)利用交錯網格有限差分求解二次磁場的雙旋度方程實現可控源音頻大地電磁模擬；殷長春等(2014)利用交錯網格有限差分求解二次電場的雙旋度方程實現海洋電磁法各向異性三維正演模擬；Haber等(2000)針對基于電或磁場的雙旋度方程在空氣電阻率無窮大引起非橢圓方程和地面強烈耦合作用造成求解不準問題，提出利用有限體積法求解Coulomb規范條件下矢勢和標勢耦合方程，以實現頻率域電磁法三維正演，其離散線性方程組是正定的，但非對稱.Hou等(2006)和張燁等(2012)利用有限體積法求解Coulomb規范條件下二次矢勢和標勢耦合方程，實現感應測井三維各向異性模擬；周建美等(2014)利用該方法實現海洋電磁法三維各向異性正演模擬.綜上所述，為了解決人工場源源項處理和場源畸變問題，通常采用數值算法求解二次電場和磁場雙旋度方程；考慮到空氣層和地下介質耦合問題，人們近年提出求解Coulomb規范條件下二次場標勢和矢勢耦合方程，但該算法具有網格節點數巨大，求解一次場非常耗時及Coulomb規范條件產生的耦合方程非對稱等缺點.

本文從麥克斯韋方程組出發，在Lorenz規范條件下建立磁矢勢和標勢耦合天然對稱方程，采用交錯采樣網格和有限體積法離散磁矢勢和標勢方程，并將不同場源類型分解成一系列短導線(電性)源組合，對離散線性方程組采用Jacobi迭代的預處理和擬最小殘差法(QMR)求解，從而實現任意場源激發下頻率域三維電磁正演.最后，我們通過與各種不同類型一維理論解進行對比驗證了該算法的有效性，并進一步對比分析不同場源激發下典型地電模型的電磁響應特征.

2　磁矢勢和標勢耦合方程及邊界條件

在電磁勘探中，對大多數地下巖石我們可以忽略位移電流.假定時間變化因子為eiω t，則準靜態條件下麥克斯韋方程式可寫為

(1a)

(1b)

(1c)

(1d)

式中，E為電場，H為磁場，Je為電流密度，μ為磁導率，ε為介電常數，σ為電導率.

(2)

假設介質磁導率μ和介電常數ε均為常數，將式(2)代入式(1b)可得

(3)

為了保證磁矢勢A和標勢φ唯一及方程的對稱性，我們引入洛倫茲(Lorenz)規范條件：

(4)

其中c=iωμ.將式(4)代入式(3)可得

(5)

對式(3)兩邊取散度,并代入Lorenz規范條件(4)可得

(7)

可以看出，該耦合方程具有天然對稱性，對離散后形成的線性方程組求解精度和穩定性有一定提升.由于地下介質為有損介質，電磁場隨傳播距離增加呈指數衰減，在正演模擬中可選擇足夠大的求解區域Ω，在求解區域Ω外的電磁場值非常小.因此可以采用簡單的截斷邊界條件.與其對應的Lorenz規范條件下磁矢勢和標勢的邊界條件可表示為

(8)

3　頻率域三維有限體積電磁正演模擬

為了求解磁矢勢和標勢耦合方程(7)，我們將系統闡述有限體積法離散過程、場源處理及線性方程組合成求解等.

3.1區域離散

首先，將求解區域Ω采用一系列平行于x、y、z坐標軸的三組平面剖分成若干個小的六面體網格單元.假設沿x軸方向被剖分成Nx段，含有l個節點(l=Nx+1)；節點編號沿x軸方向遞增i=1,2,…,l，網格間距為dxi(i=1,…,Nx)；沿y軸方向被剖分成Ny段，含有m個節點(m=Ny+1)，節點編號沿y軸方向遞增，j=1,2,…,m，網格間距為dyj(j=1,…,Ny)；沿z軸方向被剖分成Nz段，含有n個節點(n=Nz+1)，節點編號沿z軸方向遞增，k=1,2,…,n，網格間距為dzk(k=1,…,Nz);共計有Nx×Ny×Nz個長方體剖分單元.

圖1　磁矢勢和標勢交錯采樣網格Fig.1　Staggered grids of magnetic vectors and scalar potentials

3.2方程離散

為實現有限體積法離散，需要定義磁矢勢三個分量和標勢體積單元大小.我們以Ax分量為例.參見圖2a，體積單元沿x方向邊長為dxi，沿y和z方向的邊長分別為(dyj-1+dyj)/2和(dzk-1+dzk)/2.因此，該單元體積為

Vi+1/2,j,k=dxi(dzk-1+dzk)(dyj-1+dyj)/4.

(9)

按照等效電阻率原理，該體積單元電導率可以表示為

(10)

圖2　磁矢勢和標勢體積單元示意圖Fig.2　Sketch showing cells of magnetic vectors and scalar potentials

同理可得Ay,Az分量體積單元表達式.對于標勢φ體積單元(參見圖2b)，沿x，y和z方向的邊長分別為(dxi-1+dxi)/2，(dyj-1+dyj)/2，(dzk-1+dzk)/2.由此，單元體積為

(11)

而其體積單元等效電導率為

(14)

(15)

(16)

對于磁矢勢和標勢耦合方程組(2)第二式，在某個標勢φ體積單元上積分可得

(18)

仿照前面離散過程，(18)式中左端四項的離散公式為

(19)

(20)

(21)

(22)

3.3場源設置

3.4離散線性方程組合成及求解

將(14)—(17)式代入方程組(7)可以建立磁矢勢Ax分量和標勢φ的離散表達式(見圖4a)，同理將(19)—(22)式代入方程組(7)可建立標勢φ和磁矢勢Ax、Ay、Az分量的離散表達式(見圖4b).將方程組(7)離散后Ax、Ay、Az、φ四個等式，并且結合邊界條件(8)式可以合成總體線性方程組：

(23)

式中Cx，Cy，Cz，Cφ及Dx，Dy，Dz為系數項，而rx，ry，rz，rφ為源項，待求未知數個數為

N=(nx-1)×(ny-2)×(nz-2)

+(nx-2)×(ny-1)×(nz-2)

+(nx-2)×(ny-2)×(nz-1)

+(nx-2)×(ny-2)×(nz-2).

(24)

對于線性方程組求解目前主要采用直接求解和迭代求解.迭代求解通常在Krylov子空間內進行.首先采用不完全LU分解、不完全Cholesky分解等預處理方法降低系數矩陣的條件數，然后采用共軛梯度法(CG),雙共軛梯度法(BiCG)，廣義最小殘量法(GMRES),擬最小殘差法(QMR)等迭代方法進行求解.該算法需要的內存小，計算速度快，適合求解場源較少的正反演問題.直接求解方法一般直接對系數矩陣進行LU分解，然后利用各種成熟的求解庫MUMPS、PARDISO等進行求解.該方法所需內存大、計算速度慢，但特別適合多源問題.本文離散線性方程組(23)的系數矩陣為大型對稱稀疏矩陣，由于僅考慮單一發射源問題，我們選用Jacobi迭代的QMR算法進行求解.

圖3　頻率域電磁法不同類型場源類型Fig.3　Transmitting sources of different types for FDEM

圖4　標勢和磁矢勢離散關系圖Fig.4　Template for discretizing magnetic vectors and scalar potentials

4　數值模擬與結果分析

為了驗證本文的頻率域三維電磁正演算法的可行性和有效性，我們利用本文算法結果和一維正演結果進行對比，研究不同場源激發條件下頻率域視電阻率和相位響應特征.

4.1算法驗證

選取水平層狀三層H模型(如圖5a).第一層厚度500 m，電阻率100 Ωm；第二層厚度1000 m，電阻率10 Ωm；第三層電阻率為1000 Ωm.求解區域大小為10000 m×10000 m×10000 m，網格剖分為100×100×100，分別向外延拓16層，延拓指數為1.5.

首先進行平面波場X極化和Y極化模式下的電磁場三維模擬，然后計算大地電磁阻抗視電阻率和相位.計算頻率為1～1000 Hz.圖5b為水平層狀介質三維和一維正演結果響應對比.由圖可見三維模擬的視電阻率和相位與一維結果完全吻合，視電阻率相對誤差最大0.3%，相位最大相對誤差2.9%，說明本文算法對大地電磁法三維模擬具有較高精度.

為驗證本文算法的廣譜有效性，我們進一步對有限長接地線源的CSAMT進行三維正演模擬.網格參數設置和上面模型一致.假設接地長導線發射源沿x方向置于地面，長度為2000 m，中心點相對投影(0,0,0)，發射頻率為10 Hz.接收點設在位于發射源中垂線的正向y軸上，點距間隔為100 m.圖6為層狀介質CSAMT三維和一維模擬結果對比.從圖可以看出，Ex分量的實部和虛部三維模擬結果幾乎與一維結果重合，最大相對誤差不超過3%；Hz分量實部和虛部也基本與一維正演結果吻合，除了Hz分量實部受場源的畸變影響在發射源附近兩個測點相對誤差達到10%，其他測點相對誤差均小于5%.考慮到CSAMT通常在遠區或過渡區進行觀測，在場源附近發生的畸變并不會影響數據解釋精度.

圖5　三層水平層狀模型及大地電磁響應結果對比Fig.5　Comparison of 1D and 3D MT responses for a 3-layer model

圖6　接地長導線水平層狀模型CSAMT一維和三維電磁響應對比Fig.6　Comparison of 1D and 3D CSAMT responses for the 3-layer model of Fig.4

圖7　不接地回線三層水平層狀模型一維和三維電磁響應對比Fig.7　Comparison of 1D and 3D EM responses for a wire-loop transmitter over the 3-layer model in Fig.4

最后對大回線源激發條件下的三維電磁場進行正演模擬.網格參數設置和上面模型一致.設置發射源為3000 m×3000 m方形回線，置于地面，中心點相對投影(0,0,0)，發射頻率為10 Hz.接收點位于發射源中垂線的正向y軸上，點距100 m.圖6為層狀介質大回線激發下的電磁場三維和一維模擬結果對比.從圖中可以看出，Ex分量的實部和虛部三維模擬和一維模擬結果重合，最大相對誤差不超過3%；Hz分量實部和虛部也基本與一維正演結果吻合，除了Hz分量實部受場源的畸變影響在靠近發射源1.5 km附近相對誤差可達10%，但其余各測點相對誤差均小于5%.因此，本文算法對大回線激發下的電磁場三維模擬具有較高精度.

4.2典型地電模型

通過水平三層模型不同類型場源的數值模擬結果與解析解的對比，證明了本文算法適合各類場源的頻域電磁法三維正演，并且具有較高精度.下面，我們研究不同類型場源復雜模型的頻率域電磁響應特征，特別是進行視電阻率和相位響應結果的對比.如圖8，我們設計一個低阻異常體模型大小為800 m×800 m×500 m，頂部埋深為500 m，異常體電阻率為10 Ωm，圍巖電阻率為100 Ωm.將10000 m×10000 m×10000 m計算區域剖分為100×100×100個單元，沿x、y、z三個方向的擴展層數為16，空氣擴展層為12.將異常體剖分成8×8×5個單元,計算頻率為1～1000 Hz.人工源頻率域電磁法采用赤道觀測裝置，有線長接地導線長度為2000 m，磁偶源為600 m×600 m回線，收發距都為6 km.

圖8　低阻異常模型Fig.8　Low-resistivity anomaly model in a homogeneous half-space

圖9為不同場源激發下低阻模型中心點視電阻率和相位結果對比.從圖可以看出，在100 Hz以上四種不同類型場源(XY模式、YX模式、接地長導線、大回線磁偶源)的高頻段視電阻率和相位曲線基本重合；隨著頻率降低XY模式、YX模式大地電磁視電阻率響應受三維異常體的影響出現分離，而相位曲線基本重合；其中XY模式受三維異常體影響要大于YX模式.對于接地長導線、大回線源，隨頻率降低觀測場由遠區逐漸進入近區，電阻率呈現快速上升，相位快速下降趨勢.其中大回線磁偶源進入近區頻率要高于接地長導線源，而且幅度大于長導線源.因此，實際電磁勘探工作中，為便于在遠區觀測通常采用接地長導線源作為發射源.

圖10為低阻異常體模型不同類型場源激發條件下視電阻率和相位擬斷面圖.由圖可以看出，當頻率高于100 Hz時，即處于遠區條件下，CSAMT、回線源赤道裝置的視電阻率響應基本相同，在異常體上方近地表位置有一個不太明顯的小范圍高阻區，下方為低阻異常區；低阻異常區的分布范圍與模型中低阻異常體的實際寬度基本相同；相位特征也具有類似特征，呈現高值異常.隨著頻率降低，XY模式視電阻率呈現低阻視電阻率異常和高值相位異常；YX模式呈現低阻加兩個對稱高阻視電阻率異常和高值相位外加兩個對稱低值異常；然而接地長導線和大回線源的赤道裝置受近場影響，呈現高阻假異常和低值相位假異常.其中大回線源進入近場頻率要高于接地長導線，而且異常響應范圍和幅值都要小.因此開發不同場源的頻率域電磁三維反演算法可克服三維異常體和場源引起的虛假異常，有助于提高數據解釋質量.

圖9　不同場源激發下低阻模型中心點視電阻率和相位結果對比Fig.9　Comparison of apparent resistivities and phases at central point in low-resistivity model for different transmitting sources

圖10　不同類型場源低阻異常模型的視電阻率和相位響應擬斷面圖橫坐標為算術坐標，縱坐標頻率為對數坐標.Fig.10　Pseudo-cross section of apparent resistivity and phase for the model in Fig.8 for different transmitting sources

5　結論

本文從麥克斯韋方程組出發，在Lorenz規范條件下建立磁矢勢和標勢耦合對稱方程.采用交錯網格和有限體積法離散磁矢勢和標勢耦合方程，并將不同場源類型分解成一系列短導線電性源的組合；對離散線性方程組采用Jacobi迭代的預處理擬最小殘差法(QMR)進行求解，我們成功實現任意場源下頻率域三維電磁正演模擬.通過對水平層狀模型不同類型場源的三維電磁模擬可以發現，本文開發的算法能夠實現不同類型場源各種不同應用范圍頻率域三維電磁正演模擬，避免耗時的一次場計算，計算精度較高.由于采用特殊的源處理技術，使得頻率域三維電磁正反演統一數據處理和解釋平臺成為可能.對不同場源激發下典型地電模型模擬發現，不同類型場源的視電阻率和相位特征和規律不同，在實際勘探中需要考慮場源的選擇和設置，同時在電磁數據解釋中應考慮場源的畸變效應.

References

Avdeev D B.2005.Three-dimensional electromagnetic modelling and inversion from theory to application.Surveys in Geophysics,26(6):767-799.Avdeev D,Knizhnik S.2009.3D integral equation modeling with a linear dependence on dimensions.Geophysics,74(5):F89-F94.

Bakr S A,Mannseth T.2009.Feasibility of simplified integral equation modeling of low-frequency marine CSEM with a resistive target.Geophysics,74(5):F107-F117.

B?rner R U.2010.Numerical modelling in geo-electromagnetics:Advances and challenges.Surveys in Geophysics,31(2):225-245.

Cai H Z,Xiong B,Han M R,et al.2014.3D controlled-source electromagnetic modeling in anisotropic medium using edge-based finite element method.Computers &Geosciences,73:164-176.

Chen H,Deng J Z,Tan H D,et al.2011.Study on divergence correction method in three-dimensional magnetotelluric modeling with staggered-grid finite difference method.Chinese Journal of Geophysics (in Chinese),54(6):1649-1659,doi:10.3969/j.issn.0001-5733.2011.06.025.

Deng J Z,Tan H D,Chen H,et al.2011.CSAMT 3D modeling using staggered-grid finite difference method.Progress in Geophysics (in Chinese),26(6):2026-2032,doi:10.3969/j.issn.1004-2903.2011.06.017.Farquharson C G,Duckworth K,Oldenburg D W.2006.Comparison of integral equation and physical scale modeling of the electromagnetic responses of models with large conductivity contrasts.Geophysics,71(4):G169-G177.Haber E,Ascher U M,Aruliah D A,et al.2000.Fast simulation of 3D electromagnetic problems using potentials.Journal of Computational Physics,163(1):150-171.

Hou J S,Mallan R K,Torres-Verdín C.2006.Finite-difference simulation of borehole EM measurements in 3D anisotropic media using coupled scalar-vector potentials.Geophysics,71(5):G225-G233.Mackie R L,Madden T R,Wannamaker P E.1993.Three-dimensional magnetotelluric modeling using difference equations—Theory and comparisons to integral equation solutions.Geophysics,58(2):215-226.Miensopust M P,Queralt P,Jones A G.2013.Magnetotelluric 3-D inversion—a review of two successful workshops on forward and inversion code testing and comparison.Geophysical Journal International,193(3):1216-1238.Mukherjee S,Everett M.2011.3D controlled-source electromagnetic edge-based finite element modeling of conductive and permeable heterogeneities.Geophysics,76(4):F215-F226.Muoz G.2014.Exploring for geothermal resources with electromagnetic methods.Surveys in Geophysics,35(1):101-122.Nabighian M N.1988.Electromagnetic Methods in Applied Geophysics.Tulsa:Society of Exploration Geophysicists.Newman G A.2014.A review of high-performance computational strategies for modeling and imaging of electromagnetic induction data.Surveys in Geophysics,35(1):85-100.

Pardo D,Nam M J,Torres-Verdín C,et al.2011.Simulation of marine controlled source electromagnetic measurements using a parallel Fourier hp-finite element method.Computational Geosciences,15(1):53-67.Schwarzbach C,Haber E.2013.Finite element based inversion for time-harmonic electromagnetic problems.Geophysical Journal International,193(2):615-634.

Selway K.2014.On the causes of electrical conductivity anomalies in tectonically stable lithosphere.Surveys in Geophysics,35(1):219-257.

Sheard S N,Ritchie T J,Christopherson K R,et al.2005.Mining,environmental,petroleum,and engineering industry applications of electromagnetic techniques in geophysics.Surveys in Geophysics,26(5):653-669.

Shen J S.2003.Modelling of 3-D electromagnetic responses in frequency domain by using staggered-grid finite difference method.Chinese Journal of Geophysics (in Chinese),46(2):281-288.Siripunvaraporn W.2012.Three-dimensional magnetotelluric inversion:an introductory guide for developers and users.Surveys in Geophysics,33(1):5-27.

Smith J T.1996a.Conservative modeling of 3-D electromagnetic fields,Part I:Properties and error analysis.Geophysics,61(5):1308-1318.

Smith J T.1996b.Conservative modeling of 3-D electromagnetic fields,Part II:Biconjugate gradient solution and an accelerator.Geophysics,61(5):1319-1324.Smith R.2014.Electromagnetic Induction Methods in Mining Geophysics from 2008 to 2012.Surveys in Geophysics,35(1):123-156.Song Y,Kim H J,Lee K H.2002.An integral equation representation of wide-band electromagnetic scattering by thin sheets.Geophysics,67(3):746-754.Strack K M.2014.Future directions of electromagnetic methods for hydrocarbon applications.Surveys in Geophysics,35(1):157-177.

Tan H D,Tong T,Lin C H.2006.The parallel 3D magnetotelluric forward modeling algorithm.Applied Geophysics,3(4):197-202.

Tan H D,Yu Q F,Booker J,et al.2003.Magnetotelluric three-dimensional modeling using the staggered-grid finite difference method.Chinese Journal of Geophysics (in Chinese),46(5):705-711.Wannamaker P E.1991.Advances in three-dimensional magnetotelluric modeling using integral equations.Geophysics,56(11):1716-1728.

Weiss C J,Constable S.2006.Mapping thin resistors and hydrocarbons with marine EM methods,Part II—Modeling and analysis in 3D.Geophysics,71(6):G321-G332.

Yahya N,Akhtar M N,Nasir N,et al.2012.Guided and direct wave evaluation of controlled source electromagnetic survey using finite element method.Journal of Electromagnetic Analysis and Applications,4(3):135-146.

Yin C C,Ben F,Liu Y H,et al.2014.MCSEM 3D modeling for arbitrarily anisotropic media.Chinese Journal of Geophysics (in Chinese),57(12):4110-4122,doi:10.6038/cjg20141222.Yoshimura R,Oshiman N.2002.Edge-based finite element approach to the simulation of geoelectromagnetic induction in a 3-D sphere.Geophysical Research Letters,29(3):9-1-9-4.Zhang Y,Wang H N,Tao H G,et al.2012.Finite volume algorithm to simulate 3D responses of multi-component induction tools in inhomogeneous anisotropic formation based on coupled scalar-vector potentials.Chinese Journal of Geophysics (in Chinese),55(6):2141-2152,doi:10.6038/j.issn.0001-5733.2012.06.036.

Zhdanov M S,Lee S K,Yoshioka K.2006.Integral equation method for 3D modeling of electromagnetic fields in complex structures with inhomogeneous background conductivity.Geophysics,71(6):G333-G345.

Zhou J M,Zhang Y,Wang H N,et al.2014.Efficient simulation of three-dimensional marine controlled-source electromagnetic response in anisotropic formation by means of coupled potential finite volume method.Acta Physica Sinica (in Chinese),63(15):159101.

附中文參考文獻

陳輝,鄧居智,譚捍東等.2011.大地電磁三維交錯網格有限差分數值模擬中的散度校正方法研究.地球物理學報,54(6):1649-1659,doi:10.3969/j.issn.0001-5733.2011.06.025.

鄧居智,譚捍東,陳輝等.2011.CSAMT三維交錯采樣有限差分數值模擬.地球物理學進展.26(6):2026-2032,doi:10.3969/j.issn.1004-2903.2011.06.017.

沈金松.2003.用交錯網格有限差分法計算三維頻率域電磁響應.地球物理學報,46(2):281-288.

譚捍東,余欽范,Booker J等.2003.大地電磁法三維交錯采樣有限差分數值模擬.地球物理學報,46(5):705-711.

殷長春,賁放,劉云鶴等.2014.三維任意各向異性介質中海洋可控源電磁法正演研究.地球物理學報,57(12):4110-4122,doi:10.6038/cjg20141222.

張燁,汪宏年,陶宏根等.2012.基于耦合標勢與矢勢的有限體積法模擬非均勻各向異性地層中多分量感應測井三維響應.地球物理學報,55(6):2141-2152,doi:10.6038/j.issn.0001-5733.2012.06.036.

周建美,張燁,汪宏年等.2014.耦合勢有限體積法高效模擬各向異性地層中海洋可控源的三維電磁響應.物理學報,63(15):159101.

(本文編輯何燕)

A finite-volume solution to 3D frequency-domain electromagnetic modelling using Lorenz-gauged magnetic vector and scalar potentials

CHEN Hui1,2,YIN Chang-Chun1*,DENG Ju-Zhi2

1 Geo-Exploration Science and Technology Institute,Jilin University,Changchun 130026,China 2 Key Laboratory of Radioactive Geology and Exploration Technology Fundamental Science for National Defense， East China Institute of Technology，Nanchang 330013,China

To eliminate the coupling between air and the earth and to avoid calculation of primary field,we develop a universal algorithm on forward modelling of the frequency-domain electromagnetic (EM)method for different source types and applications.We first present a magnetic vector and a scalar potential with Lorenz gauge based on Maxwell equation,then divide different sources into a lot of short wires.Next,we use staggered grids and the finite volume method to discrete the potential equations and implement the quasi-minimum-residual (QMR)iteration with Jacob to solve the large sparse and symmetrical linear equations system.Finally,we accomplish frequency-domain EM modelling for different source types and application areas.Through analyses and comparison of 1D and 3D MT,CSAMT and loop source responses,we demonstrate the efficiency and accuracy of the algorithm proposed in this study.Based on this,we analyze EM responses for different source types.The algorithm and numerical results presented in this paper build a framework for 3D frequency-domain EM forward modelling and inversion.KeywordsFDEM;Finite-volume method;Forward modelling;Lorenz-gauged;Coupled potentials

陳輝,殷長春,鄧居智.2016.基于Lorenz規范條件下磁矢勢和標勢耦合方程的頻率域電磁法三維正演.地球物理學報，59(8):3087-3097,

10.6038/cjg20160831.

Chen H,Yin C C,Deng J Z.2016.A finite-volume solution to 3D frequency-domain electromagnetic modelling using Lorenz-gauged magnetic vector and scalar potentials.Chinese J.Geophys.(in Chinese),59(8):3087-3097,doi:10.6038/cjg20160831.

國家青年基金項目(41404057),國家自然科學基金項目(41164003),國家重大科研裝備研究項目(ZDYZ2012-1-03和20130523MTEM05)聯合資助.

陳輝，男，1985年8月生，博士生，主要從事電磁勘查技術正反演研究.E-mail:schoolhui@163.com

殷長春，男，1965年生，教授，國家“千人計劃”特聘專家，主要從事電磁勘探理論，特別是航空和海洋電磁方面的研究.

E-mail:yinchangchun@jlu.edu.cn

10.6038/cjg20160831

P631

2015-01-21，2015-10-08收修定稿