劉 慶,徐鳳姣,Osborne Kachaje
(長江大學 油氣資源與勘探技術教育部重點實驗室,湖北 武漢 430100)
20世紀50年代初Tikhonov(1950)[1]和Cagnird(1953)[2]分別提出利用天然電磁場進行勘探的方法。
有限單元法由于精度高,能夠模擬復雜地形,因此在國內外得到了廣泛應用。徐世浙[3-4]對有限單元法進行了深入研究,為國內學者使用有限單元法進行MT的正演奠定了理論基礎。
劉鵬茂[5]、顧關文[6]等基于CPU實現(xiàn)了大地電磁的二維與三維反演,取得了一定的加速效果。本文采用CPU(OpenMP)和GPU(CUDA)異構并行的方式進行MT的二維正演模擬,并采用全稀疏存儲的方式來存儲有限元中的剛度矩陣,有效地減少了單次MT二維正演所需的時間和對內存的占用。
根據(jù)麥克斯韋方程,角頻率為ω(時諧因子為e-iωt)的定態(tài)電磁場的方程是
▽×E=iwμH
(1)
▽×H=(σ-iωε)E
(2)
假定電性結構是二維的,取走向為y軸,x軸與y軸垂直,保持水平,z軸在垂直方向上??芍?u/?y=0,將(1)、(2)展開可得到兩個獨立的方程組,即E型波和H型波,方程組的具體形式如下。
E型:
(3)
H型:
(4)
本文選用H型波的方程,并基于有限單元法進行大地電磁二維正演。其邊值問題可以歸納為:
(5)
邊值問題(5)與下列變分問題等價:
(6)
利用4節(jié)點的等參單元來進行數(shù)值模擬,采用雙線性插值,將式(6)中的區(qū)域積分分解為各單元積分之和:
(7)
單元積分:
(8)
(9)
單元積分:
(10)
(11)
單元積分:
(12)
(13)
(14)
δF(u)=δuTKu=0
(15)
由δu得任意性,可得Ku=0。解線性方程組前,將AB線上得邊界值代入。解得線性代數(shù)方程組后,得各節(jié)點得u。
本文采用三元組的存儲方式來儲存總體剛度矩陣K。并分別采用全稀疏存儲和不做稀疏處理的存儲對不同網格大小的二維區(qū)域進行MT二維正演,所得的剛度矩陣所占用的內存對比如圖1所示。
圖1 完全存儲與稀疏存儲對比
由圖1可知,采用完全存儲的方式時,剛度矩陣中的元素將隨著網格節(jié)點的個數(shù)以二次方的量級增長,而全稀疏存儲的方式則對內存的需求較小。
本文采用CPU(OpenMP)和GPU(CUDA)異構并行的方式來進行MT二維正演計算,加速效果如表1所示。
表1 MT二維正演并行效率對比
在二維情況下,設計了如圖2所示的二維正演模型,圖3為其二維正演視電阻率斷面圖。
圖2 二維正演模型
圖3 正演結果
從圖3中可以看出正演的視電阻率斷面圖可以很好的對正演模型做出響應,進一步地驗證了MT二維并行正演程序的正確性。
1)采用全稀疏存儲的方式避免了由于存儲剛度矩陣中的零元素而造成的內存的浪費,為模擬范圍更大、網格更多、精度更高的MT二維正演提供了便利。
2)采用CPU和GPU異構并行的方式較明顯地提高了單次大地電磁二維正演的速度。