任運(yùn)通，李貅，齊彥福，曹華科

(長(zhǎng)安大學(xué) 地質(zhì)工程與測(cè)繪學(xué)院，陜西 西安 710054)

0 引言

時(shí)間域航空電磁法采用機(jī)載移動(dòng)平臺(tái)，通過(guò)線圈發(fā)射大功率一次場(chǎng)電磁波信號(hào)，利用接收線圈接收到感應(yīng)的二次磁場(chǎng)，實(shí)現(xiàn)對(duì)地下電性結(jié)構(gòu)及能源礦產(chǎn)分布的精細(xì)探測(cè)。該方法具有分辨率高、探測(cè)速度快、通行性極好等諸多優(yōu)點(diǎn),已被廣泛應(yīng)用于油氣資源、環(huán)境工程和金屬礦產(chǎn)等領(lǐng)域[1-6]。然而，傳統(tǒng)的電阻率成像和一維反演技術(shù)已經(jīng)無(wú)法滿足三維問(wèn)題解釋的要求，開展三維反演解釋技術(shù)研究是當(dāng)前的熱點(diǎn)問(wèn)題，作為反演的基礎(chǔ)，時(shí)間域航空電磁三維正演算法受到越來(lái)越多的關(guān)注。

目前，時(shí)間域航空電磁三維有限元正演方法主要分為間接法和直接法兩種。其中的間接法是基于頻率域麥克斯韋方程組，首先計(jì)算頻率域響應(yīng)，再通過(guò)時(shí)頻轉(zhuǎn)換方法轉(zhuǎn)換到時(shí)間域[7-10]。由于間接算法存在穩(wěn)定性差等諸多問(wèn)題，人們普遍將研究方向轉(zhuǎn)移到了時(shí)間域直接算法，直接求解時(shí)間域電磁響應(yīng)。2010 年，Um等[11]采用非結(jié)構(gòu)化網(wǎng)格剖分策略，對(duì)電性源的海洋三維瞬變電磁場(chǎng)進(jìn)行了時(shí)間域矢量有限元的三維正演模擬計(jì)算。2016年，李賀等[12]人研究了直接時(shí)間域矢量有限元瞬變電磁三維正演模擬方法。2017年，齊彥福等人[13-14]基于航空電磁footprint的局部網(wǎng)格策略，利用非結(jié)構(gòu)矢量有限元方法，對(duì)復(fù)雜介質(zhì)情況下的時(shí)間域三維航空電磁響應(yīng)進(jìn)行了進(jìn)一步模擬研究。然而，由于航空電磁法采用同時(shí)移動(dòng)發(fā)射源和測(cè)點(diǎn)的方式進(jìn)行觀測(cè)，且采樣密集，每次移動(dòng)測(cè)點(diǎn)均需求解一次正演方程，產(chǎn)生巨大的計(jì)算量，傳統(tǒng)的串行計(jì)算方法無(wú)法滿足未來(lái)三維反演的效率要求?？紤]到并行計(jì)算可以有效提高程序的運(yùn)算效率，解決相同數(shù)量的問(wèn)題所需執(zhí)行的時(shí)間更短，開展了基于MPI+OpenMP并行技術(shù)的時(shí)間域航空電磁快速正演算法研究。

目前廣為使用的并行加速技術(shù)有兩種，分別是MPI(Message Passing Interface)和OpenMP(Open Multi-Processing)并行技術(shù)。MPI和OpenMP均需要與計(jì)算機(jī)語(yǔ)言結(jié)合使用，因兩者較好的通信性和可移植性而受到廣泛地應(yīng)用。在尺度上，MPI多為基于多臺(tái)處理機(jī)上的跨節(jié)點(diǎn)并行方法，OpenMP主要通過(guò)一些指令集對(duì)現(xiàn)有的Fortran或C/C++等程序進(jìn)行擴(kuò)展。并行技術(shù)在提高計(jì)算效率上具有顯著優(yōu)點(diǎn)，將該方法應(yīng)用到地球物理數(shù)值模擬和數(shù)據(jù)處理中具有良好的前景。早在1997年，Newman和Alumbaugh[15]對(duì)三維大地電磁正演方法以及靈敏度矩陣進(jìn)行了深入研究，成功將并行加速技術(shù)應(yīng)用到了電磁反演當(dāng)中，實(shí)現(xiàn)了多個(gè)節(jié)點(diǎn)的并行加速計(jì)算，極大推動(dòng)了并行計(jì)算技術(shù)在地球物理領(lǐng)域的發(fā)展。國(guó)內(nèi)對(duì)于并行計(jì)算技術(shù)應(yīng)用到地球物理的研究，起步雖然稍晚于國(guó)外，發(fā)展卻并不遜色。2005年榮瑩等[16]實(shí)現(xiàn)了基于MPI的處理機(jī)集群并行計(jì)算系統(tǒng)平臺(tái)的構(gòu)建，在現(xiàn)有的硬件條件下提高計(jì)算力； 2006年，譚捍東[17]等結(jié)合MPI的優(yōu)越性，通過(guò)頻點(diǎn)并行的方式成功實(shí)現(xiàn)了大地電磁三維正演的并行計(jì)算，得到了很高的加速比，驗(yàn)證了并行算法的穩(wěn)定性以及高效性；2011年，李小康[18]研究了頻率域航空電磁法有限單元二維正演的并行計(jì)算；2015年，陳輝[19]研究了基于MPI的航空瞬變電磁一維正反演，提高了航空電磁一維模擬技術(shù)的計(jì)算效率。目前，對(duì)于并行加速技術(shù)在電磁方法數(shù)值模擬中的應(yīng)用主要集中在大地電磁法，而航空電磁正演的并行算法研究依然停留在一維正反演和二維的頻率域正演階段。

本文將并行技術(shù)與時(shí)間域有限元算法相結(jié)合，基于三維航空電磁局部網(wǎng)格之間相互獨(dú)立的特性[20]，通過(guò)采用MPI+OpenMP的并行加速策略，在保證正演結(jié)果可靠性的基礎(chǔ)上，極大提高了時(shí)間域航空電磁有限元三維正演的計(jì)算效率。首先通過(guò)與串行有限元法和有限體積法進(jìn)行對(duì)比檢驗(yàn)本文并行程序的可靠性，然后進(jìn)行計(jì)算效率分析，討論并行加速比與硬件條件的關(guān)系，最后將該并行算法應(yīng)用于復(fù)雜起伏地表模型，模擬其航空電磁響應(yīng)。

1 方法理論

1.1 時(shí)間域航空電磁有限元正演理論

時(shí)間域麥克斯韋方程組可以表示為[21]

(1)

(2)

(3)

(4)

其中，t為時(shí)間，r為位置矢量，e(r,t)、j(r,t)、d(r,t)及b(r,t)分別表示r處在t時(shí)刻電場(chǎng)強(qiáng)度、電流密度、電位移矢量以及磁感應(yīng)強(qiáng)度，q和μ分別表示累積電荷和磁導(dǎo)率。本文假定介質(zhì)的磁導(dǎo)率與自由空間磁導(dǎo)率相同，即μ=μ0=4π×10-7H/m，介電常數(shù)與自由空間介電常數(shù)ε0相同，即ε=ε0=8.85×10-12F/m。電磁場(chǎng)之間具有如下本構(gòu)關(guān)系：

j(r,t)=σe(r,t)+js(r,t) ，

(5)

d(r,t)=εe(r,t) ，

(6)

b(r,t)=μh(r,t) ，

(7)

其中：js(r,t)表示外部施加的源電流密度，h(r,t)和σ分別表示磁場(chǎng)強(qiáng)度和電導(dǎo)率。通過(guò)消去磁場(chǎng)，可以獲得電場(chǎng)擴(kuò)散方程

(8)

本文采用非結(jié)構(gòu)矢量有限元方法進(jìn)行空間離散，四面體單元中任意位置的電場(chǎng)可以表示為：

(9)

(10)

(11)

其中：M、S分別為質(zhì)量和剛度矩陣，J為電流源項(xiàng)，針對(duì)每個(gè)單元，可由下面三式給出

(12)

(13)

(14)

其中V是單元的體積。由于各個(gè)單元相互獨(dú)立，因此Mk、Sk和Jk可以采用OpenMP并行技術(shù)加速計(jì)算。針對(duì)發(fā)射源，將發(fā)射線圈分解為若干段導(dǎo)線，每段近似為一個(gè)電偶極子[24]每個(gè)電偶極子可以表示為：

js(r,t)=δ(r-rs)vI(t)dl。

(15)

上式中：I(t)表示t時(shí)刻的電流強(qiáng)度，dl表示電偶極子長(zhǎng)度，v表示電流方向，δ表示脈沖函數(shù)。然后，利用二階后推歐拉公式對(duì)有限元方程(11)進(jìn)行時(shí)間域離散：

(3A+2ΔtB)ei+2(t)=A[4ei+1(t)-ei(t)]-2ΔtSi+2，

(16)

其中：Δt為時(shí)間步長(zhǎng)，ei(t)表示第i時(shí)刻的電場(chǎng)值。式(16)亦可簡(jiǎn)寫為：

Fe=P，

(17)

其中：F表示大型稀疏系數(shù)矩陣，P是右端項(xiàng)。當(dāng)發(fā)射階躍電流波形時(shí)，在0時(shí)刻之前發(fā)射線圈中供恒定電流，在全空間產(chǎn)生穩(wěn)定的磁場(chǎng)，根據(jù)楞次定律可知空間中任意位置的電場(chǎng)均為0，故電場(chǎng)的初始條件為：

e(r,0)=0 。

(18)

采用狄利克雷邊界條件，根據(jù)電磁波在空間中的幾何指數(shù)衰減規(guī)律，可以通過(guò)加大計(jì)算區(qū)域的擴(kuò)邊范圍，使其滿足

e|Γ=0 ，

(19)

其中Γ表示計(jì)算區(qū)域的外邊界。最后，采用直接求解器Pardiso對(duì)線性方程組進(jìn)行求解，即可獲得全空間的電場(chǎng)值，再利用法拉第電磁感應(yīng)定律(式(1))計(jì)算磁場(chǎng)響應(yīng)。

1.2 局部網(wǎng)格正演計(jì)算技術(shù)

大量的數(shù)值實(shí)驗(yàn)表明三維數(shù)值模擬結(jié)果的精度受到網(wǎng)格質(zhì)量的影響非常大，因此需要對(duì)導(dǎo)電大地和發(fā)射源進(jìn)行精細(xì)的網(wǎng)格剖分來(lái)獲得高精度正演結(jié)果。傳統(tǒng)的時(shí)間域航空電磁三維正演模擬算法采用全局網(wǎng)格加密方式，即通過(guò)一套網(wǎng)格對(duì)所有測(cè)點(diǎn)的電磁響應(yīng)進(jìn)行模擬。然而，由于航空電磁法采樣密集且觀測(cè)剖面長(zhǎng)，導(dǎo)致正演網(wǎng)格單元數(shù)量巨大。如果采用全局網(wǎng)格進(jìn)行正演模擬，將產(chǎn)生巨大的計(jì)算量，嚴(yán)重降低計(jì)算效率?？紤]到航空電磁系統(tǒng)影響范圍有限(如圖1所示)，本文采用局部網(wǎng)格加密技術(shù)，針對(duì)每個(gè)測(cè)點(diǎn)或者相鄰的幾個(gè)觀測(cè)點(diǎn)分別設(shè)計(jì)獨(dú)立的網(wǎng)格，并分別在局部網(wǎng)格上進(jìn)行正演模擬，以此在保證結(jié)果可靠性的同時(shí)提高計(jì)算效率。

圖1 航空系統(tǒng)影響范圍示意Fig.1 Moving footprint for time-domian airborne EM system

局部網(wǎng)格技術(shù)根據(jù)航空電磁系統(tǒng)的影響范圍設(shè)計(jì)合理的局部精細(xì)網(wǎng)格。圖2展示了全局網(wǎng)格和局部網(wǎng)格，其中圖2a表示全局網(wǎng)格，圖2b、c、d則分別表示針對(duì)測(cè)線上不同測(cè)點(diǎn)所設(shè)計(jì)的相互獨(dú)立的局部網(wǎng)格。從圖中可知全局網(wǎng)格對(duì)測(cè)線上所有測(cè)點(diǎn)、測(cè)線下方地空分界面以及異常體部分均進(jìn)行了相對(duì)細(xì)致的網(wǎng)格剖分，因此造成網(wǎng)格數(shù)量較大，導(dǎo)致正演模擬的計(jì)算量劇增，計(jì)算效率降低。而局部網(wǎng)格僅僅對(duì)異常體、當(dāng)前測(cè)點(diǎn)位置以及其下方的地空分界面進(jìn)行局部加密，對(duì)單個(gè)測(cè)點(diǎn)進(jìn)行單獨(dú)求解，因此所使用的網(wǎng)格數(shù)量大大降低，計(jì)算量減小。

2 并行加速

2.1 基于OpenMP的并行加速

OpenMP采用分叉—合并(Fork-Join)執(zhí)行模式。一個(gè)OpenMP程序開始于一個(gè)單一的線程，該進(jìn)程又稱作主線程，通過(guò)并行指令對(duì)程序的并行區(qū)間進(jìn)行定義，在這個(gè)區(qū)間中程序塊由多個(gè)線程自動(dòng)分配任務(wù)及并行執(zhí)行，線程數(shù)由主線程決定，在一個(gè)大型程序中可以嵌套多個(gè)OpenMP并行區(qū)間，其中并行區(qū)間由一組$OMP指令進(jìn)行開啟(!OMP PARALLEL)和關(guān)閉(!OMP END PARALLEL)。圖3是OpenMP并行模式的示意圖。

2.2 基于MPI的并行加速

MPI需要用頭文件use mpi進(jìn)行聲明，然后在主進(jìn)程(即節(jié)點(diǎn))中進(jìn)行初始化，開啟MPI并行環(huán)境，獲取進(jìn)程個(gè)數(shù)及編號(hào)，再將多個(gè)任務(wù)合理地分配到各個(gè)進(jìn)程并傳輸對(duì)應(yīng)的數(shù)據(jù)，然后分別進(jìn)行計(jì)算，最后匯集到主進(jìn)程并關(guān)閉并行環(huán)境，如圖4所示。

圖2 全局網(wǎng)格與局部網(wǎng)格對(duì)比 (a為全局網(wǎng)格；b、c、d為局部網(wǎng)格)Fig.2 Global grid vs. local grid ((a) is a Global grid; (b) (c) (d) is a local grid)

圖3 OpenMP并行模式示意Fig.3 OpenMP parallel mode schematic

2.3 并行程序設(shè)計(jì)及實(shí)現(xiàn)

由于航空電磁系統(tǒng)在每個(gè)觀測(cè)點(diǎn)處的靈敏區(qū)域相對(duì)于總計(jì)算區(qū)域來(lái)說(shuō)遠(yuǎn)遠(yuǎn)不及，因此僅僅在觀測(cè)點(diǎn)的一定位置范圍內(nèi)對(duì)網(wǎng)格進(jìn)行加密即可滿足計(jì)算的精度要求，針對(duì)每個(gè)觀測(cè)點(diǎn)分別單獨(dú)設(shè)計(jì)獨(dú)立的網(wǎng)格來(lái)提高計(jì)算效率。

時(shí)間域航空電磁三維正演各個(gè)測(cè)點(diǎn)計(jì)算相互獨(dú)立，不存在數(shù)據(jù)依賴關(guān)系，具有非常好的并行性，在三維正演的基礎(chǔ)上實(shí)現(xiàn)了基于MPI的多測(cè)點(diǎn)并行計(jì)算，在獨(dú)立網(wǎng)格內(nèi)實(shí)現(xiàn)了基于OpenMP的單元矩陣并行計(jì)算。程序采用主從模式，整個(gè)程序的基本結(jié)構(gòu)及運(yùn)行流程由主進(jìn)程把控，主要負(fù)責(zé)讀取模型參數(shù)、數(shù)據(jù)傳輸以及對(duì)并行任務(wù)的分配，子進(jìn)程負(fù)責(zé)接收來(lái)自主進(jìn)程傳遞的消息、對(duì)分配的任務(wù)進(jìn)行計(jì)算以及輸出本進(jìn)程的計(jì)算結(jié)果。為了充分利用計(jì)算資源提高計(jì)算效率，主進(jìn)程在完成任務(wù)分配和數(shù)據(jù)傳輸?shù)墓ぷ骱笠矃⑴c到測(cè)點(diǎn)的并行計(jì)算當(dāng)中。

具體實(shí)現(xiàn)的流程如圖5所示，主進(jìn)程讀取模型的整套網(wǎng)格參數(shù)，隨后將其傳遞給其他子進(jìn)程，之后各進(jìn)程根據(jù)分配的網(wǎng)格參數(shù)計(jì)算模型的響應(yīng)結(jié)果，為了防止數(shù)據(jù)傳輸時(shí)發(fā)生通信沖突，采用各進(jìn)程分別輸出各自的計(jì)算結(jié)果。

圖4 MPI并行模式示意Fig.4 MPI parallel mode schematic

圖5 航空瞬變電磁三維正演并行計(jì)算流程Fig.5 Airborne transient electromagnetic three-dimensional Forward modeling parallel computing flow chart

3 模型算例

3.1 精度驗(yàn)證

為檢驗(yàn)本文所開發(fā)的時(shí)間域航空電磁多粒度并行正演方法的正確性，設(shè)計(jì)了如圖6所示的水平板狀體模型，發(fā)射波形為階躍波，發(fā)射電流強(qiáng)度為435 A，通過(guò)全局網(wǎng)格和局部網(wǎng)格兩種剖分方式，對(duì)異常體上方y(tǒng)=0 m剖面上-400～400 m范圍內(nèi)的32個(gè)測(cè)點(diǎn)，分別進(jìn)行正演模擬，且將正演結(jié)果與有限體積方法的正演結(jié)果進(jìn)行比對(duì)，進(jìn)行數(shù)值精度驗(yàn)證。

圖6 水平板狀體模型示意Fig.6 Schematic diagram of the horizontal plate model

使用的集群服務(wù)器有16個(gè)型號(hào)為Intel(R) Xeon(R) CPU E5-2609 v4 @ 1.70 GHz的邏輯CPU，每個(gè)CPU擁有2G內(nèi)存和4個(gè)線程。采用時(shí)間域有限元全局網(wǎng)格法、基于多粒度并行加速策略的時(shí)間域有限元局部網(wǎng)格法以及有限體積方法分別計(jì)算的電磁響應(yīng)結(jié)果如圖7所示，可以清楚地看出，上述3種不同方法所計(jì)算的結(jié)果吻合得非常好，且均對(duì)模型中的地下板狀體有明顯反映。這一結(jié)果有效驗(yàn)證了本文并行算法的精度。

3.2 并行性能分析

為了分析基于MPI+OpenMP并行加速策略的時(shí)間域航空電磁三維正演程序的并行加速效果，采用上述水平塊狀體模型進(jìn)行測(cè)算，根據(jù)控制變量法，控制計(jì)算測(cè)點(diǎn)數(shù)為32和線程數(shù)為4保持不變，依次對(duì)不同數(shù)量的進(jìn)程進(jìn)行比較，通過(guò)加速比和并行效率兩個(gè)參數(shù)來(lái)評(píng)估本文設(shè)計(jì)的三維正演并行程序。

令原有的串行程序在集群服務(wù)器單進(jìn)程上的運(yùn)行時(shí)間為Ts，經(jīng)過(guò)多粒度并行優(yōu)化后，運(yùn)行所需時(shí)間為Tp，其中P表示開啟的進(jìn)程數(shù)，則加速比Sp可表示為：

Sp=Ts/Tp，

(20)

并行效率可定義為：

圖7 全局網(wǎng)格串行計(jì)算和局部網(wǎng)格并行計(jì)算電磁響應(yīng)結(jié)果對(duì)比Fig.7 Comparison of electromagnetic response results between global grid serial computing and local grid parallel computing

ep=Sp/P。

(21)

測(cè)試結(jié)果見表1。在計(jì)算測(cè)點(diǎn)數(shù)一定時(shí)，隨著進(jìn)程數(shù)的增加，計(jì)算時(shí)間呈非線性下降趨勢(shì)，加速比逐漸增大，相反并行效率逐漸降低；當(dāng)開啟16個(gè)進(jìn)程時(shí)，采用多粒度并行加速后執(zhí)行耗費(fèi)的時(shí)間僅僅是串行程序耗時(shí)的1/10，同時(shí)并行效率達(dá)到最低。

表1 基于MPI+OpenMP并行加速策略的行正演計(jì)算統(tǒng)計(jì)Table 1 Line forward calculation statistics based on MPI+OpenMP parallel acceleration strategy

由圖8可以發(fā)現(xiàn)以下特點(diǎn)：并行化后的正演程序加速比與開啟的進(jìn)程數(shù)并不是整數(shù)倍關(guān)系，而是呈非線性增長(zhǎng)，且并行效率整體呈降低趨勢(shì)。這是因?yàn)樵诩悍?wù)器上開啟并行環(huán)境后，會(huì)占用一定內(nèi)存，且進(jìn)程間進(jìn)行數(shù)據(jù)傳遞需要占用時(shí)間，而且任務(wù)量與開啟的進(jìn)程數(shù)不匹配時(shí)，計(jì)算任務(wù)少的進(jìn)程會(huì)率先完成計(jì)算量，但該進(jìn)程并不會(huì)結(jié)束，而是要等待未完成任務(wù)的進(jìn)程，從而導(dǎo)致時(shí)間的損耗。所以開啟P個(gè)進(jìn)程，加速比Sp并不能達(dá)到P，且隨著進(jìn)程數(shù)的增加Sp與P的差值越來(lái)越大。

圖8 不同進(jìn)程數(shù)的加速比(a)和并行效率(b)Fig.8 Acceleration ratio (a) and parallel efficiency (b) for different node numbers

3.3 復(fù)雜模型計(jì)算

考慮到地形起伏不平對(duì)航空電磁正演結(jié)果的影響較為嚴(yán)重，采取四面體網(wǎng)格進(jìn)行加密剖分后的計(jì)算量較大，從而導(dǎo)致計(jì)算時(shí)間較長(zhǎng)，因此利用多粒度并行加速技術(shù)進(jìn)行加速計(jì)算。導(dǎo)入地形文件數(shù)據(jù)模擬起伏地表，并設(shè)計(jì)地下埋藏有塊狀及傾斜板狀良導(dǎo)體(圖9)。設(shè)高阻圍巖的電導(dǎo)率為0.01 S/m，塊狀體和板狀體的電導(dǎo)率為1 S/m；傾斜版狀體頂部埋深60 m，垂直深度250 m，塊狀體頂部埋深50 m，邊長(zhǎng)200 m。采用中心回線裝置，飛行高度30 m，發(fā)射線框半徑15 m，線圈匝數(shù)為1，采用階躍波激發(fā)。全區(qū)設(shè)置11條測(cè)線，每條測(cè)線布設(shè)171個(gè)測(cè)點(diǎn)，共計(jì)1 881個(gè)測(cè)點(diǎn)，每個(gè)網(wǎng)格單元數(shù)約90 000個(gè)，正演過(guò)程共需進(jìn)行9 350次矩陣分解，327 250次回代。若采用傳統(tǒng)串行計(jì)算時(shí)，總耗時(shí)約1 968.623 min(32.8 h)，采用并行加速開啟16個(gè)進(jìn)程后總計(jì)算時(shí)間約194.336 min(3.24 h)，消耗內(nèi)存約5.2 G，效率提高了10.13倍,如表2所示。

表2 計(jì)算情況統(tǒng)計(jì)Table 2 Calculation statistics table

圖10是主剖面(x=0測(cè)線)的多測(cè)道圖，可以看出在早期表現(xiàn)出地形的響應(yīng)，到了晚期，地形影響逐漸變小，表現(xiàn)出異常體的響應(yīng)。圖11呈現(xiàn)出整個(gè)測(cè)區(qū)在四個(gè)時(shí)間點(diǎn)時(shí)的航空電磁響應(yīng)結(jié)果，從圖中可以明顯看出在時(shí)間早期，僅存在地形的響應(yīng)，由于地形凸起，電磁響應(yīng)呈現(xiàn)出相對(duì)低異常，在x=0，y=-500附近，相對(duì)低異常達(dá)到極大值，由于塊狀異常體處于凸起地形下方，相對(duì)于傾斜板埋藏較深，故隨著時(shí)間延長(zhǎng)，傾斜板的異常響應(yīng)率先顯示出來(lái)，隨后出現(xiàn)塊狀體異常響應(yīng)，且由于地下埋藏有異常體，導(dǎo)致地形的影響減弱，到晚期時(shí)，電磁波穿過(guò)異常體，僅剩地形的電磁響應(yīng)。

圖9 復(fù)雜模型計(jì)算Fig.9 Calculation of complex models

圖10 復(fù)雜模型主剖面(x=0測(cè)線)多測(cè)道曲線Fig.10 Multiple trajectory map of complex model main section (x=0 line)

4 結(jié)論

利用MPI實(shí)現(xiàn)了三維正演模擬中的多測(cè)點(diǎn)并行計(jì)算和實(shí)現(xiàn)了獨(dú)立網(wǎng)格內(nèi)單元矩陣計(jì)算的OpenMP并行加速，最終實(shí)現(xiàn)了基于MPI+OpenMP并行加速策略對(duì)時(shí)間域航空電磁三維正演方法的并行化，并得出并行化后的正演程序，其加速比隨著開啟的進(jìn)程 數(shù)增加呈非線性增長(zhǎng)。

圖11 復(fù)雜地形響應(yīng)曲面Fig.11 Complex terrain anomaly response surface map

由于“并行開銷”的存在，開啟的進(jìn)程數(shù)越多，并行效率整體呈現(xiàn)出下降趨勢(shì)，即使只開啟一個(gè)進(jìn)程，效率也比串行效率低。地形起伏模型的數(shù)值模擬體現(xiàn)了本文并行優(yōu)化后的正演方法的高效性，正確高效的三維正演為三維反演提供了可能。