基于卷积完全匹配层的旋转交错网格高阶差分法模拟弹性波传播
基于交错网格有限差分弹性波正演模拟及波场特征分析
基于交错网格有限差分弹性波正演模拟及波场特征分析【摘要】为研究和认识多种储层中弹性波的波场特征，以利于多波地震资料解释，高精度数值模拟是有效的方法之一。
本文在弹性波方程基础上，采用高阶交错网格有限差分技术模拟地震波在各向同性介质和各向异性介质中的传播，可得到不同类型介质的弹性波场。
同时，文中也分析了各向异性系数对多波波场特征的影响。
通过对高精度数值模拟得到的波场快照对比研究表明，该方法可有效获得高精度弹性波正演结果，为研究各种复杂介质中弹性波的波场特征和传播规律奠定了基础。
【关键词】多波多分量波场特征各向异性弹性波正演1 引言随着油气田勘探技术的不断发展[1][2]，人们对地震资料的认识也不断加深，纵波地震资料在含油气的显示上存在一些不确定性，单一纵波资料解释的多解性问题尤为突出。
在地震勘探领域中，过去一直把各向同性弹性体理论作为研究地下介质的前提，但是在实际地层中普遍存在各向异性，地下介质的各向异性（如周期薄互层引起的各向异性、以及裂隙引起的各向异性）产生的弹性波场与各向同性介质产生的弹性波场存在着不可忽略的差异。
由此，多波地震勘探作为油储地球物理的主要方法之一应运而生。
在多波资料解释过程中，要求搞清楚储层的岩性与多波的波场特征之间的关系，因此，多波波场数值模拟技术显得非常重要。
高精度数值模拟技术是联系地震、地质、测井以及油藏工程的纽带，其作用主要体现在提高人们对各种复杂介质中地震波传播规律的认知，并可为新技术、新方法提供试验数据，以满足方法技术研究的需要，同时也可以检验解释结果的正确性。
弹性波波动方程高精度数值模拟可以得到全波场信息，包含了地震波的动力学和运动学特点，为准确描述地震波场特征和波的传播规律奠定基础，本文在弹性波方程基础上，采用高阶交错网格有限差分技术模拟地震波在各向同性介质和各向异性介质中的传播，比较地震波在各向同性介质和各向异性介质中的波场响应异同，并分析了各向异性系数对多波波场特征的影响，这对研究各种复杂介质中弹性波的波场特征和传播规律有着重要的意义。
完全匹配层吸收边界条件在弹性波波场分离数值模拟中的应用
ｗ ｗ ｐ＋ ｗ ５
Ｌｔ ｕ；＝ Ｌ
＋ Ｌ，
．
些 方 法容 易 实现 ， 吸 收边界 只需较 少 的空 间网格 且 点， 但它们 只在一定 的入射角度 和频 率 内有效 地衰减
外 行波。当前研 究比较广泛 的， 吸收效果最 佳的吸 且
Ｌｔ ｕｐ＝ Ｌｘ ． 】 。 ｃ
Ｌｔ Ｐ＝ Ｌｚ ｗ
１算法原理
考虑二 维情 况 ， 振春等 构建 了等 效的一 阶双 李
曲型速度一应力弹性波动方程 ：
１ ． １ ｕＰ ＋ ｕｆ
其 他一些 作者对此方法进行 了扩 展 ， 出了类 似的吸 提
收边 界 条 件 ， Ｈｉｄｎ２， ｅｎ ｌｓ］ Ｌａ 等 ， 如 ｇ ｏ ［ Ｒ ｙ ｏｄ［ ， ｉ ］ ３ ｏ 这
波 Ｐ ＭＩ吸收条件能够比较方便 地应用到一 阶双 曲型
作 者 简介 陈 可洋（ ３ ）男． １ ． 硕士． 主要从事地震资料数字处理方法研究工作。
石 油 工 业计 算 机 应 用 总第６期 ２１年 １ ５ ００ 第 期
能的虚假反射 。根据完全 匹配层构建 的基本 原理 ， 通 过频率 的复数 变换 ， 得到引人衰 减吸收 因子 的波动方
Ｌ ＝ ；（ “ Ｌｗ Ｖ Ｌ ＋ ，） Ｌ。 ＝；（ Ｌ’ Ｖ， ＋ Ｌ
．
Ｌｔ
＝Ｖ （ ＋Ｌｘ ２ ， Ｌ ｗ）
式 （ ）中 Ｌ 、 Ｌ 为时间和空间一阶导数 ， １ ｆ 、： Ｌ Ｕ和 Ｗ 分 别 为 方 向和 方 向的质点振动速度 ， 和 叫 Ｕ Ｕ 及 和 Ｗ 以此类 推 。 ‘＝ Ｌｘ Ⅺ ＋Ｌｊ ｗ
Ⅱ
和｛ ＾ ＝ ２ Ｌｗ Ｌ 一ｖ
Ⅱ
（） １
Ｌｔ
．
基于旋转交错网格的双相各向异性介质二维三分量波场模拟
基于旋转交错网格的双相各向异性介质二维三分量波场模拟林朋;彭苏萍;卢勇旭;王攀【摘要】基于Biot双相介质模型,推导了双相TTI介质二维三分量一阶速度-应力弹性波方程,采用旋转交错网格(RSG)技术建立了各向异性孔隙介质波动方程的二维三分量高精度有限差分格式(FDTD),并引入不分裂卷积完全匹配层(CPML)作为吸收边界条件.为了验证算法可行性,对均匀双相TTI介质中的弹性波场进行了模拟.结果表明:使用旋转交错网格有限差分技术能够模拟出双相TTI介质中存在的快横波、慢横波、快纵波和慢纵波;双相各向异性介质中存在明显的横波分裂、波前面尖角和三分叉现象;不分裂卷积完全匹配层对边界反射的吸收效果较好.通过对比传统交错网格(SSG)和旋转交错网格有限差分技术,证明了旋转交错网格有限差分算法稳定性较强,精度较高,是一种实用的地震波场数值模拟方法.【期刊名称】《煤炭学报》【年(卷),期】2016(041)005【总页数】9页(P1203-1211)【关键词】双相各向异性介质;旋转交错网格;传统交错网格;不分裂卷积完全匹配层【作者】林朋;彭苏萍;卢勇旭;王攀【作者单位】中国矿业大学(北京)煤炭资源与安全开采国家重点实验室,北京100083;中国矿业大学(北京)煤炭资源与安全开采国家重点实验室,北京100083;中国矿业大学(北京)煤炭资源与安全开采国家重点实验室,北京100083;中国矿业大学(北京)煤炭资源与安全开采国家重点实验室,北京100083【正文语种】中文【中图分类】P631林朋，彭苏萍，卢勇旭，等.基于旋转交错网格的双相各向异性介质二维三分量波场模拟［J］.煤炭学报，2016，41（5）：1203-1211. doi：10. 13225/ j. cnki. jccs. 2015. 1561Lin Peng，Peng Suping，Lu Yongxu，et al. Study on 2D/3C wave propagation in two-phase anisotropic media using the rotated staggered-grid method［J］. Journal of China Coal Society，2016，41（5）：1203-1211. doi：10. 13225/ j. cnki. jccs. 2015. 1561地球内部分布不均的岩石裂缝和孔隙在地震尺度上往往表现为各向异性性质，而传统的地震波传播理论仅仅适用于单相介质。
交错网格有限差分正演模拟的联合吸收边界
交错网格有限差分正演模拟的联合吸收边界胡建林;宋维琪;张建坤;邢文军;徐文会【摘要】三维声波方程交错网格有限差分正演模拟中的边界问题一直是热点问题.完全匹配层吸收边界(PML)具有较强且稳定的吸收效果,但必须具有一定的边界厚度才能吸收干净,这就增大了三维正演模拟的模型空间,即增加了运算量;Higdon边界能消除任意角度入射波的边界反射,也具有较强稳定性,但该高阶吸收边界离散化后过于复杂,而低阶时吸收效果不如PML边界.因此,基于对PML吸收层中的平面波传播规律的研究,重新推导PML最外层的Higdon吸收边界条件,得到含PML吸收系数的新的Higdon吸收边界条件.联合吸收边界不仅可使用较小厚度(相对于单纯PML边界)的PML层对分量进行衰减,而且在PML边界外层,能应用新推导的Higdon吸收边界条件对反射波进行匹配吸收.在相同吸收效果下,联合吸收边界大幅度降低了PML厚度,减小了运算量,得到精确的模拟结果.【期刊名称】《石油地球物理勘探》【年(卷),期】2018(053)005【总页数】7页(P914-920)【关键词】三维声波方程;交错网格有限差分;正演;PML边界;Higdon边界;联合吸收边界【作者】胡建林;宋维琪;张建坤;邢文军;徐文会【作者单位】中国石油大学(华东)地球科学与技术学院,山东青岛266555;中国石油大学(华东)地球科学与技术学院,山东青岛266555;中国石油冀东油田公司勘探开发研究院,河北唐山063004;中国石油冀东油田公司勘探开发研究院,河北唐山063004;中国石油冀东油田公司勘探开发研究院,河北唐山063004【正文语种】中文【中图分类】P6311 引言复杂地下介质中，地震波的传播过程繁冗，难以得到解析解，因此，一般是通过正演模拟探究地下地震波的传播。
在地震波正演模拟中，利用波动方程的正演模拟比用运动学的射线追踪法可获得更丰富的动力学信息，因此地震波场的数值模拟是地震波场传播研究中的重要手段之一[1-8]。
完全匹配层吸收条件下的交错网格弹性波模拟
完全匹配层吸收条件下的交错网格弹性波模拟
张雷;赵丽娜
【期刊名称】《内蒙古石油化工》
【年(卷),期】2015(000)004
【摘要】弹性波既包含纵波,又包含横波,同声波相比其与地震波更为相近,所以利用弹性波进行数值模拟更接近实际.本文基于速度—应力弹性波动方程进行有限差分数值模拟,采用交错网格的差分格式,将速度和应力分布于两套网格上.根据需要采用不同阶数的时间和空间差分精度,在保证结果精度的前提上,提高计算效率.为了排除虚假反射的影响,采用完全匹配层吸收条件,建立垂直阶梯模型进行数值模拟,得到对应的波场快照和单炮记录.数值结果分析表明,完全匹配层也可以很好地吸收人工边界产生的虚反射现象.综合利用这两者可以更好地探究弹性波的传播规律,有利于进行精确的数值模拟.
【总页数】4页(P29-32)
【作者】张雷;赵丽娜
【作者单位】塔里木油田分公司开发事业部轮南作业区,841000;中石化西北油田分公司塔河采油二厂,新疆塔里木841604
【正文语种】中文
【中图分类】P631.4
【相关文献】
1.完全匹配层吸收边界条件在弹性波波场分离数值模拟中的应用 [J], 陈可洋
2.完全匹配层吸收边界在孔隙介质弹性波模拟中的应用 [J], 赵海波;王秀明;王东;陈浩
3.基于卷积完全匹配层的旋转交错网格高阶差分法模拟弹性波传播 [J], 冯德山;王向宇
4.吸收边界条件下弹性波散射的数值模拟 [J], 常建梅;冯文杰
5.时域有限元弹性波模拟中的位移格式完全匹配层吸收边界(英文) [J], 赵建国;史瑞其
因版权原因，仅展示原文概要，查看原文内容请购买。
高阶交错网格有限差分弹性波场模拟的精度分析
高阶交错网格有限差分弹性波场模拟的精度分析
岳晓鹏;白超英;岳崇旺
【期刊名称】《煤田地质与勘探》
【年(卷),期】2017(045)001
【摘要】交错网格波场数值模拟是目前地震正演中广泛使用的方法,为对比分析不同阶数的差分格式下产生的计算效率和精度差异,重新推导了弹性波方程的4种时间4阶、空间2N阶的差分公式及系数,并计算了他们的稳定性条件.利用这4种差分格式进行弹性波场数值模拟,对比分析了波场快照、合成地震记录及CPU时间.结果表明:时间4阶、空间6+6阶精度的交错网格有限差分方法在进行地震波场数值模拟时具有较高的计算精度和计算效率.
【总页数】6页(P125-130)
【作者】岳晓鹏;白超英;岳崇旺
【作者单位】长安大学地质工程与测绘学院地球物理系,陕西西安 710064;许昌学院数学与统计学院,河南许昌 461000;长安大学地质工程与测绘学院地球物理系,陕西西安 710064;长安大学地质工程与测绘学院地球物理系,陕西西安 710064【正文语种】中文
【中图分类】P315.69
【相关文献】
1.各向异性介质弹性波高阶交错网格有限差分模拟 [J], 霍凤斌;李振鹏;徐发;张涛;
2.高阶交错网格有限差分法弹性波叠前逆时深度偏移 [J], 陈可洋
3.二阶各向异性弹性波动方程高阶交错网格有限差分法 [J], 姜宇飞
4.弹性波交错网格高阶有限差分法波场分离数值模拟 [J], 李振春;张华;刘庆敏;韩文功
5.双相TI介质中弹性波交错网格高阶有限差分法数值模拟 [J], 尹学爱;邱光辉因版权原因，仅展示原文概要，查看原文内容请购买。
Kelvin-Voigt黏弹性介质地震波场数值模拟与衰减特征
Kelvin-Voigt黏弹性介质地震波场数值模拟与衰减特征严红勇;刘洋【摘要】利用高阶交错网格有限差分模拟Kelvin-Voigt黏弹性介质中传播的地震波,同时将完全匹配层吸收边界条件引入到其边界处理中.数值模拟结果表明,完全匹配层吸收边界效果好,高阶有限差分能模拟得到的黏弹性介质波场精度较高.对模拟的黏弹性波场进行分析,表明介质的粘滞性使地震反射波的能量变弱,高频衰减明显,并比低频衰减得快,主频向低频方向移动,有效频带变窄,即降低了地震波的分辨率；并且反射转换波比反射纵波要衰减得快；而且还随着传播距离的增加,其峰值频率也逐渐降低.通过数值模拟分析具有不同的粘滞系数介质对地震波的吸收和衰减,结果表明随着粘滞系数的增大,地下介质对地震波的吸收衰减更明显.%This paper uses finite difference algorithm of high-order staggered-grid simulate Kelvin-Voigt viscoelaslic media of seismic waves and meanwhile introduces the perfectly matched layer(PML) absorbing boundary condition into its boundary. Numerical simulation demonstrates that the effect of this algorithm of absorbing boundary is very good and the wavefield of viscoelastic media obtained from high-order finite difference is relatively accurate. An analysis of viscoelastic wavefield simulation shows that the energy of the reflected wave becomes weaker,the attenuation of the high frequency wave is much more apparent in comparison with that of the low frequency wave,the main frequency becomes closer to the low frequency,and the effective bandwidth is narrower,which all induce low resolution of seismic wave according to the simulation of viscoelastic wavefields. Besides,the attenuation of PS-wave is much more rapid incomparison with that of PP-wave and the peak frequency becomes lower gradually with the increase of the propagating distance. It is also shown that the absorption and attenuation are more apparent with the increaseof viscosity coefficient by analysis of the absorption and attenuation of seismic wave in different viscosity coefficient media.【期刊名称】《物探与化探》【年(卷),期】2012(036)005【总页数】7页(P806-812)【关键词】黏弹性;交错网格;有限差分;完全匹配层;衰减【作者】严红勇;刘洋【作者单位】中国石油大学油气资源与探测国家重点实验室,北京102249;中国石油大学CNPC物探重点实验室,北京102249;中国石油大学油气资源与探测国家重点实验室,北京102249;中国石油大学CNPC物探重点实验室,北京102249【正文语种】中文【中图分类】P631.4目前，在讨论地震波传播理论时，绝大部分情况下是把地震波看作弹性波。
一种基于全局优化的交错网格有限差分法
一种基于全局优化的交错网格有限差分法印兴耀;刘博;杨凤英【摘要】在地震波场数值模拟中，交错网格有限差分技术得到了广泛的应用，但是在弹性模量变化较大时，通常会因插值而导致模拟误差增大。
旋转交错网格可以很好地克服这个缺点，因而适合于各向异性介质正演模拟。
但是对于同样大小的网格单元，旋转交错网格需要的步长比常规交错网格要大，这会使梯度和散度算子的误差增大因而更易产生空间数值频散。
针对这些问题，本文提出了旋转交错网格与紧致有限差分相结合的方法，并基于模拟退火算法进行全局优化，压制数值频散，拓宽波数范围。
数值模拟结果表明，此方法可以有效地压制数值频散，且具有较高的模拟精度。
%Staggered grid finite difference techniques have been widely used in numerical simulation of seismic wave field,but the interpolation might enlarge the simulation errors when the elastic moduli change dramatically.Rotated staggered grid can overcome this drawback thus be more applicable for modeling the anisotropic media.However,for the grid cells with the same size,rotated staggered grid needs larger step than conventional staggered grid does,conse-quently,the gradient operator and divergence operator generate larger errors which are more likely to cause spatial numerical dispersion.Aiming at the above problems,a combination of rotated staggered grid and compact finite difference method is proposed,simultaneously,global optimization is preformed based on simulated annealing algorithm to suppress numerical dispersion and to broaden the range of wavenumber.Numerical simulation results show thatthis method can effectively suppress the numerical dispersion and exhibits high simulation accuracy at the same time.【期刊名称】《地震学报》【年(卷),期】2015(000)002【总页数】11页(P278-288)【关键词】全局优化;旋转交错网格;紧致有限差分;数值频散;正演模拟【作者】印兴耀;刘博;杨凤英【作者单位】中国山东青岛 266580 中国石油大学华东地球科学与技术学院;中国山东青岛 266580 中国石油大学华东地球科学与技术学院;中国山东青岛 266580 中国石油大学华东地球科学与技术学院【正文语种】中文【中图分类】P315.3+1引言地震波数值模拟是在地震波传播理论的基础上，通过数值计算来模拟地震波在地下介质中的传播（董良国，2003），是研究地震波传播特性与地球介质参数关系的重要手段.交错网格有限差分技术已经广泛应用于各向异性介质（Igel et al，1995；Collino，Tsogka，2001；裴正林，王尚旭，2005；殷文等，2006；何燕，2008）、孔隙介质（Dai et al，1995；Zeng et al，2001；孙林洁，2012）、黏弹性介质（白晓寅，2008；孙成禹等，2010）的波动方程模拟中，但是常规交错网格在模拟非均匀性较强的介质和各向异性介质时需对介质参数进行平均或内插（陈浩等，2006），低阶插值可能会导致精度降低，高阶插值则会导致计算量增大.Gold等（1997）提出了旋转交错网格有限差分算法，Saenger等进一步研究了旋转交错网格的稳定性与频散的关系（Saenger et al，2000；Saenger，Bohlen，2004）.在旋转交错网格中采用沿网格对角线差分的方式，可以避免弹性模量的平均与插值，因而更准确.随后一些研究人员通过旋转交错网格的完全匹配层（perfectly matched layer，简写为PML）边界条件的实现，利用旋转交错网格技术对各向异性介质（李敏，刘洋，2012）和黏弹介质（严红勇，刘洋，2012）等进行了数值模拟.但由于旋转交错网格中差分是沿对角线进行，步长是同等条件下常规交错网格的■2倍（网格是正方形的情况下），若用常规的差分格式会导致离散的梯度算子和散度算子在计算时误差增大，因而容易出现数值频散.而如果使用高阶有限差分则会因需要更多的网格点而增大计算量与内存.紧致有限差分是一种隐式差分格式，在同等网格数的条件下，具有比常规差分格式更高的精度以及更低的数值频散，故成为目前关注的热点.Lele（1992）分析了紧致差分格式的分辨率特性；王书强等（2002）对弹性波方程的紧致差分格式及其与中心差分的误差进行了研究；Du等（2009）基于交错网格紧致有限差分进行了横向各向同性介质的正演模拟；Boersma（2011）利用6阶紧致有限差分求解了Navier-Stokes方程.实现有限差分时，差分系数的确定会影响数值模拟的精度.有限差分系数的求取一般包括泰勒展开和最优化算法两种方法（Liu，2013）.基于泰勒展开得到的差分系数在波数较小时精度高，但随波数增大其精度会降低；基于最优化算法的差分系数求取方法是在给定精度与差分算子长度的情况下尽可能拓宽波数范围（Holberg，1987）.Shan（2009），Kosloff等（2010），Zhou和Zhang（2011）以及Liu（2013）等对基于最小二乘算法的差分算子优化方法进行了深入研究；Zhang 和Yao（2012）则基于模拟退火算法对中心差分算子进行了优化.本文拟结合旋转交错网格与紧致有限差分技术，基于模拟退火全局优化算法对紧致差分算子进行优化，拓宽数值模拟的波数范围，在此基础上进行频散分析，最后通过数值模拟验证该方法的可行性和有效性.1 旋转交错网格紧致有限差分方法1.1 一阶速度－应力方程在弹性波正演模拟中，除了使用二阶方程外，还常常采用一阶速度－应力弹性波方程，其主要优点是无需对弹性常数进行空间差分（Virieux，1984）.这里根据运动平衡方程和本构关系给出二维情况下体力为零时各向异性介质的一阶速度－应力方程：式中：v，σ分别代表速度和应力，1表示x方向，3表示z方向；cij为介质弹性张量矩阵的元素，当极端各向异性介质时，该矩阵有21个独立参数；当二维时，不考虑y分量，故张量矩阵中与y分量有关参数为0，只有5个独立参数（c11，c13＝c31，c15＝c51，c33，c35＝c53，c55）.1.2 紧致有限差分在地震波正演模拟中，采用常规的中心有限差分往往需要较多的网格点才能比较有效地压制数值频散，不利于边界的处理（王书强等，2002）.紧致有限差分是一种隐式差分格式，只需要较少的网格点即可有效地压制数值频散，且同等阶次的精度要比常规中心有限差分格式高，因此本文采用紧致有限差分格式，弥补旋转交错网格因差分步长大导致梯度和散度算子在计算时误差增大因而更容易产生数值频散的不足.常规的2（M－1）点2 M阶紧致有限差分格式（Kim，Lee，1996）为式中，Δx为网格的空间步长，a和bn分别为差分的权系数，可以通过在第i个网格点的泰勒展开，然后对比对应系数来求解即可得到2（M－1）点2 M阶紧致有限差分格式的差分系数（Liu，Sen，2009），例如当M＝5时，可以得到8点10阶的紧致差分格式系数为a＝0.257 894 74，b1＝0.889 871 16，b2＝0.216 121 16，b3＝－0.004 701 2，b4＝0.000 151 55.1.3 旋转交错网格及完全匹配层吸收边界条件在二维旋转交错网格上，密度和速度的各个分量定义在相同位置，弹性模量和应力的各个分量定义在另外的相同位置上，如图1所示.旋转交错网格在计算时分为两步（陈浩等，2006）：首先，计算沿着对角线方向即˜x，˜z方向进行差分，得到相关物理量的对角线方向的空间一阶导数；然后，通过坐标旋转的换算关系，将得到的两个对角线方向的差分进行线性组合从而获得坐标轴方向即x，z方向的一阶空间导数，其换算关系为（Saenger et al，2000）：图1 旋转交错网格及完全匹配层吸收边界示意图x，z为坐标轴方向，x˜，z˜为对角线方向；灰色区域表示沿着x和z方向均进行衰减，白色区域表示只沿着x方向或者z方向进行衰减Fig.1 Schematic diagram of rotated staggeredgridand perfectly matched layer（PML）absorbing boundary xand zare the coordinate directions，andx˜andz˜ are the diagonal directions.Gray areas indicate that waves are absorbed along both xand z directions and white areas indicate that waves are absorbed only along the xor zdirection式中，∂／∂x˜和∂／∂z˜是沿对角线方向的导数；∂／∂x和∂／∂z是沿坐标轴方向的导数.利用波动方程进行数值模拟，一个关键问题就是边界条件.为了解决边界反射问题，引入弹性波方程的PML边界条件.陈浩等（2006）指出，在旋转交错网格中，PML吸收边界条件的处理方式及吸收效果与常规交错网格中几乎是一样的.为此，引入图1所示的PML吸收层.PML边界的基本做法是在研究区域四周引入PML，波在PML中传播时不会产生反射，并且随传播距离按一定规律衰减.当波传播到PML边界时，波场近似为零，也不会产生反射（王守东，2003）.以式（1）中的vx为例来说明旋转交错网格紧致有限差分算法PML的实现方法：式中表示n＋1／2时刻vx 在网格点（i，j）处沿坐标轴x方向的分量，表示n－1／2时刻的x方向分量，di为衰减因子，表示n 时刻σxx在网格点（i，j）处沿坐标轴x方向的导数.根据式（5）可推导出式中沿着所在对角线˜x与斜对角线˜z方向的导数值，根据求出.其中类似地可以推导出其它分量的旋转交错网格紧致有限差分下的PML控制方程. 1.4 紧致有限差分的全局优化根据平面波理论，令式中，β＝kΔx／2，0≤β≤π／2.为了使式（2）左右两边的误差尽可能小，须满足β－β＊≈0的β取值范围应尽可能地大.本文中，采用模拟退火全局最优算法（Kirkpatrick et al，1983）来优化旋转交错网格紧致有限差分算子.为此，建立以下目标函数：Kim和Lee（1996）指出，当θ接近最大值π／2时，会出现很多难以控制的误差，从而导致优化效果不佳，因此本文选择θ＝rπ／2，其中0＜r＜1.表1列出了优化的8点10阶，10点12阶，12点14阶以及14点16阶的紧致有限差分系数.表1 优化的10—16阶紧致有限差分系数Table 1 Optimized coefficients of 10—16-order compact finite-difference阶次 r a b1 b2 b3 b4／10－2 b5／10－5 b6／10－3 b7／10－5 10 0.85 0.384 0.728 0.369 －0.015 0.152 12 0.87 0.390 0.720 0.376 －0.016 0.171 －6.49 14 0.89 0.434 0.658 0.434 －0.022 0.361 －74 0.138 16 0.85 0.452 0.633 0.458 －0.0252 0.452 －114 0.322 －7.792 频散分析由式（12）可以看出，β接近于β＊的程度表征了优化的紧致有限差分算子的精度，故定义使α≈1的β取值范围越大，数值频散越小.图2给出了常规的基于泰勒展开得到的差分算子与本文提出的全局优化算子的精度对比.可以看出，对于常规的基于泰勒展开得到的紧致有限差分算子，随着空间阶次的增大，使α≈1的波数范围变得越宽，即数值模拟的精度随阶次的增大而增大.显然，本文得到的全局优化的紧致有限差分算子在相同最大误差下使α≈1的波数范围比常规紧致有限差分算子还要大.从图2中还可以看出，优化的8点10阶紧致有限差分算子的频散比常规14点16阶紧致有限差分算子小，即优化的10阶精度要高于常规的16阶精度，从而可以使用较少的网格点来达到较高的精度，节省计算内存.图2 10阶常规有限差分，10—16阶常规紧致有限差分及对应的优化紧致有限差分的频散对比（k代表波数，Δx表示空间网格间距）Fig.2 Dispersion comparison by 10-order conventional finite-difference，10—16-order conventional compact finite-difference and the corresponding optimized compact finitedifference.k denotes the wavenumber，andΔxis spatial grid spacing由图2中的10阶常规有限差分、10阶常规紧致有限差分与10阶优化紧致有限差分的频散曲线可以看出，10阶常规紧致有限差分格式的波数范围较10阶常规有限差分要宽，而经过优化的10阶紧致有限差分格式的波数范围则更进一步拓宽.设最大允许误差为εmax＝｜α－1｜0≤β≤βmax，取εmax＝0.2%，分别计算得到10阶常规有限差分、10阶常规紧致有限差分及10阶优化紧致有限差分的最大波数范围βmax分别为0.847 5，1.037 5，1.458 0，即10阶常规紧致有限差分的波数范围比常规差分格式拓宽了1.22倍，而经过全局优化后，波数范围提升至1.72倍.考虑到正演模拟中最常采用Δx＝Δz的情况，则对于旋转交错网格，Δr比常规交错网格差分步长Δx及Δz增大了2倍，因而可以采用优化系数的紧致有限差分来弥补这一不足，从而有效地压制数值频散.需要指出的是，优化系数的紧致有限差分格式也可以很方便地应用于PML边界条件，而不需要特殊的处理.3 模型试算3.1 模型1为了验证全局优化的旋转交错网格紧致有限差分的精度比常规旋转交错网格有限差分方法有所提升，首先采用一个简单的各向同性模型进行测试（模型大小为3 000m×3 000m，在测试中并未添加任何边界条件），震源采用主频30Hz的雷克子波，z方向集中力源激发，网格大小为Δx＝Δz＝10m，时间采样为Δt＝1ms.纵横波速度分别为4 000m／s和2 600m／s，模型密度为2 300kg／m3.图3给出了3种有限差分条件下400ms时刻x分量波场快照，均采用旋转交错网格系统.图3a采用常规10阶常规有限差分格式；图3b采用8点10阶，10点12阶，12点14阶，14点16阶常规紧致有限差分格式；图3c采用优化的8点10阶紧致有限差分格式.可以看出，10阶常规有限差分与10阶紧致有限差分对比，后者频散要小；而对于图3b中的紧致有限差分格式，随着空间阶次的增大，其精度越高，频散越小；对比图3b与图3c可以看出，优化的8点10阶紧致有限差分的数值频散甚至比常规14点16阶紧致有限差分还要小，这与图2所示是一致的.图4给出了从图3中提取的z＝500m处（即50个网格点处）的波形曲线，以此来对比以上3种差分格式的数值频散.可以看出，10阶常规有限差分算法在波形外出现了数值抖动，说明数值频散比较高，模拟精度相对较低，而10阶紧致有限差分和优化的10阶紧致差分精度较高，从虚线框中可以进一步看出优化的10阶紧致有限差分算法的精度更高一些.图3 10阶常规有限差分（a），10—16阶常规紧致有限差分（b）与优化的10阶紧致有限差分（c）得到的波场快照Fig.3 Snapshots of wave field by using 10-order conventional finite-difference（a），10－16-order conventional compact finite-difference（b）and optimized 10-order compact finite-difference（c）3.2 模型2图4 10阶常规有限差分，10阶常规紧致有限差分及10阶优化紧致有限差分波形对比图Fig.4 Comparison of waveforms by using 10-order conventional finite-difference（FD）（dark line），10-order conventional compact FD （blue line）and optimized 10-order compact FD （red line）模型2的几何构造如图5所示，其中凹陷构造上部及下部均为各向同性介质，最下层为横向各向同性（VTI）介质.模型纵向和横向长为3 000m，具体参数见表2.采用优化的旋转交错网格紧致有限差分方法进行正演模拟，网格大小为Δx＝Δz＝10m，时间采样为Δt＝1ms.震源采用30Hz雷克子波，在模型的（1 150m，1 500m）处以纵波源形式激发.模拟过程中，采用PML边界条件（PML层部分同样采用优化的旋转交错网格紧致有限差分算法），厚度层为30个网格点.为体现吸收效果包含了PML层，图6给出了600ms时刻x及z分量的波场快照.可以看出，震源激发的纵波传播至分界面时，发生了反射与透射，可以见到PP波、PS波等，并可在凹陷处看到断面波，波场复杂，但波前面清晰，数值频散小.从图7所示的单炮记录中也可以看出优化的旋转交错网格紧致有限差分算法具有较高的模拟精度，各类波均可得到清晰体现，并且能够证明所加PML边界条件的吸收效果很好，无明显人为边界反射产生.图5 凹陷模型示意图Fig.5 Schematic diagram of the depression model图6 600ms时刻优化的旋转交错网格紧致有限差分正演x分量（左）和z分量（右）的波场快照（虚线框外围代表PML层）Fig.6 Wavefield snapshots of x-（left）and z-component（right）at 600ms using optimized rotated staggered-grid compact finite-difference method（PML layers are outside the dashed box）表2 凹陷模型参数Table 2 Parameters of the depression model序号介质类型 c11／GPa c13／GPa c33／GPa c55／GPa 密度／（kg·m－3）31.77 42.68 31.77 13.75 2200第二层各向同性 42.33 47.04 42.34 18.82 2400第三层横向各向同性第一层各向同性87.88 21.22 67.60 22.56 2500图7 凹陷模型优化的旋转交错网格紧致有限差分正演x分量（左）和z分量（右）单炮记录Fig.7 Records of x-（left）and z-component（right）of the depression model using optimized rotated staggered-grid compact finite-difference method4 讨论与结论旋转交错网格在地震波正演模拟中由于其可以避免在模拟中因弹性模量插值或平均而产生的误差和其对边界条件处理的便捷性而被广泛应用.但由于其差分方向是沿网格的对角线方向，因此相对于常规的交错网格而言，其差分步长较大，梯度算子、散度算子在计算时更容易产生误差，这意味着旋转交错网格更容易产生数值频散，从而影响模拟的精度.若采用小网格，虽然可以压制数值频散，却又大大增加了运算量及内存.因此本文采用旋转交错网格与紧致有限差分技术相结合的方法来提高模拟精度，并在此基础上，基于模拟退火算法对紧致有限差分进行全局优化，进一步提高计算精度，压制数值频散.数值频散分析结果表明，优化的旋转交错网格紧致有限差分算法相对于普通旋转交错网格有限差分算法具有更宽的波数范围，这意味着优化的旋转交错网格有限差分算法可以在采用更大的空间采样间隔或者更高的震源主频时依然有较高的精度和较低的数值频散；且经过全局优化的10阶紧致有限差分算子比常规的16阶紧致有限差分算子具有更高的精度，即可以采用较少的网格点来达到较高的模拟精度，因此可以节省计算内存.数值模拟实验进一步证实了本文提出的优化的旋转交错网格紧致有限差分算法的正确性与可行性.凹陷模型中，各类波均具有比较清晰的响应，揭示了波场的传播规律.同时波场快照与单炮记录中无明显的人为边界反射产生，也证明了模拟采用的PML吸收边界条件的有效性.但是本文提出的方法亦存在不足之处.在求解紧致差分格式时需要解三角矩阵线性方程，较常规显式有限差分格式增加额外的计算负担，且不利于在GPU等平台直接并行，例如图3a和图3c中的400ms时长的计算时间分别为46.5s和85.6s，后者计算耗时约为前者的1.8倍.在科学计算中，效率与计算精度常常不可兼得，随着计算机科学计算能力的提高，可以采用GPU加速的LU分解算法或者使用LAPACK等优化程序包来克服求解大型矩阵的计算效率问题.在实际应用中需权衡利弊，选择适合的算法.衷心感谢审稿专家提出的宝贵意见和建议.参考文献白晓寅.2008.基于地震波衰减理论的地层吸收参数提取方法研究［D］.东营：中国石油大学（华东）地球科学与技术学院：24-54.Bai X Y.2008.The Study on Method of Layer Absorbing Parameters Extraction Based on the Theory of Seismic Wave Attenuation［D］.Dongying：School of Geosciences，China University of Petroleum：24-54（in Chinese）.陈浩，王秀明，赵海波.2006.旋转交错网格有限差分及其完全匹配层吸收边界条件［J］.科学通报，51（17）：1985-1994.Chen H，Wang X M，Zhao H B.2006.Rotated staggered grid finite-difference and the PML boundary［J］.Chinese Science Bulletin，51（17）：1985-1994（in Chinese）.董良国.2003.地震波数值模拟与反演中几个关键问题研究［D］.上海：同济大学海洋与地球科学学院：1-2.Dong L G.2003.Several Key Problems in Seismic Wave Numerical Simulation and Inversion［D］.Shanghai：School of Ocean and Earth Science，Tongji University：1-2（in Chinese）.何燕.2008.正交各向异性弹性波高阶有限差分正演模拟研究［D］.东营：中国石油大学（华东）地球科学与技术学院：45-80.He Y.2008.High-Order Finite-Difference Forward Modeling of Elastic-Wavein Orthorhombic Anisotropic Media［D］.Dongying：School of Geosciences，China University of Petroleum：45-80（in Chinese）.李敏，刘洋.2012.高阶旋转交错网格有限差分方法模拟TTI介质中横波分裂［J］.物探与化探，36（6）：934-940.Li M，Liu Y.2012.Modeling of the S-wave splitting in TTI media using high-order rotated staggered grid scheme［J］.Geophysical and Geochemical Exploration，36（6）：934-940（in Chinese）.裴正林，王尚旭.2005.任意倾斜各向异性介质中弹性波波场交错网格高阶有限差分法模拟［J］.地震学报，27（4）：441-451.Pei Z L，Wang S X.2005.A staggered-grid high-order finite-difference modeling for elastic wave field in arbitrary tilt anisotropic media［J］.Acta Seismologica Sinica，27（4）：441-451（in Chinese）.孙成禹，肖云飞，印兴耀，彭洪超.2010.黏弹介质波动方程有限差分解的稳定性研究［J］.地震学报，32（2）：147-156.Sun C Y，Xiao Y F，Yin X Y，Peng H C.2010.Study on the stability of finite difference solution of visco-elastic wave equations［J］.Acta Seismologica Sinica，32（2）：147-156（in Chinese）.孙林洁.2012.非均匀孔隙介质有限差分正演模拟方法研究［D］.青岛：中国石油大学（华东）地球科学与技术学院：64-108.Sun L J.2012.Finite-Difference Modeling Research in Heterogeneous Porous Media［D］.Qingdao：School of Geosciences，China University of Petroleum：64-108（in Chinese）.王守东.2003.声波方程完全匹配层吸收边界［J］.石油地球物理勘探，38（1）：31-34.Wang S D.2003.Absorbing boundary condition for acoustic wave equation by perfectly matched layer［J］.Oil Geophysi-cal Prospecting，38（1）：31-34（in Chinese）.王书强，杨顶辉，杨宽德.2002.弹性波方程的紧致差分方法［J］.清华大学学报：自然科学版，42（8）：1128-1131.Wang S Q，Yang D H，Yang K pact finite difference scheme for elastic equations［J］.Journal of Tsinghua University：Science and Technology，42（8）：1128-1131（in Chinese）.严红勇，刘洋.2012.黏弹TTI介质中旋转交错网格高阶有限差分数值模拟［J］.地球物理学报，55（4）：1354-1365.Yan H Y，Liu Y.2012.Rotated staggered grid high-order finite-difference numerical modeling for wave propagation in viscoelastic TTI media ［J］.Chinese Journal of Geophysics，55（4）：1354-1365（in Chinese）. 殷文，印兴耀，吴国忱，梁锴.2006.高精度频率域弹性波方程有限差分方法及波场模拟［J］.地球物理学报，49（2）：561-568.Yin W，Yin X Y，Wu G C，Liang K.2006.The method of finite difference of high precision elastic wave equations in the frequency domain and wave field simulation［J］.Chinese Journal of Geophysics，49（2）：561-568（in Chinese）.Boersma B J.2011.A 6th order staggered compact finite difference method for the incompressible Navier-Stokes and scalar transport equations［J］.J Comput Phys，230（12）：4940-4954.Collino F，Tsogka C.2001.Application of the perfectly matched absorbing layer model to the linear elastodynamic problem in anisotropicheterogeneous media［J］.Geophysics，66（1）：294-307.Dai N，Vafidis A，Kanasewich E R.1995.Wave propagation in heterogeneous，porous media：A velocity-stress，finitedifference method ［J］.Geophysics，60（2）：327-340.Du Q Z，Li B，Hou B.2009.Numerical modeling of seismic wavefields in transversely isotropic media with a compact staggered-grid finite difference scheme［J］.Appl Geophys，6（1）：42-49.Gold N，Shapiro S A，Burr E.1997.Modelling of high contrasts in elastic media using a modified finite difference scheme［C］∥SEG Technical Program Expanded Abstracts.Dallas，Texas：SEG：1850-1853.Holberg putational aspects of the choice of operator and sampling interval for numerical differentiation in large-scale simulation of wave phenomena［J］.Geophys Prospect，35（6）：629-655.Igel H，Mora P，Riollet B.1995.Anisotropic wave propagation through finite-difference grids［J］.Geophysics，60（4）：1203-1216.Kim J W，Lee D J.1996.Optimized compact finite difference schemes with maximum resolution［J］.AIAA J，34（5）：887-893.Kirkpatrick S，Gelatt C D Jr，Vecchi M P.1983.Optimization by simulated annealing［J］.Science，220（4598）：671-680.Kosloff D，Pestana R C，Tal-Ezer H.2010.Acoustic and elastic numerical wave simulations by recursive spatial derivative operators［J］.Geophysics，75（6）：T167-T174.Lele S pact finite difference schemes with spectral-like resolution［J］.J Comput Phys，103（1）：16-42.Liu Y，Sen M K.2009.An implicit staggered-grid finite-difference method for seismic modelling［J］.Geophys J Int，179（1）：459-474.Liu Y.2013.Globally optimal finite-difference schemes based on least squares［J］.Geophysics，78（4）：T113-T132.Saenger E H，Bohlen T.2004.Finite-difference modeling of viscoelastic and anisotropic wave propagation using the rotated staggered grid［J］.Geophysics，69（2）：583-591.Saenger E H，Gold N，Shapiro S A.2000.Modeling the propagation of elastic waves using a modified finite-difference grid［J］.Wave Motion，31（1）：77-92.Shan G.2009.Optimized implicit finite-difference and Fourier finite-difference migration for VTI media［J］.Geophysics，74（6）：WCA189-WCA197.Virieux J.1984.SH-wave propagation in heterogeneous media：Velocity-stress finite-difference method［J］.Geophysics，49（11）：1933-1957. Zeng Y Q，He J Q，Liu Q H.2001.The application of the perfectly matched layer in numerical modeling of wave propagation in poroelastic media ［J］.Geophysics，66（4）：1258-1266.Zhang J H，Yao Z X.2012.Optimized finite-difference operator for broadband seismic wave modeling［J］.Geophysics，78（1）：A13-A18. Zhou H B，Zhang G Q.2011.Prefactored optimized compact finite-difference schemes for second spatial derivatives［J］.Geophysics，76（5）：WB87-WB95.。
弹性波交错网格高阶有限差分法波场分离数值模拟
© 1994-2010 China Academic Journal Electronic Publishing House. All rights reserved.
本页已使用福昕阅读器进行编辑。
福昕软件（Ｃ）２００５－２００７，版权所有，
5 12
石油地球物理勘探
Pi-
1 2
,j
k
Pi+
1 2
, j-
1
k
Pi-
1 2
, j- 1
+ [ T - T - 2 T + k
Si+
1 2
,
j +1
k
Si-
1 2
, j +1
k
Si+
1 2
,j
+ 2 T + T - T ]} k
Si-
1 2
,j
k
Si+
1 2
,
j-
1
k
Si-
1 2
, j- 1
(5)
其他方程的精度为 O (Δt4 +Δx2 N ) 的差分方程同理 可得 ,这里从略 。
3 山东省东营市中国石油大学 (华东) 地球资源与信息学院 ,257061 本文于 2007 年 1 月 30 日收到 ,修改稿于同年 5 月 9 日收到 。 本研究得到国家自然科学基金 (40474041) 、国家 863 专题 (2006AA06Z206) 、CN PC 中青年创新基金 (04 E7040) 、中原油田博士后科研工作站 和 CN PC 物探重点实验室中国石油大学 (华东) 研究室资助 。
+
5vz 5z
5σPzz 5t
纵横波波场分离效果对比
油气勘探化 工 设 计 通 讯Petroleum ExplorationChemical Engineering Design Communications·171·第47卷第4期2021年4月弹性波正演是认识地下复杂介质的重要手段，在弹性波有限差分正演模拟中，除了二阶质点位移方程外，常用一阶速度-应力方程，其优点是不需要对弹性参数进行空间微分，就可以得到完整的地震波场响应[1]。
除了常规网格外，Madariaga 提出的交错网格，它在不增加计算量的前提下，和常规网格相比局部精度提高了4倍，收敛速度也很快，很好地压制了数值频散[2]。
弹性波纵横波场分离是进行弹性波逆时偏移的必要条件，如若不能正确分离纵横波场，则在偏移成像过程中由于纵横波的互相干扰，导致偏移成像失败。
目前大多数学者都采用传统的Helmholtz 原理分离纵横波场。
Helmholtz 分解基于弹性波混合波场是由一个无旋场和无散场组成的，这个无旋场就是纵波的势函数，无散场就是横波的势函数。
因此可以通过对混合波场做散度运算提取纵波，做旋度运算。
但是由于旋度运算和散度运算对空间做了一阶偏导数，分离后的波场的振幅和相位都会改变，且分离后的纵横波波场与原波场的物理意义不一致。
假设原波场为位移，那么经过Helmholtz 分离后的波场为质点速度相关的一个量。
而且对其做相位矫正和振幅恢复时计算量过大，对后期的逆时偏移来算不易实现。
本文采用纵横波解耦方程有限差分交错网格形式下的正演模拟，吸收边界采用CPML 。
在提高计算效率上利用MPI 和OpenMP 进行并行加速。
1 基本原理与方法1.1 弹性波各向同性一阶速度-应力方程密度非均匀，二维弹性波各向同性一阶速度-应力方程：（1）式（1）中x z 垂直分量（在多波多分量地震勘探中，也被称之为伪横波和伪纵波）；τxx 、τzz 、τxz 为质点振动x 轴向，z 轴向的正应力，和y 平面法向的剪切应力；λ、μ则为弹性波的模量，其中μ是剪切模量；ρ表示为密度参数。
弹性波理论
地震波交错网格高阶差分数值模拟研究摘要: 地震波数值模拟技术是勘探地球物理学中的重要组成部分，研究通过弹性波一阶速度——应力方程，采用交错网格高阶有限差分法实现了地震波在各向同性介质中的高精度的数值模拟，并采用完全匹配层( PML) 吸收边界来消除边界反射，可取得较好的效果。
通过模型的正演计算和复杂模型的处理结果表明，交错网格高阶有限差分法数值模拟是一种快速有效的地震波数值模拟方法。
关键词: 地震勘探; 交错网格; 有限差分; 数值模拟引言地震数值模拟是模拟地震波在介质中传播的一种数值模拟技术，随着地震波理论在天然地震和地震勘探中的应用，地震模拟技术便应运而生，并随着地震波理论和计算机技术的发展，地震数值模拟技术自20世纪60年代以来也得到了飞速发展，形成了目前具有有限差分法、有限元法、虚谱法和积分方程法等各种数值模拟方法的现代地震数值模拟技术。
有限差分法是偏微分方程的主要数值解法之一。
在各种地震数值模拟方法中，最早出现的数值模拟方法是有限差分法。
Ａlterman和Karal(1968)首先将有限差分法应用于层状介质弹性波传播的数值模拟中。
此后，Boore(1972)又将有限差分法用于非均匀介质地震波传播的模拟。
Alford等（1974）研究了声波方程有限差分法模拟的精确性。
Kelly等（1976）研究了用有限差分法制作人工合成地震记录的方法。
Virieux(1986)提出了应用速度——应力一阶方程交错网格有限差分法模拟P——SV波在非均匀介质中的传播。
交错网格方法提高了地震模拟的精度和稳定性，并消除了部分假想。
有限元法也是偏微分方程的数值解法之一。
Lysmer和Drake（1972）最早将有限元法应用于地震数值模拟。
Marfurt(1984)研究对比了模拟弹性波传播的有限差分法和有限元法的精度。
Seron等（1990,1996）给出了弹性波传播有限元模拟方法。
Padovani等（1994）研究了地震波模拟的低阶和高阶有限元法。
基于散度和旋度的弹性波波场分离数值模拟方法
Ｏ ｚ
：
＋
其 中 ， 和 ： 分别 为质 点振 动 速 度 的水 平分 量 和垂 直分 量 （ 即二维 矢 量 ＝（ ， ） ，… ｔ 和 ７ 为应 ） ｒ 力 的三个 分量 ， 和 分 别 为介 质 的纵 波 速 度 和横 Ｏ ｌ
波速 度 ， 由于密度 变化 通常较 小 ， 以公 式 中略去 了 所
度） 求导 数 。 由于 该算 法 需 要 进 行 算 子 分 裂 ， 不 这 但 增加 了计算 量 ， 而且还 增加 了计算 过程 的复 杂性 。
笔 者提 出了基 于散度 和旋 度 的弹性 波波 场分离
ａ Ｔｘ
一
ｄ ｘ Ｔ Ｔｘ
～
８ Ｔｘ ｚ
＝
－ 一 ４ －
Ｏ ｔ
Ｏ ｘ
第３ ４卷第 １期
２１ ０ ０年 ２月
．
物
探
与
化
探
Ｖ０ ． Ｉ３４． Ｎｏ．１
ＧＥＯＰ ＹＳＩ Ｈ ＣＡＬ ＆ ＧＥＯＣＨＥＭ Ｉ ＣＡＬ ＥＸＰＬ０ＲＡＴ１ ０Ｎ
Ｆ ｂ ２ １ ｅ ．。 ０ ０
基 于 散度 和旋 度 的 弹性 波 波 场分 离 数值 模 拟方 法
陈 可 洋 吴 清岭 杨 微 李来 林 范兴 才 ， ， ， ，
（． １ 中国石 油大庆 油 田有 限责任 公 司 勘 探 开 发研 究 院 ， 黑龙 江 大 庆 １ ３ １ ； ． ６ ７ ２ ２ 中国石 油 大庆 油 田有 限责任公 司 采 油六厂 ， 龙江 大庆 １ ３ １ ） 黑 ６ ７ ２
摘 要： 弹性波波场分离数值模拟方法是研究混合波场 中纯纵波和纯横波波场传播规律的一种重要手段 ， 常采用 通
获得较 好 的数值模 拟 结 果 ， 是 伪谱 法 边 界 处 理 困 但 难， 计算 量 大 ， 使得 该方 法 的推广应 用受 到 限制 。李 振春等 和 Ｚ ａｇＪ ｎｅ 等 提 出基 于一 阶 速度～ ｈｎ ｉ ｌ ａ ｉ 应 力各 向 同性 弹性 波 波 场 分 离 等价 的波 动 方程 ， 并
多轴卷积完全匹配层吸收边界条件
多轴卷积完全匹配层吸收边界条件田坤;黄建平;李振春;曹晓莉;李庆洋;路萍【摘要】在通常情况下,完全匹配层(PML)吸收边界条件较其他吸收边界条件具有更优越的吸收性能,已广泛用于地震波正演模拟.但传统的PML吸收边界条件也存在一定缺陷,如不分裂卷积完全匹配层(C-PML)吸收边界在某些介质中不稳定,分裂的多轴完全匹配层(M-PML)吸收边界的吸收效果不好.为此,本文提出了多轴卷积完全匹配层(MC-PML)吸收边界条件,既可提高C-PML吸收边界的稳定性,又可改善M-PML吸收边界的吸收效率,且采用不分裂算法实现,可有效降低计算成本.正演模拟结果表明,MC-PML吸收边界条件具有更高的计算稳定性和较好的吸收效果.【期刊名称】《石油地球物理勘探》【年(卷),期】2014(049)001【总页数】10页(P143-152)【关键词】PML;C-PML;M-PML;MC-PML;吸收边界条件;数值模拟;稳定性;吸收效果【作者】田坤;黄建平;李振春;曹晓莉;李庆洋;路萍【作者单位】中国石油大学(华东)地球科学与技术学院,山东青岛266580;中国石油大学(华东)地球科学与技术学院,山东青岛266580;中国石油大学(华东)地球科学与技术学院,山东青岛266580;中国石油大学(华东)地球科学与技术学院,山东青岛266580;中国石油大学(华东)地球科学与技术学院,山东青岛266580;中石化西北油田分公司工程监督中心,新疆轮台841600【正文语种】中文【中图分类】P6311 引言随着计算机技术的发展，对复杂介质中地震波传播规律的数值模拟研究也越来越深入，其中有限差分算法应用较为广泛［1～3］。
在有限差分算法中，计算区域通常是有限的，这就需要在人工边界处对波场进行吸收衰减处理。
到目前为止，基于模型的复杂性及研究的需要，有多种吸收边界条件，其中较有效的是Bérenger ［4］在电磁波模拟中提出的完全匹配层（PML）吸收边界条件。
基于无分裂复频移卷积完全匹配层边界的黏弹介质勒夫波模拟
基于无分裂复频移卷积完全匹配层边界的黏弹介质勒夫波模拟谢俊法;孙成禹;伍敦仕;乔志浩【摘要】本文建立了无分裂复频移卷积完全匹配层(CFS-CPML)吸收边界条件,利用交错网格下的高精度有限差分格式对黏弹性介质中的勒夫波场进行了数值模拟;分析了松弛机制个数对品质因子拟合精度的影响,验证了CFS-CPML边界条件对大角度掠射波的吸收效果.数值结果表明:本文方法所使用的5个松弛机制和空间4阶差分精度,即可在保证计算效率的前提下满足目前理论研究的需要;随着品质因子的减小,频散特征曲线的相速度逐渐向增高的方向偏离理论频散特征曲线的相速度,且各模式的高频能量也随之减弱.本文结果可为发展高精度的面波反演方法提供必要的理论依据.【期刊名称】《地震学报》【年(卷),期】2016(038)002【总页数】15页(P244-258)【关键词】勒夫波;正演模拟;复频移卷积完全匹配层(CFS-CPML);吸收边界条件;频散特征【作者】谢俊法;孙成禹;伍敦仕;乔志浩【作者单位】中国山东青岛 266580 中国石油大学(华东)地球科学与技术学院;中国兰州 730020 中国石油勘探开发研究院西北分院;中国山东青岛 266580 中国石油大学(华东)地球科学与技术学院;中国山东青岛 266580 中国石油大学(华东)地球科学与技术学院;中国山东青岛 266580 中国石油大学(华东)地球科学与技术学院【正文语种】中文【中图分类】P315.3+1勒夫波是在低速层分界面附近传播的一种SH型不均匀平面波，包含了近地表的构造信息，因此可以利用勒夫波获取近地表的横波速度、品质因子等地层信息(Lee， Solomon， 1975， 1978； Xia et al， 2002， 2013； Winsborrow et al， 2003， 2005； Safani et al， 2005， 2006； Dong et al， 2013)．地震波在实际地层的传播过程中，伴随着能量的吸收与衰减，因而采用线性黏弹性的假设能够更准确地表征地震波尤其是面波在近地表衰减作用下的传播特征．一般利用品质因子Q来表征介质的黏弹性，且其在地震频带范围内可近似为常数(Kanamori，Anderson，1977)．时间域常Q黏弹性介质中的地震波模拟，通常采用广义标准线性体(Genera-lized Zener body，简写为GZB)模型(Carcione et al， 1988； Carcione， 1993)或广义麦克斯韦尔(generalized Maxwell body define by Emmerich and Korn，简写为GMB-EK)模型(Emmerich， Korn， 1987； Kristek， Moczo， 2003)实现常Q近似，Moczo和Kristek(2005)已证明这两种模型是等价的．关于常Q黏弹性介质中体波模拟的研究较多(Liu et al， 1976；孙成禹等， 2010)，而针对面波的研究则相对较少. Carcione(1992)曾采用伪谱法对黏弹性半空间情形下的瑞雷波进行正演，并采取广义标准线性体模型进行常Q拟合，但其松弛函数遗漏了1/L(L为标准线性体个数或松弛机制个数)因子，该疏漏后来在其专著中得以修正(Carcione，2001)． Zhang等(2011)同样采用伪谱法结合广义标准线性体模型进一步分析了两层介质情况下瑞雷面波的频散特征变化，但是由于其采用的模型较为简单，高阶模式的变化规律不明显．高静怀等(2012)基于位移-应力旋转交错网格有限差分法研究了黏弹性介质中瑞雷波的正演，并与完全弹性情形进行了对比，但其采用的松弛函数同样也没有1/L因子．就勒夫波而言， Boore(1970)采用常规网格有限差分法实现了完全弹性横向非均匀介质中的勒夫波数值模拟； Luo等(2010)基于标准交错网格有限差分法研究了完全弹性介质中3种典型地层的勒夫波频散特征．而关于黏弹性常Q介质中勒夫波数值模拟的研究，特别是Q值对勒夫波频散特征的影响程度，尚不多见．波场数值模拟需要高效的边界反射处理技术．Bérenger(1994)在进行电磁波模拟时提出完全匹配层技术(perfectly matched layer，简写为PML). 由于PML具有良好的吸收效果，目前已成为地震波数值模拟中最流行的边界反射处理技术．早期的PML多采用波场分裂方式，称为分裂式完全匹配层(Split-PML)(Collino，Tsogka， 2011)，相对于旁轴近似、指数衰减等吸收边界处理技术(Clayton，Engquist， 1977； Cerjan et al， 1985)，其吸收性能有了很大的提高，但仍然存在某些局限性，例如，当震源离PML层很近或接收排列过长时，不能有效地吸收以掠射角入射的波．为了克服该问题，无需波场分裂的PML技术得以发展，其中比较典型的为复频移卷积完全匹配层技术(complex frequency shifted convolutional perfectly matched layer，简写为CFS-CPML)． Komatitsch和Martin(2007)比较了完全弹性介质中CFS-CPML与Split-PML在吸收掠射角入射波时的特性， Martin和Komatitsch(2009)进一步比较了黏弹性介质中两种边界的吸收特性．国内多位研究人员对吸收边界问题也进行了大量的研究工作(张显文等， 2009；张鲁新等， 2010；田坤等， 2013； Wang et al， 2013)，结果表明，采用CFS-CPML能够有效解决掠射情况下PML吸收效果差的问题，而且实现过程中无需分裂变量，提高了计算效率．鉴于实际地层并非完全弹性，尤其是近地表介质，其黏弹性特征较为明显. 为了模拟黏弹性介质中的勒夫波，本文将首先利用5个松弛机制的广义麦克斯韦尔模型来刻画黏弹性特征，采取无需分裂的复频移卷积完全匹配层技术处理人工边界反射，结合时间2阶、空间4阶精度的标准交错网格有限差分法，以实现黏弹性常Q介质中的勒夫波正演模拟．其次，分析不同Q值对勒夫面波频散特征的影响，以及速度递增、夹高速层和夹低速层的复杂地层的勒夫波频散特征，为利用勒夫波获取近地表信息提供必要的理论依据．在完全弹性条件下，应力与应变之间是一种瞬时对应关系，不随时间变化；而在黏弹性条件下，应力与应变均为时间的函数，三维各向同性线性黏弹性介质的本构方程可表示为(Robertsson et al， 1994)：式中，“*”表示时域卷积，ψπ(t)和ψμ(t)为应力松弛函数，σ为应力，ε为应变， t为时间，δ为克罗内克函数，下标遵循爱因斯坦求和约定．通常基于广义麦克斯韦尔模型(Emmerich， Korn， 1987； Kristek， Moczo，2003)来实现常Q近似，基于该模型得到的松弛函数为式中：μr和πr为已松弛模量； L为设置的松弛机制个数；θl为对应第l个松弛机制的松弛频率； H(t)为单位阶跃函数；和为进行常Q拟合时的待定系数，可通过最小二乘法求得(Emmerich， Korn， 1987)．将式(2)代入式(1)，最终得到黏弹性介质中的一阶速度-应力-记忆变量方程为式中: μu为未松弛模量; ξxyl和ξzyl为对应第l个松弛机制的记忆变量； v为速度分量，单位为m/s; ρ为密度，单位为kg/m3．高效的吸收边界是进行波场数值模拟必不可少的条件，而采用无分裂波场的CFS-CPML吸收边界，避免了非物理的波场分裂所引起的数值误差，同时引入了复频移技术，提高了对高角度掠射波场的吸收效果．对式(3)进行傅里叶变换，得到其频域形式：式中“^”表示对应的频域物理量．经典完全匹配层(perfectly matched layer，简写为PML)技术实际上是对PML层内的参数进行复坐标拉伸，即式中， x为原坐标，为变换后的坐标， dx为衰减系数，ω为圆频率．上式等价于空间导数的坐标变换，即式中sx和sz为复拉伸函数，具体表示为经典PML技术通常采用波场分裂的方式，对质点的速度分量和应力分量沿坐标轴方向进行分裂(Qin et al， 2009) ．当震源离PML层很近或接收排列过长时，不能有效地吸收以掠射角入射的波．为解决该问题，CFS-CPML技术(Komatitsch，Martin， 2007； Martin， Komatitsch， 2009；张显文等， 2009； Zhang，Yang， 2010；张鲁新等， 2010；田坤等， 2013；王华等， 2013)引入两个频移因子αx，αz和两个收缩因子kx， kz，将式(7)推广为将式(9)代入式(4)，经傅里叶反变换，得到基于CFS-CPML技术实现的一阶速度-应力-记忆变量时域方程如下：其中为4个内部变量. 沿x方向的衰减系数dx，收缩因子kx，频移因子αx分别为(Zhang， Yang， 2010)式中：x表示网格节点到PML层内界面的距离；L为PML层的厚度；d0，k0，α0为常数，通常取d0=-3cSlogRc/2L， cS为横波速度， Rc为反射系数．当PML层的厚度为10个网格节点时，可取Rc=0.1%， k0=2，α0=πf0， f0为震源子波的主频(Collino， Tsogka， 2001)．沿z方向的处理办法与x方向类似．由于CFS-CPML技术处理边界不需要进行波场分裂，降低了存储和求解的复杂程度，其计算效率比使用Split-PML技术高．根据图1所示的网格交错方式，对式(10)进行离散，即可实现黏弹性SH波的数值模拟．限于篇幅，这里仅给出速度分量和内部变量的时间2阶、空间4阶离散形式：其中式(10)中的其它变量，如σxy，σzy，ξxyl，ξzyl，可类似导出．利用显式有限差分进行数值模拟时，须考虑计算过程的稳定性．以黏弹介质中的未松弛(高频极限)速度代替弹性速度(田坤， 2014)，得到可以将此作为黏弹性介质稳定性的判断依据．式中Δx和Δz为横向和纵向的空间网格采样步长，Δt为时间采样间隔， cu为高频极限速度， an为相应的交错网格差分算子系数．利用式(2)的松弛函数进行常Q近似时，松弛机制个数对拟合精度具有较大的影响，松弛机制个数越多，拟合精度越高，但相应的计算量也越大；因此，在保证拟合精度足够高的前提下，需要尽量减少松弛机制的个数，从而减小计算量．图2给出了分别采用3个和5个松弛机制在1—10 Hz范围内的常Q拟合结果. 可以看出，采用3个松弛机制的拟合曲线与Q=20的曲线相差较大，而采用5个松弛机制的拟合曲线则与Q=20的曲线吻合得较好．图3给出了1—10 Hz范围内，采用5个松弛机制进行常Q拟合的数值相速度与使用Futterman(1962)公式计算所得理论相速度的对比结果，其中相速度的变化范围为2560—2680 m/s，其最大相对误差为0.15%，这充分说明采用5个松弛机制得到的拟合结果具有较高的精度，故本文的正演模拟中选取5个松弛机制．为比较地震波以大角度掠射到PML层时Split-PML与CFS-CPML的吸收特性，建立一个均匀各向同性线性黏弹性介质模型，如图4所示．该模型长2000 m，宽400 m；网格单元为Δx=Δz=2 m，模型四周为PML层．介质的密度为2300 kg/m3; 横波速度通过来确定，此处取2700 m/s，横波品质因子为QS=20； PML层的厚度为80 m；接收线的深度为84 m；炮点位于(800 m，82 m)，距离上面PML层仅2 m，会出现大角度掠射的情况．震源子波取主频f0为40 Hz的雷克子波，选取反射系数Rc=0.0001%，收缩因子k0=1，频移因子α0=πf0. 采用时间2阶、空间4阶的差分阶数，分别利用Split-PML和CFS-CPML边界条件进行两次模拟，模拟时长为0.54 s，计算步长为0.2ms．图5和图6分别给出了利用Split-PML和CFS-CPML处理吸收边界的情况下， 0.14 s， 0.3 s和0.54 s时刻的垂直分量波场快照，其大小与图4相同，且包含了吸收边界的波场信息．可以看出：图5中的上边界存在反射能量，说明采用Split-PML吸收边界对掠射波的吸收效果有限；图6中几乎没有边界反射，说明采用CFS-CPML吸收边界对掠射波具有较好的吸收效果．此外，对图4所示的模型在计算机(处理器为i5-3470，主频为3.2 GHz)进行正演时，采用Split-PML和CFS-CPML吸收边界获取0.54 s时长的地震记录所需时间分别为263.68 s和167.41 s，说明使用CFS-CPML吸收边界的计算效率更高．图7给出了采用Split-PML和CFS-CPML两种吸收边界获得的炮记录对比，两侧各40道为吸收边界，考虑到边界反射较弱，制图时对图7a和图7b均加了4%的增益．可以看到：图7a中炮检距较大的地震道存在顶部边界的反射(箭头所示区域)，而炮检距较小的地震道则没有，这是由于大炮检距处波入射到匹配层的入射角过大，匹配层吸收不完全所致；图7b中没有顶部边界的反射，说明CFS-CPML吸收边界对大角度掠射到边界的地震波具有较好的吸收效果．将图4所示的介质模型扩大成6000 m×6000 m，即网格点数为3000×3000，激发点设在(3000 m， 3000 m)点处，接收线的深度设为3002 m．在这种条件下，当计算时间tcal≤0.54 s时，直达波未到达边界，属于无边界反射记录，本文将该记录作为评价边界吸收效果的标准．选取图7a和图7b中炮检距为零的地震数据进行波形显示，如图8a所示，可以看出，采用Split-PML和CFS-CPML 吸收边界的波形与无边界反射记录非常一致．图8a下部放大图中所显示的0.05—0.1 s内振幅为10-3数量级，存在边界反射的波形与无边界反射的波形几乎一致，说明地震波以小角度入射到边界时，采用Split-PML和CFS-CPML吸收边界均有较好的吸收效果．选取图7a和图7b中炮检距为400 m的地震数据进行波形显示，如图8b所示. 从图8b下图的直达波尾部放大图可以看出，采用CFS-CPML吸收边界的波形更接近无边界反射的波形，说明对于大角度掠射的地震波，采用CFS-CPML吸收边界的处理效果较好．由于勒夫波不存在于均匀半空间中，其正演模拟需采用多层介质模型．建立最简单的两层介质模型，参数见表1. 由于密度对勒夫波的影响较小，所以模型中两层介质的密度均设为2000 kg/m3．整个模型由801×200个网格点构成，网格大小为Δx=Δz=0.5 m，时间步长为Δt=0.5 ms， CFS-CPML层由10个网格点构成；震源函数采用主频为f0=40 Hz的雷克子波，反射系数Rc=0.1%，收缩因子k0=2，频移因子α0=πf0，采用5个松弛机制在1—100 Hz范围内进行常Q 拟合．鉴于Q>150时，黏弹性衰减效应很小，因此当地层Q≥200时，该地层可视为完全弹性介质．保持底层半空间介质为完全弹性，依次使第一层的横波品质因子QS为200， 40， 20和10，采用时间2阶、空间4阶的有限差分法进行正演得到4个勒夫波记录，通过倾斜叠加算法(Park et al， 1999)计算其所对应的频散特征(图9a--d)，并将其与理论频散特征进行对比以考察表层Q值变化对勒夫波频散特征的影响．当第一层介质的横波品质因子QS值为200时，两层介质均为完全弹性，其勒夫波地震记录及相应的频散特征如图9a所示．从勒夫波记录中可以看出CFS-CPML边界的吸收效果较好．从图9a给出的频散特征(彩色区域)可以看到：该记录中包含基阶、第一高阶、第二高阶和第三高阶模式的勒夫波；对每个模式而言，随着频率的增大，相速度逐渐减小并趋于表层横波速度300 m/s；随着频率的减小，相速度趋于底层半空间横波速度600 m/s．图9a中的白色曲线为根据矩阵传递算法(Haskell， 1953)正演得到的理论勒夫波频散曲线，很容易看到，各模式的频散曲线与其吻合得较好，说明正演算法能够准确地模拟该地层条件下的勒夫波．此外，理论上已证明，在两层完全弹性介质的情况下，勒夫波高阶模式存在截止频率，即当频率低于截止频率时不存在高阶模式(Aki， Richards，2002)，这在图9a中也得以体现，即当频率低于20 Hz时只存在基阶模式，不存在高阶模式．图9b给出了当第一层介质QS值为40时的勒夫波地震记录及其相应频散特征．可以看出， 4个模式的高频成分损耗十分明显，从约80 Hz降低到近似70 Hz，这是由于吸收衰减作用所致．此外，由于黏弹性所导致的固有频散，各模式频散曲线的最低相速度均高于完全弹性情况下的理论相速度．当第一层的QS值进一步减小到20(图9c)甚至10(图9d)时，从勒夫波记录中可以看出高频成分较完全弹性情况下的地震记录损失较严重，特别是第601道之后的远偏移距部分更为严重．从图9c频散图中也可以看到，此时第三高阶模式高频能量已近乎完全衰减，第二高阶模式的能量也仅集中在40—60 Hz之间的狭窄频带内，基阶和第一高阶模式的最高频率降低至60 Hz附近．若表层介质QS值减小到10(图9d)，第三高阶模式高频能量则完全衰减，第二高阶模式的最高频率降低至50Hz左右，基阶和第一高阶模式的最高频率也降低，并且相速度在高频时趋于340 m/s，与完全弹性时的相速度300 m/s的相对偏差可达13.3%．综上可见，勒夫波频散曲线中各模式的最低相速度随第一层介质QS值的减小而增大．然而，实际地层并不是完全弹性的，近地表介质尤其如此，因此，在利用面波获取近地表速度时，第一层介质的速度受频散曲线最低相速度的影响较大，导致利用面波所获取的第一层速度与真实速度存在一定的误差，且QS值越小，该误差越大．近地表介质复杂多变，垂向上的速度并非总是递增，且近地表压实程度较低，速度和Q值通常较小．为便于对比，研究多层黏弹介质的勒夫波特征时，保持底层半空间介质为完全弹性，对3种黏弹介质进行组合，构成4层介质模型，具体参数如表2所示．横波速度为400， 600， 800 m/s和1000 m/s时，横波品质因子QS分别为10， 20， 50， 200．除模型参数设置不同外， 4层介质模型的数值模拟参数与表1所示的双层介质模型一致．从图9a可知弹性介质的勒夫波频散曲线与理论曲线吻合，因此仅分析黏弹介质的勒夫波频散曲线特征，弹性介质的频散曲线特征可根据理论曲线作出判断．下面分速度递增模型、夹高速层模型和夹低速层模型等3种情况逐一讨论．图10a--c分别给出了速度递增模型、夹高速层模型和夹低速层模型的勒夫波记录及其相应的频散曲线，其中，彩色区域是从勒夫波记录中采用倾斜叠加算法所提取的频散曲线，白色线条是理论频散曲线．可以看出：速度递增模型和夹高速层模型的频散曲线都只存在基阶模式和第一高阶模式，基阶模式的能量最强，且随着频率的增加，每个模式的频散曲线与理论频散曲线的吻合程度逐渐变差；与速度递增模型的频散曲线(图10a)相比，夹高速层模型频散曲线(图10b)的基阶模式弯曲程度更小，第一高阶模式的能量更弱；夹低速层模型频散曲线(图10c)尚存在第二高阶模式，但其能量较弱，总体频率范围和各模式的能量较速度递增模型增加，且基阶模式的频率范围更小，与理论频散曲线的吻合程度更高，此外，第一高阶模式的频率范围更大．由于黏弹性介质对高频能量的吸收较其对低频能量的吸收强，因此频散曲线中各模式高频部分与理论曲线的吻合程度相对较差．勒夫波是由一次波和多次反射的SH波相长干涉产生的，基于几种介质不同组合形成的模型得到的勒夫波频散曲线，其各阶模式存在差异，即速度递增模型、夹高速层模型和夹低速层模型的基阶模式和第一高阶在能量分布、频率范围和弯曲程度上都存在差异．反之，勒夫波包含了地下构造信息，其频散曲线中各模式的能量分布、频率范围和弯曲程度是地下构造信息的一种反映，采用适当的方法能够从勒夫波中获取地下的构造信息，目前常用的做法是拾取频散特征曲线，再结合最小二乘或者遗传算法等优化方法来反演得到地下横波速度等底层信息(Xia et al， 2002， 2013； Safani et al，2005， 2006)．本文以无分裂复频移卷积完全匹配层为吸收边界, 利用高精度交错网格有限差分，推导了黏弹性介质波动方程的一阶速度-应力-记忆变量方程，并给出了数值解法，实现了黏弹性勒夫波的数值模拟. 通过模型试算，得到结论如下：1) 相对于传统的分列式完全匹配层边界条件， CFS-CPML边界条件可以更加有效地吸收掠射到边界的地震波，且无需进行波场分裂，从而提高了计算效率．2) 本文试验表明，使用5个松弛机制和4阶差分精度，能够在保证计算效率的前提下满足目前理论研究的需要．3) 介质的黏弹性使频散特征曲线的相速度高于理论频散特征曲线的相速度，各模式的高频成分能量减弱，且Q值越小，频散特征曲线与理论频散特征曲线的偏差越大，各模式高频成分的能量越弱．实际的近地表地层具有较强的黏弹性，本文发展了基于无分裂复频移卷积完全匹配层边界的黏弹性勒夫波数值模拟方法，模拟的频散特征更加接近真实地层情况，可为面波反演提供理论依据．。
基于时间高精度-空间隐式矩形交错网格有限差分的弹性波数值模拟
基于时间高精度-空间隐式矩形交错网格有限差分的弹性波数
值模拟
徐世刚;包乾宗;任志明;刘洋
【期刊名称】《地球物理学报》
【年(卷),期】2022(65)4
【摘要】针对弹性波数值模拟,现有的时间高阶有限差分方法主要基于显式差分和正方形网格单元建立的,因此容易引起较强数值频散且缺乏灵活性.基于上述问题,本文通过结合十字形和菱形模板设计了一种改进的矩形组合差分模板,在此基础上发展了隐式交错网格差分方法求解弹性波方程.该方案采用改进差分模板和二阶离散格式共同近似时间导数,采用隐式差分格式求解空间导数.进一步给出了泰勒级数展开和最小二乘优化两种算法求取的高阶差分系数.联合应用高阶差分方案和波场分离技术产生高精度的弹性波场.将本文方法与几种现有差分方法进行了对比,包括:显式和隐式时间二阶差分法,显式时间高阶差分法.数值分析和模型算例表明,本文方法比其他方法具有更高的精度和灵活性.
【总页数】13页(P1389-1401)
【作者】徐世刚;包乾宗;任志明;刘洋
【作者单位】长安大学地质工程与测绘学院地球物理系;中国石油大学(北京)油气资源与探测国家重点实验室
【正文语种】中文
【中图分类】P631
【相关文献】
1.弹性波交错网格高阶有限差分法波场分离数值模拟
2.基于横向各向同性BISQ模型的弹性波高阶交错网格有限差分数值模拟
3.弹性介质瑞雷面波交错网格有限差分法数值模拟
4.一种全局优化的隐式交错网格有限差分算法及其在弹性波数值模拟中的应用
5.基于BISQ模型的各向同性孔隙介质弹性波三维交错网格高阶有限差分数值模拟
因版权原因，仅展示原文概要，查看原文内容请购买。
完全匹配层吸收条件下的交错网格弹性波模拟
均在 这方 面做 出 了大量 研究 。 目前 常 采 用 的 数 值 模 拟 方 法 主 要 有 有 限差 分
法 ］ 、 边界 元 法［ ６ ］ 和虚 谱 法［ ７ ］ 。有 限差 分法 是 最早
应力 表示 的一 阶 弹性波 方 程为
鲁＝ 警＋ 警 、
蛆：
ａ 出 击
、
＾ ＋
进 一步 的探究 。 。 采 用有 限差分 方法 进行 数值 模拟 的 过程 中 ， 由于 受到计 算 机 内存和 计算 时 间韵 限制 ， 有
限差分 法 只能 得 到有 限数 量 网格 点 上 的波 场 值 ， 所 以就必须 截 断计算 空 间并设 置边 界 条件 广泛 应用 的 边界 条件 有很 多 。 目前 吸 收效果 最 好 的是完 全 匹配 层（ Ｐ ＭＬ） 吸收条件 ， 许 多 文 献 都 验 证 了其 吸 收 效 果， 但是 很少 有 文 献 提及 到 将 完 全 匹 配层 吸收 条件
警＝ 警 ｚ 警
誓 警＋ 警
交错 网格 法将 速 度 和 应 力分 别分 布 于 两套 网格 上 ， 和常 规 网格 相 比提 高 了局 部 精度 和 收 敛 速度 ， 由于 其优点 ， 许 多 学 者在 交 错 网 格方 面进 行 了大 量 的研 究， 王 春燕 ［ ８ ］ 主要 研究 高 阶交错 网格 有 限差 分方 法 ， Ｔｈ ｏ ｍａ ｓ Ｂ ｏ ｈ ｌ ｅ ｎ ｍ在瑞 雷 波 交错 网格 数 值 模 拟方 面
摘
要： 弹性 波既 包含 纵 波 ， 又 包含 横 波 ， 同声 波相 比其 与地 震 波更 为相 近 ， 所 以利 用弹性 波进 行数
值 模拟 更接近 实际。 本 文基 于速 度 一应 力 弹性 波动 方程进行 有 限差 分数值 模 拟 ， 采用 交错 网格 的差分格
弹性波理论
地震波交错网格高阶差分数值模拟研究摘要: 地震波数值模拟技术是勘探地球物理学中的重要组成部分，研究通过弹性波一阶速度——应力方程，采用交错网格高阶有限差分法实现了地震波在各向同性介质中的高精度的数值模拟，并采用完全匹配层( PML) 吸收边界来消除边界反射，可取得较好的效果。
通过模型的正演计算和复杂模型的处理结果表明，交错网格高阶有限差分法数值模拟是一种快速有效的地震波数值模拟方法。
关键词: 地震勘探; 交错网格; 有限差分; 数值模拟引言地震数值模拟是模拟地震波在介质中传播的一种数值模拟技术，随着地震波理论在天然地震和地震勘探中的应用，地震模拟技术便应运而生，并随着地震波理论和计算机技术的发展，地震数值模拟技术自20世纪60年代以来也得到了飞速发展，形成了目前具有有限差分法、有限元法、虚谱法和积分方程法等各种数值模拟方法的现代地震数值模拟技术。
有限差分法是偏微分方程的主要数值解法之一。
在各种地震数值模拟方法中，最早出现的数值模拟方法是有限差分法。
Ａlterman和Karal(1968)首先将有限差分法应用于层状介质弹性波传播的数值模拟中。
此后，Boore(1972)又将有限差分法用于非均匀介质地震波传播的模拟。
Alford等（1974）研究了声波方程有限差分法模拟的精确性。
Kelly等（1976）研究了用有限差分法制作人工合成地震记录的方法。
Virieux(1986)提出了应用速度——应力一阶方程交错网格有限差分法模拟P——SV波在非均匀介质中的传播。
交错网格方法提高了地震模拟的精度和稳定性，并消除了部分假想。
有限元法也是偏微分方程的数值解法之一。
Lysmer和Drake（1972）最早将有限元法应用于地震数值模拟。
Marfurt(1984)研究对比了模拟弹性波传播的有限差分法和有限元法的精度。
Seron等（1990,1996）给出了弹性波传播有限元模拟方法。
Padovani等（1994）研究了地震波模拟的低阶和高阶有限元法。
基于高阶时域有限差分算法的电磁波传播计算
基于高阶时域有限差分算法的电磁波传播计算苏卓;谭峻东;张俊;龙云亮【期刊名称】《电波科学学报》【年(卷),期】2014(0)3【摘要】利用FDTD(2,4)高阶时域有限差分(Finite-Difference Time-Domain,FDTD)算法并结合滑动窗口的思想,对电磁波传播特性进行了仿真计算.采用的高阶FDTD算法在空间上达到四阶精度,与二阶精度的传统FDTD算法相比,在相同每波长采样数的条件下,数值色散误差能得到进一步的减少.在源脉冲传播较长距离时,数值色散的减少使得时域下脉冲扩展现象得到改善,滑动子窗口仍然能包含着激励源脉冲的全部信息,从而可更加准确地计算长距离电波传播特性.另外,在相同的数值色散误差容限下,每波长采样数比传统二阶FDTD方法有所减少,从而节省存储空间,加快计算速度.【总页数】6页(P431-436)【作者】苏卓;谭峻东;张俊;龙云亮【作者单位】中山大学电子与通信工程系,广东广州510006;中山大学电子与通信工程系,广东广州510006;中山大学电子与通信工程系,广东广州510006;中山大学电子与通信工程系,广东广州510006【正文语种】中文【中图分类】TN011【相关文献】1.二维周期钢管阵列带隙特性计算中的时域有限差分算法 [J], 王刚;赵宏刚;温激鸿;韩小云2.时域有限差分算法的计算机实现可行性研究 [J], 潘东旭;朱江3.计算电离层VLF波传输特性的时域有限差分算法 [J], 张驰;芮国胜;王瑞;薛鹏4.电磁计算中辛时域有限差分算法研究进展 [J], 查晓民; 谢国大; 沙威; 黄志祥; 吴先良5.一种高效计算各向异性磁化等离子体的时域有限差分算法 [J], 钟双英;刘崧因版权原因，仅展示原文概要，查看原文内容请购买。