第十章 偏微分方程数值解法

第十章 偏微分方程数值解法
偏微分方程问题，其求解十分困难。

除少数特殊情况外，绝
大多数情况均难以求出精确解。

因此，近似解法就显得更为重要。

本章仅介绍求解各类典型偏微分方程定解问题的差分方法。

§1 差分方法的基本概念
1.1 几类偏微分方程的定解问题
椭圆型方程：其最典型、最简单的形式是泊松（Poisson ）方程
),(22
2
2y x f y
u x u u =∂∂+∂∂=∆ 特别地，当
0),(≡y x f 时，即为拉普拉斯（Laplace ）方程，又称
为调和方程
22
22
=∂∂+∂∂=∆y
u
x u u Poisson 方程的第一边值问题为
⎪⎩
⎪⎨
⎧Ω
∂=Γ=Ω∈=∂∂+∂∂Γ∈),(),(),(),(),(22
22y x y x u y x y x f y u
x u y x ϕ 其中
Ω为以Γ为边界的有界区域，Γ为分段光滑曲线，
ΓΩ
称为定解区域，),(y x f ，),(y x ϕ分别为Ω，Γ上的已知连
续函数。

第二类和第三类边界条件可统一表示为
),(),(y x u u y x ϕα=⎪⎪⎭
⎫ ⎝⎛+∂∂Γ∈n 其中n 为边界Γ的外法线方向。

当0=α时为第二类边界条件， 0≠α时为第三类边界条件。

抛物型方程：其最简单的形式为一维热传导方程
2
20(0)u u
a a t x
∂∂-=>∂∂ 方程可以有两种不同类型的定解问题：
初值问题
⎪⎩
⎪⎨⎧+∞
<<∞-=+∞<<-∞>=∂∂-∂∂x x x u x t x u a t
u )()0,(,00
22
ϕ
初边值问题
2
212
00,0(,0)()0(0,)(),(,)()0u u
a t T x l t x u x x x l
u t g t u l t g t t T
ϕ⎧∂∂-=<<<<⎪∂∂⎪⎪
=≤≤⎨⎪==≤≤⎪⎪⎩
其中
)(x ϕ，)(1t g ，)(2t g 为已知函数，且满足连接条件
)0()(),
0()0(21g l g ==ϕϕ
边界条件)(),(),(),0(21t g t l u t g t u ==称为第一类边界条
件。

第二类和第三类边界条件为
)()()()(22101t g u t x u t g u t x u l
x x =⎥⎦
⎤⎢⎣⎡+∂∂=⎥⎦⎤
⎢⎣⎡-∂∂==λλ
T t ≤≤0
其中0)(1≥t λ，0
)(2≥t λ。

当
0)()(21≡=t t λλ时，为第二
类边界条件，
否则称为第三类边界条件。

双曲型方程:
最简单形式为一阶双曲型方程
=∂∂+∂∂x
u
a t u 物理中常见的一维振动与波动问题可用二阶波动方程
2
2
2
22x
u a t u
∂∂=∂∂
描述，它是双曲型方程的典型形式。

方程的初值问题为
⎪⎪⎪⎩⎪
⎪⎪⎨⎧+∞
<<∞-=∂∂=+∞
<<∞->∂∂=∂∂=x x t
u x x u x t x u a t u t )()()0,(,0022
222ψϕ
边界条件一般也有三类，最简单的初边值问题为
22
2
220
12
00,0(,0)(),
()0(0,)(),(,)()0t u u a t T x l t x u u x x x x l t u t g t u l t g t t T
ϕψ=⎧∂∂==<<<<⎪∂∂⎪⎪∂⎪
==≤≤⎨∂⎪
⎪==≤≤⎪⎪⎩
1.2 差分方法的基本概念
差分方法又称为有限差分方法或网格法，是求偏微分方程定 解问题的数值解中应用最广泛的方法之一。

它的基本思想是：先对求解区域作网格剖分，将自变量的连 续变化区域用有限离散点（网格点）集代替；将问题中出现的连 续变量的函数用定义在网格点上离散变量的函数代替；通过用网 格点上函数的差商代替导数，将含连续变量的偏微分方程定解问 题化成只含有限个未知数的代数方程组（称为差分格式）。

如果 差分格式有解，且当网格无限变小时其解收敛于原微分方程定解 问题的解，则差分格式的解就作为原问题的近似解（数值解）。

因此，用差分方法求偏微分方程定解问题一般需要解决以下问题： （1）选取网格；
（2）对微分方程及定解条件选择差分近似，列出差分格式； （3）求解差分格式；
（4）讨论差分格式解对于微分方程解的收敛性及误差估计。

下面，用一个简单的例子来说明用差分方法求解偏微分方程 问题的一般过程及差分方法的基本概念。

设有一阶双曲型方程初值问题。

⎪⎩
⎪⎨⎧=+∞
<<∞->=∂∂+∂∂)()0,(,00x x u x t x u a t
u
ϕ
（1） 选取网格：
分，最简单
kh
x x k ==(0,1,2
,0,1,2,)j t t j k j τ===±±=将
D 分成许
多小矩形
区域。

这些直线称为网格线，其交点称为网格点，也称为节点，
h 和τ
分别称作
x 方向和t
方向的步长。

这种网格称为矩形网格。

（2） 对微分方程及定解条件选择差分近似，列出差分格式： 如果用向前差商表示一阶偏导数，即
2
211(,)12(,)(,)(,)(,)2(,)(,)(,)
2
k j k j k j k j k j x x t k j k j k j t x t u x t u x t u
h u x h t x h u x t u x t u
u x t t θτθττ++-∂''=-+∂-∂''=-+∂
其中
1,021<<θθ。

方程 0
u u
a t x
∂∂+=∂∂ 在节点
),(j k t x 处可表示为
h
t x u t x u a
t x u t x u j k j k j k j k )
,(),()
,(),(11-+-++τ
),(2
),(2122
2
j k x
j k t t h x u ah
t x u θτθτ+''++''=
),2,1,0,,2,1,0()
,( =±±==j k t x R j k
其中
(,0)()(0,1,2,)k k u x x k ϕ==±±。

由于当
τ
,h 足够小时，在式
中略去
),(j k t x R ，就得到一个与方程相近似的差分方程
1,1,,,0
k j k k k j
j j
u u u u a
h
τ
++--+=
此处，
j k u ,可看作是问题的解在节点),(j k t x 处的近似值。

同初值
条件
),2,1,0()
(0, ±±==k x u k k ϕ
结合，就得到求问题的数值解的差分格式。

式
)
()
,(2
),(2),(1222h O t h x u ah
t x u t x R j k x j k t j k +=+''++''=τθτθτ
称为差分方程的截断误差。

如果一个差分方程的截断误差为
)(p
q h O R +=τ，则称差分方 程对
t
是
q 阶精度，对x 是p 阶精度的。

显然，截断误差的阶数
越大，差分方程对微分方程的逼近越好。

若网格步长趋于0时，差分方程的截断误差也趋于0，则称 差分方程与相应的微分方程是相容的。

这是用差分方法求解偏微 分方程问题的必要条件。

如果当网格步长趋于0时，差分格式的解收敛到相应微分方 程定解问题的解，则称这种差分格式是收敛的。

§2 椭圆型方程第一边值问题的差分解法
本节以Poisson 方程为基本模型讨论第一边值问题的差分方法。

2.1 差分格式的建立
考虑Poisson 方程的第一边值问题
⎪⎩
⎪⎨⎧Ω
∂=Γ=Ω∈=∂∂+∂∂Γ∈),(),(),()
,(),(2222y x y x u y x y x f y u
x
u y x ϕ
取τ
,h 分别为
x 方向和y 方向的步长，如图所示，以两族平行
线
kh x x k ==，(,0,1,2,)j y y j k j τ===±±将定
解区域剖分成矩形网 格。

节
点
的
全
体
记
为
{(,),,,}k j k j R x y x kh y j k j τ===为整数。

定解区
域内部的节点称为内点，记内点集
Ω R 为τ
h Ω。

边界
Γ与网格
线的交点称为边界点，边界点全体记为
τ
h Γ。

与节点
),(j k y x 沿x 方
向或
y
方向只差一个步长的点
),(1j k y x ±和),(1±j k y x 称为节点
),(j k y x 的
相邻节点。

如果一个内点的四个相邻节点均属于ΓΩ ，称为正
则内点，正内点的全体记为
)
1(Ω
，至少有一个相邻节点不属于
ΓΩ
的内点称为非正则内点，非正则内点的全体记为)
2(Ω。

问题是要
求出第一边值问题在全体内点上的数值解。

为简便，记
)
,(),(j k y x j k =，
(,)(,)
k j u k j u x y =，
),(,j k j k y x f f =。

对正则
内点)
1(),(Ω∈j k ，由二阶中心差商公式
442
2
(2
(,)
2
(4)2
22(4)2
2(,)
(1,)(,)(,)(1,)
12
(1,)2(,)(1,)(12
(,1)2(,)(,1)(,12
x
k j x k j
y k j u k j u k j u k j u k j u h h h u x h u k j u k j u k j h u x h u
u k j u k j u k j u x y y ττ+----∂=-∂+-+-=-∂+-+-=-∂
Poisson 方程
22
22
(,)u u
f x y x y
∂∂+=∂∂ 在点
),(j k 处可表示为
22
(1,)2(,)(1,)(,1)2(,)(,
u k j u k j u k j u k j u k j u k h τ
+-+-+-++
其中
(),(12
),(12),(2
2)4(21)4(244+=+++=τ
τθτθh O y x u y h x u h j k R j k x j k x
为其截断误差表示式，略去
),(j k R ，即得与方程相近似的差分方程
j
k j k j k j k j
k j k j k f u u u h
u u u ,2
1
,,1,2
,1,,122=+-+
+--+-+τ
式中方程的个数等于正则内点的个数，而未知数j k u ,则除了包含
正则内点处解
u 的近似值外，还包含一些非正则内点处u 的近
似值，因而方程个数少于未知数个数。

在非正则内点处Poisson 方程的差分近似不能按上式给出，需要利用边界条件得到。

边界条件的处理可以有各种方案，下面介绍较简单的两种。

（1）直接转移
用最接近非正则内点的边界点上的u 值作为该点上u
值的
近似，这就是边界条件的直接转移。

例如，点),(j k P 为非正则内点，其最接近的边界点为Q 点，则有
)
2(,),()
()(Ω∈==j k Q Q u u j k ϕ
上式可以看作是用零次插值得到非正则内点处u 的近似值，容易
求出，其截断误差为)(τ+h O 。

将上式代入，方程个数即与未知数
个数相等。

（2）线性插值
这种方案是通过用同一条网格线上与点P 相邻的边界点与
内点作线性插值得到非正则内点),
(j k P 处u
值的近似。

由点
R 与
T
的线性插值确定
)(P u 的近似值j k u ,，得
)()(,T d
h d
R d h h u j
k ϕϕ+++=
其中RP
d =，其截断误差为
)(2
h O 。

将其与方程相近似的差分方
程联立，得到方程个数与未知数个数相等的方程组，求解此方程 组可得Poisson 方程第一边值问题的数值解。

上面所给出的差分格式称为五点菱形格式，
j
k j k j k j k j
k j k j k f u u u h
u u u ,2
1
,,1,2
,1,,122=+-+
+--+-+τ
实际计算时经常取τ=h ，此时五点菱形格式可化为
1,1,,1,1,,21
(4)k j k j k j k j k j k j u u u u u f h
+-+-+++-= 简记为
21h
◇
j k j k f u ,,= 其中◇
j k j k j k j k j k j k u u u u u u ,1,1,,1,1,4-+++=-+-+。

例1 用五点菱形格式求解拉普拉斯（Laplace ）方程第一边值问题
⎪⎩
⎪⎨⎧Ω∂=Γ++=Ω∈=∂∂+∂∂Γ])1lg[(),(),()
,(222222y x y x u y x y x f y u
x
u 其中}1,0),{(≤≤=Ωy x y x 。

取3
1
==τh 。

（0,0）
[解]程组
⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧=-+++=-+++=-+++=-+++0
)4(10)4(10)4(10)4(12,21,23,22,12,322,11,13,12,02,221,20,22,21,11,321,10,12,11,01,22u u u u u h
u u u u u h u u u u u h u u u u u h 代入边界条件
1,02,00,10,21,32,33,13,21625
lg ,lg
991013
lg ,lg
992534
lg ,lg
993740
lg ,lg
99
u u u u u u u u ⎧
==⎪⎪
⎪==⎪⎨
⎪==⎪⎪⎪==⎪⎩
其解为 2756919
.01
,1=u ，
4603488
.01,2=u ，
3467842
.02,1=u ，5080467
.02
,2=u
当
τ
=h 时
，
对
1,1,,1,1,,21
(4)k j k j k j k j k j k j u u u u u f h
+-+-+++-= 利用点),(j k ，)1,1(-±j k ，)1,1(+±j k 构造的差分格
式，称为五点 矩形格式，简记为
2
21h □
j k j k f u ,,=
其
中
□
j k j k j k j k j k j k u u u u u u ,1,11,11,11,1,4-+++=--+--+++，其
截断误差为
2444
4
4244(,)
(,)6()12k j h u u u R k j O h O x x y y ⎛⎫∂∂∂=+++= ⎪∂∂∂∂⎝⎭
五点菱形格式与矩形格式的截断误差均为
)(2
h O ，称它们具有 二阶精度。

如果用更多的点构造差分格式，其截断误差的阶数可
以提高，如利用菱形格式及矩形格式所涉及的所有节点构造出的 九点格式就是具有四阶精度的差分格式。

§3 抛物型方程的差分解法
以一维热传导方程
)0(022
>=∂∂-∂∂a x
u a t u
为基本模型讨论适用于抛物型方程定解问题的几种差分格式。

3.1 差分格式的建立
首先对
xt 平面进行网格剖分。

分别取τ,h 为x
方向与
t 方向
的步长，用两族平行直线
(0,1,2,)k x x kh k ===±±，
)2,1,0( ===j j t t j τ，将xt
平面剖分成矩形网格，节点
为(,)
(0,1,2,
,0,1,2)k j x t k j =±±=。

为
简便，记
)
,(),(j k t x j k =，
)
,(),(j k t x u j k u =，)
(k k x ϕϕ=，
)
(11j j t g g =，
221122(),(),()j j j j j j g g t t t λλλλ===。

（一）微分方程的差分近似
在网格内点
),(j k 处，对
t
u ∂∂分别采用向前、向后及中心差商公式
(,)
(,)
2(,)
(,1)(,)
()(,)(,1)
()
(,1)(,1)()
2k j k j k j u u k j u k j O t u u k j u k j O t u u k j u k j O t
ττ
ττ
ττ
∂+-=
+∂∂--=
+∂∂+--=+∂
一维热传导方程
)0(022>=∂∂-∂∂a x
u
a t u
可分别表示为
()
,1(),(2),1()
,()1,(2
O h
j k u j k u j k u a j k u j k u +=-+-+--+ττ
()
,1(),(2),1()1,()1,(2
O h
j k u j k u j k u a j k u j k u +=-+-+---+ττ
22
2
(,1)(,1)(1,)2(,)(1,)(2u k j u k j u k j u k j u k j a O h h
ττ+--+-+--=+
由此得到一维热传导方程的不同差分近似
022
,1,,1,1,=+----++h
u u u a
u u j
k j k j k j
k j k τ
022
,1,,11
,,=+----+-h
u u u a
u u j
k j k j k j k j k τ
,1,1
1,,1,2
202k j k j k j k j k j
u u u u u a
h
τ
+-+---+-=
上述差分方程所用到的节点各不相同。

其截断误差分别为
)(2
h O +τ，)(2h O +τ和)(22h O +τ。

因此，它们都与一维热传导方程相容。

如果将式
2
2
(,1)(,1)(1,)2(,)(1,)(2u k j u k j u k j u k j u k j a O h τ
τ+--+-+--=
中的j k u ,用
)(2
1
1,1,-++j k j k u u 代替，则可得到又一种差分近似
022
,11,1,,11
,1,=+------++-+h
u u u u a
u u j
k j k j k j k j k j k τ
差分方程用到四个节点。

由Taylor 公式容易得出
)()(2
1
21,1,,τO u u u j k j k j
k ++=-+
故其的截断误差为
⎪⎪⎭
⎫ ⎝
⎛++2
222
)(h O h O τ
τ。

因而不是对任意的
0,→τh ，
此差分方程都能逼近热传导方程
)0(022>=∂∂-∂∂a x
u
a t u ，
仅当
()o h τ=时，才成立。

综上可知，用不同的差商公式可以得到微分方程的不同的差
分近似。

构造差分格式的关键在于使其具有相容性、收敛性和稳 定性。

前面三个方程都具有相容性，而此方程则要在一定条件下 才具有相容性。

（二）初、边值条件的处理
为用差分方法求解定解问题初值问题
⎪⎩
⎪⎨⎧+∞
<<∞-=+∞<<-∞>=∂∂-∂∂x x x u x t x u a t
u )()0,(,00
22
ϕ
初边值问题
221200,0(,0)()0(0,)(),(,)()0u u
a t T x l t x u x x x l
u t g t u l t g t t T ϕ⎧∂∂-=<<<<⎪∂∂⎪⎪
=≤≤⎨⎪==≤≤⎪⎪⎩
还需对定解条件进行离散化。

对初始条件及第一类边界条件，可直接得到 k
k k x u u ϕ==)0,(0,
(0,1, 0,1,
,)k k n =±=或
j
j j n j
j j g t l u u g t u u 2,1,0),(),0(====
)1,,1,0(-=m j
其中
τ
T
m h l n ==,
对第二、三类边界条件
)
()()()(2210
1t g u t x u t g u t x u l
x x =⎥⎦
⎤
⎢⎣⎡+∂∂=⎥⎦⎤⎢⎣⎡-∂∂==λλ
T
t ≤≤0
需用差分近似。

下面介绍两种较简单的处理方法。

（1）
在左边界)0(=x 处用向前差商近似偏导数
x
u ∂∂，
在右边界
)(l x =处用向后差商近似
x
u
∂∂，即
)
(),1(),()
(),0(),1()
,()
,0(h O h
j n u j n u x
u h O h j u j u x u j n j +--=∂∂+-=∂∂
),,1,0(m j =
则得边界条件的差分近似为
⎪⎪⎩
⎪⎪⎨⎧=+-=---j
j n j j n j n j
j j j
j g u h
u u g u h
u u 2,2,1,1,01,0,1λλ
),,1,0(m j =
其截断误差为
)(h O 。

（2）用中心差商近似
x
u ∂∂,即
2(0,)2(,)(1,)(1,)
()
2(1,)(1,)
()
2j n j u u j u j O h x h u u n j u n j O h x h
∂--=+∂∂+--=+∂
),,1,0(m j =
则得边界条件的差分近似为
⎪⎪⎩
⎪⎪⎨
⎧=+-=---+-j j n j j n j n j j j j
j g u
h u u g u h
u u 2,2,1,11,01,1,12λλ
),,1,0(m j =
其截断误差为
)(2
h O 。

误差的阶数提高了，但出现定解区域外的 节点),1(j -和),1(j n +，这就需要将解拓展到定解区域外。

可以通
过用内节点上的
u
值插值求出
j u ,1-和j n u ,1+，也可以假定热传导方
程在边界上也成立，将差分方程扩展到边界节点上，由此消去j u ,1-
和
j
n u ,1+。

（三）几种常用的差分格式
以热传导方程的初边值问题
221200,0(,0)()0(0,)(),(,)()0u u
a t T x l t x u x x x l
u t g t u l t g t t T ϕ⎧∂∂-=<<<<⎪∂∂⎪⎪
=≤≤⎨⎪==≤≤⎪⎪⎩
为例给出几种常用的差分格式。

（1）古典显式格式
令
2
h
a r τ=，则
022
,1,,1,1,=+----++h
u u u a
u u j
k j k j k j
k j k τ
可
改
写
成
j k j k j k j k ru u r ru u ,1,,11,)21(-+++-+=
将其与初始条件及第一类边界条件 k
k k x u u ϕ==)0,(0,
(0,1, 0,1,
,)k k n =±=或
j
j j n j
j j g t l u u g t u u 2,1,0),(),0(====
)1,,1,0(-=m j
结合，我们得到求解此问题的一种差分格式
⎪⎩⎪
⎨⎧=====-=-=+-+=-++),,1,0(,),,1,0(,,1,0,1,,2,1()21(2,1,00,,1,,11,m j g u g u n k u m j n k ru u r ru u j j n j j
k
k j k j k j k j k ϕ
由于第0层
)0(=j 上节点处的u 值已知)(0,k k u ϕ=，由此即可算出
u 在
第一层
)1(=j 上节点处的近似值1,k u 。

重复使用此式，可以逐层计
算出所有的
j k u ,，因此此差分格式称为古典显式格式。

又因式中
只出现相邻两个时间层的节点，故此式是二层显式格式。

（2）古典隐式格式
将
式
22
,1,,11
,,=+----+-h
u u u a
u u j
k j k j k j k j k τ
整理并与初始条件及第一类边界条件式联立，得差分格式如下
,1,1,1,11,1,00,1,2(2)(1,2,,1,0,(0,1,,),(0,1,,)k j k j k j k j k j k k
j j n j j
u u r u u u k n j u k n u g u g j m ϕ++++-+⎧=+-+=-=⎪
==⎨⎪===⎩
其中
2h
a r τ
=。

虽然第0层上的
u
值仍为已知，但不能由上式直
接计算以上各层节点上的值
j k u ,，必须通过解下列线性方程组
,11,1,11,11,1,11,(2)(12)k j k j k j k j k j k j k j u r u u u ru r u ru ++++-+-++++
--+=-++-
1,11,2,12,2,12,
1,11,120
1
2001212j j j j n j n n j n j r
r u u rg
r r r u u u u u u rg r r r r r ++-+--+-+-⎡⎤⎢⎥+-+-⎡⎤⎡⎢⎥⎢⎥⎢⎢⎥⎢⎥⎢⎢⎥⎢⎥⎢=⎢⎥⎢⎥⎢
⎢⎥⎢⎥⎢
⎢⎥
⎢⎥⎢+-+-⎢⎥
⎣⎦⎣
⎢⎥-+⎣
⎦
)1,,1,0(-=m j
才能由
j k u ,计算1,+j k u ，故此差分格式称为古典隐式格式。

此方程
组是三对角方程组，且系数矩阵严格对角占优，故解存在唯一。

（3）Richardson 格式
Richardson 格式是将式
,1,1
1,,1,2
20
2k j k j k j k j k j
u u u u u a
h
τ
+-+---+-=
整理后与初始条件及第一类边界条件式联立。

其计算公式为
⎪⎩
⎪
⎨⎧=====-=+-+=-+-+),,1,0(,),,1,0(,1,,2,1()2(22,1,00,,1,,11,1,m j g u g u n k u n k u u u r u u j j n j j k
k j k j k j k j k j k ϕ
这种差分格式中所涉及的节点出现在
1,,1+-j j j 三层上，故为三层
显式格式。

Richardson 格式是一种完全不稳定的差分格式，因此它在实际计算中是不能采用的。

（4）杜福特-弗兰克尔（DoFort-Frankel ）格式
DoFort-Frankel 格式也是三层显式格式，它是由式
22
,11,1,,11
,1,=
+------++-+h
u u u u a
u u j
k j k j k j k j k j k τ
与初始条件及第一类边界条件式结合得到的。

具体形式如下：
⎪⎪⎩
⎪
⎪⎨⎧
=====-=+-+++=--++),,1,0(,),,1,0(,,2,1(2121)(2122,1,00,1,,1,11,m j g u g u n k u n k u r r u u r r u j j n j j k
k j k j k j k j k ϕ
用这种格式求解时，除了第0层上的值
0,k u 由初值条件得到，必
须先用二层格式求出第1层上的值
1,k u ，然后再按上式逐层计
算
),,3,2(,m j u j k =。

（5）六点隐式格式
对二阶中心差商公式
22
2
(,)
(1,)2(,)(1,)k j u
u k j u k j u k j x
h
∂+-+-=∂
如果用
2
2
x u ∂∂在点
)1,(+j k 与点),(j k 处的二阶中心差商的平
均值
近似
2
2x u ∂∂在
)2
1
,(+j k 处的值，即 22,22
,,11,11,1,1)
2
1,(2
2
h u
u u u u u x u j k j k j k j k j k j k +-++-=
∂∂++-++++
同时
t
u ∂∂在点
)2
1
,(+j k 处的值也用中心差商近似，即 ,1,1
(,)
2
k j k j
k j u u u x
τ
++-∂=
∂
这样又得到热传导方程的一种差分近似
2,2(2,,11,11,1,12,1,-++---++-++++k j k j k j k j k j
k j k u u u u u h
a
u u τ
其截断误差为
)(22h O +τ，将上式与初始条件及第一类边界条件式 联立并整理，得差分格式
,1,1,1,11,11,,,00,1,2(2)(222(1,2,,1,0,1,,1)(0,1,,),(0,1,,)k j k j k j k j k j k j k j k k k j j n j j
r r u u u u u u u u k n j m u k n u g u g j m ϕ++++-++⎧
=+-++-+⎪⎪⎪=-=-⎨
⎪==⎪
===⎪⎩
此格式涉及到六个节点，它又是隐式格式，故称为六点隐式格式。

与古典隐式格式类似，用六点格式由第
j 层的值j k u ,计算第1+j
层的值1,+j k u 时，需求解三对角方程组
1,12,12,11,11/20
/21/200
/21/2/21j j n j n j r r u r r r u u u r r r r r ++-+-++-⎡⎤
⎢⎥-+-⎡⎤
⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥
⎢⎥-+-⎢⎥⎣⎦
⎢⎥-+⎣⎦
)1,,1,0(-=m j
此方程组的系数矩阵严格对角占优，故仍可用追赶法求解。

例2 用古典显式格式求初边值问题
⎪⎪⎪⎩
⎪⎪⎪⎨
⎧≤≤==≤≤=<<<<=∂∂-∂∂3
09),3(,0),0(30)0,(30,30022
2t t u t u x x x u x t x u
t u 的数值解，取
5.0,1==τh 。

[解] 这里
21,0.5
a a r h
τ
===，
2
)(x
x =ϕ，
9)(,0)(21==t g t g 。

由格式
⎪⎩⎪
⎨⎧=====-=+-+=-++),,1,0(,),,1,0(,1,,2,1()21(2,1,00,,1,,11,m j g u g u n k u n k ru u r ru u j j n j j
k
k j k j k j k j k ϕ
可得到
⎪⎩⎪⎨⎧========+=-++)
6,,1,0(9
,0)3,2,1,0()5,,1,0,2,1()(5.0,3,0220,,1,11, j u u k k
x u j k u u u j j
k k j k j k j k 将初值
0,k u 代入上式，即可算出
2)04(5.0)(5.00,00,21,1=+⨯=+⨯=u u u 5)19(5.0)(5.00,10,31,2=+⨯=+⨯=u u u
将边界条件
01,0=u ，91,3=u 及上述结果代入又可求得
9,5.5,5.2,02,32,22,12,0====u u u u
如此逐层计算，得全部节点上的数值解为
0,31,32,33,30,40, 2.75, 5.75,9,0,
u u u u u ====
=
§4 双曲型方程的差分解法
对二阶波动方程
2
2
2
22x
u a t u ∂∂=∂∂
如果令1u
v t
∂=∂，x u
v ∂∂=2，则方程可化成一阶线性双曲型方程
组
⎪⎪⎩⎪⎪⎨
⎧∂∂=
∂∂∂∂=∂∂x v t
v x v a t v 1
22
21
记
T
v )
,(21v v =，则方程组可表成矩阵形式
x A x a t ∂∂=∂∂⎥⎦
⎤⎢⎣⎡=∂∂v v v 0102
矩阵
A 有两个不同的特征值a ±=λ，故存在非奇异矩阵P ，使得
Λ=⎥
⎦
⎤⎢⎣⎡-=-a a PAP 001
作变换
T
w w P )
,(21==v w ，方程组可化为
x
t ∂∂Λ=∂∂w
w
方程组由二个独立的一阶双曲型方程联立而成。

因此本节主要讨 论一阶双曲型方程的差分解法。

4.1几种简单的差分格式
考虑一阶双曲型方程的初值问题
⎪⎩
⎪
⎨⎧+∞
<<∞-=+∞<<∞->=∂∂+∂∂x x x u x t x u a t
u )()0,(,00ϕ 将xt 平面剖分成矩形网格，取x 方向步长为t h ,方向步长为τ，
网
格
线
为
(0,1,2,)
k x x kh k ===±±，
)2,1,0( ===j j t t j τ。

为简便，记
),(),(j k t x j k =，)(),,(),(k k j k x t x u j k u ϕϕ==。

以不同的差商近似偏导数，可以得到方程的不同的差分近似
0,,1,1,=-+-++h u u a
u u j
k j k j
k j k τ
0,1,,1,=-+--+h
u u a u u j
k j k j k j k τ
02,1,1,1,=-+--++h
u u a
u u j
k j k j
k j k τ
截断误差分别为)(h O +τ，)(h O +τ与
)(2
h O +τ。

结合离散化的初始条件，可以得到几种简单的差分格式
)
,2,1,0,,2,1,0()(0,,,1,1, =±±=⎪⎩⎪⎨⎧=--=++j k u u u ar u u k
k j k j k j k j k ϕ
)
,2,1,0,,2,1,0()(0,,1,,1, =±±=⎪⎩⎪⎨⎧=--=-+j k u u u ar u u k
k j k j k j k j k ϕ
,1,1,1,,0()
(0,1,2,,0,1,2,)2k j k j k j k j k k ar u u u u k j u ϕ++-⎧
=--⎪=±±=⎨
⎪=⎩
其中
r h
τ
=。

如果已知第
j
层节点上的值
j k u ,，按上面三种格式就
可求出第1+j 层上的值1,+j k u 。

因此，这三种格式都是显式格式。

如果对t
u ∂∂采用向后差商，x
u
∂∂采用向前差商，则方程可化成
)()
,(),1()1,(),(=++-++--h O h
j k u j k u a j k u j k u ττ
相应的差分格式为
)
,2,1,0,,2,1,0()(0,1,1,1,1, =±±=⎪⎩⎪⎨⎧=--=++++j k u u u ar u u k
k j k j k j k j k ϕ
此差分格式是一种隐式格式，必须通过解方程组才能由第
j 层节
点上的值j k u ,求出第1+j 层节点上的值1,+j k u 。

例对初值问题
⎪⎩
⎪⎨⎧+∞
<<∞-=+∞<<∞->=∂∂+
∂∂x x x u x t x u
t u )()0,(,00ϕ 其中
⎪⎪⎩⎪
⎪⎨⎧>=<=0
00
2
101)(x x x x ϕ 用差分格式求其数值解
)
4,3,2,1(,=j u j k ，取
2
1
=
=h r τ。

[解] 记
),2,1,0( ±±==k kh x k ，由初始条件
⎪⎪⎩⎪
⎪⎨⎧==--===
,2,100
2
1,2,11)(k k k x k k
ϕϕ
按差分格式
2
,1,0,,2,1,0()(0,,,1,1, =±±=⎪⎩⎪⎨⎧=--=++j k u u u ar u u k
k j k j k j k j k ϕ
计算公式为
j
k j k j k u u u ,1,1
,2
1
23++-=
计算结果略。

如果用差分格式
2
,1,0,,2,1,0()(0,,1,,1, =±±=⎪⎩⎪⎨⎧=--=-+j k u u u ar u u k
k j k j k j k j k ϕ
求解，计算公式为
)(2
1
,1,1
,j k j k j k u u u -++=
计算结果见表略。

与准确解 11(,)200
x t u x t x t
x <⎧⎪⎪
==⎨⎪>⎪⎩ 比较知，按前一个差分格式所求得的数值解不收敛到初值问题
的解，而后一个差分格式的解收敛到原问题的解。