《2022年第九章偏微分方程差分方法汇总》由会员分享，可在线阅读，更多相关《2022年第九章偏微分方程差分方法汇总（14页珍藏版）》请在金锄头文库上搜索。
1、170 第 9 章偏微分方程的差分方法含有偏导数的微分方程称为偏微分方程。由于变量的增多和区域的复杂性，求偏微分方程的精确解一般是不可能的，经常采用数值方法求方程的近似解。偏微分方程的数值方法种类较多，最常用的方法是差分方法 。差分方法具有格式简单，程序易于实现，计算量小等优点，特别适合于规则区域上偏微分方程的近似求解。本章将以一些典型的偏微分方程为例，介绍差分方法的基本原理和具体实现方法。9.1 椭圆型方程边值问题的差分方法9.1.1 差分方程的建立最典型的椭圆型方程是Poisson （泊松）方程Gyxyxfyuxuu),(),()(2222（9.1 ）G 是 x, y 平面上的有界区域，其
2、边界为分段光滑的闭曲线。当f(x, y) 0 时，方程（9.1 ）称为Laplace( 拉普拉斯）方程。椭圆型方程的定解条件主要有如下三种边界条件第一边值条件),(yxu（9.2 ）第二边值条件),(yxnu（9.3 ）第三边值条件),()(yxkunu（9.4 ）这里， n 表示上单位外法向，(x,y),(x,y),(x,y)和 k( x,y) 都是已知的函数，k( x,y) 0。满足方程（ 9.1 ）和上述三种边值条件之一的光滑函数u( x,y) 称为椭圆型方程边值问题的解。用差分方法求解偏微分方程，就是要求出精确解u(x,y) 在区域 G 的一些离散节点（ xi，yi）上的近似值ui,j
3、( xi，yi） 。差分方法的基本思想是，对求解区域G 做网格剖分，将偏微分方程在网格节点上离散化，导出精确解在网格节点上近似值所满足的差分方程，最终通过求解差分方程，通常为一个线性方程组，得到精确解在离散节点上的近似值。设 G=0xa, 0 y0, B( x, y) Bmin 0, E(x, y) 0。引进半节点, 12121hxxii,22121hyyii利用一阶中心差商公式，在节点（i, j）处可有)(2), 1(), 1(),()(), 1(),(),(), 1(1)(),21)(),21)(1),)(211211,211,211211hOhjiujiujixuhOhjiujiuAhj
4、iujiuAhhOjixuAjixuAhjixuAxjiji精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 3 页，共 14 页173 对yuyuBy),(类似处理，就可推得求解方程（9.9 ）的差分方程hjijijijijijijijijijiGjijifuauauauaua),(),(,1,1,1,1, 1, 1, 1,1（9.10 ）其中jijijijijijijijijijijijijijijijijijiEBBhAAhaDhBhaDhBhaChAhaChAha,21,21,22,21,2121,221,221,221,221,1,2121
5、,1,1,2121,1)()()2()2()2()2(（9.11 ）显然， 当系数函数A( x, y)=B( x, y)=1, C( x, y)=D( x, y)= E( x, y)=0 时， 椭圆型方程（9.9 ）就成为Poisson 方程（ 9.1 ） ，而差分方程（9.10 ）就成为差分方程（9.6 ） 。容易看出，差分方程（9.10 ）的截断误差为)(2221hhO阶。9.1.2 一般区域的边界条件处理前面已假设G 为矩形区域， 现在考虑G 为一般区域情形，这里主要涉及边界条件的处理。考虑 Poisson 方程第一边值问题),(),(),(),(yxyxuGyxyxfu（9.12 ）其
6、 中G可 为 平 面 上 一 般 区 域 ， 例 如 为 曲 边 区 域 。 仍 然 用 两 组 平 行 直 线 ：x=x0+ih1, y=y0+jh2, i, j=0, 1, , 对区域 G 进行矩形网格剖分，见图9-3 。如果一个内节点（i, j）的四个相邻节点（i+1, j） ， （ i-1, j） ， （i,j+1）和（ i, j-1 ）属于GG，则称其为 正则内点 ，见图 9-3 中打“。 ”号者；如果一个节点 （i, j）属于G且不为正则内点，则称其为非正则内点 ，见图9-3中打“ . ”号者。记正则内 点 集 合 为hG， 非 正 则 内 点 集 合 为h。 显 然 ， 当G 为
7、 矩 形 区 域 时 ，精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 4 页，共 14 页174 hhhhGG,成立。在正则内点（ i, j）处，完全同矩形区域情形，可建立五点差分格式hjijijijijijijiGjifuuuhuuuh),(,2121,1,1,22, 1, 121（9.13 ）在方程（ 9.13 ）中，当（i, j）点临近边界时，将出现非正则内点上的未知量，因此必须补充非正则内点处的方程。若非正则内点恰好是边界点，如图9-4 中D点，则利用边界条件可取uD=(D) 对于不是边界点的非正则内点，如图 9-4 中 B点，一般可采用
8、如下两种处理方法。a. 直接转移法 . 取与点 B 距离最近的边界点 （如图 9-4 中 E点）上的 u 的值作为u( B) 的近似值 uB，即 uB=u(E)=(E) 直接转移法的优点是简单易行，但精度较低，只为一阶近似。b. 线性插值法 . 取 B 点的两个相邻点（如图9-4 中边界点A 和正则内点C 作为插值节点对u(B) 进行线性插值)()()()(21hOCuxxxxAuxxxxBuACABACBC则得到点B 处的方程ABCBxxuhAhhu,)(111线性插值法精度较高，为二阶近似。对每一个非正则内点进行上述处理，将所得到的方程与（9.13 ）式联立，就组成了方程个数与未知量个数相
9、一致的线性代数方程组。求解此方程组就可得到一般精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 5 页，共 14 页175 区域上边值问题（9.12 ）的差分近似解。对于一般区域上二阶椭圆型方程（9.9 ）的第一边值问题，可完全类似处理。第二、三边值条件的处理较为复杂，这里不再讨论。9.2 抛物型方程的差分方法本节介绍抛物型方程的差分方法，重点讨论差分格式的构造和稳定性分析。9.2.1 一维问题作为模型，考虑一维热传导的初边值问题Ttlxtxfxuatu0,0),(22（9.14 ）lxxxu0),()0,(（9.15 ）Tttgtlutgtu0),
10、(),(),(),0(21（9.16 ）其中 a 是正常数，)()(),(),(21tgtgxtxf和都是已知的连续的函数。现在讨论求解问题（9.14 ）-(9.18)的差分方法。首先对求解区域G=0 x l, 0tT 进行网格剖分。取空间步长h=l/ N, 时间步长 =T/M, 其中 N,M 是正整数，作两族平行直线MkkttNjjhxxkj, 1 ,0, 1 , 0,将区域 G 剖分成矩形网格，见图9-5 ，网格交点（ xj，tk）称为节点。用差分方法求解初边值问题（9.14 ）-（ 9.16 ）就是要求出精确解u( x, t) 在每个节点 （xj，tk）处的近似值),(kjkjtxuu。
11、为简化记号， 简记节点 （xj，tk）=u( j, k) 。利用一元函数的Taylor展开公式，可推出下列差商表达式精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 6 页，共 14 页176 )(),() 1,(),(Okjukjukjtu（9.17 ）)() 1,(),(),(Okjukjukjtu（9.18 ）)(2) 1,()1,(),(2Okjukjukjtu（9.19 ）)(), 1(),(2), 1(),(2222hOhkjukjukjukjxu（9.20 ）1. 古典显格式在区域 G 的内节点 ( j, k) 处，利用公式（9.17
12、） 和 （9.20 ） ， 可将偏微分方程 （9.14 ）离散为)(), 1(),(2), 1(),() 1,(22hOfhkjukjukjuakjukjukj其中),(kikjtxff。舍去高阶小项)(2hO，就得到节点近似值（差分解）kju所满足的差分方程kjkjkjkjkjkjfhuuuauu21112（9.21 ）显然，在节点 ( j, k) 处，差分方程（9.21 ） 逼近偏微分方程 （9.14 ） 的误差为)(2hO，这个误差称为截断误差 ，它反映了差分方程逼近偏微分方程的精度。现将（9.21 ）式改写为便于计算的形式，并利用初边值条件（9.15 ）与（ 9.16 ）补充上初始值和
13、边界点方程，则得到MktgutguNjxuMkNjfruurruukkNkkjjkjkjkjkjkj, 1 , 0),(),(1, 2, 1),(1, 1 , 0, 1,2, 1)21(2100111（9.22 ）其中2har称为 网比 。与时间相关问题差分方程的求解通常是按时间方向逐层进行的。对于差分方程（9.22 ） ，当第k 层节点值kju已知时，可直接计算出第k+1 层节点值1kju。这样，从第0 层已知值)(0ijxu开始，就可逐层求出各时间层的节点值。差分方程（9.22 ）的求解计算是显式的，无须求解方程组，故称为古典显格式 。此外，在式（9.22 ）中，每个内节点处方程仅涉及k
14、和 k+1 两层节点值，称这样的差分格式为精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 7 页，共 14 页177 双层格式 。差分方程（ 9.22 ）可表示为矩阵形式011, 1 , 0,uMkFAuukkk（ 9.23 ）其中rrrrrrrrA212121TNTkkNkNkkkkTkNkkxxtrgffftrgfFuuu)(,),()(,),(),(11212211112. 古典隐格式在区域 G 的内节点 ( j, k) 处， 利用公式（9.18 ） 和 （9.20 ） ， 可将偏微分方程 （9.14 ）离散为)(), 1(),(2), 1(
15、) 1,(),(22hOfhkjukjukjuakjukjukj舍去高阶小项)(2hO，则得到如下差分方程kjkjkjkjkjkjfhuuuauu21112（9.24 ）它的截断误差为)(2hO，逼近精度与古典显格式相同。改写（9.24 ）式为便于计算的形式，并补充上初始值与边界点方程，则得到MktgutguNjxuMkNjfuruurrukkNkkjjkjkjkjkjkj, 1 ,0),(),(1,2, 1),(, 1 ,0, 1,2, 1)21 (2100111（9.25 ）精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 8 页，共 14 页1
16、78 与古典显格式不同，在差分方程（9.25 ）的求解中，当第k-1 层值1kju已知时，必须通过求解一个线性方程组才能求出第k 层值kju，所以称（ 9.25 ）式为 古典隐格式 ，它也是双层格式。差分方程（ 9.25 ）的矩阵形式为01,2, 1,uMkFuBukkk（9.26 ）其中rrrrrrrrB212121向量,kkFu同（ 9.23 ）式中定义。从（9.26 ）式看到，古典隐格式在每一层计算时，都需求解一个三对角形线性方程组，可采用追赶法求解。3.Crank-Nicolson格式（六点对称格式）利用一元函数Taylor展开公式可得到如下等式)()1,(),(21)21,()(),
17、()1,()21,(22222222OkjxukjxukjxuOkjukjukjtu使用这两个公式，在)21,(kj点离散偏微分方程（9.14 ） ，然后利用（9.20 ）式进一步离散二阶偏导数，则可导出差分方程2121121111112221kjkjkjkjkjkjkjkjkjfhuuuhuuuauu（9.27 ）其截断误差为)(22hO，在时间方向的逼近阶较显格式和隐格式高出一阶。这个差分格式称为Crank-Nicolson格式 ，有时也称为 六点对称格式，它显然是双层隐式格式。改写（9.27 ）式，并补充初始值和边界点方程得到精选学习资料 - - - - - - - - - 名师归纳总结
18、 - - - - - - -第 9 页，共 14 页179 MktgutguNjxuMkNjfruurruruurrukkNkkjjkjkjkjkjkjkjkj, 1 ,0),(),(1,2, 1),(1, 1 ,0, 1,2 , 12)1(2)1(22100211111111（9.28 ）它的矩阵形式为1, 1 ,0,2, 1,)()(0211MkuMkFuAIuBIkkk（9.29 ）在每层计算时，仍需求解一个三对角形方程组。4. Richardson格式利用公式（9.19 ）和（ 9.20 ） ，可导出另一个截断误差为)(22hO阶的差分方程kjkjkjkjkjkjfhuuuauu211
19、1122称之为 Richardson格式 。可改写为kjkjkjkjkjkjfuuuruu2)2(21111 (9.32) 这是一个三层显式差分格式。在逐层计算时，需用到1kju和kju两层值才能得到k+1层值1kju。这样，从第0 层已知值)(0jjxu开始，还须补充上第一层值1ju，才能逐层计算下去。可采用前述的双层格式计算1ju。除上述四种差分格式外，还可构造出许多逼近偏微分方程（9.14 ）的差分格式，但并不是每个差分格式都是可用的。一个有实用价值的差分格式应具有如下性质：（ 1）收敛性 。对任意固定的节点（xj，tk）, 当剖分步长0,h时，差分解kju应收敛到精确解u（xj，tk）
20、 。（ 2）稳定性 。当某一时间层计算产生误差时，在以后各层的计算中，这些误差的传播积累是可控制而不是无限增长的。理论上可以证明，在一定条件下，稳定的差分格式必然是收敛的。因此，这里主要研究差分格式的稳定性。精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 10 页，共 14 页180 作为例子，先考查Richardson格式的稳定性。设kju是当计算过程中带有误差时， 按 Richardson格式 （9.30 ） 得到的实际计算值。kju是理论值， 误差kjkjkjuue。假定右端项kjf的计算是精确的，网比21r，则kje满足)2(1111kjk
21、jkjkjkjeeeee（9.31 ）设前 k-1 层计算时精确的，误差只是在第k 层0j点发生，即01,0, 00jjeeekjkjkj。则利用（ 9.31 ）式可得到误差的传播情况，见表9-1。表 9-1 r=1/2时 Richardson格式的误差传播jkj0-4 j0-3 j0-2 j0-1 j0j0+1 j0+2 j0+3 j0+4 k0 0 0 0 0 0 0 0 k+1 0 0 0 -2 0 0 0 k+2 0 0 -4 7-4 0 0 k+3 0 -6 17-24 17-6 0 k+4 -831-68 89-68 31-8k+5 -10 49-144 273-388 273-1
22、44 49-10 k+6 71-260 641-1096 1311-1096 641-26071从表中看出，误差是逐层无限增长的。表中的计算虽然是就网比21r进行的，实际上对任何r0 都会产生类似现象，所以Richardson格式是不稳定的。利用误差传播图表方法考查差分格式的稳定性虽然直观明了，但只能就具体取定的 r 值进行，并且也不适用于隐式差分格式。9.2.2 差分格式的稳定性前节构造的几种双层差分格式都可以表示为如下的矩阵方程形式01uFHuukkk（9.32 ）其中 H 称为 传播矩阵 。对于显格式H=A, 隐格式H=B-1, 六点对称格式H=( I+B) -1(I+ A) 。一般的三
23、层格式也可以转化为双层格式。为了讨论方便，设在初始层产生误差0，且假定右端项kF的计算是精确的。用ku精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 11 页，共 14 页181 表示当初始层存在误差0时，由差分格式（9.32 ）得到的计算解，则ku满足方程001uFuHukkk（9.33 ）记误差向量kkkuu，则k满足方程为初始误差01, 2, 1,kHkk（ 9.34 ）定义 9.1称差分格式 （9.32 ） 是稳定的， 如果对任意初始误差0， 误差向量k在某种范数下满足000,0,kCk（9.35 ）其中 C 为与 h, 无关的常数。这个定
24、义表明，当差分格式稳定时，它的误差传播是可控制的。从（ 9.34 ）式递推得到TkHkk0,0因此，差分格式稳定的充分必要条件是TkCHk0,（9.36 ）定理 9.3 （稳定性必要条件）差分格式稳定的必要条件是存在与h,无关的常数 M,使谱半径MH1)(（9.37 ）定理 9.4 （稳定性充分条件）设 H 为正规矩阵，即HHHH*, 则（ 9.37 ）式也是差分格式稳定的充分条件。下面讨论几种差分格式的稳定性。为便于讨论，引进N-1 阶矩阵0110110110S这个特殊矩阵的特征值为精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 12 页，共 1
25、4 页182 1, 2, 1,cos2NjNjSj（ 9.38 ）例 9-1 古典显格式此时 H=A=(1-2r)I+rS。 利用（9.38 ）式和三角函数公式，可求得 H 的特征值为1,2, 1),2(sin412NjNjrj为使稳定性条件式（9.39 ）成立，必须且只须21r。由于 H =A 为实对称矩阵，所以古典显格式稳定的充分必要条件是网比21r例 9-2 古典隐格式此时 H=B-1,B=(1+2r)I-rS 。利用（ 9.38 ）式可求得H 的特征值为1,2, 1,)2(sin4112NjNjrj显然，对任意r0 ，条件（ 9.37 ）成立。注意，H=B-1仍为实对称矩阵，所以古典隐
26、格式对任何网比r0 都是稳定的，称为绝对稳定 。例 9-3 六点对称格式此时 H=(I+B)-1(I+A)，利用矩阵A 和 B 的特征值可得到矩阵 H 的特征值为1,2, 1,)2(sin42)2(sin4222NjNjrNjrj则对任意r0, 条件（ 9.37 ）成立。由于A 和 B 均为实对称矩阵，且AB=BA，则可验证 H 也是实对称矩阵。所以六点对称格式是绝对稳定的。习题九9-1 试用五点差分格式求解Poisson 方程的边值问题),(,0|),(,162222yxuGyxyuxu其中 G=-1 x, y1。取步长h=0.5 求解。9-2 试写出求解Laplace 方程边值问题精选学习
27、资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 13 页，共 14 页183 40, 0)3,(,4sin)0 ,(30,0), 4(),3(), 0(30,40,162222xxuxxuyyuyyyuyxyuxu的五点差分格式，取步长h=1, 并写出差分方程的矩阵形式。9-3 试用古典显格式求解热传导方程定解问题TttutuxxxxuTtxxutu0, 0),1 (),0(10),1(4)0 ,(0, 10,22只计算 k=1,2 两层上的差分解，取网比r=1/6, h=0.2 。9-4 用古典隐格式求解热传导方程定解问题3.00,0), 1(),0(10,sin)0,(3 .00, 10,22ttutuxxxutxxutu取2.0, 1.0h精确解为xetxutsin),(2。9-5 导出求解热传导方程的差分格式)2(1)1(11211kjkjkjkjkjkjkjuuuhuuuu的截断误差，并选取使其达到二阶。精选学习资料 - - - - - - - - - 名师归纳总结 - - - - - - -第 14 页，共 14 页
