李袁傲, 任政勇,2,3,4*, 湯井田,2,3,4, 吳啟紅

1 中南大學(xué)地球科學(xué)與信息物理學(xué)院, 長沙 410083 2 中南大學(xué)有色金屬成礦預(yù)測與地質(zhì)環(huán)境監(jiān)測教育部重點實驗室, 長沙 410083 3 有色資源與地質(zhì)災(zāi)害探查湖南省重點實驗室, 長沙 410083 4 自然資源部覆蓋區(qū)深部資源勘查工程技術(shù)創(chuàng)新中心， 合肥 230001 5 成都大學(xué)建筑與土木工程學(xué)院, 成都 610106

0 引言

大地電磁法(magnetotellurics，MT)是一種利用天然交變電磁場研究地球內(nèi)部電性結(jié)構(gòu)的地球物理方法(Cagniard,1953)，被廣泛運用于礦產(chǎn)資源勘察、巖石圈電性結(jié)構(gòu)研究以及石油天然氣勘探等領(lǐng)域(Chave and Jones, 2012).當前，大地電磁數(shù)據(jù)采集趨向于三維方式，采集區(qū)域地質(zhì)情況復(fù)雜且常具有起伏地形，因此需要研發(fā)高效的大地電磁數(shù)據(jù)三維反演方法.反演的迭代過程中需要多次調(diào)用正演算法來計算模型更新量，正演響應(yīng)的準確性保證了更新模型的正確性.因此，研發(fā)快速高精度的大地電磁三維正演算法是大地電磁數(shù)據(jù)解釋的研究熱點.

大地電磁正演算法主要基于積分方程法，有限差分法和有限單元法.積分方程法只離散化異常體，適合于較小異常體模型的計算.由于受限于簡單模型才有格林函數(shù)，積分方程法一般只能用于簡單模型電磁場響應(yīng)的計算與分析(Hohmann，1975; Lajoie and West，1976; Zhdanov et al.，2006; 王若等，2009).有限差分法采用差分來近似微分，具有方法簡單、編程容易等優(yōu)點，在積分方程法之后迅速受到了普及，被廣泛用于三維地電磁感應(yīng)問題的求解 (Mackie et al.，1994; Newman and Alumbaugh，2002; Siripunvaraporn et al.，2002; 沈金松，2003; 譚捍東等，2003; Streich，2009; 李展輝等，2009; Dong and Egbert，2019) .但是，由于差分近似一般適用于直線、不適合于復(fù)雜折線和曲線，有限元差分方法在處理任意起伏地形時存在一定的困難.有限元法采用在小區(qū)域內(nèi)利用簡單多項式來近似復(fù)雜未知場的思想，在小區(qū)域幾何形狀的選取上具有任意性，可以精確逼近三維復(fù)雜地質(zhì)情況以及起伏地形，從而計算出高精度的電磁場響應(yīng)，因此有限單元法逐漸成為了當前地電磁感應(yīng)正演問題的主流算法.

Coggon(1971)率先基于有限單元法計算了地球電磁場響應(yīng).隨后，地電磁感應(yīng)問題的有限元計算方法在網(wǎng)格離散化技術(shù)、線性方程組求解技術(shù)、并行加速技術(shù)、后處理方法與技術(shù)等方面取得了長足進步(Raiche，1994; Avdeev，2005; B?rner，2010; Newman，2014).與此相對的是，在微分控制方程最優(yōu)化設(shè)計方向的研究成果較少.目前，地電磁感應(yīng)問題的有限元計算一般采用電場方程、磁場方程和勢場方程.利用電場方程，Nam等(2007)和Farquharson和Miensopust (2011)進行了大地電磁問題的有限元計算；Ren等(2013)將面向目標自適應(yīng)有限元法應(yīng)用到大地電磁三維正演計算中.利用磁場方程，Siripunvaraporn等(2002)和Franke等(2007)進行了大地電磁三維場分布的計算.勢場方程方面，Haber(2000)基于A-φ勢實現(xiàn)了大地電磁三維有限元模擬；Mitsuhata和Uchida(2004)實現(xiàn)了基于T-Ω勢的大地電磁三維有限元模擬.

空解現(xiàn)象產(chǎn)生的本質(zhì)原因為在空氣和地球中采用了同樣的控制方程，由于空氣和地下具有完全不同的電性特征(電導(dǎo)率對比度為無窮大)，采用同樣的控制方程必然導(dǎo)致病態(tài)特征.為解決此問題，本文提出了一種新的空氣與地球分離混合公式.解決思路為：首先，將求解區(qū)域分解為空氣、地下兩個獨立的區(qū)域，分別采用不同的控制方程；然后，利用地空分界面上電磁場連續(xù)性條件，進行界面兩側(cè)的耦合與信息交換，保證電磁場解的唯一性.基于該解決思路，空氣層中磁場H分解為電矢量勢T和磁標量勢φ，利用Coulomb規(guī)范條件，使得矢量勢T和標量勢φ滿足具有Laplace結(jié)構(gòu)的勢場方程.地下區(qū)域中，采用雙旋度結(jié)構(gòu)的電場方程.在地下區(qū)域采用電場雙旋度方程的原因在于：地下介質(zhì)的電導(dǎo)率大于零，這使得雙旋度系統(tǒng)在地下區(qū)域避免了零空間問題且有良好的穩(wěn)定性.

最后，我們計算了國際標準測試模型COMMEMI 3D-1、DTM-1與梯形山模型的大地電磁響應(yīng)，將計算結(jié)果與其他算法進行了對比，來驗證本文算法的正確性，以及處理起伏地表問題的能力.我們還進行了系統(tǒng)條件數(shù)與特征值的計算和對比分析，來檢驗本文算法對電場雙旋度系統(tǒng)病態(tài)性的改進.

1 節(jié)點矢量混合有限元系統(tǒng)

利用大地電磁的場源位于求解區(qū)域外的一般性假設(shè)，在求解區(qū)域內(nèi)(見圖1)，大地電磁問題的電磁場響應(yīng)滿足頻率域Maxwell方程組(選取時諧因子為e-iω t)：

(1)

其中，E為電場矢量，H為磁場矢量，i為純虛數(shù)單位，ω為角頻率，μ為磁導(dǎo)率，

(2)

其中，J為電流密度，σ為電導(dǎo)率，ε為介電常數(shù)，

(3)

(4)

圖1 MT地電模型示意圖Ω0表示空氣區(qū)域，Ω1表示地下區(qū)域，Γ表示地空邊界，Γ0表示空氣外部邊界，Γ1表示地球外部邊界，n0為空氣區(qū)域外邊界的外法向向量，n1為地球區(qū)域外邊界的外法向向量，n01為地空邊界上空氣一側(cè)的外法向向量，方向由空氣指向地球，n10為地空邊界上地球一側(cè)的外法向向量，方向由地球指向空氣，n為約定的地空邊界上法向向量，方向由地球指向空氣.Fig.1 Schematic diagram of MT geoelectric modelΩ0 represents the air area, Ω1 represents the underground area, Γ represents the air-earth interface，Γ0 represents the outer boundary of air，Γ1 represents the outer boundary of earth，n0 is the outer normal vector of the outer boundary of the air region, n1 is the outer normal vector of the outer boundary of the earth region, n01 is the outer normal vector of the air side on the air-earth interface, and the direction is from the air to the earth, n10 is the outer normal vector on the side of the earth on the air-earth interface, and the direction is from the earth to the air, n is the normal vector on the agreed air-earth interface, with the direction from the earth to the air.

1.1 空氣層電磁場積分弱形式

(5)

將公式(5)代入公式(2)中得到：

(6)

其中H0表示空氣中的磁場.引入標量函數(shù)φ的負梯度來表示公式(6)中的矢量部分，可得到：

(7)

將公式(7)代入公式(2)并取旋度得到：

即：

(8)

將公式(8)中的矢量勢T的雙旋度結(jié)構(gòu)進行化簡，強加Coulomb規(guī)范條件來唯一化矢量勢T：

(9)

(10)

(11)

綜上，空氣區(qū)域的控制方程為

(12)

為了求解公式(12)，采用Galerkin節(jié)點有限元法(Jin，2002)，先假設(shè)余量rT為

(13)

余量rφ為

(14)

采用線性插值基函數(shù)，四面體單元中任意點處的矢量勢T以及標量勢φ可由節(jié)點基函數(shù)表示為

(15)

其中Txi，Tyi，Tzi分別為矢量勢T的三個分量，Li為節(jié)點基函數(shù)，在空氣區(qū)域中滿足以下條件：

?Ω0Li·rTdv=0,

?Ω0Li·rφdv=0.

(16)

利用Green恒等式：

(17)

(18)

綜合公式(15)—(18)，給出公式(12)的積分弱形式：

高溫環(huán)境下GH3536基體軟化，導(dǎo)致其粘著趨勢上升，存在外部載荷時易發(fā)生塑性變形，在磨球滑動過程中以剝層的方式不斷脫落的磨屑被膠合在一起，形成層狀物(見圖5a)并逐漸被撕裂。另外，基體磨痕的局部區(qū)域可見沿滑動方向的犁溝，這可能是由于部分磨屑在磨損過程中嵌入到基體中，在滑動中推擠基體使之塑性流動并犁出溝槽，從而引起磨粒磨損。由圖5b可見，NiAlW涂層的磨痕比基體更為光滑和平整，沒有明顯的單個顆粒沉積物脫落和涂層整塊粘著或撕裂現(xiàn)象。涂層的磨損機制以微切削為主，其表面的單個顆粒沉積物凸起部分首先被切削，在外界的反復(fù)撞擊與擠壓下發(fā)生破碎，最后脫落并造成涂層的磨損[15-16]。

(19)

寫成分量形式為

(20)

公式(19)與公式(20)即為空氣中的控制方程積分弱形式.

1.2 地球電磁場積分弱形式

地下區(qū)域中，電場滿足雙旋度方程：

(21)

采用Galerkin矢量有限元法(Jin，2002)，求解公式(21).先假設(shè)余量RE為

(22)

采用棱邊基函數(shù)，四面體單元中的電場可表示為

(23)

其中E1i為對應(yīng)棱邊上的電場切向分量，Ni為棱邊基函數(shù)，計算區(qū)域內(nèi)滿足以下條件：

?Ω0Ni·REdv=0.

(24)

(25)

其中a,b分別為矢量函數(shù)，則

(26)

(27)

(28)

公式(28)即為地下介質(zhì)區(qū)域控制方程的積分弱形式.

1.3 空氣-地球耦合系統(tǒng)

結(jié)合公式(19)和公式(28)以及圖1中所定義的邊界外法向向量，空氣地球耦合系統(tǒng)的積分弱形式為

(29a)

(29b)

(29c)

根據(jù)矢量公式：

(30)

公式(29a)可以進一步改寫為

(31)

將公式(29)中與地空界面有關(guān)的項移到等式左邊，且在地空界面上有n=n10=-n01，可以得到：

(32)

(33)

在地空分界面上施加如下電磁場連續(xù)性條件:

(34)

(35a)

(35b)

(35c)

公式(35a)—(35c)中的右端項均可由1D層狀介質(zhì)的解析解給出，將公式(35)寫成矩陣形式：

(36)

其中：

公式(35)和公式(36)即為空氣地球分離策略的有限元離散方程.

其中i,j=x,y(x≠y)，Zij為阻抗，ρij為視電阻率，φij為相位，E1i和H1j分別為電場E1、磁場H1在對應(yīng)方向上的分量.

2 數(shù)值試驗

本小節(jié)所有計算均在CPU主頻為2.30 GHz，RAM超過200 GB的小型AMAX工作站上完成.

2.1 COMMIME 3D-1模型

3D-1模型是國際標準測試COMMEMI項目(Zhdanov et al., 1997)的一個模型.它是均勻半空間中嵌套一個導(dǎo)電棱柱，其頂埋深為250 m，背景電阻率為100 Ωm，棱柱電阻率為0.5 Ωm，尺寸為1000 m×2000 m×2000 m，空氣電阻率為1016Ωm,計算區(qū)域尺寸為70 km×70 km×140 km.

計算網(wǎng)格包含148576個四面體單元.傳統(tǒng)雙旋度方程共173335個自由度，本文提出的地空分離公式共168978個自由度，頻率為0.1 Hz.我們將測試結(jié)果(Air-Earth decomposition)與高精度自適應(yīng)有限元計算結(jié)果(Curl-Curl system)(Ren et al., 2013)進行對比，圖2展示了視電阻率和相位的對比，圖3展示了這兩種方法的相對誤差.空氣地球分離公式計算用時 19.0 s，電場雙旋度公式計算用時 21.8 s.

圖2和圖3顯示混合公式系統(tǒng)與雙旋度公式系統(tǒng)能夠獲得一致的視電阻率與相位曲線.視電阻率值從測線端點到y(tǒng)=0處逐漸降低，過y=0之后逐漸增大，曲線具有對稱性.視電阻率最大相對誤差0.34%，相位最大相對誤差0.18%，驗證了算法的正確性和精度.為了進一步研究混合系統(tǒng)的數(shù)值特征，我們計算并對比了特征值分布(如圖4和圖5).

圖2 3D-1模型視電阻率與相位，測試頻率為0.1 Hz, Curl-Curl system代表了自適應(yīng)有限元結(jié)果(Ren et al., 2013),Air-Earth decomposition 代表了本文計算的結(jié)果Fig.2 Comparison of calculated apparent resistivity and phase on the 3D-1 model among our solution (denoted by Air-Earth decomposition) and the reference solution (Curl-Curl system, Ren et al., 2013), at a frequency of 0.1 Hz

圖3 3D-1模型視電阻率與相位相對誤差(高精度自適應(yīng)有限元計算結(jié)果為參考值(Ren et al., 2013))，測試頻率為0.1 HzFig.3 Relative error of apparent resistivity and phase by comparing to the reference solution (Curl-Curl system,Ren et al., 2013) for the 3D-1 model at a frequency of 0.1 Hz

圖4表明兩個系統(tǒng)的最大特征值分布具有相似性.圖5表明地空分離公式系統(tǒng)有更多的非零或模長更大的最小特征值.當一個系統(tǒng)的特征值全部位于復(fù)數(shù)平面的左半部分時，系統(tǒng)是收斂的.與此相反，當特征值全部位于復(fù)數(shù)平面的右半部分時，系統(tǒng)則是發(fā)散的.從圖4與圖5可以直觀的看出，空氣地球分離公式的特征值更多地集中在復(fù)數(shù)平面的左半部分且距離原點更遠，因而空氣地球分離公式相對于電場雙旋度公式有更好的穩(wěn)定性.根據(jù)條件數(shù)的定義(最大特征值與最小特征值之比)，地空分離系統(tǒng)擁有更小的條件數(shù).為了定量研究地空分離系統(tǒng)的條件數(shù)，計算了不同頻率下的條件數(shù)(從0.001 Hz，0.01 Hz，…，10 Hz共六個對數(shù)等間隔分布的頻點).測試結(jié)果如圖6所示，結(jié)果表明混合系統(tǒng)的條件數(shù)相對于雙旋度系統(tǒng)均降低，下降可達3個數(shù)量級.

圖4 3D-1模型上地空分離公式和雙旋度電場公式對應(yīng)的前1000個最大特征值分布的對比分析，測試頻率為0.1 HzFig.4 The distribution of 1000 maximum eigenvalues for two methods (one using the Air-Earth decomposition,and one using the Curl-Curl system) on the 3D-1 model, at a frequency of 0.1 Hz

圖5 3D-1模型上地空分離公式和雙旋度電場公式對應(yīng)的后1000個最小特征值分布的對比分析，測試頻率為0.1 HzFig.5 The distribution of 1000 minimum eigenvalues for two methods (one using the Air-Earth decomposition,and one using the Curl-Curl system) on the 3D-1 model at a frequency of 0.1 Hz

圖6 3D-1模型地空分離公式和雙旋度電場公式的條件數(shù)對比Fig.6 Comparison of condition numbers for two methods (one using the Air-Earth decomposition, and one using the Curl-Curl system) on the 3D-1 model at different frequencies

2.2 梯形山模型

為了測試對地形的適應(yīng)性，計算了梯形山模型(Nam et al., 2007).如圖7所示，梯形山頂部為450 m×450 m的正方形，底部為2000 m×2000 m的正方形且高為450 m.地形電阻率為100 Ωm，空氣的電阻率為1016Ωm，測線方向為x=0方向，計算區(qū)域尺寸為25 km×25 km×50 km.計算區(qū)域包含186168個四面體單元，雙旋度系統(tǒng)共227241個自由度，地空分離系統(tǒng)共219957個自由度，測試頻率為2 Hz.地空分離公式和雙旋度電場公式的模型視電阻率與相位曲線對比如圖8所示，相應(yīng)的相對誤差見圖9.空氣地球分離公式計算用時 54.26 s，電場雙旋度公式計算用時 59.66 s.

圖7 梯形山模型Fig.7 The trapezoidal hill model

圖8 梯形山模型視電阻率與相位，測試頻率為2 Hz，Curl-Curl system代表了自適應(yīng)有限元結(jié)果(Ren et al., 2013),Air-Earth decomposition代表了本文計算的結(jié)果Fig.8 Comparison of calculated apparent resistivity and phase on the trapezoidal hill model among our solution (denoted by Air-Earth decomposition) and the reference solution (Curl-Curl system, Ren et al., 2013) at a frequency of 2 Hz

測試結(jié)果顯示：視電阻率最大相對誤差0.25%，相位最大相對誤差0.50%.在視電阻率曲線中，ρxy在梯形山基底部位的電阻率稍高于圍巖的電阻率，而在起伏地形以及山頂部分則低于圍巖電阻率，ρyx則呈現(xiàn)兩端高中間低的電阻率變化趨勢，這是由于電流的地形效應(yīng)引起的.xy模式中的電流沿y方向(垂直測線方向)傳播并在山體的側(cè)面上向山頂彎曲，使得沿測線方向山腳處的電流發(fā)散而山脊和山頂處電流聚集，導(dǎo)致了山腳的電阻率高于圍巖而山脊和山頂?shù)碾娮杪实陀趪鷰r.在yx模式中，電流沿x軸方向(沿測線方向)傳播，電流的傳播路徑發(fā)生彎曲且都在x=0方向上匯聚，因而在yx模式中，ρyx低于圍巖電阻率.

為了研究混合系統(tǒng)與雙旋度系統(tǒng)的數(shù)值特征，我們計算了特征值分布，測試結(jié)果如圖10和圖11.與3D-1模型類似，圖10表明兩個系統(tǒng)具有相似的最大特征值分布，圖11表明地空分離系統(tǒng)擁有模長更大的最小特征值，且地空分離系統(tǒng)的特征值更加集中在復(fù)數(shù)平面的左半部分且距離原點更遠.

圖10 梯形山模型上地空分離公式和雙旋度電場公式對應(yīng)的前1000個最大特征值分布的對比分析，測試頻率為1 HzFig.10 The distribution of 1000 maximum eigenvalues for two methods (one using the Air-Earth decomposition,and one using the Curl-Curl system) on the trapezoidal hill model at a frequency of 1 Hz

圖11 梯形山模型上地空分離公式和雙旋度電場公式對應(yīng)的后1000個最小特征值分布的對比分析，測試頻率為1 HzFig.11 The distribution of 1000 minimum eigenvalues for two methods (one using the Air-Earth decomposition,and one using the Curl-Curl system) on the trapezoidal hill model at a frequency of 1 Hz

我們在這個模型中選取0.01 Hz，0.1 Hz，…，100 Hz共五個頻點測試了系統(tǒng)的條件數(shù).如圖12所示，模型的條件數(shù)分布與COMMEMI 3D-1模型的條件數(shù)分布特征一致，條件數(shù)有1～3個數(shù)量級的改善，在10 Hz頻率下兩系統(tǒng)條件數(shù)差異最大.

圖12 梯形山模型地空分離公式和雙旋度電場公式的條件數(shù)對比Fig.12 Comparison of condition numbers of the Air-Earth decomposition and the Curl-Curl system on the trapezoidal hill model at different frequencies

2.3 DTM-1模型

該模型由嵌套在均勻半空間中的三個異常體組成，其模型示意圖如圖13，其中ρ1=10 Ωm,ρ2=1 Ωm,ρ3=10000 Ωm，介質(zhì)電阻率最大對比度可達10000∶1，模型外邊界距離地面250 km.計算區(qū)域尺寸為250 km×250 km×500 km，測線沿y=0方向分布.

圖13 DTM-1模型結(jié)構(gòu)Fig.13 Structure of DTM-1 model

計算區(qū)域包含289178個四面體單元，雙旋度系統(tǒng)共343978個自由度，地空分離系統(tǒng)共343615個自由度.選取測試頻率為0.01 Hz進行計算，空氣地球分離公式計算用時 70.47 s，電場雙旋度公式計算用時 71.86 s.模型的電磁響應(yīng)曲線如圖14所示，相對誤差如圖15所示，結(jié)果顯示兩種系統(tǒng)的響應(yīng)曲線吻合度十分高，視電阻率最大相對誤差0.20%，相位最大相對誤差0.43%.與梯形山模型視電阻率曲線分析類似，該模型的視電阻率的變化趨勢依舊可以根據(jù)電流的傳播方向來分析.由于ρ3=10000 Ωm高于圍巖電阻率，使得通過異常體的電流會向周圍發(fā)散，則在x=[-25000 m,0]的范圍內(nèi)電流向測點聚集，反應(yīng)在視電阻率曲線上則是這一段范圍內(nèi)的視電阻率呈現(xiàn)下降趨勢.相反，ρ1=10 Ωm,ρ2=1 Ωm低于圍巖電阻率，電流向異常體聚集，所以在x=[0,25000 m]范圍內(nèi)測點處的電流發(fā)散，視電阻率升高.

圖14 DTM-1模型視電阻率與相位，測試頻率為0.01 Hz，Curl-Curl system代表了自適應(yīng)有限元結(jié)果(Ren et al., 2013),Air-Earth decomposition代表了本文計算的結(jié)果Fig.14 Comparison of calculated apparent resistivity and phase on the DTM-1 model among our solution (denoted by Air-Earth decomposition) and the reference solution (Curl-Curl system, Ren et al., 2013) at a frequency of 0.01 Hz

圖15 DTM-1模型電磁響應(yīng)相對誤差(高精度自適應(yīng)有限元計算結(jié)果為參考值(Ren et al., 2013))，測試頻率為0.01 HzFig.15 Relative error of our apparent resistivity and phase by comparing to the reference solution (Curl-Curl system,Ren et al., 2013) for the DTM-1 model at a frequency of 0.01 Hz

同樣，我們也在該工作頻率下測試了兩個系統(tǒng)的特征值分布，如圖16和圖17所示,DMT-1模型的特征值分布特征與3D-1模型以及梯形山模型的測試結(jié)果相似，空氣地球分離公式系統(tǒng)矩陣的特征值更多地分布在復(fù)數(shù)平面左半部分距離原點更遠處.我們在這個基礎(chǔ)上選擇了從0.001 Hz，0.01 Hz，…，100 Hz共六個對數(shù)等間隔分布的頻點進行了條件數(shù)的測試，結(jié)果如圖18所示.根據(jù)圖18可以得知，混合系統(tǒng)對雙旋度系統(tǒng)的病態(tài)性做出了改進，在1 Hz的工作頻率下，條件數(shù)的改善最明顯，可達3個數(shù)量級.

圖16 DTM-1模型上地空分離公式和雙旋度電場公式對應(yīng)的前1000個最大特征值分布的對比分析，測試頻率為0.01 HzFig.16 The distribution of 1000 maximum eigenvalues for two methods (one using the Air-Earth decomposition,and one using the Curl-Curl system) on the DTM-1 model at a frequency of 0.01 Hz

圖17 DTM-1模型上地空分離公式和雙旋度電場公式對應(yīng)的后1000個最小特征值分布的對比分析，測試頻率為0.01 HzFig.17 The distribution of 1000 minimum eigenvalues for two methods (one using the Air-Earth decomposition,and one using the Curl-Curl system) on the DTM-1 model at a frequency of 0.01 Hz

圖18 DTM-1模型地空分離公式和雙旋度電場公式的條件數(shù)對比Fig.18 Comparison of condition numbers of the Air-Earth decomposition and the Curl-Curl system on the DTM-1 model at different frequencies

3 結(jié)論

我們推導(dǎo)了一種新的基于空氣地球分離策略的大地電磁混合計算公式，在空氣中采用Laplace算子取代雙旋度算子，將雙旋度系統(tǒng)分解為Helmholtz-Curl耦合系統(tǒng)，有效改善了系統(tǒng)矩陣在低頻條件下的病態(tài)性特征.

模型測試驗證了新系統(tǒng)的正確性與精度.在相同網(wǎng)格的條件下，本文算法的自由度數(shù)低于電場雙旋度方程的自由度數(shù).在計算速度上，空氣地球分離公式稍高于電場雙旋度公式.特征值分布與條件數(shù)分布表明新系統(tǒng)具有更好的特征值分布與更低的條件數(shù)，在數(shù)值上更加穩(wěn)定，對雙旋度系統(tǒng)的病態(tài)性做出了有效的改進，因而新系統(tǒng)應(yīng)更適用于低頻條件下的大地電磁正演工作.

本研究基于空氣地球分離思路推導(dǎo)了一種具有良好條件數(shù)的混合公式，為尋求大地電磁最優(yōu)化控制方程開啟了新的求解思路.基于新思路，我們期待更多、性能更佳的大地電磁混合公式.