基于快速推进法的三维随机介质地震波走时计算

合集下载

地震波初至旅行时计算策略

收稿日期:2007-12-18基金项目:新世纪优秀人才支持计划(NCET -07-0845)、国家高技术研究发展计划(2006AA06A202)、中国石油天然气集团公司石油科技中青年创新基金(05E7028)和中国石油天然气集团公司物探重点实验室开放基金(GPK L0802)联合资助作者简介:杜启振(1969-),男(汉族),山东单县人,教授,博士,博士生导师,主要从事地震波全弹性波传播理论与偏移成像方法研究。

文章编号:167325005(2009)022*******地震波初至旅行时计算策略杜启振,侯　波(中国石油大学地球资源与信息学院,山东东营257061)摘要:以波前扩展法为基础,提出了适合强速度纵横向变化的有限差分法地震波初至旅行时外推策略,即在计算时采用二阶精度的有限差分法对程函方程进行了近似计算。

在没有出现对负数开平方的情况下只计算三个可能的体波,当出现已知条件不满足或算法不稳定(如对负数开平方)的情况时就计算可能的3个首波和两个散射波,两种情况下都取其最小值作为地震初至旅行时。

采用二分法的旅行时查找方法,完成了初至旅行时算法的改进。

复杂模型试算结果表明,该方法提高了初至旅行时计算的精度和稳定性。

关键词:有限差分法;程函方程;初至旅行时;二分法;地震波中图分类号:P 631.443 文献标识码:AStra teg i es of se is m i c wave f i rst break travel 2ti m e ca lcul a ti onDU Q i 2zhen,HOU Bo(College of Geo 2R esources and Infor m ation in China U niversity of Petroleum ,D ongying 257061,China )Abstract :A strategy of seis m ic wave first break travel 2ti m e calculati on with finite difference method was p r oposed on the basis of wavefr ont s p read method .The strategy can adap t t o the str ongly lateral and /or vertical vel ocity variati on media .The first break travel 2ti m es were calculated with the second order finite difference app r oxi m ati on in s olving the eikonal equati on .W ithout the circu m stances of extracting the square r oots of the negative nu mber,only three p r obable body waves were calculated .W hether the given conditi ons are not satisfied or the alg orith m is not stabilized (extracting the square r oot of the negative nu m 2ber ),three p r obable body waves and t w o scattering waves must be calculated si m ultaneously .The m ini m izati on of travel 2ti m es was selected as the seis m ic wave first break travel 2ti m e in both cases listed above .The dichot o my was e mp l oyed in l ooking for the l ocati on of wavefr ont point .Thus,the alg orith m of first break travel 2ti m e calculati on was i m p r oved .The method was app lied successfully with a co mp licated model in the seis m ic wave first break travel 2ti m e calculati on,which sho ws its high p recisi on and stability .Key words :finite difference method;eikonal equati on;first break travel 2ti m e;dichot omy;seis m ic wave 波动方程法模拟波场[124]虽然精确,但计算量太大。

反射波法隧道超前地质预报地震波走时模拟研究

反射波法隧道超前地质预报地震波走时模拟研究徐若格;孙辉;陈瑞;任洪涌;张健;李猛;高福柳;赵晓彦【期刊名称】《工程地质学报》【年(卷),期】2024(32)2【摘要】隧道施工安全常常由于复杂的地质情况面临着巨大的挑战,反射波法地震勘探具有探测深度大、抗干扰能力强等优点,是确保隧道施工安全的重要超前地质预报手段。

反射地震波走时计算与地震数据处理中的纵横波分离、速度建模、偏移成像等多个步骤息息相关,直接影响着隧道地质预报的结果。

本文将仿照反射波法隧道超前地质预报的实际施工区域进行速度建模,以求解程函方程为基础,结合快速推进法(FMM)中的迎风有限差分方法以及窄带扩展技术实现针对隧道模型的初至地震波走时模拟。

并通过对隧道模型中反射界面位置信息的提取,构建反射波初始窄带以及结合窄带扩展技术实现针对隧道模型反射波走时的计算。

通过数值分析以及复杂隧道模型计算表明,文中给出的反射波法隧道超前地质预报地震波走时模拟具有良好的计算精度以及较强的适应能力。

【总页数】9页(P623-631)【作者】徐若格;孙辉;陈瑞;任洪涌;张健;李猛;高福柳;赵晓彦【作者单位】西南交通大学地球科学与工程学院【正文语种】中文【中图分类】TU45【相关文献】1.基于TRT6000弹性波反射法、地质雷达法、地质编录的隧道综合超前预报技术研究2.基于TRT6000弹性波反射法、地质雷达法、地质编录的隧道综合超前预报技术研究3.基于TRT6000弹性波反射法、地质雷达法、地质编录的隧道综合超前预报技术研究4.隧道三维地震反射波法超前地质预报偏移成像应用研究5.地震反射波法技术及其在隧道超前地质预报中的应用研究因版权原因，仅展示原文概要，查看原文内容请购买。

地震波初至走时的计算方法综述_赵烽帆

地震波初至走时的计算方法综述
ｗｏｆｔｈｅｔｒａｖｅｌｌｃｕｌａｔｉｏｎｍｅｔｈｏｄｓｏｆｓｅｉｓｍｉｃｆｉｒｓｔｂｒｅａｋＡｒｅｖｉｅｔｉｍｅｃａ－
ＣＥＮＣＢＣ１．ａｋｅｒｋｓｒ，５，ｈｉｎａａｒｔｈｕｅｔｗｏｅｎｔｅｅｉｉｎ０００４ｈｉｎａｑｊｇ１
，ａｔｏｒｓｈｅｒｉｃｔｉｏｎ，ｔｕｔｅｓ，ｓｅｍｃｅｓ，９２．ＳｔａｔｅＫｅＬａｂｏｒｏＬｉｔｈｏＥｖｏｌｕＩｎｓｔｉｏＧｅｏｌｏａｎｄＧｅｏｈｓｉｃＣｈｉｎｅＡｃａｄｅｏＳｃｉｅｎＢｅｉｉｎ０００２Ｃｈｉｎａｙｐｙｙｙｆｆｇｐｙｆｊｇ１３．ＵｎｉｖｅｏＣｈｉｎｅＡｃａｄｅｏＳｃｉｅｎＢｅｉｉｎ０００３Ｃｈｉｎａｒｓｉｔｓｅｍｃｅｓ，９，ｆｆｊｇ１ｙｙ
，ＡｂｓｔｒａｃｔＩｎｔｈｅｓｅｉｓｍｉｃｗａｖｅｆｉｅｌｄｔｈｅｆｉｒｓｔａｒｒｉｖａｌｔｒａｖｅｌｔｉｍｅ－，ｎｉｍｏｒｔａｎｔｒｏｌｅｉｎｔｈｅｆｉｅｌｄｏｆｓｅｉｓｍｏｌｏｗｈｉｃｈｉｓｄｕｅｌａｓａｐｇｙｐｙｔｏｔｈｅｆｉｒｓｔａｒｒｉｖａｌｈａｓｅｓｃａｎｂｅｅａｓｉｌｔｒａｃｅｄａｎｄｉｄｅｎｔｉｆｉｅｄ．ｐｙＴｈｅｓｅｉｓｍｉｃｆｉｒｓｔａｒｒｉｖａｌｔｒａｖｅｌｅｌｄｉｓｗｉｄｅｌｕｓｅｄｉｎｔｉｍｅｆｉ－ｙ，ｖ，ｓｓｅｉｓｍｉｃｒｅｓｔａｃｋｍｉｒａｔｉｏｎｅｌｏｃｉｔａｎａｌｓｉｓｅｉｓｍｉｃｔｒａｖｅｌ－ｐｇｙｙａｎｄｅａｒｔｈｕａｋｅｌｏｃａｔｉｏｎ．Ｔｈｉｓａｅｒｍａｉｎｌｍｏｒａｈｔｉｍｅｔｏｑｐｐｙｇｐｙｉｎｔｒｏｄｕｃｅｓｏｕｒｉｃａｌｅｉｓｍｉｃｉｒｓｔｒｒｉｖａｌｒａｖｅｌｔｉｍｅｆｔｓｆａｔ－ｙｐ，，（ｍｅｔｈｏｄｂａｓｅｄｏｎｌｃｕｌａｔｉｏｎｍｅｔｈｏｄｓｔｈａｔｉｓ１）ｓｅｉｓｍｉｃｒａｃａｙ，ｕｓｅｄｍｅｔｈｏｄｓａｒｏｘｉｍａｔｉｏｎｔｈｅｃｏｍｍｏｎｌｔｈｅｈｉｈｃｆｒｅｕｅｎ－ｐｐｙｇｙｑａｒｅＳｈｏｒｔｅｓｔＰａｔｈＭｅｔｈｏｄ（ＳＰＭ）ａｎｄＭｏｄｉｆｉｅｄＳｈｏｒｔｅｓｔＰａｔｈ；（Ｍｅｔｈｏｄ（ＭＳＰＭ）２）ｔｈｅｍｈｏｄｂａｓｅｄｏｎｔｈｅｎｕｍｅｒｉｃａｌｅｔ，ｗｉｎｃｌｕｄｅＦｉｎｉｔｅｓｏｌｕｔｉｏｎｏｆｅｉｋｏｎａｌｅｕａｔｉｏｎｈｉｃｈｍａｉｎｌｑｙ，Ｍｅｔｈｏｄ（ＦＭＭ）ａｎｄＤｉｆｆｅｒｅｎｃｅＭｅｔｈｏｄ（ＦＤ）ＦａｓｔＭａｒｃｈｉｎｇ；（Ｍｅｔｈｏｄ（ＦＳＭ）３）ＷａｖＦａｓｔＳｗｅｅｉｎｒｏｎｔＣｏｎｓｔｒｕｃｔｉｏｎｅｆｐｇ（；ＷＦＣ）ｍｅｔｈｏｄｂａｓｅｄｏｎｔｈｅＨｕｅｎｓｒｉｎｃｉｌｅａｎｄ（４）ｔｈｅｙｇｐｐｄｏｍａｉｎｗａｖｅｅｕａｔｉｏｎ（ＦＷＱ）．ｍｅｔｈｏｄｂａｓｅｄｏｎｔｈｅｆｒｅｕｅｎｃｑｙｑａｎｄｂｅｔｔｅｒＴｈｅｆｉｒｓｔｏｎｅｈａｓｈｉｈｅｒｃｏｍｕｔａｔｉｏｎａｌａｃｃｕｒａｃｇｐｙ，ｗ，ｗｒｅｓｕｌｔｉｎｈｉｌｅｉｔｎｅｅｄｓｍｏｒｅｒｉｄｏｉｎｔｓｈｉｃｈｍａｓｔａｂｉｌｉｔｇｐｙｙｌｏｗｃｏｍｕｔａｔｉｏｎａｌｅｆｆｉｃｉｅｎｃ．Ｔｈｅｃａｌｃｕｌａｔｉｏｎｏｆｒａａｔｈｓｉｓｐｙｙｐ，ｓｎｏｔｒｅｕｉｒｅｄｉｎｔｈｅｓｅｃｏｎｄｏｎｅｏｉｔｈａｓａｄｖａｎｔａｅｓｉｎｑｇ，ｂ，ｓｃｏｍｕｔａｔｉｏｎａｌｅｆｆｉｃｉｅｎｃｔａｂｉｌｉｔａｎｄｒｅａｌｉｚａｔｉｏｎｕｔｔｈｅｐｙｙ，ｗｉｓｌｏｗｅｒｈｉｃｈｃａｎｂｅｉｍｒｏｖｅｄｂｃｏｍｕｔａｔｉｏｎａｌａｃｃｕｒａｃｐｙｐｙｔｈｅｈｉｈｉｎｔｒｏｄｕｃｉｎｉｎｉｔｅｄｉｆｆｅｒｅｎｃｅｓｃｈｅｍｅ．Ｔｈｅｔｈｉｒｄｏｒｄｅｒｆ－ｇｇ，ｏｎｅｃａｎｃａｌｃｕｌａｔｅｔｒａｖｅｌｔｉｍｅｗｉｔｈｈｉｈａｃｃｕｒａｃａｎｄｓｔａｂｉｌｉｔｙｇｙｗｈｉｌｅｉｔｒｅｕｉｒｅｓｔｈｅｒｉｄｔｒａｎｓｉｔｉｏｎｂｅｔｗｅｅｎｒａｒｉｄａｎｄｒｅｕｌａｒｑｇｙｇｇ，ｗｒｅｓｕｌｔｉｎｌｏｗｃｏｍｕｔａｔｉｏｎａｌｅｆｆｉｃｉｅｎｃ．Ｔｈｅｒｉｄｈｉｃｈｍａｐｙｇｙｌａｓｔｏｎｅｃａｎａｄａｔａｎｃｏｍｌｅｘｍｅｄｉｕｍ，ｂｕｔｔｈｅｃｏｍｕｔａｔｉｏｎａｌｐｙｐｐａｃｃｕｒａｃａｎｄｅｆｆｉｃｉｅｎｃａｒｅｌｏｗｅｒ．ｙｙ；；Ｋｅｗｏｒｄｓｆｉｒｓｔｂｒｅａｋ；ｓｈｏｒｔｅｓｔａｔｈｅｉｋｏｎａｌｅｕａｔｉｏｎｆｉｎｉｔｅｐｑｙ；ｗ；ｆ；ｆｄｉｆｆｅｒｅｎｃｅａｓｔａｒｃｈｉｎａｓｔｗｅｅｉｎａｖｅｆｒｏｎｔｍｓｇｐｇ；ｄｏｍａｉｎｗａｖｅｅｃｏｎｓｔｒｕｃｔｉｏｎｆｒｅｕｅｎｃｎｕａｔｉｏｑｑｙ

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

基于快速推进迎风双线性插值法的三维地震波走时计算孙章庆;孙建国;岳玉波;江兆南【摘要】三维地震波走时计算技术是三维地震反演、层析成像、偏移成像等诸多地震数据处理技术中非常重要的正演计算工具.为了获得精度高且兼顾效率的三维走时计算方法:首先,在常规双线性插值公式推导过程中,充分利用平面波双线性假设的结论,获得了二元极小值超越方程的解析解,进而推导出了准确的局部走时计算公式,同时构造性地证明了该计算公式满足地震波的传播规律和Eikonal方程;其次,引入迎风差分的基本思想,提出迎风双线性插值的局部走时计算策略,该计算策略能简化算法、提高效率且保证无条件稳定性;然后,将上述计算公式和迎风双线性插值策略与常规快速推进法中的窄带技术结合,获得了一种新的基于快速推进迎风双线性插值法的三维地震波走时计算方法;最后,通过精度和效率分析检验了新算法的精度、效率和正确性,并通过计算实例验证了算法在面对复杂介质时的稳定性和有效性.【期刊名称】《地球物理学报》【年(卷),期】2015(058)006【总页数】13页(P2011-2023)【关键词】三维;平面波双线性假设;迎风双线性插值;窄带技术;地震波走时计算【作者】孙章庆;孙建国;岳玉波;江兆南【作者单位】吉林大学地球探测科学与技术学院,长春 130026;国土资源部应用地球物理综合解释理论开放实验室—波动理论与成像技术实验室,长春130026;吉林大学地球探测科学与技术学院,长春 130026;国土资源部应用地球物理综合解释理论开放实验室—波动理论与成像技术实验室,长春130026;东方地球物理公司物探技术研究中心,河北涿州 072751;重庆市地质矿产勘查开发局208水文地质工程地质队,重庆408300【正文语种】中文【中图分类】P6311 引言地震波走时是描述地震波运动学特征的一个很重要的物理参数，其刻画了地震波由震源点激发后到达地下介质空间中各个位置坐标点所需耗消的走时，而该走时和地下介质的速度参数的分布情况密切相关，所以地震波走时计算方法常常被应用于走时反演、层析成像、偏移成像、速度分析等一些重要的地震数据处理技术中.地震波走时计算的精度、效率、稳定性等性质直接影响着这些地震数据处理的效果和效率（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 computationformulas of algorithm：the bilinear interpolation formulas为了推导计算tF的公式，首先以点A为坐标原点，向量AB方向为x轴正方向，向量AD方向为y轴正方向，向量AF方向为z轴正方向，建立局部笛卡尔坐标系，则有点A、B、C、D、F的坐标分别为A（0，0，0），B（h，0，0），C（h，h，0），D（0，h，0），F（0，0，h）；其次，假定到达F点的射线经过了走时值分布为已知的平面ABCD内的任意一点E，设E点的坐标为 E（x，y，0），因为E 点被限定位于平面ABCD内，则有0≤x≤h，0≤y≤h.过E点作平行于y轴的直线，该直线与线段AB、CD分别相交于E1 和E2 点，则点E1、E2 的坐标分别为E1（x，0，0）、E2（x，h，0）.根据平面ABCD上的走时值为双线性分布的假设，线段AB、CD及E1E2上的走时值均为线性分布，则有对上述（1）式进行整理可得其中k1 ＝（tB－tA）/h，k2 ＝（tD －tA）/h，k3 ＝（tA＋tC－tB－tD）/h2.设在局部正方体网格单元内地震波传播速度为均匀，并基于上述到达F点的射线经过平面ABCD内的E（x，y，0）点的假定，可得其中s为当前正方体网格单元内的地震波传播速度的倒数（即慢度）.将式（2）代入（3）中，则有根据Fermat原理，只要确定能使得（4）式取最小值的E（x，y，0）点的位置，即可获得计算tF的公式，所以对（4）式tF分别关于变量x，y求偏导数，并令该偏导数为零，即分析方程组（5）可知：只要求出该方程组的解，即可获得E（x，y，0）点的位置坐标，然后再将其带回到（4）式即可获得计算tF的公式.但是，一些研究已表明方程组（5）实际上是一个超越方程组，无法求出相应的解析解（其可能解分布在一条双曲线上），所以一些学者分别采用网格界面剖分法（张东等，2009）、快速插值法（李培明等，2013；梅胜全等，2010）、最速下降法（刘锋等，2012）来获得一个约束条件下接近于真实解的一个近似解.与以往的双线插值法的近似处理不同，在此笔者将充分利用双线性假设的结论，来获取方程组（5）的精确解析解，并构造性地证明计算tF的公式的正确性以及其背后满足的地震波传播波的物理事实.根据走时插值法的提出者Asakawa的线性假设（Asakawa and Kawanaka，1993），实际上即是平面波的假设，其阐述的是如果当前传播的为平面波，则该波传播经过的二维空间中的走时分布在二维空间中的任意一条网格线上均为线性，而在三维空间中的一个平面上的走时分布则均为双线性.同样基于平面波的假设，并假定在局部的相对较小的一个网格单元内地震波的传播速度为匀速.如图3所示，则无论从那个方向（假定箭头方向为平面波的任意一个传播方向）途经平面ABCD均有tD－tA＝tC－tB，因为如图3所示，设该平面波与平面ABCD相交于如图3所示的一系列平行的虚线，其中经过点A、B、C、D的虚线分别位于平面波tA、tB、tC、tD 时刻的波前上，因为这些波前是相互平行的，并且网格四边形ABCD为一正方形，所以必有tA时刻的波前传到tD时刻花费的走时等于tB时刻的波前传到tC时刻花费的走时，也即是tD－tA＝tC－tB. 图3 平面波假设分析Fig.3 The analysis of plane wave assumption基于上述双线性假设，即平面波假设的进一步结论tD－tA＝tC－tB，可以得出k3＝0.将该结论带回到（5）式中，则（5）可简化为再分析（6）式，若在k1、k2均不为零的情况下，则用下式除以上式再整理有：y ＝k2x/k1，即方程组（6）的解实际上是位于一条直线上（而并非一条双曲线上），并且将该直线代回方程组（6）中的任意一个方程，均可以获得解析的x、y的解.不过要获得各种情况下的方程组（6）的解，还得分情况讨论k1、k2取不同值时，方程组（6）的解的情况：①k1＝0且k2＝0时此时，将k1、k2代回方程组（6）可得：x＝y＝0，再将该结论代回（4）式中可得计算tF的公式为②k1＝0且k2＜0时在此种情况下，将k1＝0代入方程组（6）可得x＝0，再将该结论代回方程组（6）中可将其简化为融入0≤y≤h的限定条件，求解（8）式可得其中Δt＝tD－tA，且需满足：（此是基于0≤y≤h的限定条件）.综合x＝0的结论，可知 E点的坐标为0），再将其代入（4）式中即可获得此种情况下计算tF的公式：③k2＝0且k1＜0时这种情况与情况②相似，不同之处仅在于交换了变量x与y的角色，所以同理情况②的推导过程可得：其中：Δt＝tB－tA，且需满足：（此是基于0≤x≤h的限定条件）.同时E点的坐标为计算tF 的公式为：④k1＜0且k2＜0这种情况下，直接求解方程组（6），即将上述结论代入方程组（6）的上式，并融入0≤x≤h、0≤y≤h的限定条件，则可得方程（6）的解：其中，Δt1＝tB－tA、Δt2＝tD－tA，且需满足：0≤（此是基于0≤x≤h、0≤y≤h 的限定条件）.将（14）式代回（5）式中，即可获得计算tF的公式：同时，此时E点的坐标为：⑤k1＞0或k2＞0这种情况下，因为上述将E点限定于平面ABCD内的限定条件有0≤x≤h、0≤y≤h，所以很容易分析出，此时方程组（6）无解.这种无解的情况表面上表明：上述计算公式可能出现不稳定的情况.不过，这种所谓的不稳定的情况，在将上述算法融入一定的局部实现策略和算法的整体实现策略后，是不可能出现的，关于该问题将在后续内容中阐述.经过以上的推导，实际上已经获得了当k1、k2取不同值的情况时，点E（x，y，0）的位置坐标和tF的计算公式.如图4所示，下面将构造性地证明这些公式的正确性，以及它们对应的不同情况下的地震波传播的物理事实.图4 双线性插值公式的分析Fig.4 The analysis of bilinear interpolation formulas上述情况①：当k1＝0且k2＝0时，即有tB＝tA，tD＝tA，也即tA＝tB＝tD，该式对应的物理事实为点A、B、D位于地震波传播过程中的同一个等时面上.此外，根据不在同一条直线上的三点确定一个平面的数学事实和上述基于的平面波假设，则有与A、B、D点位于同一平面上的C点也位于该等时面上，所以可以进一步得出tA＝tB＝tC＝tD，也即平面ABCD为平面波的一个等时面.同时，根据局部笛卡尔坐标系可知，直线AF与平面ABCD垂直，再根据射线路径和等时面相互垂直的物理规律，可以得出在这种情况下平面波沿着z轴正方向传播，而经过平面ABCD到达F点的射线与平面ABCD相交于A点，所以此时自然有：x＝y＝0，tF＝tA＋sh，也就是说上述情况①得的方程组（6）的解和tF的计算公式符合在这种情况下平面波的传播规律，即是正确的.上述情况②：当k1＝0且k2＜0时，即有tB＝tA，tD＜tA，此时再根据双线性假设，因为tB＝tA，所以平面ABCD上的走时分布沿x方向没有变化，而仅沿着y方向有变化，所以此时的插值问题实际上已简化为二维插值问题，而E点则被固定到距F更近的线段AD上，即有x＝0.此外，基于二维空间中Asakawa提出的走时线性插值法同样可以获得其中Δt＝tD－tA.再分析tF的计算公式，实际上其可以变形为该方程实际上是二维Eikonal方程（∂t/∂z）2＋（∂t/∂y）2＝s2在网格上以A点为中心的差分离散方程，也就是说此处得到的tF的计算公式，在局部上满足Eikonal方程，所以其是正确的.同时上述0的限定条件还能保证tF的计算公式不会出现负数开平方的不稳定情况.上述情况③：当k2＝0且k1＜0时，此时的情况与情况②是类似的，也可以简化为一个二维插值问题，不同之处仅在于坐标方向x与y交换了角色，所以在此不再重复阐述.上述情况④：当k1＜0且k2＜0时，即有tB＜tA，tD＜tA，在这种情况下插值问题才真真变为一个三维双线性插值问题，在此条件下方程组（6）有解.求解该方程组即可获得E点的坐标为x＝－hΔt1/同时将它们代回（5）式得到tF的计算公式：tF＝tA＋其中Δt1＝tB－tA、Δt2＝tD－tA. 同样其可以变形为该方程实际上是三维 Eikonal方程（∂t/∂x）2 ＋（∂t/∂y）2 ＋（∂t/∂z）2＝s2在网格上以A点为中心的差分离散方程，也就是说此处得出的tF的计算公式在局部上是满足Eikonal方程的，所以其是正确的.此外，再仔细分析E点的坐标可得：y＝kx，其中k＝Δt2/Δt1，也就是说，实际上E点是位于坐标平面xAy内的一条直线上，该直线的斜率为k.k的取值实际上和已知条件tA、tB、tD的取值情况有关：k＝Δt2/Δt1 ＝（tD－tA）/（tB－tA）.如图4所示，tD ＜tB时有k＞1，此时直线途径ACD半区域，也即E点位于走时分布较小（因为双线假设和tD＜tB）的区域.这实际上与Fermat原理和地震波传播的基本规律是相符合的，即地震波总是从走时值相对更小的区域传向下一个区域.而当tB＜tD时有k＜1，此时有如上同样的结论.最后当tB＝tD时，则k＝1，此时E点位于正方形ABCD的对角线AC 上.综上所述，在情况④下得出的E点的坐标位置是满足地震波的传播规律和Fermat原理的，tF的计算公式在局部上满足Eikonal方程.同时0≤2s2 h2/3的限定条件，还能保证tF的计算公式不会出现负数开平方的不稳定情况.上述情况⑤：当k1＞0或k2＞0时，即有tB＞tA，tD＞tA，此时上述的结论为无解.实际上，再回顾方程组（6），此时并非方程组（6）无解，只是解的情况为x＜0，y＜0，即E点不在平面ABCD内.实际上，为了公式推导的方便和局部网格计算的限定，将E点限定于走时值已知的平面ABCD内，而在该限定条件下x＜0，y＜0的解不能满足该限定条件，所以致使方程组（6）无解.事实上，x＜0，y＜0的解的情况说明：在这种情况下E点位于其他的邻近的网格单元内.同时tA＜tB、tA＜tD，再加上双线性假设，暗含了其他相邻的网格单元上的分布的走时值会比网格单元ABCD上的更小.而根据Fermat原理，此时计算tF的公式则应该通过其他网格单元上推导获得，也即在计算时应该采用其他网格单元进行计算.该问题实际上引伸出了局部走时计算策略的问题：如图5所示，在进行三维局部走时计算时，包围被算点F周围的类似于ABCD的网格单元实际上有24个，从几何的角度讲，到达F点的射线必与用于封闭F点的类似于ABCD的24个网格单元中的一个相交，而该交点正是我们上述需要求解的E点.当然要是在24个网格单元中逐次搜寻计算E点，必定是一种繁琐且计算效率低下的算法.为了解决该问题，此处借鉴迎风差分的思想，提出了一种非常简洁的迎风双线性插值的局部实现策略，下面将详细阐述该策略.4 算法的局部实现策略：迎风双线性插值法如图5所示，在三维走时计算的某一时刻，当前需要计算的是点F的走时值tF.如图2所示，上述推导了网格单元ABCD走时分布为已知的情况下点F的走时值tF 的计算公式.但实际上，在被算点F周围的24个类似于ABCD的网格单元中，走时分布为已知的网格单元可能有很多个，那么究竟应该采用那个网格单元作为已知条件是合理的和稳定的？在此，引入迎风差分的基本思想（Sethian and Popovici，1999），该思想是指地震波总是从地震波走时值分布更小的区域向走时值未知的区域传播，而在局部数值计算实现时就应该迎着走时值分布更小的区域去做差分代替微分的近似.迎风差分格式的基本思想实际上是Huygens原理和Fermat原理的综合利用，主要体现在其隐含的将被算点周围的走时值最小的网格节点当作子震源点，也正因为此迎风差分格式具有无条件稳定性.同样借鉴该思想，在此提出迎风双线性插值法，该方法的核心借用迎风差分格式选取最小邻近走时点的方法，选取被算点周围的最小走时分布网格单元，然后再采用类似于图2所示的推导的公式进行局部走时计算.总体上讲，实际就是要构建相应的迎风双线性插值公式.如图5a所示的局部正方体网格，以点F为中心的大正方体的六个面分别由以A1、A2、A3、A4、A5、A6为中心的4个类似于图2中的ABCD的网格单元组成，这样就总共有24个类似于ABCD的网格单元.在当前计算点F的走时值tF的时刻，这24个网格单元上的26个网格节点的走时值可能为已知，可能为未知.由后续将阐述的算法的整体实现策略的初始化步骤（详见后续阐述内容的“算法的整体实现步骤：快速推进法”部分）可知，若网格节点的走时值为未知（即没有被计算出来），其值为固定的一个非常大的数（这个数为人为规定的，其值远大于所有网格节点最终可能算出的走时值），若网格节点的走时值被计算完成，其值即为最终的走时计算结果.在一般情况下，当前被算点F周围的26个网格节点的走时值为一部分被计算完成，一部分未开始计算（如图6a所示）.而在被计算完成的部分里面，先完成的走时值一定是小于后完成的走时值（该结论由后续内容的窄带技术的阐述可知），也就是说与点F邻近的26个网格节点均有一个走时值，而用于计算tF的网格单元的走时分布应为最小.具体挑选最小走时分布区域的方法为：①确定与点F直接相邻的6个网格节点的走时值最小的网格节点Amin，其中有 tAmin ＝min（tA1，tA2，tA3，tA4，tA5，tA6）；② 在以Amin为中心的平面上分别确定其他两个方向上的最小走时网格节点Bmin，Dmin，若设Amin＝A5则分别有tBmin ＝min（tB51，tB52），tDmin ＝ min（tD51，tD52）；③Amin，Bmin，Dmin所在的网格单元即为最终的类似于公式推导过程中的网格单元ABCD，与该网格单元相对应，tF的计算公式相应地变为图5 算法的局部实现策略：迎风双线性插值Fig.5 The local implementation strategy of algorithm：the upwind bilinear interpolation其中Δt1＝tBmin－tAmin，Δt2＝tDmin－tAmin.公式（16）即为最终的迎风双线性插值公式，与常规的双线性插值公式相比，该公式实际上是一种优化的公式.该公式隐含着一个优化的局部实现策略，而该优化的局部实现策略将Huygens原理和Fermat原理融入到算法的局部计算过程中，使得算法是在满足地震波传播的基本规律的基础上实现的.与迎风差分格式类似，此处的迎风双线性插值法也具有无条件稳定性.5 算法的整体实现策略：窄带技术上述内容分别建立了算法的局部走时计算公式和局部实现策略，实际上要想完成整个计算区域的地震波走时计算，还需要拟定一个算法的整体实现步骤，即算法的整体实现策略.算法的整体实现策略往往决定着算法的稳定性和计算效率.基于对计算效率、无条件稳定性、对迎风双线性插值的适应性等因素的综合考虑，在此引入快速推进法中的窄带技术（Sethian and Popovici，1999）作为算法的整体实现策略.如图6b所示，窄带技术的核心思想是采用一个包含当前所有波前点的窄带来近似模拟地震波前，通过窄带的扩展演化来模拟地震波前的传播过程，其充分考虑了地震波前传播过程中的熵守恒理论，同时此处因为迎风差分思想的融入，所以窄带技术具有无条件稳定性.窄带技术在实现过程中将计算空间中的网格节点分为三种类型：完成点，其走时计算已经完成的点，在图6b中用黑色填充的点表示；窄带点，当前窄带内的点，也即近似波前上的网格节点，在图6b中用灰色填充的点表示；远离点，窄带扩展还未到达的点，也即还未进行走时计算的网格节点，在图6b中用白色填充的点表示.图6 算法的整体实现策略：窄带技术Fig.6 The global implementation strategy of algorithm：the narrow band technique窄带技术的实现过程可以分为初始化和扩展计算两个部分.在初始化时，震源点为唯一的完成点，其走时值为零；与震源点直接相邻的正方体网格单元内的网格节点为初始窄带点，它们构成初始窄带，即初始波前，它们的走时值通过局部的解析公式求取；其它所有剩下的计算空间中的网格节点为远离点，其初始值设置为一个相对很大的数（远大于最终可能计算出来的最大的走时值即可，这主要是为了实现上。

基于最大瞬时输入能的时程分析地震波选择与调幅方法[发明专利]

专利名称：基于最大瞬时输入能的时程分析地震波选择与调幅方法
专利类型：发明专利
发明人：吴琛,项洪,陈思坚,杨超
申请号：CN201810184143.0
申请日：20180306
公开号：CN108416140A
公开日：
20180817
专利内容由知识产权出版社提供
摘要：基于最大瞬时输入能的时程分析地震波选择与调幅方法，包括：将一列地震波输入到结构中，由振型分解反应谱法计算三维有限元模型的顶点位移d、基底剪力Q和第一自振周期T；将Q、d 代入公式可反推得到输入三维有限元模型的统计意义上的最大瞬时输入能ΔE；输入的地震动按《抗规》的最大加速度值A进行调幅，输入到结构自振周期为T的单自由度体系中，计算得到最大瞬时输入能ΔE；代入公式可得到该条地震动基于统计意义上的最大瞬时输入能的最大加速度值A；地震动按最大加速度值A进行调幅，输入三维有限元结构中，进行弹性时程分析。

本发明针对地震波的选择与调幅，能较精确、高效地从地震波数据库选择与建筑场地条件相符合的地震波。

申请人：福建工程学院
地址：350000 福建省福州市闽侯县上街镇福州地区大学新校区学园路
国籍：CN
代理机构：福州市鼓楼区京华专利事务所(普通合伙)
代理人：宋连梅
更多信息请下载全文后查看。

复杂介质情况下的三维地震旅行时计算的开题报告

复杂介质情况下的三维地震旅行时计算的开题报告一、问题描述在地震勘探中，对于许多复杂的地质构造，需要使用三维地震旅行时计算方法来获得合理的成像结果。

对于复杂介质情况下的三维地震旅行时计算，需要考虑多种因素，包括波场方程、介质参数、计算方法等。

因此，本文将针对复杂介质情况下的三维地震旅行时计算进行研究。

二、研究目的本文旨在探究复杂介质情况下的三维地震旅行时计算方法，深入研究波场方程、介质参数和计算方法等关键因素，并通过数值模拟验证研究结果的可行性。

三、研究内容1.三维地震旅行时计算的波场方程：首先，本文将介绍在复杂介质情况下进行三维地震旅行时计算所需的波场方程，并对其进行深入研究。

针对波场方程中的一些关键参数，例如介质密度、纵横波速度等，将从理论和实验两个方面进行探究。

2.三维地震旅行时计算的介质参数：介质参数是影响三维地震旅行时计算精度的另一个重要因素。

本文将在探究波场方程后，进一步研究如何确定适合复杂介质的介质参数，并分析不同介质参数对三维地震旅行时计算结果的影响。

3.三维地震旅行时计算的计算方法：对于复杂介质情况下的三维地震旅行时计算，需要使用高效的计算方法来避免计算结果出现误差。

本文将探究不同的计算方法，并评估其在复杂介质情况下的适用性和精度。

4.数值模拟：最后，本文将通过数值模拟的方法，验证所得研究结果的可行性，并对比不同方法在模拟数据上的表现，以进一步确定采用何种方法获得最佳计算结果。

四、研究意义本文所研究的复杂介质情况下的三维地震旅行时计算方法，能够为地震勘探领域提供更加准确的成像结果，提高勘探效率。

同时，研究结果也对于深入理解地质构造和介质物性有所帮助。

五、研究方法本文将采用文献调研、数值模拟等方法进行研究。

六、进度安排本项目预计用时 6 个月，具体进度安排如下：第一周：确定研究课题并进行背景调研第二周-第六周：深入调研波场方程和介质物性相关文献第七周-第十周：研究三维地震旅行时计算的计算方法第十一周-第十四周：探究复杂介质情况下的三维地震旅行时计算的方法和模型第十五周-第十七周：进行数值模拟并整理结果第十八周-第二十周：撰写研究报告并进行修改和完善七、预期结果本文预期可以得出适用于复杂介质情况下的三维地震旅行时计算的方法，并通过数值模拟证明其可行性和高精度。

波前快速推进法三维走时计算技术

波前快速推进法三维走时计算技术
朱丹
【期刊名称】《大庆石油地质与开发》
【年(卷),期】2001(020)005
【摘要】三维叠前深度偏移简称3DPSDM,是解决复杂构造和强速变化条件下地
质体成像问题的有效工具,走时计算是其配套技术.波前快速推进算法是一种非常快速、准确稳定的计算走时的方法,不仅计算速度极快(大约比有限差分法快4～5倍),而且对任何复杂速度场都具有无条件稳定性.该算法可用于三维克希霍夫叠前偏移、三维克希霍夫叠后偏移、三维速度分析和三维克希霍夫正演模型计算等.根据此算法,在SGI并行机上研制了并行计算软件,并对其进行了基准点测试.同时,还给出了SEG/EAGE模型的计算实例.
【总页数】2页(P65-66)
【作者】朱丹
【作者单位】大庆油田有限责任公司勘探开发研究院
【正文语种】中文
【中图分类】P631.1+2
【相关文献】
1.波前快速推进法起伏地表地震波走时计算 [J], 孙章庆;孙建国;韩复兴;杨昊
2.基于快速推进迎风双线性插值法的三维地震波走时计算 [J], 孙章庆;孙建国;岳玉波;江兆南
3.基于完全三叉树的快速推进法地震波走时计算 [J], 王乾龙;孙建国;孙辉;黄兴国
4.基于快速推进与波前构建的联合3D走时计算（英文） [J], 孙辉;孙建国;孙章庆;韩复兴;刘明忱;刘志强;高正辉;石秀林
5.地震波走时场模拟的快速推进法和快速扫描法比较研究 [J], 兰海强;张智;徐涛;白志明;梁锴
因版权原因，仅展示原文概要，查看原文内容请购买。

地震波走时表及体波量规函数在软件设计中的使用方法

地震波走时表及体波量规函数在软件设计中的使用方法和跃时;孙文斌;孟宪森
【期刊名称】《防灾减灾学报》
【年(卷),期】2005(021)001
【摘要】本文介绍了通过计算机软件设计实现地震波走时表的自动查找和体波量规函数的插值计算,在台站处理地震和地震速报中取得了满意的效果.
【总页数】6页(P37-42)
【作者】和跃时;孙文斌;孟宪森
【作者单位】牡丹江地震台,黑龙江,牡丹江,157009;黑龙江省地震局,黑龙江,哈尔滨,150090;黑龙江省地震局,黑龙江,哈尔滨,150090
【正文语种】中文
【中图分类】P315.61;P315.69
【相关文献】
1.编制区域地震波走时表方法探讨 [J], 王秀娟
2.气枪震源地震波走时变化影响因素分析及地震相关的走时变化 [J], 周青云;刘自凤;贺素歌
3.三维多值走时地震波场重建及格林函数计算 [J], 张钋;杨长春;李幼铭
4.二滩遥测台网地震波速度模型的建立和走时表的编制 [J], 赵珠
5.甘表区域地震波走时表 [J], 张诚;高世垒
因版权原因，仅展示原文概要，查看原文内容请购买。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

基于快速推进法的三维随机介质地震波走时计算地球内部介质的自组织性广泛存在着,由于岩石矿物成分的非均匀性及成岩过程中诸如温度、压力的差异和构造运动的影响,地层介质的物性参数必然表现出非均匀性,即自组织性。

传统的建模方法很难描述这种非均匀性,随机介质建模的出现为描述小尺度非均匀性开辟思路。

波场模拟作为研究随机介质中地震波波场传播的工具,已然发挥很大作用。

然而随机介质的波场信息绕射较严重,不利于分析其规律。

地震波走时作为地震波运动学信息的一个重要参量,通过走时的计算研究波在随机介质中的传播规律较直观。

地震波走时的计算方法主要分为基于程函方程的有限差分法和基于射线理论的射线追踪法,通过有限差分直接求解由波动方程高频近似得到的程函方程获得地震波走时是一种高效准确的方法。

针对常规有限差分法的不足,快速推进法(FMM)采用迎风差分解决算法的不稳定问题。

窄带扩展和堆排序是实现快速推进法的两个关键技术。

窄带扩展用于近似模拟波前的传播,堆排序用于选取窄带中的最小走时点,通过两者有机地结合精确高效完成整个空间的计算。

由于窄带扩展每次都需要选取最小走时点,因此堆排序的效率直接影响着FMM的计算效率。

一般地,堆排序算法是基于完全二叉树的。

事实上,这种堆排序并不是最优的,而且堆排序的方式仍可改进以契合迎风差分计算走时的特点。

针对随机复杂介质建模和快速推进法地震波走时计算两方面问题,本文对三维随机介质速度建模和快速推进法中的问题开展了相关研究,针对性地提出了有效的改进方案。

论文的具体研究内容是:(1)将基于统计学随机过程的谱分解理论“叠加型”速度建模方法由二维扩展到三维,重点研究三维随机复杂介质建模过
程。

对三维随机介质建模作了补充,详细深入地分析自相关函数、三个方向上的自相关长度和均方差对三维(平稳)随机介质建模结果的影响;(2)对FMM中的堆
排序做了改进,将传统基于完全二叉树的堆排序扩展为基于完全三叉树的堆排序,并在完全三叉树堆排序的基础上仍采用空位下沉法,有效地结合迎风差分法的特点。

然后将改进后的算法应用到均匀模型中,对比分析改进前后的计算效率。

(3)将改进的FMM应用到三维复杂随机介质的地震波走时计算中,详细阐述
了三维迎风双线性快速推进法地震波走时计算的实现方法,主要包括三维局部算法和窄带扩展的实现策略。

通过三维均匀介质模型试算对算法的计算精度和计算效率进行分析,然后对三维复杂介质模型(含随机介质)进行计算并给出具体分析。

通过模型计算,“叠加型”随机复杂介质建模的实现方案是可行的,丰富了现有随机介质建模方法,在一定程度上提高了速度模型的复杂度,更加贴切地描述
实际地层的真实信息。

基于完全三叉树堆排序的快速推进法相比于基于完全二叉树堆排序快速推进法更高效,且无论在二维还是在三维随机复杂介质中地震波走时计算都是可行的、稳定的。