差分方程_有限差分法(FDM)实现声波方程数值模拟

0. 前言

实现声波方程的数值模拟，常见的方法包括有限差分法（Finite Difference Method-FDM）、伪谱法（Pseudo-spectral Method-PSM）、有限单元法（Finite Element Method-FEM），这里先介绍如何使用FDM-2、FDM-4（即三点有限差分法和五点有限差分法）求解声波方程，在之后的文章中在逐一介绍其他方法。

---

1. 模型建立

声波方程的一般形式如公式（1）所示，可以发现方程中存在对波动u的二阶时间偏导，以及对空间坐标x，y的二阶偏导。要直接求解（1）的解析解是比较困难的，我们常常会使用各种近似的数值解来逼近解析解。

有限差分法便是其中常用的一种数值计算方法，FDM的基本原理类似于我们求解斜率的方法，其本质上是使用了泰勒级数展开并略去高阶项，这里就不展开讲了。

• 1.1 偏微分项的处理

对于时间项的二阶偏导，使用三点差分格式：

我们在这里约定u的下角标代表空间上的离散点索引（index），上标代表时间上的差分索引，k时刻代表当前时刻，k+1代表下一时刻，k-1代表前一时刻。

对于空间上，在x和y方向上进行离散，我们假设在x和z方向上离散点的间距相同，即

。

对空间离散并网格化

关于空间的二阶偏导，我们可以使用三点计算，也可以使用五点计算。对于三点差分格式（称之为FDM-2），其形式与时间差分格式相同，对于五点差分格式(称之为FDM-4)，其形式为：

其中的系数是通过泰勒展开式求解出来的。

• 1.2 速度项以及源的处理

对于均匀介质，声波的传播速度相同，因此公式（1）中的

是一个常数。

代表波动源，我们可以使用雷克子波（零相位子波）或者高斯函数作为震源函数。这里使用雷克子波，其表达式为:

其中的

代表波的主频。它的差分形式也很简单，只需要将

换成

代入计算就好了。需要注意的是，为了避免使用突变信号，需要将上面的震源公式进行一定量的时移，即将

设置成