郭雪豹，王建民，王維紅，王云專(zhuān)，柯 璇，劉詩(shī)竹

（1.東北石油大學(xué) 地球科學(xué)學(xué)院，黑龍江 大慶 163318； 2.大慶油田有限責(zé)任公司 勘探開(kāi)發(fā)研究院，黑龍江 大慶163712）

0 引言

隨著我國(guó)各大油田勘探的不斷深入，探區(qū)內(nèi)有利儲(chǔ)藏大多已勘探完畢，而剩余的油氣藏多具有油水分布復(fù)雜、薄儲(chǔ)層變化快等特點(diǎn).對(duì)于日趨復(fù)雜多樣化的油氣藏類(lèi)型，常規(guī)的偏移方法已無(wú)法對(duì)地下深層次的復(fù)雜構(gòu)造準(zhǔn)確成像.基于波動(dòng)理論的逆時(shí)偏移與同樣基于波動(dòng)理論的克希霍夫偏移和單程波偏移相比，波動(dòng)方程不存在任何假設(shè)，既不會(huì)受到速度橫向劇烈變化的影響，也不存在成像傾角的限制，理論上能夠?qū)Χ啻尾?、回轉(zhuǎn)波等準(zhǔn)確成像，是目前成像精度最高的方法.

逆時(shí)偏移方法最早是在1983年被提出的，Baysal E等［1］提出不同的逆時(shí)偏移概念；Whitmore N D、Mc Mechan G A、Loeuenthal D等［2-4］將它應(yīng)用于疊后偏移.隨著油氣資源的復(fù)雜多樣化，常規(guī)的地面地震逆時(shí)偏移方法已經(jīng)無(wú)法滿(mǎn)足精細(xì)化勘探的需求，尤其對(duì)地下深層一些微幅構(gòu)造及井旁的小斷層等信息.VSP技術(shù)是與地面地震觀(guān)測(cè)方法相對(duì)應(yīng)的，它將檢波器置于井中垂向布置.因此，比地面檢波器更加接近目的層，減少地下淺層對(duì)地震信號(hào)的干擾，直接接收來(lái)自目的層的反射，同樣還能夠接收到一些地面檢波器所接收不到的陡傾角信息，具有比地面地震記錄更高的分辨率.郭建［5］分析VSP技術(shù)的應(yīng)用現(xiàn)狀.朱金明等［6］通過(guò)按時(shí)間逆的次序求解雙程無(wú)反射波動(dòng)方程，實(shí)現(xiàn)VSP逆時(shí)偏移.孫文博等［7］采用偽譜法實(shí)現(xiàn)逆時(shí)偏移方法，分別對(duì)分離波場(chǎng)偏移處理，得到較好的效果.Alok K S等［8］采用包含多次波在內(nèi)的全波場(chǎng)成像，提高照明范圍.Sun Wenbo等［9］在角度域?qū)崿F(xiàn)VSP逆時(shí)偏移.Mark E W等［10］利用walkaway VSP資料在鹽丘側(cè)翼上取得較好的成像效果.VSP資料在對(duì)地下深層構(gòu)造成像方面具有優(yōu)勢(shì)，并且其內(nèi)部豐富的波場(chǎng)信息也有助于更高精度成像.

采用VSP數(shù)據(jù)進(jìn)行逆時(shí)偏移，利用VSP資料的高分辨率、高信噪比的特點(diǎn)，實(shí)現(xiàn)地下微幅構(gòu)造的準(zhǔn)確成像.根據(jù)設(shè)計(jì)的VSP觀(guān)測(cè)系統(tǒng)，筆者正演VSP地震記錄，利用VSP資料的豐富波場(chǎng)信息進(jìn)行逆時(shí)偏移，以彌補(bǔ)常規(guī)地面地震在深層構(gòu)造及井旁小構(gòu)造方面的成像不足；逆時(shí)偏移算法中波場(chǎng)計(jì)算部分由GPU執(zhí)行，其余邏輯判斷部分由CPU執(zhí)行，充分挖掘兩者的計(jì)算特性來(lái)獲得最高加速比.算法經(jīng)GPU加速后，計(jì)算效率得到了顯著提升.

1 VSP觀(guān)測(cè)系統(tǒng)

VSP觀(guān)測(cè)系統(tǒng)是一種由地表震源激發(fā)，井中檢波器接收的觀(guān)測(cè)方式（見(jiàn)圖1）.首先，與傳統(tǒng)地面檢波器接收到的反射信號(hào)路徑相比，VSP檢波器接收的反射波路徑明顯更短，避免了反射波回到地表時(shí)所要經(jīng)過(guò)的淺層部分，減少地表低速帶對(duì)反射波信息中高頻成分的吸收，振幅畸變相對(duì)更小.其次，由于檢波器布置在地下，與地表的檢波器相比更加貼近目的層，因此能夠直接接收到來(lái)自目的層的反射信息，有助于研究來(lái)自目的層的單一反射信息；在地下存在陡傾角構(gòu)造時(shí)，地表檢波器范圍有限，可能損失部分反射信息，而VSP的檢波器則能接收到.根據(jù)波場(chǎng)信息，VSP波場(chǎng)較地面地震波場(chǎng)更為豐富，不僅存在上行波（圖中以上行一次反射波為例），還有下行波（圖中以下行直達(dá)波為例），將有利于逆時(shí)偏移方法對(duì)VSP波場(chǎng)信息的充分利用.因此，與地面地震記錄相比，VSP數(shù)據(jù)的分辨率和信噪比更高，將給地下深層成像中帶來(lái)明顯優(yōu)勢(shì).

2 逆時(shí)偏移原理

逆時(shí)偏移過(guò)程可以分為震源波場(chǎng)正傳和檢波點(diǎn)處波場(chǎng)反傳.首先進(jìn)行震源波場(chǎng)正傳，在正傳過(guò)程中為削弱用計(jì)算機(jī)有限存儲(chǔ)模擬地下無(wú)限介質(zhì)的波場(chǎng)傳播所帶來(lái)的邊界反射，需要添加邊界條件，削弱邊界反射對(duì)波場(chǎng)的影響，并且保存每一時(shí)刻的波場(chǎng)信息；然后進(jìn)行檢波點(diǎn)處波場(chǎng)反傳，與對(duì)應(yīng)時(shí)刻的震源波場(chǎng)應(yīng)用成像條件得到逆時(shí)偏移結(jié)果，對(duì)逆時(shí)偏移結(jié)果去噪處理得到最后的成像剖面.

二維聲波方程：

式中：p為波場(chǎng)值；v為介質(zhì)速度.

通過(guò)有限差分法對(duì)波場(chǎng)進(jìn)行數(shù)值模擬，時(shí)間二階、空間2L階精度的有限差分格式［12］為

式中：Δx、Δz分別為x、z方向的網(wǎng)格大小；Δt為時(shí)間步長(zhǎng)；i、j分別為x、z方向網(wǎng)格點(diǎn)的位置；n為時(shí)間點(diǎn)；al為有限差分系數(shù).

3 GPU加速技術(shù)

盡管逆時(shí)偏移方法成像精度極高，但需解決計(jì)算量大和存儲(chǔ)需求大的問(wèn)題.GPU即圖形處理單元，經(jīng)NVIDIA公司發(fā)布CUDA后，使它可以被開(kāi)發(fā)人員編寫(xiě)指令程序，便于讓GPU執(zhí)行CPU所分布的計(jì)算任務(wù).GPU具有比CPU更多的計(jì)算核心，更易處理大量的并行計(jì)算任務(wù)，對(duì)于逆時(shí)偏移過(guò)程中的大量波場(chǎng)值計(jì)算更加具有優(yōu)勢(shì)［11］，因此可以顯著的提高算法的計(jì)算效率.

在逆時(shí)偏移過(guò)程中，計(jì)算量主要體現(xiàn)在模擬波場(chǎng)傳播計(jì)算中.對(duì)于一個(gè)二維的剖面，將它按縱橫向進(jìn)行離散化，得到用網(wǎng)格點(diǎn)表示的剖面.在偏移過(guò)程中，不斷利用式（2）計(jì)算每一個(gè)網(wǎng)格點(diǎn)上不同時(shí)刻的波場(chǎng)值.對(duì)于同一時(shí)刻的波場(chǎng)剖面，每個(gè)網(wǎng)格點(diǎn)上的波場(chǎng)值計(jì)算并無(wú)先后順序，即同一時(shí)刻的所有網(wǎng)格點(diǎn)的計(jì)算是并行的.首先在GPU上開(kāi)辟同一時(shí)刻所需要的線(xiàn)程數(shù)，然后將每一個(gè)網(wǎng)格點(diǎn)上的波場(chǎng)值計(jì)算放到對(duì)應(yīng)的線(xiàn)程中，GPU上的計(jì)算核心多（計(jì)算核心即流多處理器，其本身包含標(biāo)量流處理器，每個(gè)線(xiàn)程都是被發(fā)射到一個(gè)標(biāo)量流處理器上執(zhí)行），因此比CPU更加適合并行計(jì)算［13-15］.

由式（2）可知，在計(jì)算一點(diǎn)的波場(chǎng)值時(shí)需要多個(gè)點(diǎn)的波場(chǎng)值參與運(yùn)算，因此需要在GPU上開(kāi)辟部分存儲(chǔ)以存放波場(chǎng)值.在GPU的存儲(chǔ)器中，全局存儲(chǔ)器是存儲(chǔ)空間最大的，但訪(fǎng)問(wèn)速度遠(yuǎn)不如其內(nèi)部的共享存儲(chǔ)器.由于共享存儲(chǔ)器存儲(chǔ)有限，首先將計(jì)算的數(shù)據(jù)先讀入全局存儲(chǔ)器；然后在每個(gè)線(xiàn)程運(yùn)算時(shí)，將所需要的數(shù)據(jù)再?gòu)娜执鎯?chǔ)器復(fù)制到共享存儲(chǔ)器中，在計(jì)算時(shí)可以直接訪(fǎng)問(wèn)共享存儲(chǔ)器，有效加快計(jì)算式的訪(fǎng)問(wèn)速度，從而進(jìn)一步提高計(jì)算效率.GPU加速計(jì)算流程（見(jiàn)圖2）.

4 模型測(cè)試

4.1 鹽丘模型

采用鹽丘模型測(cè)試算法（見(jiàn)圖3），網(wǎng)格大小及時(shí)間步長(zhǎng)根據(jù)穩(wěn)定性條件選取.模型橫、縱向網(wǎng)格點(diǎn)數(shù)分別為800、600，橫、縱向網(wǎng)格大小分別為5、5m.采用地面放炮、井中接收的觀(guān)測(cè)方式，震源從地表左端50m處開(kāi)始，向右每隔50m放一炮，共80炮；在分別距模型左端0、4 000m處布置2口井，井中每隔5m布置一個(gè)檢波器，每口井放置600個(gè)；震源處雷克子波頻率為40Hz，時(shí)間步長(zhǎng)為0.5ms，時(shí)間采樣點(diǎn)為7 500.

在數(shù)值模擬過(guò)程中，采用保存部分波場(chǎng)信息的完全匹配層邊界條件［16-19］，吸收由計(jì)算機(jī)有限存儲(chǔ)所帶來(lái)的邊界反射衰減，成像條件為互相關(guān)成像條件［20］，對(duì)逆時(shí)偏移結(jié)果采用拉普拉斯去噪［21-23］處理.鹽丘模型運(yùn)算采用GPU加速節(jié)點(diǎn)運(yùn)行，硬件配置：CPU：E5-2630（2.3GHz／15M／6核心）＊2；內(nèi)存：64 GB；系統(tǒng)硬盤(pán)：1TB；GPU：Nvidia Tesla K10＊3（3 072核心／4.57T單精度／0.19T雙精度／8G顯存）.

VSP逆時(shí)偏移的成像剖面見(jiàn)圖4.由圖3和圖4可知，地下構(gòu)造基本成像清晰，淺層部分成像略顯不足，井旁地層成像準(zhǔn)確且清晰，鹽丘輪廓清楚，證實(shí)VSP逆時(shí)偏移對(duì)井旁高分辨率成像的優(yōu)勢(shì)，能夠有效利用來(lái)自地下深層的反射信息.

4.2 Marmous模型

利用Marmous模型進(jìn)行測(cè)試，網(wǎng)格大小及時(shí)間采樣由模型發(fā)布方參數(shù)決定.模型縱、橫向網(wǎng)格點(diǎn)數(shù)分別為750、993，縱、橫向網(wǎng)格大小分別為6.25、4.00m（見(jiàn)圖5）；采用地面放炮、井中接收的觀(guān)測(cè)方式，震源從地表左端62.50m處開(kāi)始，向右每隔62.50m放一炮，共99炮；在分別位于模型的62.50、3 125.00、6 187.50m處3口井，井中每隔4.00m布置一個(gè)檢波器，每口井放置750個(gè)；震源處雷克子波頻率為40Hz，時(shí)間步長(zhǎng)為0.4ms，時(shí)間采樣點(diǎn)為10 000.

VSP逆時(shí)偏移的成像剖面見(jiàn)圖6，地面地震的偏移結(jié)果見(jiàn)圖7.由圖5可知，Marmous模型地下地層較多，地勢(shì)起伏大，給成像帶來(lái)一定難度.由圖7可知，地面地震逆時(shí)偏移可以較清晰地對(duì)它進(jìn)行成像，尤其在淺層部分，地層清晰連續(xù)，但在2 000m以下的部位開(kāi)始模糊.在圖6中2 000m以下的部位，構(gòu)造形態(tài)基本正確，成像清晰，突出VSP逆時(shí)偏移對(duì)于深層構(gòu)造的成像優(yōu)勢(shì).在計(jì)算過(guò)程中，采用CPU偏移一炮需要2 844s，在相同參數(shù)下采用GPU偏移一炮僅需36 s，有效提高偏移的速度70倍以上.

5 結(jié)論

（1）逆時(shí)偏移方法能夠充分利用VSP數(shù)據(jù)的豐富波場(chǎng)信息及高分辨率、高信噪比的優(yōu)勢(shì)，在井旁及深層構(gòu)造成像方面具有比地面地震逆時(shí)偏移更好的成像效果.但由于對(duì)淺層反射信號(hào)接收不足，在淺層部分成像效果不如地面地震逆時(shí)偏移成像清晰.

（2）逆時(shí)偏移算法對(duì)計(jì)算機(jī)的計(jì)算能力要求較高，充分利用GPU上的共享存儲(chǔ)器，使其達(dá)到顯著的加速比，增強(qiáng)算法的實(shí)用性.

（3）二維VSP逆時(shí)偏移在一定程度上彌補(bǔ)地面地震逆時(shí)偏移對(duì)深層成像不足的問(wèn)題，但隨著勘探目標(biāo)的日趨復(fù)雜，二維VSP逆時(shí)偏移也難以達(dá)到復(fù)雜構(gòu)造精細(xì)勘探的需求，開(kāi)展三維VSP逆時(shí)偏移技術(shù)是未來(lái)的研究方向.

［1］Baysal E，Kosloff D D，Sherwood W C.Reverse time migration［J］.Geophysics，1983，48（11）：1514-1524.

［2］Whitmore N D.Iterative depth migration by backward time propagation［C］∥SEG Technical Program Expanded Abstracts，1983：382-385.

［3］McMechan G A.Migration by extrapolation of time-dependent boundary values［J］.Geophysical Prospecting，1983，31（3）：413-420.

［4］Loewenthal D，Mulfti I R.Reversed time migration in spatial frequency domain［J］.Geophysics，1983，48（5）：627-635.

［5］郭建.VSP技術(shù)應(yīng)用現(xiàn)狀及發(fā)展趨勢(shì)［J］.勘探地球物理進(jìn)展，2004，27（1）：1-8.Guo Jian.The application status and development trends of VSP technology［J］.Progress in Exploration Geophysics，2004，27（1）：1-8.

［6］朱金明，顏俊華.VSP逆時(shí)偏移［J］.石油地球物理勘探，1991，26（5）：564-570.Zhu Jinming，Yan Junhua.VSP reverse-time migration［J］.Oil Geophysical Prospecting，1991，26（5）：564-570.

［7］孫文博，孫贊東.基于偽譜法的 VSP逆時(shí)偏移及其應(yīng)用研究［J］.地球物理學(xué)報(bào)，2010，53（9）：2196-2203.Sun Wenbo，Sun Zandong.VSP reverse time migration based on the pseudo-spectral method and its applications［J］.Chinese J.Geophys，2010，53（9）：2196-2203.

［8］Alok K S，Xander S，Eric V.VSP imaging using all multiples：Full wavefield migration approach［C］∥SEG Technical Program Expanded Abstracts，2012：1-6.

［9］Sun Wenbo，Sun Zandong.Angle-domain VSP reverse-time migration and its application［C］∥SEG Technical Program Expanded Abstracts，2010：4339-4343.

［10］Mark E W，Lu Rongrong，Xander C，et al.A novel application of time-reversed acoustic：Salt-dome flank imaging using walkaway VSP survey［J］.Geophysics，2006，71（2）：A7-A11.

［11］劉紅偉，李博，劉洪.地震疊前逆時(shí)偏移高階有限差分算法及 GPU實(shí)現(xiàn)［J］.地球物理學(xué)報(bào)，2010，53（7）：1725-1733.Liu Hongwei，Li Bo，Liu Hong.The algorithm of high order finite difference pre-stack reverse time migration and GPU implementation［J］.Chinese J.Geophys，2010，53（7）：1725-1733.

［12］劉洋，李承楚，牟永光.任意偶數(shù)階精度的有限差分?jǐn)?shù)值模擬方法［J］.石油地球物理勘探，1998，33（1）：1-10.Liu Yang，Li Chengchu，Mu Yongguo.Finite-difference numerical modeling of any even-order accuracy［J］.Oil Geophysical Prospecting，1998，33（1）：1-10.

［13］李博，劉紅偉，劉國(guó)峰，等.地震疊前逆時(shí)偏移算法的CPU／GPU實(shí)施對(duì)策［J］.地球物理學(xué)報(bào)，2010，53（12）：2938-2943.Li Bo，Liu Hongwei，Liu Guofeng，et al.Computational strategy of seismic pre-stack reverse time migration on CPU／GPU［J］.Chinese J.Geophys，2010，53（12）：2938-2943.

［14］石穎，陸加敏，柯璇，等.基于 GPU 并行加速的疊前逆時(shí)偏移方法研究［J］.東北石油大學(xué)學(xué)報(bào)，2012，36（4）：111-115.Shi Ying，Lu Jiamin，Ke Xuan，et al.Prestack reverse time migration based on GPU parallel accelerating algorithm［J］.Journal of Northeast Petroleum University，2012，36（4）：111-115.

［15］柯璇，石穎，劉詩(shī)竹.地震資料逆時(shí)偏移中的圖像處理器加速算法［J］.計(jì)算機(jī)系統(tǒng)應(yīng)用，2013，22（11）：115-118.Ke Xuan，Shi Ying，Liu Shizhu.GPU acceleration algorithm of the reverse time migration with seismic data［J］.Computer Systems＆ Applications，2013，22（11）：115-118.

［16］郭念民，吳國(guó)忱.基于PML邊界的變網(wǎng)格高階有限差分聲波方程逆時(shí)偏移［J］.石油地球物理勘探，2012，47（2）：256-265.Guo Nianmin，Wu Guochen.Highorder finite difference method in reverse-time migration with variable grids based on PML boundary condition［J］.Oil Geophysical Prospecting，2012，47（2）：256-265.

［17］王守東.聲波方程完全匹配層吸收邊界［J］.石油地球物理勘探，2003，38（1）：31-34.Wang Shoudong.Absorbing boundary condition for acoustic wave equation by perfectly matched layer［J］.Oil Geophysical Prospecting，2003，38（1）：31-34.

［18］王永剛，邢文軍，謝萬(wàn)學(xué)，等.完全匹配層吸收邊界條件的研究［J］.中國(guó)石油大學(xué)學(xué)報(bào)：自然科學(xué)版，2007，31（1）：19-24.Wang Yonggang，Xing Wenjun，Xie Wanxue，et al.Study of absorbing boundary condition by perfectly matched layer［J］.Journal of China University of Petroleum：Natural Science Edition，2007，31（1）：19-24.

［19］胡昊，劉伊克，常旭，等.逆時(shí)偏移計(jì)算中的邊界處理分析及應(yīng)用［J］.地球物理學(xué)報(bào)，2013，56（6）：2033-2042.Hu Hao，Liu Yike，Chang Xu，et al.Analysis and application of boundary treatment for the computation of reverse-time migration［J］.Chinese J.Geophys，2013，56（6）：2033-2042.

［20］薛東川，張?jiān)迄i，朱振宇，等.幾種疊前逆時(shí)偏移成像條件的比較［C］.中國(guó)地球物理，2012：581.Xue Dongchuan，Zhang Yunpeng，Zhu Zhenyu，et al.Comparison of several imaging condition for prestack reverse-time migration［C］.The Chinese Geophysical，2012：581.

［21］劉紅偉，劉洪，鄒振，等.地震疊前逆時(shí)偏移中的去噪與存儲(chǔ)［J］.地球物理學(xué)報(bào)，2010，53（9）：2171-2180.Liu Hongwei，Liu Hong，Zou Zhen，et al.The problems of denoise storage in seismic reverse time migration［J］.Chinese J.Geophys.，2010，53（9）：2171-2180.

［22］Zhang Y，Sun J.Practical issues in reverse time migration：true amplitude gathers，noise removal and harmonic source encoding［J］.First Break，2009，26（1）：29-35.

［23］陳可洋.基于拉普拉斯算子的疊前逆時(shí)噪聲壓制方法［J］.巖性油氣藏，2011，23（5）：87-95.Chen Keyang.Pre-stack reverse-time noise suppressing method based on Laplacian operator［J］.Lithologic Reservoirs，2011，23（5）：87-95.