２００７年ｌＯ月 石油墟球幽 参探 第４２卷第５期
处理方法·
弹性波交错网格高阶有限差分法
波场分离数值模拟
李振春 ① 张 华① 刘庆敏① 韩文功②
（①中国石油大学（华东）地球资源与信息学院；②中石化胜利油田分公司）
李振春，张华，刘庆敏，韩文功．弹性波交错网格高阶有限差分法波场分离数值模拟．石油地球物理勘探，２００７， ４２（５）：５１０４５１５
摘要 地震波中的弹性波传播会产生纵波和横波，所以用完全弹性波波动方程进行弹性波波场数值模拟时，只 能得到纵横波耦合的混合波场。本文从Ｐ波波场为无旋场，而ｓ波波场为无散场的思路出发，推出满足此条件 的一阶速度一应力弹性波波场分离方程，并利用交错网格高阶有限差分法对波场分离方程进行数值模拟。模 拟实例表明，此方法不但成功地将Ｐ波波场和ｓ波波场从混合波场中分离出来，而且这种方法的稳定性好、模 拟精度高，可用于弹性波传播规律研究及地震资料处理。
关键词 弹性波 交错网格 高阶有限差分 波场分离
１ 引言
多分量地震记录的每个分量包含着不同的波
型。尤其在垂直分量和径向分量的记录中既有纵波
反射波，也有横波反射波。当在各向同性介质中用
Ｐ波震源激发时，在垂直分量Ｚ和径向水平分量Ｘ
的记录上都存在ＰＰ反射波和ＰＳ转换波。因此分
解出ＰＰ波和ＰＳ波是进行偏移成像等地震资料数
据处理的前提。
弹性波波场分离最初由Ｄａｎｋｂａａｒｌ１ 提出；Ｄｅｖ—
ａｎｅｙ等 利用极化波分解对弹性波ＶＳＰ数据中的
ＰＳ波场进行了分离；Ｄｅｌｌｉｎｇｅｒ等人 在空间一频
率域中利用散度和旋度因子对各向异性介质中的
ＰＳ波场进行了分离。本文根据马德堂等人ｌ４ 提出
的满足Ｐ波为无旋场、ｓ波为无散场的等价方程思
路，推导出了既满足完全弹性波动方程又满足Ｐ
波、ｓ波方程的一阶速度一应力等价方程，并利用交
错网格高阶有限差分法对其波场进行模拟，从而分
别获得混合波场、Ｐ波波场和ｓ波波场，实现了弹性
波波场分离数值模拟。 ２ 一阶速度一应力弹性波波场分离
方程
在二维ＶＴＩ介质中，假定体力为零，用速度一
应力表示的弹性波方程为
Ａ ＋Ｂ —Ａ＿十 ＿ （１） 口 ０ ０Ｚ 其中：ｗ一（ ， ，ｄ ，ａｚ．ｚ， ） 是由速度分量 （ｉ
— ，ｚ）和应力分量ｄ （ ，Ｊ— ，ｚ）构成的列向量，Ａ，
Ｂ为系数矩阵，分别为
Ａ＝
Ｂ一 ０ ０ Ｐ－ ０ ｐ－ ０ ０ ０ ０ ０ ０ ０ ０ ０ ０
０ ｐ－ ０ ｐ－ ０ ０ ０ ０ ０ ０ ０ ０ Ｏ Ｏ Ｏ
矩阵中Ｐ为介质密度；Ｃ 为介质的弹性常数。在各
＊山东省东营市中国石油大学（华东）地球资源与信息学院，２５７０６１ 本文于２００７年１月３Ｏ日收到，修改稿于同年５月９日收到。 本研究得到国家自然科学基金（４０４７４０４１）、国家８６３专题（２００６ＡＡ０６Ｚ２０６）、ＣＮＰＣ中青年创新基金（０４Ｅ７０４０）、中原油田博士后科研工作站 和ＣＮＰＣ物探重点实验室中国石油大学（华东）研究室资助。 Ｏ Ｏ ＯＯ
Ｏ 第４２卷第５期 李振春等：弹性波交错网格高阶有限差分法波场分离数值模拟
向同性介质中，Ｃ 一Ｃ。。一 ＋２／１一ｐＶ；，Ｃ 。一 —
ｌＤ（ 一２Ｖ｝），Ｃ ： 一ＩＤＶ｝； 和 为介质的拉梅常
数； 和Ｖｓ分别为纵、横波速度。为构造等价方程
在式（１）中引入混合波场新变量Ｕ一｛ ， ｝，Ｐ波
波场新变量Ｕ 一｛ ， ｝和Ｓ波波场新变量Ｕ。一 ｛伙，伙 ｝，经变换后得如下方程组
＝ ＂Ｕｐ ￣－ＶＳ
３Ｖｐ ｄｔＴｐ Ｄ— 一— ‘ ３ｔ ３ｘ ３ＶＰ ａ Ｄ— ＝— ‘ ３ｔ ３ｚ
（ ＋篝）
３ｔ一 （凳＋箦）
ａ Ｓ ａ ｌ ３ａｓ ＩＤ 一言＋
３ｖ ａ ３ａ ｌＤ 一 ＋
２ 篝
２ 凳
１／２／３一ｖ￣３ ｖ，）
式中： （ｉ—ｚ，ｚ）、 （ｉ—ｚ，ｚ）分别为Ｐ波、Ｓ波波
场的速度分量； （ｉ，Ｊ—ｚ，ｚ）、 ｓ（ｉ， —ｚ，ｚ）分别
为Ｐ波、Ｓ波波场的应力分量； （ｉ＝ｚ，ｚ）为混合波 波场的速度分量。
在均匀各向同性介质中，全弹性波波场（或矢量
弹性波波场）可以分解成纯纵波和纯横波两部分＿６］。
通过对波场分别取散度和旋度可以提取纯纵波和纯
横波。纵波是无旋场，即 ×Ｕ 一０；横波为无散
场，即 ·Ｕ。一０。下面就对式（２）中横波的散度和
纵波的旋度进行分析。在式（２）中，由于
３ （ ·Ｕｓ） ３ ｒａ仉 ．３Ｖｓ：］ ａ ａ Ｌ ３ｚ‘３ｚ
３ｘ３ｔ（ ３ｔ）＋ ３ｚ３ｔ（ ３ｔ）＼ ，’ ＼ ／
３ｘ３ｔ［等＋鲁］＋ ３ｚ３ｔ［鲁＋鲁］ Ｌ ａｚ ’ ａｚ Ｊ’ Ｌ ａｚ ’ ａｚ ２咐＿＿３＿３ ｖ￣ａ ＋ ＋ 一
由此可知， ·Ｕ。关于时间的二阶偏导数为零，即 ·Ｕｓ关于时间的一阶偏导数为常数，故 ·Ｕ。要
么为常数，要么为线性函数；由波动特性知： ·Ｕ。 关于时间只能为常数，同时由解波动方程的零初始
条件知： ·Ｕｓ应等于零，因此式（２）中的横波波场
Ｕ。是无散场；同理可推得 一ｏ，即纵波波
场Ｕ 是无旋场。所以理论上，通过求解式（２）可得
到弹性波混合波场Ｕ一｛ ， ｝，Ｐ波波场Ｕ 一 ｛ ， ｝和Ｓ波波场Ｕ。一｛ 。， 。｝，且有Ｕ—ＵＰ＋ Ｕｓ，即实现弹性波波场分离。
３一阶速度一应力交错网格任意偶数
阶精度差分格式
对方程（２）进行差分时，用有限差分计算ｔ时刻 的应力，并将该时间差分进行泰勒展开，得到在时间
域的２Ｍ阶精度差分递推公式，即
（件 ）一（卜 ）＋２ ×
×（ ） －Ｉ－Ｏ（ＡｔＴＭ）（３）
在交错网格技术中，变量的导数是在相应变量网格
点之间的半程上计算的。为此，用式（４）计算式（２） 中的一阶空间导数，即
ｃ ｛ ［ｚ＋等ｃ２ 一 ］一
Ｉ ｚ一等（２ 一１）ｌ｝＋Ｏ（Ａｘ ） （４）
式中：系数ｃ 由在ｚ处将 Ｉ ｚ＋ （２ｎ－－１）ｌ和
Ｉ ｚ一等（２ｎ－－１）ｌ进行Ｔａｙｌｏｒ展开而得到 。
在式（３）的运算中，可利用式（４）把应力的空间
导数代替速度的时间导数，用速度的空间导数代替 应力的时间导数。这样，对时间的高阶差分转化为
对空间的高阶差分，不再涉及到除了ｔ—Ａｔ／２和ｔ 以外的时间层，可得到时间为２Ｍ阶，空间为２Ｎ阶
精度的有限差分格式。
令ｕ 、 卅丢分别为混合波场Ｕ中速度分
量 、 的离散值；Ｕ 、 ．Ｊ＋毒、 粤
分别为Ｐ波波场Ｕ 中速度分量 、
，应力分量 石油地球物理勘探
（ＹＰ的离散值；己，皇ｉ ｊ、Ｖ 麓．，十｛、Ｒ§ １、 ＋ １、
Ｈｉ ＋ １分别为Ｓ波波场Ｕｓ中速度分量ＶＳ应
力分量 、 、０＂ｓ 的离散值，假设Ａｘ—Ａｚ，对式
（２）中关于纵波速度分量的方程，其精度取ｏ（Ａｔ ＋ Ａｘ ）时的差分格式为
ｕ簪一ｕ譬＋ ｛ ｃ Ｘ
Ｘ … Ｒｋ ］｝＋ （ ） ×
Ｘ｛［Ｒ 。 一３Ｒ ． ＋３Ｒ 一Ｒ 。］＋
＋［Ｒ ｛３Ｒ ＋３Ｒ Ｒ 寻］＋
＋２［Ｈｋｓ １一Ｈ １．Ｊ｛一２Ｈｋｓ ＋ １＋
＋２Ｈｋｓ １＋Ｈ ｛一Ｈ ｛］＋
＋［ ｛．，＋ 一 ．，＋ 一２Ｔ； ＋｛＋
＋２ ｛＋ 丢．ｒ．１］＋
＋［ ＋｛ 一 ｛ 一２ 告＋
＋２ ＋ ｛｛．ｒ ］） （５）
其他方程的精度为０（Ａｔ ＋Ａｘ ）的差分方程同理
可得，这卑从略。
４边界条件及稳定性条件说明
波动方程法的人为边界条件构造方法主要有四
种：衰减区法、傍轴近似法、吸收边界法、多次入射法 等。为了构造较好的人为反射边界条件，这里采用
吸收边界和衰减边界相结合的组合边界条件 ］。利 用特征分析法吸收边界条件 ］，根据不同边界区域
有选择性的压制边界反射波；再使用最佳匹配层法 （ＰＭＬ）［１０，１１］来构造衰减边界条件，在匹配层内得到
： 主ｌｏｌ 墨 ／网格点数 ｌＯｌ
《 萎ｌｏｌ 个衰减的解，使得传播到边界的入射波能量衰减
至很弱，其反射波能量更弱，几乎可以忽略。但是为 了有较好效果，最佳匹配层法必须要求较宽宽度，从
而大大增加了计算量。为了结合两种边界各自优
点，这里将吸收边界与最佳匹配层的衰减边界有机 地结合构造了组合边界条件：先采用宽度较小的最
佳匹配层衰减边界，使得传播到边界的入射波能量 衰减；然后再采用吸收边界条件，对反射波进行吸
收。这样做的优点在于计算量增加不大，吸收也较
彻底，较好地消除了人为反射边界。 稳定性条件是有限差分数值模拟中十分重要的
问题，这里给出一个简单的、比通用稳定性条件更为 严格的交错网格高阶有限差分法的稳定性条件。假
设 一△ —ｚ，最大速度值为Ｖ册 ，则时间为２Ｍ阶、 空间为２Ｎ阶精度的有限差分法稳定性条件为【１ ］
二 ：二 ｆ Ｖ △￡ Ｉ（２ｍ一１）！＼ｌｐ ／一
Ｎ Ｘ『∑（一１） Ｃ ｝≤１ ｌ 其中：ｌ０为密度；Ｃ 为式（５）中作为２Ｍ阶空间展开 时的系数。
５数值模型试算
５．１均匀介质波场分离模拟
图１为均匀介质模型模拟结果。水平、垂直方 向采样间隔Ａｘ、Ａｚ分别为５ｍ，时间采样间隔△￡为
０．５ｍｓ；纵波速度为３５００ｍ／ｓ；泊松比为０．２５；在此
模型中心用集中力源在ｚ方向激发波场。精度为
０（Ａｔ ＋Ａｘ ）的交错网格有限差分获得的地震波 波场快照和单炮记录如图１所示。
由图１的模拟效果分析可知：在均匀介质体中
网格点数 ｌＯｌ
难ｌ０ｌ 墨
２０ｌ ／网格点数 ｌＯｌ
图１ａ混合波场水平分量（左）、Ｐ波波场水平分量（中）和Ｓ波波场水平分量（右）快照（ｆ一３５０ｍｓ）