应用PDE讲义13_扩散方程差分解

2 ∆
整理成关于 维未知向量
， 的常微分方程组
， 1，2， ，
，，
T
0 其中 Jacobi 矩阵
21 12 1
12 1 ∆
12 1 12
如果把问题的解 ， 看成 ， 平面上方的曲面，上述过程相当
于沿一系列直线计算曲面的横切面，如下图所示。由于 Jacobi 矩阵的
特征值 介于
4 0
∆
4
致使随着空间网格∆ 加密，常微分方程组变得“刚性”，必须特别给 予处理。
等因素有时需要选用隐式方法，但使用显式算法远比隐式方便．
在隐式算法中，不能直截了当地求出 ，需要借助于迭代方法，
而迭代过程的实质是逐步显式化．例如，在第 步迭代，可以有
，
，
1
，
这是最简单的预报－校正算法．
根据积分中值定理可知，存在
， ，使得
，
，
但是一般说来 是不可能准确地知道．很自然地想到如果在 ， 中取几个点，然后用这些点处函数值的加权平均来近似 ， ，这样可得到更有效的显式 Runge‐Kutta 类 级算法的一
，
，
则可以得到 或
，
Euler 格式
10
，
向后 Euler 格式
对 Euler 方法还可以给出其它两种不同的解析解释，以 Euler 格式为 例。其一是数值微分，用最简单的向前差商来近似代替导数，即
，
其二是幂级数展开，将 的项，有
在 r 作 Taylor 展开，并略上 以上
，
，
2!
，
上述的每一种解释，都可对 Euler 方法进行改进，导出更为有效的办 法的途径，用数值微分来近似似乎是最自然的，但其实它是三种推广
应用偏微分方程与科学计算 讲义（十三）
Lecture Notes on Applied Partial Differential Equations and
Scientific Computing No. 13
马石庄
2010．11．01．北京
1
第 13 讲 发展问题差分近似：扩散方程
教学目的：研究一维扩散方程的有限差分的完全离散和半离散格式， 讨论差分格式的经济实用性，建立稳定性分析方法，探讨半离散格式 的刚性问题解法。 主要内容： §1 线上方法 ............................................................................................ 3
， 1，2， ， 相应的特征向量记成 ，解有形式
假定 Re 0 ， ∞时，有
1，2， ， ，即系统是渐近稳定的。当
0
则称此为系统的暂态解，而称 为稳态解．各个 称作系统的齐
次方程组的解分量（简称解分量）．暂态解可表示成为解分量
的线性组合．实部 确定
， 1，2， ， 的衰减特性，而虚部 确定振荡特性．
般形式
，
其中 是在配置点 处的加权因子。例如经典的四级 Runge‐Kutta 算 法为
12
，
，
2
2
，
2
2
，
22 6 无论是那种—步算法，都可以用下面的公式统一地描述．
，， ，
其中 ， ， 称为增量函数．显式单步法为
，，
对于单次方法，设 是方程的精确解，把
，，
称为在 的局部截断误差．所谓“局部”，是假定在 前各步没有
如上图示，可以进一步表明常微分方程（组）初值问题的解法在科学 计算中的重要地位．如果能有效地解决了这一中心地位的问题，对其 建立有效的算法和软件，并把各类计算问题都成功地转化为这个问题 求解，有可能建立起科学计算的新途径。 1.1 “刚性”问题
既然任何一个高阶常微分方程都可通过变量代换化为一阶常微 分方程组初值问题，不失一般性可以考虑线性系统
最基本常用的数值方法是单步方法，即计算 公式中只明显地 包含 ．设计单步算法的出发点是将微分方程两端对 在 ， 上 积分，得到
，
对上式中的定积分采用不同的数值积分方法，也就得到了不同的单步 方法．
科学计算中常碰到的一个问题是一维数值积分，即计算
，； 这个问题可以直截了当地化成求解如下常微分方程问题
0 用任何一种常微分方程数值求解方法计算出 就得到了所要的积 分值。如果在 ， 上用矩形公式计算积分
不但如此，对于求解常微分方程两点边值问题，例如
，， ，
，
，
已经建立了差分方法．其实，可以用所谓“打靶法”来求解。其基本
思想是设法将边值问题转化为初值问题来求解，即依 据边界条件
0，
0寻求与它等价的初始条件
，， ，
，
，
也就是说，反复调整初始时刻的斜率
值
，使初值问题的积分曲
线
能“命中”
．
计算过程很简单，若凭经验能够
6
其中 是待求的的 维向量函数，而 是已知的向量函数， 是独 立变量，可看成是时间， 是 矩阵．
当 1时， ，解为
受 的 Maclaurin 展开式的启发，定义 维方阵 的指数
1
exp
!
其中
，则有
1
1
!
1!
因此
非奇异矩阵 存在分解
Λ 当矩阵 对称时，有正交分解
其中，
，，
ΛT 是列为特征向量的正交矩阵，且
对于一个稳定系统， 的实部 一定是负的，称量 1/ *
为时间常数，用来表征 的衰减速度。 的振荡频率为 /2 ，
8
越大，振荡越快．一般 是衰减的或振荡衰减的，但各个 的 衰减速度之间的差异可能是很大的．
可以看到，刚性性质是数学问题本身的性质．它不依赖于求解 这个问题的数值方法．但是正是由于这个性质，使得传统的常微分方 程的数值积分方法遇到极大的困难．为了克服这个困难，刚性常微分 方程数值积分方法的研究成为数值方法中最活跃的方向之一。
就是此方法是相容的。
上面讨论的问题均没有考虑到计算时舍入误差的影响．而实际上
不可能得到 的准确值．如果考虑舍入误差，得到的是 而不是
。如果在某一点舍入误差的影响在后续计算中恶性积累以至于无
2
上一讲主要讨论椭圆型方程的有限差分法，位势（椭圆型）方程 描写的状态（如温度、电位等）不随时间 改变，称为驻定问题；现 在开始讨论与时间 有关的发展问题：扩散（抛物型）方程和波动（双 曲型）方程．驻定问题可看成是某一发展问题当 ∞的渐近状态， 所以，当用渐近方法（例如迭代法）求解验定问题时，只关心最终状 态，而不管瞬时状态或中间过程；相反发展问题的瞬时状态有物理意 义，需要在考虑偏微分方程的数值解法时，要注意到这两类问题的这 些联系和区别。
提供， 的两个预测值 ， ，分 别按这两个斜率值“试射”——求解
5
相应的初值问题，从而获得 的两个结果 ， ．如果均不满足预 定的精度，就用线性插值方法校正 ， 得新的斜率值
然后再按斜率值 试射，求解相应的初值问题，又得新的结果 ．继续这一过程，直到计算结果 与 相当符合为止。显
然这样方法法过分依赖于经验，局限性大．但是，从一个方面是说明 了常微分方程的初值问题的重要性。
§1 线上方法
扩散方程属于发展方程的一种，介乎位势方程和波动方程之间， 比较适合发展一种半离散化方法，即先把偏微分方程转化为大规模常 微分方程，然后再对常微分方程组的初值问题求解。以齐次扩散方程 的混合问题
3
0， ，0
0， 0，1 ， 0 1， 0 ， 0
， 0，1
为例。以步长 取均匀空间网格
， 0，1，2， ， 1 对二阶空间导数作中央差分近似，得到
解．就是寻求解 在一系列离散节点
上的近似值为 ， ， ， ，
． 如不特别说明，相邻两个
节点的间距
，总是假定
为定数，
，
0，1， 。
9
初值问题的数值解法都有“步进式”的基本特点，即求解过程顺 着节点排列的次序一步一步地向前推进．描述这类算法，只要给出用 已知信息从 ， ， ， 计算从 的逐推公式即可，从 出发可以依次算出 的近似值。这种计算公式称为 步差分格式．
1.1 “刚性”问题 ............................................................................. 6 1.2 从 Euler 方法出发 ....................................................................... 9 1.3 稳定性 ....................................................................................... 14 §2 完全差分格式 .................................................................................. 19 2.1 区域网格 .................................................................................... 19 2.2 连续算子离散化 ....................................................................... 21 2.3 格式的经济实用性..................................................................... 24 §3 Fourier 稳定性分析........................................................................... 28 2.1 离散方程的特解 ....................................................................... 28 3.2 解析解与离散解的比较............................................................ 32 3.3 Von Neumann 分析 .................................................................... 34 附录：发展问题，电子计算机和计算方法 ......................................... 37 习题 13 ................................................................................................... 39
发展方程是常微分方程初值问题的自然的推广，常微分方程的确 也能看成是没有空间变化的发展方程。常微分方程与发展方程在数值 处理方面有许多相似之处，事实上计算后者最有效的方法之一就是把 发展问题的偏微分方程近似转化为一个常微分方程组。但是必须注意， 这种相似不总是靠得住的，数值求解发展偏微分方程要在空间和时间 上离散，一个成功的算法中这两个过程不是独立发展的。发展偏微分 方程的那些数值分析比起常微分方程来说，更复杂也更微妙。
1.2 从 E源自文库ler 方法出发
考察的一阶方程的初值问题 ，
只要函数 ， 适当光滑，譬如关于 满足 Lipschitz 条件
，
，
在理论上就可以保证初值问题的解
存在并且唯一，其中 ·
是范数。
虽然求解常微分方程有各种各样的解析方法，但只能用来求解殊
类型的方程．实际问题中归结出来的微分方程主要靠数值解法求
途径中收效最差的。
如果数值积分改用梯形公式近似，则得到
， 2 或更进一步的加权方法
，
梯形 格式
，
1
，
格式
11
其中0 关于
1。Euler格式与其它几个格式有着本质的区别，前者是 的一个直接的计算格式，这类格式是显式的；而后者几个格
式的右端含有未知的 ，它实际上是关于 的一个函数方程，这
类格式是隐式的．显式与隐式两类方法各有特点，考虑到数值稳定性
误差，如果
，
则称这个方法为 阶精度的．而称
为整体截断误差。对于固定的
，如果当 0时，
，即 0，则称该方法是收敛的。
假设显式单步法有 阶精度，且增量函数满足 Lipschitz 条件
，，
，，
|
|
又设初值是准确的，
，那么其整体截断误差
13
进一步地，如果单步法的增量函数 满足
， ，0
，
条件，就称此方法与原始初值问题相容。显式单步法收敛的充要条件
diag ， ， ，故有
1
ΛT
T
!
矩阵 关于向量 的 2－范数 的矩阵范数为
7
sup | | max| | 也有
max 如果矩阵 是正定的，线性系统的解是渐近稳定的
尽管如此，在发展方程的半离散化得到的线性系统初值问题的系 数矩阵 经常是“病态”的。假定矩阵 的 Jordan 标准型是对角矩阵， 其特征值