2. 山东科技大学山东省沉积成矿作用与沉积矿产重点实验室, 山东青岛 266590;
3. 河北工程大学地球科学与工程学院, 河北邯郸 056038
2. Shandong Provincial Key Laboratory of Depositional Mineralization & Sedimentary Mineral, Shandong University of Science and Technology, Qingdao Shandong 266590, China;
3. School of Earth Science and Engineering, Hebei University of Engineering, Handan Hebei 056038, China
线性Born近似是一种常用的地震波场正演模拟方法(Tarantola, 1984; Beylkin and Burridege, 1990; Chapman, 2004),该方法描述了地下介质扰动同地震散射波场之间的线性映射关系,可以准确的模拟地下的一次反射波和绕射波波场信息,具有重要的研究意义和应用价值.
在实际应用中,Born正演通常结合常规的地震射线方法进行实现.通过追踪背景速度场中地震射线,我们可以计算上行和下行格林函数的走时和振幅信息,然后通过求取关于地下散射点的体积分即可得到最终的射线Born正演模拟结果. Červeny和Coppoli (1992)利用射线Born正演对包含单个绕射点的速度模型进行地震记录的合成,克服了普通的射线正演无法模拟地下绕射波场的难题;Gibson等(1993)利用射线Born正演进行了裂隙薄层介质中VSP(Vertical Seismic Profiling,垂直地震剖面)数据的合成;Moser(2012)验证了声波射线Born正演在速度模型包含突变边界和尖点时的绕射波模拟效果;Šachl(2013)探讨了三维射线Born近似实现流程,并且提出在模型边界处应用余弦窗来消除模型截断所造成的虚假边界绕射;Tengesdal(2013)和Minakov等(2017)对比了射线Born正演同有限差分法正演的应用效果,探讨了参数的选择对结果的影响,并且将该方法应用到全波形反演之中.射线Born正演虽然具备常规射线方法灵活、高效的优点,但是也保留了射线类方法的弊端,例如射线的焦散区、阴影区等波场奇异性区域会对射线Born正演的精度造成影响,更为重要的是,由于多次走时表的存取和编程实现的困难,射线Born正演往往使用的是单次走时表,这样也就使得该方法无法准确的模拟多次走时波场,严重影响了该方法在复杂模型中的应用效果.
为了克服射线Born正演的不足,Huang等(2016)提出了一种基于高斯束的Born正演模拟方法.该方法使用高斯束取代常规的地震射线进行地下散射点走时和振幅信息的计算(Červeny, 1983; 岳玉波等, 2012;Gao et al., 2017),在理论上克服了常规射线方法的弊端,并且具备模拟多次走时波场的能力.然而,他们所提出基于最小虚值走时的频率域等时线叠加算法,无法保证全部波场信息的利用,波场模拟精度有限且计算效率也不存在优势.岳玉波等(2019)提出了一种兼具模拟精度和计算效率的高斯束Born正演的实现方法,该方法利用高斯束所携带的走时和振幅信息,将地下的反射率映射为地表稀疏束中心位置处的局部平面波,然后通过逆倾斜叠加将局部平面波累加到临近束中心的接收点波场中.
在岳玉波等(2019)研究成果的基础上,我们将高斯束Born正演应用到弹性波地震波场的模拟中.在弹性各向同性介质中,弹性波的主分量波场可以等效为标量波场,相关的正反演问题可以得到极大的简化(Schleicher et al., 2007).以此为理论基础,我们推导了地下不同波型的反射率同弹性波主分量波场之间的线性映射关系,进而建立了弹性波多分量波场Born正演公式;我们利用高斯束的走时、振幅和极性信息将地下反射率映射为不同波型的局部平面波,然后通过逆倾斜叠加将所合成的局部平面波转化为时空域的多分量地震记录.本文首先对理论公式进行了推导,接下来给出了相关算法的实现过程,最后通过两个数值模型对本文方法的应用效果进行了验证.
1 方法原理在三维弹性各向同性介质中,假设震源xs处的P波全方位点源所激发的地震波,经由地下构造的反射和转换后产生反射PP波、PS1波和PS2波,被位于自由空间地表的接收点xr接收形成三分量地震记录u(xr, xs, ω).
1.1 弹性波主分量波场根据零阶射线理论(Červeny, 2005),我们可以使用射线中心坐标系中的三个满足右手系且相互正交的基矢量e1, e2, e3来表示S1波、S2波和P波主分量的极化方向.其中,e3沿着射线的切线方向,是固定的不变的,而e1, e2则位于垂直于射线的平面内,可以在满足上述基矢量准则下任意选取.如果我们特别的选取e1使其位于地下散射点处的入射平面内(图 1中的êSV),那么e2将会垂直于该平面(图 1中的êSH),此时的S1波、S2波就成为了通常所说的SV波和SH波,而且两者之间是完全解耦的.基于该选择,弹性波的传播过程可以由两个独立的P-SV传播系统和SH传播系统来描述,对应的反射系数矩阵也由9个分量简化为了5个,从而有效的简化了相关的正反演问题.
在确定地下散射点处x的横波极化矢量êSV和êSH后,我们可以相应的计算接收点xr处的极化矢量eSV, eSH, eP(原理见图 1),然后将该点接收到的多分量地震波场u(xr, xs, ω)通过式(1)转化为弹性波的主分量uSV, uSH, uP
(1) |
其中
(2) |
为由射线中心坐标系到笛卡尔坐标系的旋转矩阵.经过式(1)转换得到的主分量波场,可以视为对应不同波型的标量波场,从而可以使用相对简单的标量传播算子进行波场的正向或者反向延拓.Jackson等(1991),Takahashi (1995),Goertz(2002),Lüth等(2005)均使用了上述方式来解决弹性波的偏移成像问题.需要注意的是,由于我们所选取的接收点处的极化矢量是由震源出射且经过地下散射点后到达接收点的射线路径所决定.因此,式(1)所描述的波场转换过程只能在波场的模拟过程中动态的进行.
1.2 弹性波Born正演公式根据线性散射理论(Aki and Richards, 2002; Stolt and Weglein, 2012),我们可以得到各向同性介质中某个弹性波主分量波场uwv(xr, xs, ω)的Born正演公式:
(3) |
其中,上标w, v分别代表了上行波和下行波的波型,F(ω)是震源子波,Mwv(x)是地下散射点x处的散射势,也可以视作该点的反射率(Stolt and Weglein, 2012),Gw(x, xr, ω)是上行波格林函数,Gv(x, xs, ω)是下行波格林函数.在考虑到所有可能的上行和下行波波型组合后,我们得到如下以算子形式表示的弹性多分量波场Born正演公式:
(4) |
其中,M (x)为地下散射点x处的反射率矩阵,含有5个非零分量:
(5) |
G (x, xr, ω)是上行波3×3格林函数矩阵:
(6) |
G (x, xs, ω)是下行波3×1格林函数矩阵:
(7) |
某个波型u的格林函数Gu(x, x ′, ω),可以表示为一系列由地表x ′出射且对地下散射点x有贡献的高斯束的叠加(Hill,2001; Bai et al., 2016):
(8) |
其中,pu=(px, py, pz)T为高斯束中心射线在地表处的初始慢度矢量,Au(x, x ′),τu(x, x ′)分别为高斯束的复值振幅和走时(Červeny, 2005):
(9) |
(10) |
式中,下标R, I分别代表复数的实部和虚部.将式(8)代入式(4)后,我们得到以高斯束表征的弹性波Born正演公式:
(11) |
其中,m=1, 2, 3分别代表三分量波场的x, y, z分量.ARwv(xr, xs),AIwv(xr, xs)分别为高斯束振幅乘积项的实部和虚部,τRwv(xr, xs),τIwv(xr, xs)分别为高斯束走时之和的实部和虚部,具有如下的形式:
(12) |
直接计算式(11)需要在每个接收点进行高斯束的计算,为了提高计算效率,我们在稀疏的束中心位置L计算高斯束(Hill, 2001; 岳玉波, 2011;Yuan et al., 2017),并且引入相位校正因子exp[-iωpw·(xr- L)]来近似表示束中心附近的接收点的格林函数
(13) |
为了减少上述近似所引入的误差,我们在式(11)中应用下述的单位分解公式(Hill, 2001; Gray and Bleistein, 2009)
(14) |
从而得到最终的弹性波Born正演公式:
(15) |
其中,
(16) |
式(15)的作用是通过局部的逆倾斜叠加将束中心处多分量局部平面波累加到接收点波场中,进而形成模拟的多分量地震记录.为保证计算效率,我们可以在频率波数域内进行该式的计算.式(16)的作用是利用高斯束所包含的走时、振幅和极性信息,将地下的反射率映射转化为束中心处的单分量局部平面波.由于该式的计算位于高斯束Born正演的最内层,其不但决定了波场模拟精度,还决定了最终的计算效率.岳玉波等(2019)对该式的计算进行了详细的探讨,并且给出了一套高效时间域算法.
1.4 具体实现过程弹性波高斯束Born正演的具体实现过程概括如下:
(1) 根据所选择的模拟参数信息,确定地表的束中心位置、高斯束的角度范围和角度间隔;
(2) 根据平滑后的P波和S波速度场,分别在震源和束中心位置进行P波和S波高斯束的试射计算,并存储高斯束走时、振幅和极化矢量等信息的参数表;
(3) 对于每一个束中心位置,根据存储的高斯束参数表以及地下散射点处PP和PS波的反射率,通过计算式(16)来求取束中心处对应不同波型的单分量局部平面波;
(4) 循环步骤(2)和(3)直到所有的束中心处理完成.计算式(15),将对应不同的束中心位置、不同波型的局部平面波通过逆倾斜叠加累加到接收点波场中,得到最终的多分量地震记录.
2 模型试算在本节中,我们通过两个二维模型数据进行弹性波高斯束Born正演的试算.首先使用一个层状介质模型来展示该方法的具体实现过程.接下来,使用Marmousi2模型测试该方法在复杂模型中的应用效果,并同弹性波动方程有限差分正演结果进行对比分析.
2.1 层状介质模型首先,使用一个层状模型进行测试.图 2a,2b分别展示了该模型的纵波和横波速度场,密度假设为常数.模型横向采样点数为1000,采样间隔为10 m,纵向采样点数为550,采样间隔为5 m.图 3a,3b分别为对速度场进行平滑后计算得到的PP和PS波反射率剖面.在进行单炮正演时,我们将炮点放在模型地表横向位置2.4 km处,并使用主频为20 Hz的雷克子波作为震源;240个接收点均匀的分布在炮点两侧,接收点间隔为20 m;每个接收道道长为2.7 s,时间采样间隔为2 ms.
在正演的过程中我们使用了高斯束偏移所使用的参数选取准则(Hill, 2001),该炮共有19个束中心,每个束中心计算32个高斯束.我们以x分量为例,展示高斯束Born正演的过程.首先,我们根据高斯束的走时和振幅以及所处理波型的极化矢量信息计算式(16),将PP波和PS波反射率映射到tau-p域,得到的结果如图 4a和4b所示;接下来,我们将其分别同地震子波进行褶积,得到如图 5a和5b所示的对应波型分别为PP和PS的x分量局部平面波;最后,我们利用式(15)所示的逆倾斜叠加,将PP和PS局部平面波累加到接收点波场中,从而得到图 6a所示的x分量单炮记录.应用相似的高斯束Born正演过程,我们得到了图 7a所示的z分量单炮记录.为了对上述结果的正确性进行验证,我们使用相同观测和模拟参数进行了波动方程有限差分法正演,得到分别如图 6b、图 7b所示的x、z分量记录.可以看到,高斯束Born正演得到的多分量记录同有限差分结果非常相似,从而证明了本文方法的正确性.不过,由于高斯束Born正演无法模拟直达波、广角反射波和多次波,因此在箭头所标识的位置同波动方程正演结果存在一定的差异.
计算效率是高斯束Born正演的重要优势,我们在单炮正演的过程中统计了两种正演方法的运算时间.其中,高斯束Born正演的计算时间是3.1 s,而有限差分法是162.5 s,高斯束Born正演的计算效率是有限差分法的近50倍.
2.2 Marmousi2模型接下来,使用Marmousi2模型来测试高斯束Born正演在复杂模型中的应用效果.图 8展示了Marmousi2模型的P波速度场,假设S波速度是P波速度的0.6倍,密度为常数.该模型横向采样点数为2701,采样间隔为10 m,纵向采样点数为701,采样间隔为5 m.在正演时,我们使用30 Hz的雷克子波作为震源函数,并将炮点放置在横向位置为15 km处的地表位置;该炮共有301个接收道,道间隔为20 m,且沿炮点对称分布.由图 8可以看到,该炮对地下照明的区域恰好为模型中构造最复杂的部分.
首先使用高斯束Born正演模拟了该炮的x和z分量单炮记录,得到的结果分别如图 9a和10a所示,然后使用有限差分法计算了具有相同观测和模拟参数的x和z分量单炮记录,得到的结果分别如图 9b和10b所示.为了便于两种方法的对比,我们切除了有限差分正演结果中的直达波,并且在高斯束Born正演结果中应用了相同的切除函数.通过对比可以看到,两者整体上具有非常相似的波场分布特征,由于高斯束Born正演无法模拟多次波等波场信息,因此两者在某些局部细节上还存在一定的差异.
我们同样对比了两种方法的计算时间,其中高斯束Born正演为7.8 s,而有限差分正演为403.1 s.高斯束Born正演具有50倍以上的计算效率优势.虽然高斯束Born正演只能模拟一次散射波(一次反射和绕射),无法像有限差分法那样得到全部的波场信息,然而在地震勘探中,我们感兴趣且能有效利用的也正是一次散射波场.此外,当模型中横波速度变低时,有限差分法需要相应的减小计算网格以保证计算过程的稳定性,计算时间会大大增加,此时高斯束Born正演计算效率方面的优势会更加的明显.因此可以说,高斯束Born正演为复杂构造模型弹性波的正演模拟提供了一种精确、高效的实现方式.
3 结论本文提出了一种弹性波高斯束Born正演方法,可以用于模拟各向同性介质一次散射弹性波地震波场.由于使用高斯束作为波场的传播算子,该方法不但保留了射线类方法高效的优点,还具备了模拟多次走时波场的能力,保证了波场模拟的精度.文中的两个数值模型的例子证明,弹性波高斯束Born正演方法不但可以准确的模拟复杂构造模型中的弹性波一次散射波场,而且具有比有限差分法弹性波正演高得多的计算效率.接下来,我们将以该方法为基础构建相应的共轭算子,开展弹性波最小二乘高斯束偏移技术的相关研究.
Aki K, Richards P G. 2002. Quantitative Seismology. 2nd ed.. Sausalito, CA: University Science Books.
|
Bai M, Chen X H, Wu J, et al. 2016. Q-compensated migration by Gaussian beam summation method. Journal of Geophysics and Engineering, 13(1): 35-48. DOI:10.1088/1742-2132/13/1/35 |
Beylkin G, Burridge R. 1990. Linearized inverse scattering problems in acoustics and elasticity. Wave Motion, 12(1): 15-52. DOI:10.1016/0165-2125(90)90017-X |
Červeny V. 1983. Synthetic body wave seismograms for laterally varying layered structures by the Gaussian beam method. Geophysical Journal International, 73(2): 389-426. |
Červeny V, Coppoli A D M. 1992. Ray-Born synthetic seismograms for complex structures containing scatterers. Journal of Seismic Exploration, 1: 191-206. |
Červeny V. 2005. Seismic Ray Theory. Cambridge, UK: Cambridge University Press.
|
Chapman C. 2004. Fundamentals of Seismic Wave Propagation. Cambridge ,UK: Cambridge University Press.
|
Gao Z H, Sun J G, Sun X, et al. 2017. A fast algorithm for depth migration by the Gaussian beam summation method. Journal of Geophysics and Engineering, 14(1): 173-183. DOI:10.1088/1742-2140/aa52dd |
Gibson R L Jr, Toksöz M N, Batini F. 1993. Ray-Born modelling of fracture-zone reflections in the Larderello geothermal field. Geophysical Journal International, 114(1): 81-90. DOI:10.1111/gji.1993.114.issue-1 |
Goertz A. 2002. True-Amplitude Multicomponent Migration of Elastic Wavefields. Logos Verlag Berlin.
|
Gray S H, Bleistein N. 2009. True-amplitude Gaussian-beam migration. Geophysics, 74(2): S11-S23. DOI:10.1190/1.3052116 |
Hill N R. 2001. Prestack Gaussian-beam depth migration. Geophysics, 66(4): 1240-1250. DOI:10.1190/1.1487071 |
Huang X G, Sun H, Sun J G. 2016. Born modeling for heterogeneous media using the Gaussian beam summation based Green's function. Journal of Applied Geophysics, 131: 191-201. DOI:10.1016/j.jappgeo.2016.06.004 |
Jackson G M, Mason I M, Lee D. 1991. Multicomponent commonreceiver gather migration of single-level walk-away seismic profiles. Geophysical Prospecting, 39(8): 1015-1029. DOI:10.1111/gpr.1991.39.issue-8 |
Lüth S, Buske R, Giese R, et al. 2005. Fresnel volume migration of multicomponent data. Geophysics, 70(6): S121-S129. DOI:10.1190/1.2127114 |
Minakov A, Keers H, Kolyukhin D, et al. 2017. Acoustic waveform inversion for ocean turbulence. Journal of Physical Oceanography, 47(6): 1473-1491. DOI:10.1175/JPO-D-16-0236.1 |
Moser T J. 2012. Review of ray-Born forward modeling for migration and diffraction analysis. Studia Geophysica et Geodaetica, 56(2): 411-432. DOI:10.1007/s11200-011-9046-0 |
Šachl L. 2013. Validation of 3D synthetic seismograms based on the ray-Born approximation. Studia Geophysica et Geodaetica, 57(1): 84-102. |
Schleicher J, Tygel M, Hubral P. 2007. Seismic True-Amplitude Imaging. SEG Geophysical Developments Series No. 12. SEG.
|
Stolt R H, Weglein A B. 2012. Seismic Imaging and Inversion:Volume 1:Application of Linear Inverse Theory. Cambridge: Cambridge University Press.
|
Takahashi T. 1995. Prestack migration using arrival angle information. Geophysics, 60(1): 154-163. DOI:10.1190/1.1443742 |
Tarantola A. 1984. Inversion of seismic reflection data in the acoustic approximation. Geophysics, 49(8): 1259-1266. DOI:10.1190/1.1441754 |
Tengesdal H C B. 2013. Ray-Born modelling and full waveform inversion[Master's thesis]. Bergen, Norway: The University of Bergen.
|
Yuan M L, Huang J P, Liao W Y, et al. 2017. Least-squares Gaussian beam migration. Journal of Geophysics and Engineering, 14(1): 184-196. DOI:10.1088/1742-2140/14/1/184 |
Yue Y B. 2011. Study on Gaussian beam migration in complex medium[Ph. D. thesis] (in Chinese). Qingdao: China University of Petroleum.
|
Yue Y B, Li Z C, Qian Z P, et al. 2012. Amplitude-preserved Gaussian beam migration under complex topographic conditions. Chinese Journal of Geophysics (in Chinese), 55(4): 1376-1383. DOI:10.6038/j.issn.0001-5733.2012.04.033 |
Yue Y B, Qian Z P, Zhang X D, et al. 2019. Gaussian beam based Born modeling method for single-scattering waves in acoustic medium. Chinese Journal of Geophysics (in Chinese), 62(2): 648-656. DOI:10.6038/cjg2019M0367 |
岳玉波. 2011.复杂介质高斯束偏移成像方法研究[博士论文].青岛: 中国石油大学. http://cdmd.cnki.com.cn/article/cdmd-10425-1011286385.htm
|
岳玉波, 李振春, 钱忠平, 等. 2012. 复杂地表条件下保幅高斯束偏移. 地球物理学报, 55(4): 1376-1383. DOI:10.6038/j.issn.0001-5733.2012.04.033 |
岳玉波, 钱忠平, 张旭东, 等. 2019. 声波介质一次散射波场高斯束Born正演. 地球物理学报, 62(2): 648-656. DOI:10.6038/cjg2019M0367 |