三維逆時(shí)偏移中的隨機(jī)邊界條件

三維逆時(shí)偏移中的隨機(jī)邊界條件

0多波地震資料偏移處理的意涵基于標(biāo)量波理論的反射縱波勘探技術(shù)在過(guò)去的能源開(kāi)發(fā)過(guò)程中發(fā)揮了重要作用。隨著油氣開(kāi)發(fā)對(duì)地震勘探精度的要求的不斷提高，勘探目標(biāo)結(jié)構(gòu)和巖性復(fù)雜性的不斷增加，單因素?cái)U(kuò)散的縱波采集精度的提高越來(lái)越受到其理論上的局限性和波場(chǎng)信息的個(gè)體化。碳酸鹽巖裂縫儲(chǔ)層、煤層和巖漿巖儲(chǔ)層的勘探難以以德容忍。近年來(lái),以矢量波動(dòng)理論為基礎(chǔ)的多波地震勘探技術(shù)引起了業(yè)界的廣泛重視,同常規(guī)縱波勘探技術(shù)相比,多波地震技術(shù)的理論假設(shè)與實(shí)際情況更為接近,且橫波在復(fù)雜介質(zhì)中具有獨(dú)特的傳播機(jī)理,故理論上多波地震技術(shù)解決實(shí)際問(wèn)題的能力更強(qiáng),更有利于提高勘探精度或?qū)崿F(xiàn)某些特殊目標(biāo)的精確勘探。疊前深度偏移是多波地震技術(shù)的重要研究?jī)?nèi)容,對(duì)多波地震資料進(jìn)行疊前深度偏移處理的意義在于:一方面,地下構(gòu)造的精確成像需要借助高精度的疊前偏移技術(shù)來(lái)實(shí)現(xiàn);另一方面,多波偏移剖面中所包含的地震波動(dòng)力學(xué)特征還可用于巖性和流體識(shí)別;同時(shí)偏移得到的共成像點(diǎn)道集還可為振幅隨偏移距變化(AVO)、振幅隨入射角變化(AVA)或其它疊前反演工作提供輸入數(shù)據(jù)。因此,多波地震資料偏移應(yīng)當(dāng)實(shí)現(xiàn)兩個(gè)目標(biāo):(1)地下構(gòu)造的縱、橫波準(zhǔn)確成像;(2)獲取保真的縱橫波疊前反演道集和疊加剖面?；趶椥圆ɡ碚摰亩喾至康卣鹳Y料聯(lián)合逆時(shí)偏移技術(shù)是實(shí)現(xiàn)多波地震縱橫波偏移成像的有效工具逆時(shí)偏移思想提出之初,在疊后地震資料偏移成像中獲得了良好的成像剖面目前的野外采集普遍向高密度、三維方向發(fā)展,為確保復(fù)雜構(gòu)造觀測(cè)波場(chǎng)的連續(xù)性,多采用小炮距、小接收線距、長(zhǎng)排列接收,部分探區(qū)單炮記錄過(guò)萬(wàn)道本文針對(duì)有效邊界存儲(chǔ)在三維情況下存儲(chǔ)過(guò)大的情況,提出應(yīng)用隨機(jī)邊界條件,在不保存邊界波場(chǎng)的情況下,實(shí)現(xiàn)正向延拓波場(chǎng)的重構(gòu);考慮三維逆時(shí)偏移情況下,開(kāi)辟數(shù)組空間內(nèi)存較大,提出應(yīng)用區(qū)域分解技術(shù)將計(jì)算數(shù)據(jù)分配到不同節(jié)點(diǎn),實(shí)現(xiàn)基于MPI+CUDA的協(xié)同并行。1縱波場(chǎng)疊加響應(yīng)的基本原理彈性波逆時(shí)偏移將多分量地震記錄當(dāng)作彈性波場(chǎng)的邊值問(wèn)題彈性波逆時(shí)偏移的基本原理為:假定地下介質(zhì)由一系列繞射點(diǎn)組成,地面、海底或井中接收到的三分量記錄為各繞射點(diǎn)產(chǎn)生的縱橫繞射波場(chǎng)在接收點(diǎn)處的疊加響應(yīng),其中各繞射點(diǎn)縱橫繞射波的產(chǎn)生時(shí)間為炮點(diǎn)波場(chǎng)的主能量到達(dá)該點(diǎn)的時(shí)間1.1時(shí)間上1階差分根據(jù)各向同性介質(zhì)中,位移-應(yīng)變關(guān)系、應(yīng)力-應(yīng)變關(guān)系、應(yīng)力-位移關(guān)系其中,υ在交錯(cuò)網(wǎng)格空間中對(duì)式(1)進(jìn)行差分離散可得到時(shí)間2階,空間任意偶數(shù)階的差分格式,以式(1)中第一個(gè)方程為例:其中,Δx、Δy、Δz分別表示坐標(biāo)系中x、y、z方向上的網(wǎng)格間距;Δt表示時(shí)間采樣間隔;i,j,k分別代表x、y、z方向上離散點(diǎn)的序號(hào);n表示時(shí)間上的離散點(diǎn)的序號(hào)。且ρ、μ半網(wǎng)格點(diǎn)上的值由同一方向上相鄰整網(wǎng)格點(diǎn)求取算術(shù)平均得到。2N為差分階數(shù),D1.2無(wú)分裂pml邊界彈性波方程交錯(cuò)網(wǎng)格有限差分格式屬于顯示差分型遞推格式,計(jì)算當(dāng)前時(shí)刻波場(chǎng)值時(shí),需利用上一時(shí)刻波場(chǎng)值。由于數(shù)值計(jì)算截?cái)嗾`差的影響,隨著時(shí)間遞推誤差會(huì)累積,破壞數(shù)值穩(wěn)定性甚至引起計(jì)算溢出通過(guò)對(duì)三維各向同性介質(zhì)一階雙曲型彈性波方程進(jìn)行平面諧波分析式(4)中,Δx、Δy、Δz分別表示坐標(biāo)系中x、y、z方向上的網(wǎng)格間距;Δt表示時(shí)間采樣間隔;υ采用無(wú)分裂式PML邊界條件解決階段邊界的偽反射問(wèn)題上式中Ω其中,d(x)、d(y)、d(z)為三個(gè)方向的衰減因子;根據(jù)公式(5)、(6)可依此類推基于無(wú)分裂PML邊界條件的彈性波控制方程,其計(jì)算步驟為:(1)對(duì)包含PML邊界的整個(gè)區(qū)域進(jìn)行常規(guī)的彈性波波場(chǎng)求解;(2)在PML邊界區(qū)域?qū)⒌玫降某Ｒ?guī)波場(chǎng)減去衰減量,然后按照時(shí)間軸方向更新輔助變量。由控制方程和計(jì)算步驟可知,無(wú)分裂的PML邊界條件避免了對(duì)彈性波波場(chǎng)進(jìn)行分裂,降低了存儲(chǔ)。由于采用同一套方程也不會(huì)出現(xiàn)在邊界處兩套方程耦合的編程難題。1.3偏移剖面的基本原理成像條件決定著逆時(shí)偏移成像剖面的質(zhì)量。目前業(yè)界廣泛應(yīng)用互相關(guān)成像條件,其基本原理為:利用炮點(diǎn)波場(chǎng)和接收點(diǎn)逆時(shí)延拓波場(chǎng)做零延時(shí)互相關(guān),在反射點(diǎn)位置的波場(chǎng)將得到加強(qiáng),非反射點(diǎn)處的波場(chǎng)將削弱,最后所有時(shí)刻波場(chǎng)相加得到偏移剖面?；ハ嚓P(guān)成像條件的實(shí)現(xiàn)算法為:其中,I2橫波場(chǎng)的重構(gòu)上述彈性波逆時(shí)偏移基本流程中,在利用互相關(guān)成像條件時(shí)需同時(shí)應(yīng)用相同時(shí)刻的檢波點(diǎn)波場(chǎng)和炮點(diǎn)波場(chǎng),上述兩個(gè)波場(chǎng)的時(shí)間序列相反,因此需要重構(gòu)炮點(diǎn)波場(chǎng)以實(shí)現(xiàn)縱橫波成像;本文采用隨機(jī)邊界條件散射人工邊界的偽反射,實(shí)現(xiàn)波場(chǎng)的重構(gòu)。為了解決三維情況下單卡GPU顯存不足的難題,利用區(qū)域分解技術(shù)將數(shù)據(jù)劃分,實(shí)現(xiàn)大規(guī)模數(shù)據(jù)的三維彈性波逆時(shí)偏移。2.1機(jī)速度無(wú)線性的影響在解決炮點(diǎn)、檢波點(diǎn)波場(chǎng)時(shí)間序列不一致的問(wèn)題中,國(guó)內(nèi)外學(xué)者先后提出Checkpointing方法、有效邊界存儲(chǔ)策略隨機(jī)邊界條件實(shí)質(zhì)為在計(jì)算區(qū)域外部隨機(jī)賦速度值,利用隨機(jī)速度的無(wú)規(guī)律性,將波場(chǎng)值散射,這樣在進(jìn)行互相關(guān)成像和疊加時(shí),無(wú)規(guī)律的波場(chǎng)值能量得到減弱,反射點(diǎn)成像值相對(duì)得到增強(qiáng)。所以,隨機(jī)邊界條件的關(guān)鍵之處在于隨機(jī)邊界層內(nèi)隨機(jī)速度的選擇,常用隨機(jī)速度賦值公式為:式中,V為了驗(yàn)證隨機(jī)邊界條件的適用性,本文采用二維均勻速度模型,驗(yàn)證隨機(jī)邊界的重構(gòu)和散射效果。二維均勻速度模型大小1500m×1500m,縱波速度為2500m/s,橫波速度1800m/s,炮點(diǎn)位置位于(750m,750m),縱波源激發(fā),記錄不同時(shí)刻v2.2區(qū)域分解和多卡并行數(shù)據(jù)傳遞GPU中有大量晶體管,是專為執(zhí)行復(fù)雜的數(shù)學(xué)和幾何計(jì)算而設(shè)計(jì)的。當(dāng)今的GPU已不再將這些計(jì)算局限于圖形渲染,其通用計(jì)算技術(shù)的發(fā)展和天生的并行處理機(jī)制在計(jì)算方面提供了高于CPU數(shù)十倍乃至上百倍于CPU的性能。業(yè)界已經(jīng)利用GPU的并行特性在二維逆時(shí)偏移取得了一系列突破在對(duì)彈性波方程的差分格式,進(jìn)行有限差分計(jì)算編程時(shí),將一個(gè)維度上(以y方向?yàn)槔?的區(qū)域進(jìn)行劃分,按照實(shí)際資料的大小和顯卡顯存的限制來(lái)決定劃分的數(shù)量。利用區(qū)域分解進(jìn)行差分計(jì)算,需要解決的核心問(wèn)題在于如何隱藏卡與卡間的數(shù)據(jù)交換。如圖5所示,為了保證差分計(jì)算的正確性,在每一個(gè)時(shí)間步需要在虛線區(qū)域進(jìn)行數(shù)據(jù)交換。以其中的兩塊卡為例,說(shuō)明分解區(qū)域后每塊卡的計(jì)算與數(shù)據(jù)傳遞。如圖6所示,每塊卡可以分為三部分(圖中分別用A、B和C表示),C區(qū)域?yàn)閮蓧K卡之間的重疊部分,厚度由差分格式來(lái)決定,B區(qū)域?yàn)閭鬟f部分與C區(qū)域保持一樣的厚度,A區(qū)域?yàn)轱@卡本身的計(jì)算區(qū)域。在正演和延拓計(jì)算的每一個(gè)時(shí)間步內(nèi),C區(qū)域的網(wǎng)格處在邊界而無(wú)法進(jìn)行差分運(yùn)算,所以需要將B區(qū)域的數(shù)據(jù)傳遞到C。這樣計(jì)算過(guò)程就分為了如下步驟:首先是每塊卡計(jì)算B區(qū)域;然后,每塊卡計(jì)算A區(qū)域的同時(shí),傳遞B區(qū)域數(shù)據(jù)到相鄰顯卡。計(jì)算能力高于2.0的NVIDIA顯卡支持在同節(jié)點(diǎn)上進(jìn)行顯卡與顯卡之間的數(shù)據(jù)交換。而本文提出的多卡并行是基于多節(jié)點(diǎn)的并行,在數(shù)據(jù)傳遞時(shí)需要將傳遞數(shù)據(jù)通過(guò)內(nèi)存中轉(zhuǎn)傳遞,由MPI的每一個(gè)進(jìn)程控制每塊顯卡,在MPI進(jìn)程組內(nèi)進(jìn)行非阻塞的數(shù)據(jù)交換。整個(gè)的實(shí)現(xiàn)流程圖如圖7所示:(1)由主進(jìn)程讀取到炮記錄和速度模型并分發(fā)于進(jìn)程;(2)各進(jìn)程控制GPU顯卡實(shí)現(xiàn)每一時(shí)間步的圖4B區(qū)域波場(chǎng)計(jì)算;(3)GPU回傳數(shù)據(jù)并在內(nèi)存中交換,同時(shí)實(shí)現(xiàn)圖4A區(qū)域的波場(chǎng)計(jì)算;(4)循環(huán)迭代波場(chǎng),最后由主進(jìn)程規(guī)約得到最終的偏移結(jié)果。3鹽丘模型三維彈性波逆時(shí)偏移試驗(yàn)SEG/EAEG鹽丘模型是一個(gè)國(guó)際上通用測(cè)試三維復(fù)雜構(gòu)造成像效果的地質(zhì)模型。本文采用SEG/EAEG鹽丘速度模型作為縱波模型,通過(guò)固定的縱橫波速度比構(gòu)造橫波模型,炮集記錄由一階彈性波方程有限差分正演獲得。模型參數(shù)如下:模型大小3380m×3380m×1050m,網(wǎng)格大小Δx=Δy=Δz=5m,時(shí)間采樣間隔0.35ms,記錄長(zhǎng)度2.8s,震源為35Hz雷克子波。觀測(cè)系統(tǒng)參數(shù)為:總炮數(shù)4690炮,67條測(cè)線,炮線間距50m,炮間距60m,接收線間距10m,接收道間距5m,每條測(cè)線70炮,每炮30400道接收;偏移算法采用本文所提出的基于MPI+CUDA的三維彈性波逆時(shí)偏移并行算法。圖8為三維鹽丘模型的縱橫波速度模型以及多分量逆時(shí)偏移剖面。另外,PP、PS的成像效果較SP、SS的成像效果更好。原因?yàn)?PP,PS波是在縱波傾斜入射到彈性界面時(shí)產(chǎn)生能量較強(qiáng),而SP、SS是由縱波傾斜入射到彈性界面產(chǎn)生的轉(zhuǎn)換橫波再次入射得到的地震記錄,該情況涉及波型的多次轉(zhuǎn)換,情況復(fù)雜。因此PP成像結(jié)果主要代表了P波的反射系數(shù),PS主要成像結(jié)果代表了PS波的反射系數(shù);而SP、SS成像結(jié)果并不對(duì)應(yīng)一種波的反射系數(shù),多次轉(zhuǎn)換后并不清楚它的能量分布。PP波和PS波均可以對(duì)鹽丘模型構(gòu)造得到準(zhǔn)確的偏移成像結(jié)果。PP波主頻比PS波主頻低,其原因在于:模型資料的縱橫波記錄是由同一縱波源激發(fā)得到,且轉(zhuǎn)換橫波在空間域中的波長(zhǎng)小于縱波波