有限差分法

我们现在考虑方程(4.2.1)中 p 为常数的二维的情况，我们可以得到下面的方程：
∂ 2φ + ∂ 2φ + f (x, y)φ = q(x, y). ∂x 2 ∂y 2
设函数φ 在区域 D 内满足方程（4.2.6）式(区域 D 的边界为 G)。
(4.2.6)
7
区域 D 的离散化：
即通过任意的网络划分方法把区域 D 离散为许许多多的小单元。原则上讲这种网格分割是可以任 意的，但是在实际应用中，常常是根据边界 G 的形状，采用最简单，最有规律，和边界的拟合程度最 佳的方法来分割。常用的有正方形分割法和矩形分割法（如图 4.2.1）。有时也用三角形分割法（见 图 4.2.2）。对圆形区域，应用图（4.2.3）所示的极网络格式也许更方便些。这些网络单元通常称 为元素，网络点称为节点。
第四章 有限差分方法
4.1 引言
有限差分法：数值求解常微分方程或偏微分方程的方法。
物理学和其他学科领域的许多问题在被分析研究之后, 往往可以归结为常微分方程或偏微分方 程的求解问题。一般说来，处理一个特定的物理问题，除了需要知道它满足的数学方程外，还应当同 时知道这个问题的定解条件，然后才能设计出行之有效的计算方法来求解。
φ |G = g(s).
(4.2.3)
（2）第二类边界条件，或称诺伊曼(Neumann)问题 (g1 ≠ 0, g2 = 0) 。
∂φ ∂n
|G
=
g ( s ).
（3）第三类边界条件，或称混合问题 (g1 ≠ 0, g2 ≠ 0) 。
对于算符 L 为斯杜-刘维尔(Sturm-Liouville)算符的特定情况，即
1
有限差分法的具体操作分为两个部分：
(1)用差分代替微分方程中的微分，将连续变化的变量离散化，从而得到差分方程组的数学形式； (2)求解差分方程组。 在第一步中，我们通过所谓的网络分割法，将函数定义域分成大量相邻而不重合的子区域。通 常采用的是规则的分割方式。这样可以便于计算机自动实现和减少计算的复杂性。网络线划分的交点 称为节点。若与某个节点 P 相邻的节点都是定义在场域内的节点，则 P 点称为正则节点；反之，若节 点 P 有处在定义域外的相邻节点，则 P 点称为非正则节点。 在第二步中，数值求解的关键就是要应 用适当的计算方法，求得特定问题在所有这些节点上的离散近似值。
−
φ
0
)
⎤ ⎥
⎦
+
f 0φ0
=
q0
（4.2.18）
如果在 x 和 y 方向的步长分别相等， 即 h1 = h3 = hx 和 h2 = h4 = hy 时，则上式化为
φ1
− 2φ 0 hx 2
+φ3
+ φ2
− 2φ 0 hy 2
+φ4
+
f 0φ 0
=
q0
，
（4.2.19）
一般可以用角标来表示节点的标记，将上式写为
..., x = x + h = b , 然
n+1
n
后求出 f(x)在这些点上的近似值。显然步长 h 越小，近似解的精度就越好。与节点 xi 相邻的节点有
xi − h 和 xi + h ，因此在 xi 点可以构造如下形式的差值：
f (xi + h) − f (xi ), f (xi ) − f (xi − h),
有限差分法以变量离散取值后对应的函数值来近似微分方程中独立变量的连续取值。
在有限差分方法中，我们放弃了微分方程中独立变量可以取连续值的特征，而关注独立变量离 散取值后对应的函数值。但是从原则上说，这种方法仍然可以达到任意满意的计算精度。因为方程的 连续数值解可以通过减小独立变量离散取值的间格，或者通过离散点上的函数值插值计算来近似得 到。这种方法是随着计算机的诞生和应用而发展起来的。其计算格式和程序的设计都比较直观和简单， 因而，它的实际应用已经构成了计算数学和计算物理的重要组成部分。
⎞ ⎟
的差分表达式：
⎝ ∂y 2 ⎠ 0
⎛ ⎜ ⎝
∂
2
φ
⎞ ⎟
∂y 2 ⎠ 0
≈
2
h4 (φ 2
− φ 0 ) + h2 (φ 4 h2h4 (h2 + h4 )
−φ0)
.
当采用等步长 h2 = h4 = hy 时， 有
⎛ ⎜
∂
2φ
⎞ ⎟
≈ φ 2 − 2φ 0 + φ 4 .
⎝ ∂y 2 ⎠ 0
hy 2
⎞ ⎟ ⎠
0
≈
2
h3 (φ 1
−φ0) + h1h3 (h1
h1 (φ 3 + h3 )
−
φ0)
.
当用等步长 h1 = h3 = hx 时，上式成为
（4.2.14）
10
⎛ ⎜
∂
2φ
⎞ ⎟
≈ φ1 − 2φ 0 + φ 3 .
⎝ ∂x2 ⎠ 0
hx 2
用完全相同的计算方法，我们可以推导出
⎛ ⎜
∂
2φ
本书中我们将略去对差分法稳定性和收敛性理论的讨论，尽管这方面的内容是相当重要的。以下 的讨论中所讲到的各种差分格式，我们均假定求解方法满足稳定性和收敛性的要求。
5
4.2 有限差分法和偏微分方程
利用上节所介绍的微分的差分表示，我们就很容易地将微分方程离散化为差分方程组的形式。但 是由差分方程所得的解完全取决于待求微分方程的特性。正如我们在物理上所知道的，边界条件的情 况变化将会引起差分方程组的不同。在求解微分方程中，我们会遇到两类问题：一类是初始值问题； 另一类是边值条件的问题。在初始值问题中，部分边界上的函数值和部分的函数偏导值是给定的。通 常在这类问题中的独立变量之一是时间 t。在边界值问题中，边界上的信息是给定的。本书中我们仅 讨论后一类问题。
f ( xi + h) − f ( xi − h).
节点 xi 的一阶向前差分 节点 xi 的一阶向后差分 节点 xi 的一阶中心差分
2Leabharlann 与 xi 点相邻两点的泰勒展开式可以写为
f (xi
− h) =
f
( xi
)
−
hf
'
( xi
)
+
h2 2
f
''
( xi
)
−
h3 3!
f
'''
( xi
)
+
h4 4!
f
⎛ ∂φ ⎞ ⎝⎜ ∂x ⎠⎟ 0
≈
φ1 −φ0 h1
.
或
（4.2.7）
⎛ ∂φ ⎞ ⎝⎜ ∂x ⎠⎟ 0
≈
φ0 −φ3 h3
.
显然这种单侧差商的误差较大。
（4.2.8）
如果要寻求更精确的差分格式，我们可以引入待定常数 α , β ，由 φ1 和 φ3 的泰勒展开，构造出如下的 关系式：
α (φ1
− φ0
y
y
D
D
0
x0
x
4.2.1 求解区域的矩形分割。 4.2.2 求解区域的三角形分割。 4.2.3 求解区域的极网络分割。
8
用节点上的函数值来表示节点上偏导的数值。
设区域内部某节点 0 附近的各节点如图 4.1.1 所示。这里我们取步长 h 不相等的最一般情况。以
φ0,φ1,φ2,φ3,φ4 分别代表在节点 0，1，2，3，4 处φ 的函数值。如前所述，0 点的一阶偏导数可以通过先 前或向后的差商，由 1 和 3 节点近似写出
假定某方程形式上可以写为：
Lφ = q.
其中 L 为含偏微商的算符．
(4.2.1)
它的边界条件一般可写为：
φ
|G
+ g1(s)
∂φ ∂n
|G =
g2 (s).
G 表示场域 D 的边界， g1(s), g2 (s) 为边界上 s 点的逐点函数。
(4.2.2)
6
三类边界条件：
（1）第一类边界条件，或称为狄利克莱(Dirichlet)问题 (g1 = 0, g2 ≠ 0) 。
二阶偏导数的差分表达式（“五点格式”或“菱形格式”）
在（4.2.9）式中，如果令
⎛ ⎝⎜
∂φ ∂x
⎞ ⎠⎟
0
的系数为零，则有
α
和
β
间存在关系式：
α = h3 β. h1
（4.2.13）
将上式（4.2.13）代入（4.2.9）式中，并忽略 h 三阶以上的高次项，则得到表达式：
⎛ ⎜ ⎝
∂ 2φ ∂x 2
有限差分法的差分格式：
一个函数在 x 点上的一阶和二阶微商，可以近似地用它所临近的两点上的函数值的差分来表示。
如对一个单变量函数 f(x)，x 为定义在区间[a,b]的连续变量。以步长 h= Δ x 将[a,b]区间离散化，我
们得到一系列节点 x = a , 1
x=x +h,
2
1
x =x
3
2
+ h = a + 2Δx ,
上述差分步骤应用于偏微分：
例如，对于 f
=
f
(x,
y)
的情况，拉普拉斯算符在
0
点作用在此函数上的值
⎛ ⎜
∇
2
⎝
f
=
⎛ ⎜ ⎝
∂ ∂
2f x2
+∂2f ∂ y2
⎞ ⎟
⎞ ⎟
，也可
⎠⎠
以用临近的点上的函数值来表示出来。(见图 4.1.1, 且 h1 = h2 = h3 = h4 = h 时)
∇2 f ≈
⎛⎝⎜
∂φ ∂x