k k 1
u t
t k *t
* t
1 u
2
2 t
2
t k *t
* t o(t )
2 2
(4-3)
27
将上两式相加，略去高阶小量，整理得(i,j)点k时刻 的二阶时间微商为：
u
2
t
2

u i , j 2u i , j u i , j
k
k 1
k 1
叠加。该方法计算速度较快，但由于射线追踪中存在着 诸如焦散、多重路径等问题，故其一般只能适合于较简 单的模型，难以模拟复杂地层的波场信息。
19
傅里叶变换法是利用空间的全部信息对波场函数
进行三角函数插值，能更加精确地模拟地震波的传播
规律，同时，利用快速傅里叶变换(FFT)进行计算，还
可以提高运算效率，其主要优点是精度高，占用内存
k
u z
2

1 z
{ 2
1 12
[u
k i , j 2
u
k i, j2
]
4 3
[u
k i , j 1
u
k i , j 1
]
5 2
ui , j }
k
31
同理可得时间二阶、空间四阶精度的声波方程差分格式为：
u
k 1 i, j
2u
k i, j
u
2
k 1 i, j
2

v t
2
p
2
x
2
2

p
2
z
2
) S(t)
P t v x
一阶压力-速度方程组：
x 1 P ( ) t x 1 P ( ) t y
C (
v x

v y y

v z z
)
v y
v z
1 P ( ) t z
12
二、波动方程类型及其局限性
26
时间二阶、空间二阶差分格式推导如下： 将
k 1
u 在(i,j)点k时刻用Taylor展式展开：
k
k 1 i, j
ui, j ui, j
u t
t k *t
* t
1 u
2
2 t 2
t k *t
* t o( t )
2 2
(4-2)
ui, j ui, j
方程等）、震源和地下几何界面、物性参数（岩层
密度、速度等）情况下，研究弹性波或声波的传播
规律。
10
内容提纲
一、地震勘探基本原理
二、波动方程类型及其 局限性
三、数值算法类型及其局限性
四、声波方程的有限差分法数值模拟
11
二、波动方程类型及其局限性
1、声波方程： 二阶标量声波方程：
p
2
t
2
v (
k
s (t ) * (i i0 ) * ( j j0 )
(4-8)
32
s (k为震源函数， )
一般使用一个理论的雷克型子波代替，即：
s (t )

e
( 2f / ) t
2
2
cos 2ft
t

为时间
f 为中心频率，一般取为20-40HZ
为控制频带宽度的参数，一般取2-5 确定震源位置
实际地下介质非常复杂，所得到的炮集记录也包含
更多的地下信息。实际的炮集记录见图1-1和1-2。
4
图1-1 陆上某区实际地震记录
5
图1-2 海上某区实际地震记录
6
广义的地震反演即是从地震炮集记录出发，经 过复杂的去噪、速度分析以及偏移成像处理等手段 得到反映地下的地质结构的地震剖面。实际的地震 剖面见图1-3和1-4。
i 2, j2 i 2, j 1 i 1, j 1 i 1, j2 i 1, j 1 i 2, j2 i 2, j 1
i, j 2
i, j 1
i 2, j
i 1, j
i 1, j 1
i 1, j2
i, j
i 1, j
i 1, j 1 i 1, j2
2
2
h
1 12
2
{
k i 2, j
1 12
[u
k i 2, j
u
4 3
k i 2, j
]
4 3
[u
k i 1, j
u
5 2
k i 1, j
]
5 2
ui , j }
k

v t h
2
{
[u
u
k i 2, j
]
[u
k i 1, j
u
k i 1, j
]
ui, j
对于二维速度-深度模型，地下介质中地震波
的传播规律可以近似地用声波方程描述：
u
2
t
2
v (
2
u
2
x
2

u
2
z
) S (t ) 2
(4-1)
v ( x, z ) 是介质在点（x , z）处的纵波速度，
u
s (t )
为描述速度位或者压力的波场， 为震源函数。
24
空间模型网格化（如图4-1所示）：
胡克定律（这种理想介质称虎克固体），但波在实
际介质中传播时，是有能量损耗的，这就是所谓的
弹性波吸收。波在传播过程中，实际介质的不同部 位之间会出现某种摩擦力，称为内摩擦力或粘滞力。 这种力导致机械能向其他形式能量转换，最终转化 为热能消耗掉。
16
在地震勘探中，地震波传播的实际介质是十分
复杂的。在一定条件下，即震源作用时间短，作用
内容提纲
一、地震勘探基本原理
二、波动方程类型及其 局限性Biblioteka 三、数值算法类型及其局限性
四、声波方程的有限差分法数值模拟
18
三、数值算法类型及其优缺点
地震波波动方程数值模拟方法主要包括克希霍夫积
分法、傅里叶变换法、有限元法和有限差分法等。
克希霍夫积分法引入射线追踪过程，本质上是波动
方程积分解的一个数值计算，在某种程度上相当于绕射
t
2
(4-4)
28
同理可得(i,j)点k时刻的二阶空间微商分别为：
u
2
x
2
2

u i 1, j 2u i , j u i 1, j
k k k
h
k k
2
(4-5)
u z
2

u i , j 1 2u i , j u i , j 1
k
h
2
(4-6)
29
这就实现了用网个点波场值的差商代替了偏微 分方程的微商，将上三个式子代入(4-1)式中得：
散，发生畸变，且随走时的增加而增加。
35
Dablain给出了一个能有效减少频散的经验公式为：
h v min /( Gf N )
f N 为Nyquist频率，一般取为 主频的两倍
(4-9)
G为每个波长所占的网格点数，时间、空间为两阶差
分的情况G取8，而时间、空间为四阶差分的情况G
取4。
36
当所给震源函数、空间网格间隔、时间采样间隔以
33
(i i0 ) * ( j j 0 )
4.1 稳定性条件
对于特定的偏微分方程只有特定的几种有限差分
格式是无条件或有条件稳定的，（4-7）、（4-8）式即 是已被证明的有条件稳定格式，其稳定性条件分别为：
v max * t / h
2/2
(4-9)
vmax * t / h 3 / 8
ui , j (v * t / h) * (ui 1, j ui 1, j ui , j 1 ui , j 1 )
2 k k k k k 1
(
2h 4v t
2 2
2
h
2
)ui , j ui , j s (k ) * (i i0 ) * ( j j0 )
边界的处理是地震波场正演数值模拟的一个重要课题。
本次作业不涉及边界条件的使用，可通过增大模型
来避免边界反射干扰，有精力的同学可通过查阅资料获
得解决边界问题的方法。
38
本次大作业的具体要求为：
1、应用声波方程作为正演模拟的波动方程；
2、将所提供震源函数离散后绘图；
声波方程数值模拟
――地球物理场论
基础Ⅰ期末作业（1） 任课教师： 宋鹏
1
内容提纲
一、地震勘探基本原理
二、波动方程类型及其 局限性
三、数值算法类型及其局限性
四、声波方程的有限差分法数值模拟
2
一、地震勘探基本原理
＊ ▽ ▽ ▽ ▽ ▽
x
t
3
同相轴为双曲线，即反射波的时距曲线为双
曲线，反射波一个同相轴可带来一个地层的信息。
i 2, j
i 2, j 1 i 2, j2
z
i 2, j 1 i 2, j2
i, j 1
i, j 2
x
图4-1 差分网格划分示意图
25
网格间隔长度 时间采样步长
，
h
t
x ih
z jh
t kt
u
k i, j
u i, j
k
表示(i,j)点k时刻的波场值
能够描述纵、横波的传播规律，包括直达波、反
射波、透射波、折射波以及转换波等。 需要的已知条件包括： 1）震源函数 2）地层速度或根据方程的类型需要提 供的地层的其它弹性参数 3）边界条件