基于CPU/GPU異構(gòu)平臺的疊前逆時偏移成像系統(tǒng)

高新成 (東北石油大學(xué)現(xiàn)代教育技術(shù)中心，黑龍江 大慶 163318)石穎 (東北石油大學(xué)地球科學(xué)學(xué)院，黑龍江 大慶 163318)

?

基于CPU/GPU異構(gòu)平臺的疊前逆時偏移成像系統(tǒng)

高新成 (東北石油大學(xué)現(xiàn)代教育技術(shù)中心，黑龍江 大慶 163318)石穎 (東北石油大學(xué)地球科學(xué)學(xué)院，黑龍江 大慶 163318)

針對目前地震資料偏移成像過程中代碼復(fù)雜分散、可移植性差和可視化操作弱等問題，結(jié)合疊前逆時偏移業(yè)務(wù)需求，研究了CPU/GPU協(xié)同并行計算技術(shù)優(yōu)化疊前逆時偏移成像算法，設(shè)計了逆時偏移系統(tǒng)流程與功能，開發(fā)了一套基于CPU/GPU異構(gòu)平臺的疊前逆時偏移成像系統(tǒng)。利用自設(shè)模型和傳統(tǒng)模型數(shù)據(jù)對系統(tǒng)進行測試，測試結(jié)果表明系統(tǒng)操作方便可行，運行效果良好，能夠有效的對復(fù)雜構(gòu)造地質(zhì)體精確成像。

逆時偏移；CPU/GPU并行加速；PML邊界條件；QT

隨著地震資料偏移成像應(yīng)用推廣，在逆時偏移成像處理中積累了很多的應(yīng)用程序和方法，如何將這些程序和方法集成整合起來，實現(xiàn)偏移成像的可視化操作和高效管理成為當前研究的重要課題[1]。目前疊前逆時偏移是獲得復(fù)雜構(gòu)造精確成像最有效的方法[2]。為此，筆者結(jié)合疊前逆時偏移業(yè)務(wù)需求研究地震波場正演模擬、成像條件、邊界吸收條件和噪音壓制算法等，形成了一套疊前逆時偏移成像方法；并研究了CPU/GPU異構(gòu)平臺并行加速技術(shù)，將其應(yīng)用于疊前逆時偏移成像方法中地震波場正演計算、成像條件和波場反傳計算等計算環(huán)節(jié)，開發(fā)了一套基于CPU/GPU異構(gòu)平臺的疊前逆時偏移成像系統(tǒng)，為油田開發(fā)資料的解釋提供技術(shù)支持。

1 CPU/GPU異構(gòu)平臺的逆時偏移流程

1.1 疊前逆時偏移成像算法

筆者采用逆時偏移算法，在波場正演過程中采用隨機邊界條件和PML完全吸收邊界條件相結(jié)合的2種邊界處理方式，成像過程中采用互相關(guān)成像條件，成像后采用相應(yīng)濾波法進行去噪[3～7]。逆時偏移算法實現(xiàn)流程如圖1所示。算法實現(xiàn)步驟如下：

步1 使用波動方程有限差分法沿時間正向延拓波場計算到Tmax時刻，保存最后2個時刻波場值。

步2 采用PML邊界條件，使用波動方程有限差分法沿時間正向延拓波場計算到Tmax時刻，保存每個時刻的檢波點數(shù)據(jù)。

步3 讀取最后2個時刻波場值，使用波動方程有限差分法沿時間反向延拓計算到0時刻，計算每個時刻波場值。

步4 讀取檢波點數(shù)據(jù)，使用波動方程有限差分法沿時間反向延拓計算到0時刻，計算每個時刻波場值。

步5 讀取步3和步4同一時刻2個波場值，利用互相關(guān)成像條件進行成像，從Tmax到0時刻。

步6 循環(huán)計算多炮數(shù)據(jù)到Nshot，完成多炮成像結(jié)果疊加。

步7 采用濾波法對成像結(jié)果進行壓制低頻噪聲。

步8 輸出成像結(jié)果。

1.2 CPU/GPU協(xié)同加速逆時偏移流程

將CUDA編程思想引入到逆時偏移波場延拓和相關(guān)成像算法中，采用基于CUDA架構(gòu)的CPU/GPU協(xié)同工作模式來處理逆時偏移計算任務(wù)。其主要思想是GPU端負責完成波場延拓和相關(guān)成像等并行計算任務(wù)，CPU端負責數(shù)據(jù)輸入和輸出等相關(guān)操作控制工作。

圖2 CPU/GPU協(xié)同加速逆時偏移流程

首先，將激發(fā)點初始波場值和檢波點最大時刻波場值由內(nèi)部存儲器(內(nèi)存)傳至設(shè)備存儲器(顯存)中[8]。其次，設(shè)置GPU計算維度，將GPU多核處理器劃分為相應(yīng)個數(shù)的計算塊Block；同時將每個Block劃分為若干個計算線程Thread；開始執(zhí)行炮點波場正傳和檢波點波場反傳計算任務(wù)。根據(jù)波場計算相關(guān)成像條件，需要同一時刻正向傳播的震源波場和反向傳播的檢波點波場，因此互相關(guān)成像計算也是由GPU來完成[9，10]，具體實現(xiàn)流程如圖2所示。

2 系統(tǒng)流程與功能設(shè)計

2.1 系統(tǒng)設(shè)計流程

根據(jù)逆時偏移處理流程，將逆時偏移系統(tǒng)數(shù)據(jù)處理分為工區(qū)創(chuàng)建、數(shù)據(jù)加載、參數(shù)設(shè)置、作業(yè)發(fā)送和任務(wù)執(zhí)行等多個步驟，并結(jié)合CPU/GPU協(xié)同并行計算模式[11]設(shè)計了逆時偏移數(shù)據(jù)處理系統(tǒng)流程(見圖3)。

圖3 逆時偏移系統(tǒng)處理流程

2.2 系統(tǒng)功能設(shè)計

根據(jù)逆時偏移的實際業(yè)務(wù)需求，設(shè)計并開發(fā)了疊前逆時偏移成像系統(tǒng)，主要完成二維和三維疊前逆時深度偏移成像功能。系統(tǒng)主要包括資源管理、參數(shù)設(shè)置、數(shù)據(jù)處理、任務(wù)管理、偏移成像和用戶管理等6個功能模塊(見圖4)。

1)資源管理。對所有資源節(jié)點進行實時有效地管理，能夠?qū)崟r監(jiān)控各個計算節(jié)點的運行狀態(tài)，獲取節(jié)點的性能參數(shù)，完成計算節(jié)點的獲取、添加和刪除等操作。

2)參數(shù)設(shè)置。完成逆時偏移成像前所需要的各類參數(shù)設(shè)置，包括工區(qū)創(chuàng)建、偏移成像前基本參數(shù)設(shè)置和高級參數(shù)設(shè)置。

圖4 疊前逆時偏移系統(tǒng)功能模塊圖

3)數(shù)據(jù)處理。完成偏移成像前數(shù)據(jù)的預(yù)處理，包括導(dǎo)入原始數(shù)據(jù)和導(dǎo)入速度模型等操作。

4)偏移成像。包括地震波場正演計算、檢波點波場反傳計算和互相關(guān)成像等操作。包括二維數(shù)據(jù)疊前逆時深度偏移成像和三維數(shù)據(jù)疊前逆時深度偏移成像功能，其中包括2D CPU RTM、2D GPU RTM、3D CPU RTM和3D GPU RTM。

5)任務(wù)管理。對逆時偏移成像任務(wù)的執(zhí)行情況進行實時監(jiān)測和管理，能夠監(jiān)管整體任務(wù)和各個計算節(jié)點任務(wù)的執(zhí)行情況。包括任務(wù)添加、任務(wù)刪除、任務(wù)重啟、任務(wù)中斷以及查看任務(wù)執(zhí)行時間等操作，同時，也能實時監(jiān)測出每個節(jié)點運行性能參數(shù)。

6)用戶管理。主要包括用戶管理、權(quán)限設(shè)置和日志查詢等功能。用戶管理包括用戶的添加、刪除和修改等操作，能夠根據(jù)用戶使用系統(tǒng)功能分配相應(yīng)權(quán)限，根據(jù)權(quán)限不同，用戶登錄后顯示的操作界面也不同，并能夠記錄用戶操作日志。

3 系統(tǒng)主要功能實現(xiàn)

根據(jù)疊前逆時偏移數(shù)據(jù)處理流程與實現(xiàn)算法，通過研究可視化平臺開發(fā)技術(shù)和疊前逆時偏移方法，在Linux環(huán)境下應(yīng)用QT圖形界面開發(fā)平臺[12]，設(shè)計并實現(xiàn)了一套疊前逆時偏移成像系統(tǒng)。針對逆時偏移實現(xiàn)過程中二維CPU、二維GPU、三維CPU和三維GPU偏移算法和處理方法不相同，分別設(shè)計了二維 CPU逆時偏移、二維GPU逆時偏移、三維 CPU逆時偏移和三維 GPU逆時偏移等功能。

3.1 疊前逆時偏移成像實現(xiàn)

疊前逆時偏移成像主要包括2D CPU逆時偏移、2D GPU逆時偏移、3D CPU逆時偏移和3D GPU逆時偏移4個模塊(見圖5)。筆者以2D GPU逆時偏移模塊為例描述系統(tǒng)的實現(xiàn)過程，其他模塊操作類似：首先創(chuàng)建工區(qū)，導(dǎo)入炮集數(shù)據(jù)和速度模型文件，設(shè)置偏移成像的相關(guān)參數(shù)；選擇2D GPU RTM處理方式和CPU計算節(jié)點，完成整體偏移成像。2D GPU逆時偏移具體操作如圖6所示。

圖5 2D CPU逆時偏移模塊

圖6 2D CPU逆時偏移成像

偏移成像界面中可以通過Edit按鈕修改參數(shù)，Inport按鈕載入?yún)?shù)卡文件，參數(shù)設(shè)定后，通過Export輸出參數(shù)卡文件。所有參數(shù)設(shè)置完成后，點擊Start按鈕發(fā)送作業(yè)開始執(zhí)行逆時偏移計算，作業(yè)執(zhí)行過程中點擊Refresh按鈕可以更新進度條，查看作業(yè)執(zhí)行進度。

3.2 多節(jié)點任務(wù)管理的實現(xiàn)

圖7 多節(jié)點作業(yè)流程圖

系統(tǒng)設(shè)計了多節(jié)點間任務(wù)監(jiān)控功能，支持異構(gòu)節(jié)點間作業(yè)的發(fā)送與終止。通過對節(jié)點文件信息的讀取，識別節(jié)點設(shè)備名稱；根據(jù)節(jié)點名與登錄用戶名信息生成唯一的shell文件，作業(yè)發(fā)送或終止命令發(fā)出時，各節(jié)點執(zhí)行自身匹配的shell文件完成作業(yè)操作。多節(jié)點作業(yè)發(fā)送流程如圖7所示。

在完成多節(jié)點作業(yè)發(fā)送功能基礎(chǔ)上，添加了多節(jié)點作業(yè)進度監(jiān)控功能。在多節(jié)點任務(wù)執(zhí)行過程中，作業(yè)流會定期更新進度文件并保存在工區(qū)指定目錄中。通過定期發(fā)送更新指令，系統(tǒng)從進度文件中讀取作業(yè)進度信息，完成進度顯示模塊的更新。

4 系統(tǒng)應(yīng)用與效果分析

4.1 系統(tǒng)應(yīng)用可行性分析

為驗證系統(tǒng)的可行性及算法集成后的有效性，試驗中選取一個斷層模型的200炮數(shù)據(jù)完成正演計算和偏移成像試驗。試驗節(jié)點為1個Intel E5-2630六核CPU處理器和Nvidia Tesla K10 GPU卡，試驗具體參數(shù)如表1所示。

表1 系統(tǒng)軟件測試

逆時偏移過程中采用時間二階、空間十階的有限差分數(shù)值算法完成地震波場正傳和檢波點波場反傳計算，采用震源歸一化互相關(guān)成像條件提取成像點，并采用拉普拉斯算子濾波法對圖像進行去噪處理，得到最終成像效果。試驗中的速度模型、逆時偏移成像后剖面和采用拉普拉斯算法去噪后剖面分別見圖8(a)、(b)和(c)所示。

圖8 逆時偏移系統(tǒng)測試

通過試驗成像結(jié)果發(fā)現(xiàn)，系統(tǒng)成像結(jié)果與源代碼方式成像結(jié)果基本一樣，沒有明顯區(qū)別。同時，試驗中源代碼計算與系統(tǒng)計算的執(zhí)行時間分別為134.28s和133.46s，系統(tǒng)執(zhí)行時間略少?？梢?，系統(tǒng)完成逆時偏移從成像效果和執(zhí)行時間等方面具有可行性，成像效果良好。

4.2 系統(tǒng)2D RTM成像效果分析

表2 試驗參數(shù)

為驗證系統(tǒng)的成像效果，選取二維鹽丘(Salt-2D)模型進行疊前逆時偏移試驗。試驗環(huán)境為4個計算節(jié)點，每個節(jié)點配備1個Intel E5-2630六核CPU處理器和Nvidia Tesla K10 GPU卡，試驗參數(shù)如表2所示。

在Salt-2D模型試驗中，采用時間二階、空間十階的有限差分數(shù)值算法完成波場數(shù)值計算，采用互相關(guān)成像條件提取成像點，成像后采用拉普拉斯算子濾波法去噪。速度模型如圖9所示，逆時偏移成像結(jié)果如圖10所示。

從試驗結(jié)果可以看出，筆者設(shè)計的逆時偏移系統(tǒng)能夠有效的對地下結(jié)構(gòu)精確成像，深層構(gòu)造成像清楚，整體趨勢與速度模型數(shù)據(jù)基本一致，地層分界明顯，微幅構(gòu)造成像清晰。

5 結(jié)語

筆者主要圍繞逆時偏移系統(tǒng)設(shè)計與應(yīng)用展開研究，通過研究逆時偏移成像算法和CPU/GPU并行加速技術(shù)，設(shè)計逆時偏移系統(tǒng)數(shù)據(jù)處理流程，并根據(jù)CPU和GPU程序代碼、處理方法不同設(shè)計了2D CPU逆時偏移、2D GPU逆時偏移、3D CPU逆時偏移和3D GPU逆時偏移等功能模塊，利用QT環(huán)境開發(fā)了基于CPU/GPU異構(gòu)平臺的疊前逆時偏移成像系統(tǒng)。試驗結(jié)果表明系統(tǒng)具有可行性，與源代碼方式相比成像結(jié)果一致，能夠?qū)?fù)雜地下構(gòu)造獲得良好的成像效果。同時，利用GPU多線程實現(xiàn)了大規(guī)模并行計算，大大提高了偏移計算效率。

圖9 Salt-2D速度模型

圖10 Salt-2D逆時偏移結(jié)果

[1]陳可洋. 地震波逆時偏移方法研究綜述[J]. 勘探地球物理進展，2010，33(3)：153～159.

[2] Liu Y，Chang X，Jin D，et al. Reverse time migration of multiples for subsalt imaging[J]. Geophysics，2011，76(5)：209～216.

[3] Costa J C，Silva Neto F A.Obliquity-correction imaging condition for reverse time migration[J]. Geophysics，2009，74(3)：57～66.

[4] 高新成，李春生.基于CPU/GPU協(xié)同加速疊前逆時偏移方法研究[J].陜西理工學(xué)院學(xué)報(自然科學(xué)版)，2016(2)：44～49.

[5] 石穎，柯璇，張瑩瑩. 逆時偏移邊界條件與存儲策略分析[J].地球物理學(xué)進展，2015(2)：581～585.

[6] Clapp R G. Reverse time migration with random boundaries[J]. 79th Annual International Meeting，SEG Expanded Abstracts，2009，28：2809～2813.

[7]Jiaze He，F(xiàn)uh-Gwo Yuan.Damage identification for composite structures using a cross-correlation reverse-time migration technique[J]. Structural Health Monitoring，2015， 146：224～226.

[8] Ying Shi，Yanghua Wang. Reverse time migration of 3D vertical seismic profile data[J]. Geophysics，2016，81(1)：31～38.

[9] 石穎，陸加敏，柯璇，等. 基于GPU并行加速的疊前逆時偏移方法[J].東北石油大學(xué)學(xué)報，2012，36(4)：111～115.

[10] 劉守偉，王華忠，陳生昌，等. 三維逆時偏移GPU/CPU機群實現(xiàn)方案研究[J]. 地球物理學(xué)報，2013，56(10)：3487～3496.

[11] 王維紅，柯璇，郭雪豹，等. 基于QT的Walkaway VSP逆時偏移軟件開發(fā)研究[J]. 石油物探，2015，54(4)：452～458.

[12] 陳楠，祝媛媛，張光德，等. 基于QT的地震勘探可擴展平臺研發(fā)與應(yīng)用[J]. 華北地震科學(xué)，2013，32(4)：31～36.

[編輯] 辛長靜

2016-09-18

國家自然科學(xué)基金項目(41574117)；黑龍江省教育科學(xué)規(guī)劃重點課題(GJB1215013)

高新成(1979-)，男，博士，副教授，現(xiàn)主要從事地震數(shù)據(jù)處理與智能計算方面的教學(xué)與研究工作；E-mail：gxc@nepu.edu.cn。

TP399

A

1673-1409(2016)34-0015-06

[引著格式]高新成，石穎.基于CPU/GPU異構(gòu)平臺的疊前逆時偏移成像系統(tǒng)[J].長江大學(xué)學(xué)報(自科版)，2016,13(34)：15～20.