縱橫波波場(chǎng)分離效果對(duì)比

郭政宏，胡葉正

（西南石油大學(xué) 地球科學(xué)與技術(shù)學(xué)院，四川成都 610500）

彈性波正演是認(rèn)識(shí)地下復(fù)雜介質(zhì)的重要手段，在彈性波有限差分正演模擬中，除了二階質(zhì)點(diǎn)位移方程外，常用一階速度-應(yīng)力方程，其優(yōu)點(diǎn)是不需要對(duì)彈性參數(shù)進(jìn)行空間微分，就可以得到完整的地震波場(chǎng)響應(yīng)[1]。除了常規(guī)網(wǎng)格外，Madariaga提出的交錯(cuò)網(wǎng)格，它在不增加計(jì)算量的前提下，和常規(guī)網(wǎng)格相比局部精度提高了4倍，收斂速度也很快，很好地壓制了數(shù)值頻散[2]。

彈性波縱橫波場(chǎng)分離是進(jìn)行彈性波逆時(shí)偏移的必要條件，如若不能正確分離縱橫波場(chǎng)，則在偏移成像過(guò)程中由于縱橫波的互相干擾，導(dǎo)致偏移成像失敗。目前大多數(shù)學(xué)者都采用傳統(tǒng)的Helmholtz原理分離縱橫波場(chǎng)。Helmholtz分解基于彈性波混合波場(chǎng)是由一個(gè)無(wú)旋場(chǎng)和無(wú)散場(chǎng)組成的，這個(gè)無(wú)旋場(chǎng)就是縱波的勢(shì)函數(shù)，無(wú)散場(chǎng)就是橫波的勢(shì)函數(shù)。因此可以通過(guò)對(duì)混合波場(chǎng)做散度運(yùn)算提取縱波，做旋度運(yùn)算。但是由于旋度運(yùn)算和散度運(yùn)算對(duì)空間做了一階偏導(dǎo)數(shù)，分離后的波場(chǎng)的振幅和相位都會(huì)改變，且分離后的縱橫波波場(chǎng)與原波場(chǎng)的物理意義不一致。假設(shè)原波場(chǎng)為位移，那么經(jīng)過(guò)Helmholtz分離后的波場(chǎng)為質(zhì)點(diǎn)速度相關(guān)的一個(gè)量。而且對(duì)其做相位矯正和振幅恢復(fù)時(shí)計(jì)算量過(guò)大，對(duì)后期的逆時(shí)偏移來(lái)算不易實(shí)現(xiàn)。本文采用縱橫波解耦方程有限差分交錯(cuò)網(wǎng)格形式下的正演模擬，吸收邊界采用CPML。在提高計(jì)算效率上利用MPI和OpenMP進(jìn)行并行加速。

1 基本原理與方法

1.1 彈性波各向同性一階速度-應(yīng)力方程

密度非均勻，二維彈性波各向同性一階速度-應(yīng)力方程：

式（1）中：vx、vz為彈性波質(zhì)點(diǎn)振動(dòng)速度的水平分量和垂直分量（在多波多分量地震勘探中，也被稱之為偽橫波和偽縱波）；τxx、τzz、τxz為質(zhì)點(diǎn)振動(dòng)x軸向，z軸向的正應(yīng)力，和y平面法向的剪切應(yīng)力；λ、μ則為彈性波的模量，其中μ是剪切模量；ρ表示為密度參數(shù)。

其中，密度，縱橫波速度參數(shù)和彈性波的模量參數(shù)關(guān)系如下：

彈性波一階速度-應(yīng)力方程交錯(cuò)網(wǎng)格各個(gè)分量的配置方式如圖1所示，其中密度位于半網(wǎng)格點(diǎn)上。

各向同性彈性波介質(zhì)中，縱波和橫波是完全解耦：縱波是一個(gè)有源無(wú)旋散度場(chǎng)，橫波是一個(gè)有旋無(wú)散的散度場(chǎng)，混合波場(chǎng)是由一個(gè)散度場(chǎng)（縱波波場(chǎng)）和一個(gè)旋度場(chǎng)（橫波波場(chǎng)）疊加而成。

基于Helmholtz彈性波縱橫波場(chǎng)分離理論，可以從質(zhì)點(diǎn)振動(dòng)速度的水平分量和垂直分量提取散度獲得縱波波場(chǎng)，從中提取旋度獲得橫波波場(chǎng)。

時(shí)空域直接采用散度算子和旋度算子如下：

從式（3）中可以看出，Helmholtz分離后的縱橫波場(chǎng)，其相位整體移動(dòng)了90°，縱橫波振幅比也發(fā)生了變化，同時(shí)得到的縱橫波波場(chǎng)和原波場(chǎng)的物理意義也不相同。

1.2 縱橫波解耦的一階速度-應(yīng)力方程

密度非均勻，各向同性縱橫波波場(chǎng)分離方程：

應(yīng)用有限差分進(jìn)行彈性波波動(dòng)方程正演時(shí)，由于正演的空間有限，需要在節(jié)段的人工邊界上引入吸收邊界條件。當(dāng)前吸收效果較好的邊界條件是完全匹配層（PML）邊界，它能夠使大角度入射的彈性波在匹配層中無(wú)反射地通過(guò)，并按照指數(shù)規(guī)律進(jìn)行衰減。由于在彈性波正演模擬中，縱波速度、橫波速度和密度的差異變化復(fù)雜，使用常規(guī)的PML吸收邊界會(huì)留下一些反射波和大量的隱失波，而采用CPML算法后則可以大大減少反射波的存在，并且對(duì)隱失波的吸收具有很好的效果。

2 模型測(cè)試

選取Mamousi模型的一塊區(qū)域（圖1（a）），該區(qū)域的網(wǎng)格大小為（1 201×401），空間采樣間隔為5m，縱橫波波速比為1.73，不考慮密度模型（密度默認(rèn)為1），縱波震源激發(fā)，零相位雷克子波，震源主頻30Hz，時(shí)間采樣間隔0.000 5s，有限差分階數(shù)16階，CPML吸收層數(shù)32。

圖1 速度模型及波場(chǎng)分離

基于彈性波各向同性一階速度-應(yīng)力方程，并采用Helmholtz原理提取的縱橫波場(chǎng)如圖1所示。使用相同正演模擬參數(shù)，基于縱橫波解耦一階速度-應(yīng)力彈性波方程分離的縱橫波如1（b）所示。

由于采用彈性波進(jìn)行正演模擬，正演記錄中出現(xiàn)了特有的面波現(xiàn)象，這是由于縱橫波或者部分分量耦合在一起導(dǎo)致的，是聲波正演模擬所不具有的特征。由于Marmousi的模型相對(duì)復(fù)雜，其反射波比較雜亂，波場(chǎng)也相對(duì)復(fù)雜。

對(duì)比圖1（b）和圖2（a），分離出的縱橫波振幅與圖1（a）相比更加符合真實(shí)的縱橫波振幅比值，其反射波的細(xì)節(jié)更加清晰。同時(shí)縱橫波解耦分離的縱橫波保持原有的矢量特征，這與通過(guò)Helmholtz分離出來(lái)的標(biāo)量波場(chǎng)有本質(zhì)不同，分離出來(lái)的波場(chǎng)的物理意義和原來(lái)的波場(chǎng)的物理意義相一致。

圖2 縱橫波解耦

3 結(jié)論

基于縱橫波解耦方程交錯(cuò)網(wǎng)格下進(jìn)行數(shù)值模擬，并借助Marmousi模型進(jìn)行測(cè)試。同時(shí)和基于Helmholtz縱橫波分解的方程相比，前者提取到縱橫波場(chǎng)振幅和相位沒(méi)有發(fā)生變換，波場(chǎng)的細(xì)節(jié)清楚，其分離出的波場(chǎng)干凈，更加適應(yīng)于后續(xù)的彈性波逆時(shí)偏移成像等工作流程。