阮 倫 程玖兵

(①同濟(jì)大學(xué)海洋地質(zhì)國家重點實驗室，上海 200092； ②同濟(jì)大學(xué)海洋與地球科學(xué)學(xué)院，上海 200092)

0 引言

地震體波在地球內(nèi)部非均勻介質(zhì)中的傳播過程中存在縱波(P波)和橫波(S波)兩種模式(即波型)之間的相互轉(zhuǎn)換與耦合。為了獲得具有明確物理含義的彈性波偏移圖像，縱橫波模式解耦被認(rèn)為是至關(guān)重要的步驟[1-4]。模式解耦預(yù)條件對于彈性全波形反演也很有意義[5-6]。

在各向同性介質(zhì)中，根據(jù)亥姆赫茲分解理論，彈性矢量波場可以分解成相應(yīng)的縱波和橫波成分[7-8]，馬德堂等[9]、李振春等[10]提出了P波和S波近似解耦的彈性波方程；吳瀟等[11]比較了一系列各向同性P波和S波解耦方法。然而，地球內(nèi)部介質(zhì)普遍存在不同程度的彈性各向異性，表現(xiàn)為傳播速度與偏振特征的方向相關(guān)性。此時，縱波的偏振方向與傳播方向不平行，橫波偏振方向與傳播方向不垂直，而橫波還會分裂成偏振方向相互正交的快橫波和慢橫波[12]，它們依次被稱為qP波、qS1波和qS2波。qS1和qS2波的偏振方向隨著傳播方向的變化通常是不連續(xù)的，而且當(dāng)兩種qS波相速度相等(即發(fā)生橫波奇異性)時，無法確定各自的偏振方向，故而難以區(qū)分兩種qS波模式。在對稱性較高的橫向各向同性(TI)介質(zhì)中，qS1和qS2的偏振方向要么在對稱軸平面以內(nèi)(qSV波)，要么與對稱軸平面垂直(qSH波)[13-16]。

各向異性彈性波模式解耦理論和算法研究有著近三十年的發(fā)展歷程。一方面，人們一直試圖從方程解耦的角度實現(xiàn)qP和qS波場的解耦，代表性的工作包括擬聲學(xué)近似[17]、純模式近似[18]和偽純模式波動方程[19-20]等，這些近似方程都存在動力學(xué)誤差。另一方面，受亥姆赫茲分解理論的啟發(fā)，人們也對各向同性介質(zhì)中縱、橫波分離的算法進(jìn)行了擴(kuò)展。最初，Dellinger等[21]提出在波數(shù)域求取偽散度和偽旋度算子，通過偏振投影分離qP與qS波。為了適應(yīng)非均勻各向異性介質(zhì)，Yan等[14]將偏振投影算子構(gòu)建成非穩(wěn)相空間濾波形式，魏石磊等[22]在此基礎(chǔ)上研究了空間域濾波算子的特征。然而非均勻各向異性介質(zhì)中基于偏振投影的模式解耦計算成本很高[16]，為了突破計算瓶頸，Yan等[23]設(shè)計了一種類似于單程波偏移成像中的相移加插值(PSPI)的混合域解耦算法。而后，Cheng等[16]將模式分解表達(dá)成空間—波數(shù)域的廣義傅里葉積分算子，并借用矩陣低秩近似技術(shù)[24]，提出了一種分離qP、qSV和qSH波標(biāo)量和矢量場的快速算法；Wang 等[25]建議對介質(zhì)模型進(jìn)行分塊，再利用低秩近似解耦算法。此外，也有學(xué)者將Yan等[14]提出的非穩(wěn)相空間濾波算子簡化成更經(jīng)濟(jì)的偽散度/偽旋度算子，實現(xiàn)qP和qS波近似解耦[26]。總之，從精度、成本和可操作性等角度發(fā)展各向異性彈性波模式解耦新算法仍然非常必要。

本文在回顧各向異性彈性波模式解耦的基本原理基礎(chǔ)上，首先改進(jìn)了Zhou等[27]提出的基于主能量傳播與偏振方向偏差角構(gòu)造的qP/qS波解耦偽散度、偽旋度算法，一定程度上考慮偏差角隨傳播方向的變化，進(jìn)而建立精度更高的偽散度、偽旋度算子。然后，探討低秩近似解耦算法及其模型分區(qū)、CPU-GPU異構(gòu)平臺實現(xiàn)策略。最后，應(yīng)用SEG Hess VTI模型檢驗和評價本文的模式解耦高效算法。

1 彈性波模式解耦基本原理

在各向同性介質(zhì)中，彈性波場P/S波分離通常采用散度和旋度運算[8]

(1)

(2)

式中：K=(Kx,Ky,Kz)T表示波矢量；U或U為u或u的傅里葉變換。在各向同性介質(zhì)中，波的傳播方向與P波的偏振方向一致，因此上面的模式解耦也可視為偏振投影的范疇。

在各向異性介質(zhì)中，采用散度和旋度運算會造成模式分離不徹底，出現(xiàn)模式泄漏或串?dāng)_。為此，Dellinger等[21]提出以qP波偏振方向為基礎(chǔ)的偽散度和偽旋度算子，從理論上解決了各向異性qP與qS波分離問題，即

(3)

式中aP(K)為歸一化的qP波偏振矢量。

為了求取偏振矢量，需要求解平面波解對應(yīng)的Christoffel方程

(4)

式中：m為qP、qS1、qS2；I為單位對角矩陣；G為Christoffel張量

Gik=Cijklnjnli,j,k,l=1,2,3

(5)

在VTI介質(zhì)中，qS波一般分解成偏振方向連續(xù)變化的qSV和qSH波，前者在對稱軸平面內(nèi)垂直qP波偏振，后者在對稱軸平面法向方向上偏振[13]。因此，TI介質(zhì)彈性波模式解耦對應(yīng)波數(shù)域的偏振投影公式為

U(m)(K)=iam(K)·U(K)

(6)

式中qSV、qSH波歸一化偏振矢量滿足

(7)

式中s為TI介質(zhì)的對稱軸方向。

不難發(fā)現(xiàn)，如同散度和旋度運算一樣，Dellinger等[21]構(gòu)造的偽散度、偽旋度運算分離得到標(biāo)量縱波和矢量橫波場(其中橫波場在二維情況下僅剩一個分量)。而式(6)表示偏振投影解耦得到了三種波模式的標(biāo)量場。為了從總的彈性波場分離出各個波模式量綱一致、振幅和相位相對保真的矢量場，可采用Zhang 等[15]提出的矢量分解算子

U(m)(K)=am(K)[am(K)·U(K)]

(8)

對于非均勻VTI介質(zhì)，偏振投影算子比較復(fù)雜，因為偏振方向在空間上隨介質(zhì)參數(shù)變化而改變。例如，Yan等[14]把式(6)改寫為空間域的非穩(wěn)相濾波(投影)算子

(9)

2 偏振投影空間域濾波近似算法

在各向異性介質(zhì)中，qP波傳播方向與偏振方向通常是不一致的，它們之間存在一個偏差角。一旦給定傳播方向，可以借助Christoffel方程計算精確的偏振方向。為了避免采用計算復(fù)雜度非常高的精確的偏振投影算子，有學(xué)者利用相對容易估算的相(傳播)方向或群方向去估計傳播與偏振方向之間的夾角，進(jìn)而獲得偏振方向的近似估計，從而構(gòu)造出比較容易計算的、近似的qP/qS波模式分離算子[27]。

三維情況下，qP波偏振矢量aP(K)的方向可用方位角α0(K)和傾角β0(K)表示，波矢量K的方向可用方位角α1(K)和傾角β1(K)表示，因此，可以用φ(K)和θ(K)表示aP(K)與K的方向偏差角，即φ(K)=α0(K)-α1(K)，θ(K)=β0(K)-β1(K)。理論上，除了極個別特殊方向(如對稱軸方向)，這兩個偏差角隨傳播方向變化，且都不為零。利用偏差角可從傳播方向估計qP波的偏振方向，即

(10)

將式(10)施加到偏振投影算子并變換到空間域，就與式(9)完全等價。對非均勻TI介質(zhì)而言，空間域非穩(wěn)相濾波算子遠(yuǎn)比散度/旋度算子復(fù)雜，計算成本極高。

2.1 基于主能量傳播與偏振方向夾角的偽散度與偽旋度算子

為了在精確的偏振投影非穩(wěn)相濾波算子和散度/旋度算子之間找到比較恰當(dāng)?shù)恼壑?，Zhou等[27]提出了qP波偏振方向的近似估算方法，進(jìn)而構(gòu)造出了一種非常經(jīng)濟(jì)的偽散度與偽旋度算子。其思路是，首先用散度運算獲得近似qP波標(biāo)量場u(P)，然后估計其能流密度矢量(坡印廷矢量)e=(e1，e2，e3)T

(11)

表征空間各點主能量傳播方向K。針對一階速度—應(yīng)力方程計算的矢量波場，也可采用

(12)

計算的群方向去估計相方向，從而避免估計qP波標(biāo)量場引起的誤差。式中：τjk表示應(yīng)力分量；vk表示振動速度分量；j，k=1，2，3。再由Christoffel方程計算主能量方向qP波的偏振方向以及它與傳播方向的夾角。為了簡化投影濾波算子，相同空間點上可采用不隨傳播方向變化的偏差角去逼近qP波的偏振方向[27]，即

(13)

上述近似雖可以由傳播方向估計qP波的偏振方向，但無法估計qSV和qSH波的偏振方向，至多只能分離qP與qS波。另外，當(dāng)qP、qS波從不同方向同時達(dá)到該點時，忽視偏差角隨傳播方向變化也會引入明顯誤差。

根據(jù)式(13)，偏振投影(濾波)算子可表示成旋轉(zhuǎn)后的散度/旋度算子(即偽散度/偽旋度算子)，即

(14)

式中

(15)

與旋度算子類似，偽旋度算子僅從矢量總波場中分離出矢量形式的qS波，無法區(qū)分qSV和qSH波。

在Zhou等[27]僅推導(dǎo)了二維情況下的偽散度/偽旋度算子

(16)

以及相應(yīng)的模式解耦非穩(wěn)相濾波算子

(17)

本質(zhì)上，這些近似的偏振投影運算用傳播與偏振方向夾角相關(guān)的系數(shù)去顯式地對波場的一階偏導(dǎo)數(shù)進(jìn)行加權(quán)求和。由于這些近似的權(quán)系數(shù)可以提前計算，偏導(dǎo)數(shù)波場在傳播模擬過程中也預(yù)先計算好，這種空間非穩(wěn)相濾波的計算復(fù)雜度大幅度降低。值得注意的是，由于傳播方向是從每個時刻波場的坡印廷矢量估計的，因此每個延拓時間步上都需要構(gòu)造并施加相應(yīng)的偽散度/偽旋度算子。

如圖1所示，在一個兩層VTI介質(zhì)模型中，震源位于(2000m，0)，通過彈性波方程數(shù)值模擬計算得到900ms時刻的質(zhì)點速度波場。然后采用主能量方向偏振投影近似投影算法，基本實現(xiàn)大部分空間qP、qSV波的分離。只是在一些波前面交叉區(qū)域(圖1d方框所示)，存在明顯的模式泄漏。其原因在于，傳播與偏振矢量夾角隨傳播方向變化，僅依據(jù)主能量方向估計的一個偏差角構(gòu)建偏振投影算子，誤差較大，不能處理不同方向同時達(dá)到多種波模式的分離問題。其中采用主能量方向偏振投影近似算法得到的是標(biāo)量qP和qSV波。同時，結(jié)合式(1)和式(2)，波場空間偏導(dǎo)數(shù)運算相當(dāng)于在波數(shù)域乘以iK，虛數(shù)單位i會使得縱、橫波相位改變90°，與非歸一化的波數(shù)K相乘會使二者的振幅發(fā)生改變(對比圖1b和圖1d)。為了保證分離得到的縱、橫波在相位和振幅上盡可能保真，可以采用以式(8)為基礎(chǔ)的矢量分解模式解耦算法，圖1f～圖1k為精確偏振投影算法分離的矢量qP和qSV波，可以看出矢量qP和qSV波的水平分量和垂直分量相位上分別與質(zhì)點速度的水平分量和垂直分量保持一致。

圖1 兩層VTI介質(zhì)模型偏振投影近似算法分離結(jié)果

2.2 基于多個主能量傳播與偏振方向夾角的偽散度與偽旋度算子

一般來講，波傳播到某一點時能量隨球面擴(kuò)散，傳播過程中qP波與qS波在同一時刻都經(jīng)過空間某點也是非常普遍的。為了提升Zhou等[27]方法的精度，本文適當(dāng)考慮傳播與偏振方向偏差角隨傳播方向的變化，通過多方向子波場的偽散度與偽旋度運算實現(xiàn)更精確的模式解耦。

首先，將原始波場各個分量變換到波數(shù)域，按角度面元拆分成若干方向扇區(qū)(Ωi)的子波場[28]

(18)

式中u(x)和U(K)分別表示空間域、波數(shù)域的子波場。N一般取4～8，具體取值視波場復(fù)雜性和精度要求而定。然后，分別在每個方向扇區(qū)計算相應(yīng)的坡印廷矢量ei，指示多個主能量傳播方向。再利用Christoffel方程估計不同帶區(qū)的偏差角φi和θi，進(jìn)而建立各個方向扇區(qū)內(nèi)偏振方向的估算公式

(19)

據(jù)此構(gòu)造出每個方向扇區(qū)內(nèi)的偏振投影算子，將其作用于相應(yīng)帶區(qū)的子波場就獲得解耦后的qP與qS波子波場。最后所有方向帶區(qū)解耦子波場疊加就得到最終分離的qP和qS波場。

如圖2所示，對8個主能量方向子波場施加相應(yīng)的偏振投影再求和，明顯提升了qP與qS波重疊區(qū)域的模式分離精度。與精確偏振投影算法處理結(jié)果(圖1f和圖1g)相比，整體上效果是可接受的。不足的是，基于多個主能量傳播與偏振方向夾角的偽散度與偽旋度算子分離的標(biāo)量qP和qSV波在相位和振幅上仍存在一定的誤差(對比圖1f和圖1g與圖2a和圖2b)。

圖2 兩層VTI介質(zhì)模型改進(jìn)的偏振投影近似算法分離結(jié)果

3 偏振投影低秩近似算法及其模型分區(qū)與GPU加速策略

在非均勻VTI介質(zhì)中，彈性波模式解耦的偏振投影對應(yīng)空間域非穩(wěn)相濾波算子。無論是Yan等[23]的PSPI算法，還是基于主能量傳播與偏振方向夾角估算的非穩(wěn)相濾波算法，都會損失偏振投影算子的精度，不足以應(yīng)對非均勻性和各向異性都很強的地質(zhì)情況。

3.1 偏振投影低秩近似算法

為了既能夠保證模式解耦的可靠性，又能夠大幅度降低計算復(fù)雜度，Cheng等[16]提出了精確表征偏振投影的空間—波數(shù)域積分算子及其低秩近似快速算法。對非均勻VTI介質(zhì)，偏振投影被表達(dá)成廣義傅里葉積分公式

(20)

式中am=(amx，amy，amz)T表示某種波模式對應(yīng)的偏振向量。假設(shè)空間或波數(shù)采樣點數(shù)為M，則這種混合域積分的計算復(fù)雜度為O(M2×lgM)，現(xiàn)有計算條件下很難承擔(dān)這樣的任務(wù)。根據(jù)式(8)，還可以構(gòu)建縱、橫波矢量分解的精確偏振投影算子，表示為廣義傅里葉積分形式(以x方向為例)

(21)

根據(jù)式(20)和式(21)可知，式中的積分核同時依賴于空間域模型參數(shù)和傅里葉域波數(shù)成分。廣義傅里葉積分算子可以通過核函數(shù)近似來降低積分計算的成本。通常來說，地球介質(zhì)具有分層或分區(qū)變化的彈性特征，物理上也認(rèn)為該混合域核函數(shù)矩陣是低秩的。因此根據(jù)低秩近似理論[24]，該混合域核函數(shù)矩陣可以表示為三個矩陣的乘積

(22)

式中：W表示混合域核函數(shù)矩陣；B和C分別表示降維之后的波數(shù)和空間維度的混合域矩陣；Dmn為L×N矩陣，其中L和N分別代表該算子矩陣在低秩近似意義下的列秩與行秩。如果給定矩陣分解的誤差門檻值為1.0×10-6，常見各向異性介質(zhì)模型偏振投影算子矩陣的秩均較小，通常為幾或十幾這樣的數(shù)量級。對于圖1a中兩層VTI模型，L和N均等于2。偏振投影模式解耦算子矩陣之所以存在上述低秩近似形式，主要是因為受地質(zhì)作用的控制，大多數(shù)地層都呈現(xiàn)一定的層狀或分塊特征，從而使該算子被那些最能代表模型非均勻性和各向異性特征的空間位置和傳播方向所主導(dǎo)[16]。

依據(jù)上述低秩分解(式(22))，式(20)中的積分式可退化成

(23)

式中j={x，y，z}。它們主要的計算量體現(xiàn)在N個矢量波場的正、反傅里葉變換，以及空間波場的N和L次加權(quán)求和。因為L和N遠(yuǎn)遠(yuǎn)小于空間或波數(shù)采樣點數(shù)M，所以計算復(fù)雜度降低至O(L×M×lgM)或O(N×M×lgM)。這樣就以較小的計算成本合理逼近了式(20)中的廣義傅里葉積分運算[16]。更重要的是，偏振投影算子矩陣的低秩分解完全是根據(jù)介質(zhì)模型本身特點自動進(jìn)行的，不存在Yan等[23]的PSPI算法中設(shè)置參考模型類似的人為因素，這種模型自適應(yīng)特點對實際應(yīng)用非常重要。

3.2 模型分區(qū)策略

Wang[25]等認(rèn)為，對于一些大規(guī)模應(yīng)用問題，可把模型空間劃分成若干區(qū)域，分別在每個區(qū)域構(gòu)造低秩近似偏振投影算子，進(jìn)一步提高算法效率和并行擴(kuò)展?jié)摿?。如圖3a所示，首先需要將模型進(jìn)行分區(qū)(黃色部分)，區(qū)域之間保留一定的重疊(綠色部分)，模型外圍有PML吸收邊界(藍(lán)色部分)，然后分別對每塊區(qū)域求取偏振矢量并構(gòu)建低秩近似投影算子，對彈性波場也進(jìn)行同樣地分區(qū)并施加窗函數(shù)。如圖3b所示，重疊區(qū)域的邊界衰減函數(shù)采用的是正弦和余弦函數(shù)，且滿足二者權(quán)重系數(shù)求和為1； 最后，對不同區(qū)域的波場分別施加低秩近似偏振投影，并將分離后的波場按區(qū)域拼接。

圖3 模型分區(qū)示意圖(a)及重疊區(qū)域邊界衰減系數(shù)設(shè)置(b)

盡管模型分區(qū)一定程度上可減小低秩近似解耦算法中FFT的總體成本，且一定程度降低了每個分區(qū)內(nèi)偏振投影算子矩陣的秩，有助于節(jié)省單個時間片波場的模式解耦計算成本。但是，模型和波場分區(qū)需要每個分區(qū)周圍的鑲邊處理，既帶來額外計算開銷，又在邊界附近引入一定的誤差，這些問題對彈性波逆時偏移或全波形反演等實際應(yīng)用中大量延拓時間步上的模式解耦來說是不可忽略的。

3.3 GPU加速策略

GPU并行加速在大規(guī)模波場數(shù)值模擬、逆時偏移成像領(lǐng)域已經(jīng)得到廣泛應(yīng)用，常采用的有OpenCL和CUDA兩種并行編程框架，現(xiàn)已開始使用OpenACC[29]。

從前文分析可知，正、反傅里葉變換在低秩近似模式解耦算法中占據(jù)大部分計算量，因此提高傅里葉變換計算效率非常關(guān)鍵。利用基于CUDA架構(gòu)的快速傅里葉變換并行算法庫cuFFT很好地解決了這個問題。其次，彈性波波模式解耦算法本身具有很好的數(shù)據(jù)獨立性，因此可利用GPU進(jìn)行細(xì)粒度的并行加速。

對于模型空間太大、單卡內(nèi)存不足的情況，需要通過多卡并行方式來解決。多GPU并行算法的基礎(chǔ)是沿著慢維方向?qū)δＰ涂臻g進(jìn)行分區(qū)，而核心問題是GPU間的通信及其效率。單節(jié)點內(nèi)GPU間的實時通信可以通過P2P技術(shù)實現(xiàn)，本文以單節(jié)點雙卡并行作為示例。如圖4所示，整個模型空間分成兩部分，分別對應(yīng)在GPU1和GPU2上計算。每個計算分區(qū)都包含計算區(qū)域、數(shù)據(jù)發(fā)送區(qū)域(S1、S2)和數(shù)據(jù)接收區(qū)域(R1、R2)，其中S1(S2)、R1(R2)分別對應(yīng)為GPU1(GPU2)上的數(shù)據(jù)交換(發(fā)送、接收)內(nèi)存區(qū)域。數(shù)據(jù)傳輸過程如圖中箭頭所示。P2P技術(shù)避免以CPU作為數(shù)據(jù)傳輸中轉(zhuǎn)站，減小了數(shù)據(jù)傳輸延遲。

圖4 雙卡并行模型分區(qū)示意圖

4 數(shù)值試驗

本文以SEG Hess VTI介質(zhì)二維模型為例，測試本文發(fā)展的模式解耦高效算法。圖5為模型右半部分，尺寸為11.25km×9.38km(對應(yīng)的空間網(wǎng)格數(shù)為1800×1500)，震源位于(4.75km，3.81km)處。圖6為1400ms時的質(zhì)點速度水平和垂直分量的波前快照。圖7分別展示了本文依據(jù)4個主能量方向傳播與偏振矢量夾角構(gòu)建的偽散度/偽旋度算子、偏振投影低秩近似解耦算法以及基于模型分區(qū)(以分成四個區(qū)域為例)的偏振投影低秩近似解耦算法得到的qP和qSV波場，可見幾種算法都取得了比較好的模式解耦效果。

圖5 SEG Hess VTI模型

圖6 SEG Hess VTI模型彈性波質(zhì)點速度場1400ms時刻的水平(a)和垂直(b)分量

圖7 SEG Hess VTI模型改進(jìn)的偏振投影近似算法和低秩近似算法波場分離結(jié)果

如果采用偽散度/偽旋度算子進(jìn)行模式解耦，需要在每個時間步計算坡印廷矢量，進(jìn)而構(gòu)建每個時間步的空間非穩(wěn)相濾波器。相反，對于偏振投影低秩近似解耦算法，僅需要提前一次性構(gòu)建出解耦濾波器，然后在各個時間步直接應(yīng)用。在單核CPU(Intel(R) Xeon(R) CPU E5-2650 v4 @ 2.20GHz)平臺上，以2000個時間步彈性qP/qSV波模式解耦為例，采用依據(jù)4個主能量方向傳播與偏振矢量夾角構(gòu)建的偽散度/偽旋度算子，平均每個時間片波場模式解耦耗用時間為6.8s，而偏振投影低秩近似算法的平均耗時為5.8s，兩種算法在計算效率方面比較接近。

對于SEG Hess VTI模型低秩近似模式解耦，基于模型分區(qū)的低秩近似解耦GPU算法相比于常規(guī)的低秩近似解耦算法，單GPU卡(Nvidia Tesla K80)上改進(jìn)后的低秩近似算法實現(xiàn)的加速比達(dá)到了4以上，同時從分離結(jié)果上可以看出前者精度也基本沒有損失(適當(dāng)選取區(qū)域劃分的個數(shù))。對于彈性波逆時偏移和全波形反演等處理而言，需要每個時間的波場都施加模式解耦，基于模型分區(qū)和GPU加速的低秩近似解耦算法在保證精度的前提下明顯縮短了處理周期。此外，對于圖1的兩層VTI模型而言，基于模型分區(qū)和GPU加速的低秩近似解耦算法反而比常規(guī)的耗時多，原因在于CPU-GPU異構(gòu)平臺涉及它們之間的數(shù)據(jù)通信開銷，只有當(dāng)模型規(guī)模達(dá)到一定程度時，基于模型分區(qū)和GPU加速算法優(yōu)勢才能凸顯。

對比圖7a、圖7b與圖7g～圖7j可知，偏振投影低秩近似矢量解耦算法可以使得分離后的縱橫波場的相位和振幅保真，同時結(jié)合模型分區(qū)和GPU加速策略，使算法的精度和計算效率達(dá)到最佳。

5 結(jié)論

針對VTI介質(zhì)qP與qS波模式解耦高昂的計算成本，本文一方面基于現(xiàn)有的由主能量傳播與偏振方向偏差角構(gòu)造的偽散度/偽旋度算法，通過波場方向分解估計的多個主能量方向逼近偏振投影算子，實現(xiàn)了更精確的空間域非穩(wěn)相濾波算法，較好地解決了原算法在qP與qS波同時經(jīng)過相同空間區(qū)域時的模式泄漏問題。由于該方法利用坡印廷矢量估計的傳播與偏振矢量夾角逼近qP波的偏振方向(無法分別估計qSV和qSH波的偏振方向)，它構(gòu)建的模式解耦濾波只能實現(xiàn)qP與qS波分離，無法區(qū)分VTI介質(zhì)中的qSV與qSH波。另一方面，基于現(xiàn)有的偏振投影低秩近似算法，通過采用模型分區(qū)和GPU加速策略進(jìn)一步提高了算法效率。將模型空間分區(qū)和線程塊分配有機(jī)結(jié)合，利用CUDA架構(gòu)下的cuFFT快速傅里葉變換算法以及P2P數(shù)據(jù)通信技術(shù)，可以縮短這種高精度解耦算法的計算周期。對于二維VTI介質(zhì)彈性波模式解耦問題，本文兩種高效算法性價比相當(dāng)。第二種算法因為能夠?qū)S波進(jìn)一步解耦，因此三維擴(kuò)展性更強。今后將結(jié)合三維情況下的彈性波逆時偏移和全波形反演實際應(yīng)用發(fā)揮該算法的優(yōu)勢。