基于快速推进迎风双线性插值法的三维地震波走时计算

基于快速推进迎风双线性插值法的三维地震波走时计算

孙章庆;孙建国;岳玉波;江兆南

【摘要】三维地震波走时计算技术是三维地震反演、层析成像、偏移成像等诸多

地震数据处理技术中非常重要的正演计算工具.为了获得精度高且兼顾效率的三维

走时计算方法:首先,在常规双线性插值公式推导过程中,充分利用平面波双线性假设的结论,获得了二元极小值超越方程的解析解,进而推导出了准确的局部走时计算公式,同时构造性地证明了该计算公式满足地震波的传播规律和Eikonal方程;其次,引入迎风差分的基本思想,提出迎风双线性插值的局部走时计算策略,该计算策略能简

化算法、提高效率且保证无条件稳定性;然后,将上述计算公式和迎风双线性插值策

略与常规快速推进法中的窄带技术结合,获得了一种新的基于快速推进迎风双线性

插值法的三维地震波走时计算方法;最后,通过精度和效率分析检验了新算法的精度、效率和正确性,并通过计算实例验证了算法在面对复杂介质时的稳定性和有效性.【期刊名称】《地球物理学报》

【年(卷),期】2015(058)006

【总页数】13页(P2011-2023)

【关键词】三维;平面波双线性假设;迎风双线性插值;窄带技术;地震波走时计算

【作者】孙章庆;孙建国;岳玉波;江兆南

【作者单位】吉林大学地球探测科学与技术学院,长春 130026;国土资源部应用地

球物理综合解释理论开放实验室—波动理论与成像技术实验室,长春130026;吉林

大学地球探测科学与技术学院,长春 130026;国土资源部应用地球物理综合解释理

论开放实验室—波动理论与成像技术实验室,长春130026;东方地球物理公司物探

技术研究中心,河北涿州 072751;重庆市地质矿产勘查开发局208水文地质工程地质队,重庆408300

【正文语种】中文

【中图分类】P631

1 引言

地震波走时是描述地震波运动学特征的一个很重要的物理参数，其刻画了地震波由震源点激发后到达地下介质空间中各个位置坐标点所需耗消的走时，而该走时和地下介质的速度参数的分布情况密切相关，所以地震波走时计算方法常常被应用于走时反演、层析成像、偏移成像、速度分析等一些重要的地震数据处理技术中.地震

波走时计算的精度、效率、稳定性等性质直接影响着这些地震数据处理的效果和效率（Sun，1998，1999，2000，2004；井西利等，2007；瞿辰等，2007），同时考虑到三维地震勘探正在全面大范围展开的大背景，所以对三维地震波走时计算进行研究具有很重要的理论意义和实际价值.

在走时计算方面，目前采用的主要方法有两点射线追踪法（Julian and Gubbins，1977）、最短路径射线追踪法（Moser，1991；张美根等，2006）、波前构建

法（Vinje et al.，1993）、走时插值法（Asakawa and Kawanaka，1993；Wang and Ma，1999）、有限差分法（Vidale，1988；Sethian and Popovici，1999；Sun et al.，2011）、辛几何算法（秦孟兆和陈景波，2000；高亮等，2000）等.其中，两点射线追踪法在早期天然地震数据处理中起着非常重要的作用，但在现代三维地震数据处理中其计算效率相对较低.最短路径射线追踪法具有不错

的走时计算精度，但由于其需要设置很多网格线上的计算节点，所以在求解三维问题时其需要耗费相对更大的计算和存储成本.波前构建法能同时计算走时和射线路

径，并且其计算精度也相对不错，不过在三维算法时，其中的网格定位、插值问题等都相对非常复杂.实际上，走时插值法和有限差分法是目前被广泛采用的三维走时算法，其中有限差分法具有计算效率高、局部走时算法简洁等优势，但其计算精度有限，而具有更高计算精度的高阶有限差分需要花费相对更多的计算成本.走时插值法具有不错的计算精度，但是其算法的整体实现策略的效率相对有限差分法低很多.因此，本文将采用走时插值法作为局部走时算法，而在整体实现策略上采用有限差分算法中的快速算法.

目前，三维走时插值法主要采用B样条插值法（张东等，2013）、双线性插值法（李培明等，2013；刘锋等，2012；梅胜全等，2010；张东等，2009）等算法.在双线插值法中，有两方面关键算法.一方面是局部走时算法，其包括计算公式和局部计算策略等技术环节；另一方面是算法的整体实现策略，也即算法的整体实现步骤.在局部走时算法方面，目前普遍认为双线性插值法求解极小值的方程为一个超越方程，无法求出解析解（其可能解分布在一条双曲线上），所以一些学者分别采用网格界面剖分法（张东等，2009）、快速插值法（李培明等，2013；梅胜全等，2010）、最速下降法（刘锋等，2012）来获得一个约束条件下接近于真实解的一个近似解.而在走时计算的整体实现步骤方面，目前普遍采用“向前处理”的按行列扫描计算（张东等，2009）、最短路径射线追踪算法的实现步骤（李培明等，2013；梅胜全等，2010）等.针对以上两个方面的关键技术，笔者充分利用平面波双线性假设的结论，推导出了双线性法求解极小值超越方程的解析计算公式，同时构造性地证明了该解析计算公式的正确性和相对应的地震波传播的物理规律.此外，笔者还借鉴快速推进法中的迎风差分格式的构建思想，提出了迎风双线性插值的局部走时计算的实现策略，并以快速推进法中的窄带技术（这是一种满足地震波传播基本规律的、计算效率高且灵活的波前扩展算法）作为算法的整体实现策略.最后的精度分析和计算实例表明了本文算法的正确性和在面对复杂介质时的稳定性

和有效性.

2 算法的描述

为了获得一种精度高且兼顾效率的三维走时计算方法，本文将提出一种快速推进迎风双线性插值法的三维地震波走时算法.如图1所示，这是一种网格算法，和以往基于网格的双线插值法、有限差分法、最短路径追踪法等类似，该算法主要包括三个方面的核心内容：局部计算公式、局部实现策略以及整体实现策略.

图1 算法的描述Fig.1 The description of the algorithm

如图1所示，局部计算公式的建立可描述为如下问题：在计算的某一个时刻一个平面上的点A、B、C、D的走时值已知，且该平面内的走时值为双线性分布，则如何构建计算与该平面邻近的网格节点F的走时值的计算公式；局部实现策略则可以描述为如下问题：在计算过程中，与走时值待求的F点相邻近的、类似于平面ABCD的、走时值已知且分布满足双线性分布的平面可能会有很多种情况，而在不同情况时应该分别采用怎样对应的不同的处理方法；最后，整体实现策略可以描述为如下问题：从源点S处给定的初始条件出发，应该采用怎样的实现步骤或过程来逐次地计算整个计算空间内所有网格节点的走时值.

与如上描述的三个问题相对应，这实际上涉及到算法的局部计算公式、算法的局部实现策略、算法的整体实现策略，所以接下来本文将分别详细阐述这三方面问题.最后，还将通过算法的精度分析和计算实例来验证算法的精度、正确性和有效性.

3 算法的局部计算公式：双线性插值解析公式

为了计算三维地震波走时，首先需要推导建立局部走时计算公式.如上所述，局部走时计算公式的推导可描述为：如图2所示，设点A、B、C、D的走时值已知，分别为：tA、tB、tC、tD；正方体网格的网格间距为h；假设平面ABCD上的走时值为双线性分布；求F点的走时值tF的计算公式.

图2 算法的局部计算公式：双线性插值公式Fig.2 The local computation