偏微分课程课件10-椭圆型方程的有限差分方法(II).ppt

1、1求解代数方程组求解代数方程组Hu=g的方法的方法直接方法直接方法:高斯消去法高斯消去法,三角分解三角分解,追赶法追赶法,QR分解等分解等迭代方法迭代方法:基本迭代法基本迭代法,预处理迭代法预处理迭代法,多重网格法等多重网格法等Jacobi迭代法迭代法Gauss-Seidel迭代法迭代法超松弛迭代法超松弛迭代法其它迭代法其它迭代法Jacobi迭代法求解代数方程迭代法求解代数方程Au=fAD L U :DdiagonalLlower triangularUupper triangular 对对角角下下三三角角 上上三三角角 将将A分解分解41104000140104001041004001140

2、004 000001101000000110000001011 00000 3Jacobi迭代法求解代数方程迭代法求解代数方程Au=fAD L U A为为n阶矩阵阶矩阵,分解为分解为12(,)Tnuu uu 11 122 2:iiinn na ua uDua ua u ,i ja:DdiagonalLlower triangularUupper triangular 对对角角下下三三角角 上上三三角角 4Jacobi迭代法求解代数方程迭代法求解代数方程Au=fAD L U A为为n阶矩阵阶矩阵,分解为分解为12(,)Tnuu uu ,i ja 21 1131 132 2,11 12 2,110

3、:ii jjjnnn nna uLua ua ua ua ua uau 5Jacobi迭代法求解代数方程迭代法求解代数方程Au=fAD L U A为为n阶矩阵阶矩阵,分解为分解为12(,)Tnuu uu ,i ja 12 213 31,23 32,11,:0n nn nni jjj inn na ua ua ua ua uUua uau 6Jacobi迭代法求解代数方程迭代法求解代数方程Au=fAD L U :DdiagonalLlower triangularUupper triangular 对对角角下下三三角角 上上三三角角 AuD L U uf 将将A分解分解 DuL U uf 11u

4、DL U uD f 11DDA uD f 11ID A uD f 1JR uD f Jacobi迭代法迭代法(1)()1kkJuR uD f 1JRID A ,),()0()0(1)0(Tnxxx,/1)()1(iinijjkjijikiaxabx(1,2,)(0 1.ink,分量形式分量形式8Guass-Seidel迭代法求解代数方程迭代法求解代数方程Au=fAD L U :DdiagonalLlower triangularUupper triangular 对对角角下下三三角角 上上三三角角 AuD L U uf 将将A分解分解 DL uUuf 11uDLUuDLf 1()GSR uDL

5、f G-S迭代法迭代法(1)()1()kkGSuR uDLf 1()GSRDLU A为为T矩阵矩阵时矩阵矩阵时,G-S迭代法收敛速度是迭代法收敛速度是Jacobi法两倍法两倍G-S迭代法矩阵形式迭代法矩阵形式分量迭代法分量迭代法计算公式为计算公式为bAx(初始向量),Tnxxx),()0()0(1)0(1(1)(1)()11/,inkkkiiijjijjiijjixba xa xa 10超松弛迭代法超松弛迭代法求解代数方程求解代数方程Au=fAD L U :DdiagonalLlower triangularUupper triangular 对对角角下下三三角角 上上三三角角 A大型稀疏矩阵

6、大型稀疏矩阵将将A分解分解1(1)(1)()11/,inkkkiiijjijjiijj ixba xa xa Guass-Seidel迭代法解迭代法解(1)()(1)(1)kkkiiixxx1(1)()(1)()11(1)inkkkkiiiiiiiijjijjjj ia xa xba xa x (1)()(1)()(1)kkkkDxDxbLxUx超松弛迭代法超松弛迭代法求解代数方程求解代数方程Au=f (1)-1()-1()(1)()kkuDLDU uDLf 1 超超松松弛弛1 G-S 即即迭迭代代法法 (1)()()(1)kkDL xDU xb (1)1()1()(1)()kkxDLDU x

7、DLb(1)()kkxL xf (1)()(1)(1)kkkiiixxx12逐次超松弛迭代法逐次超松弛迭代法(SOD)Successive Over Relaxation Method(0)(0)(0)1(,),Tnxxx 1(1)()(1)()1/,inkkkkiiiijjijjiijj ixxba xa xa (1)-1()-1()(1)()kkuDLDU uDLf SOD迭代的矩阵形式迭代的矩阵形式SOD迭代的分量形式迭代的分量形式13例：解线性方程组例：解线性方程组解解:SOR 迭代格式迭代格式 (1)()()()()()111234(1)()(1)()()()221234(1)()(

8、1)(1)()()331234(1)()(1)(1)(1)()4412341 44144144144kkkkkkkkkkkkkkkkkkkkkkkkuuuuuuuuuuuuuuuuuuuuuuuu 122341101140111041101141uuuu 1(1)()(1)(),1()inkkkkiiii jji jjjj iiiuufa ua ua 14对于矩形区域对于矩形区域possion方程第一边值问题差分格式方程第一边值问题差分格式11214 1,kIkkkkIkuuuuufhkIJ (1)()2(1)(1)()()()1144kkkkkkkkkkk Ikkkk Iuuh fuuuuu

9、 (0)kuK(1)(1)()()(),11,1,1(1)()2,44kkkkki jiji jiji jkki ji ji juuuuuuuh f （二）五点差分格式的性质（二）五点差分格式的性质极值原理极值原理1.存在唯一性存在唯一性5.1,0(0),maxmax(minmin)hhhhijhhijhijijijijijDDDDuDuuuuuu 定定理理（极极值值定定理理）设设是是上上的的函函数数 若若则则 bua在在 或或 取取到到最最大大值值0,uxa b 0,(,)ux yD uD在在 的的边边界界取取到到最最大大值值的解存在唯一。的解存在唯一。差分方程边值问题差分方程边值问题定理定

10、理 hjiijijhjiijijhDyxuDyxfu),(,),(,2.2 0,(,)Pr.0,(,)hijijhijijhuxyDoofuxyD 由极值原理只需证明只需证明齐次方程齐次方程只有零解只有零解,(,),(,)h ijijijhijijijhvfx yDvx yD uuvmaxmax00minmin0hhhhijijDDijijijDDuuuuu 2.差分方程解的收敛性差分方程解的收敛性定理定理:设设 是定义在是定义在 上的函数，那么有上的函数，那么有 其中其中a为矩形区域为矩形区域D的的x方向的边长。方向的边长。ijuhD|max2|max|max2ijhDijDijDuauuh

11、hh 收敛性：收敛性：h0,k0时，差分方程的解逼近于时，差分方程的解逼近于 微分方程的解。微分方程的解。21(,),max|2hijijijhijDx yxuuAAu 令其中证明：证明：定义定义则则由定义由定义2012ijhija且22maxmaxmax|.22hhhijijijijDDDaauuuAuA=+0hijhijuuA 2maxmax=maxmax2hhhhijijijijijDDDDauuuAuA ijijuu 定理：如果第一边值问题定理：如果第一边值问题2222(,),(,)(,)(,)uuf x yx yDxyu x yx yD 二阶收敛二阶收敛DD 的解在的解在 上有四阶连

12、续的偏导数，上有四阶连续的偏导数，则五点差分格式收敛并有估计则五点差分格式收敛并有估计22max|(,)|().hijijDuu x yK hk证明：证明：设设u(x,y)是微分方程之解，是微分方程之解，(,)(,),(,)(,)0,(,)hijijijijhijijijuu xyE u xyxyDuu xyxyD iju是差分方程之解是差分方程之解(,)(,)(,),ijhijijE u x yu x yu x y 2max|(,)|max(,),2hhijijijDDauu xyE u xy(,)(,)(,),hijijiju xyf xyE u xy(,),hijijijuff xy|m

13、ax2|max|max2ijhDijDijDuauuhhh 222334434(,)(,)(,)(,)+2!(,)+(,)3!4!ijijijijijiju xyu xyhu xh yu xyhxxu xyhhu xh yxx 22max|(,)|().hijijDuu x yK hk44221244(,)(,)(,)1(,)(,)12ijhijijijijE u xyu xyu xyhu xh yku xykxy 2max|(,)|max(,),2hhijijijDDauu xyE u xy 2222221 sin(),02,01.(0,)sin(),(2,)sin(),01.(,0)0,(

14、,1)0,02.uuyxyxyuyyuyeyyu xu xx 应用五点差分格式计算如下问题：应用五点差分格式计算如下问题：).sin(),(yeyxux10)()1(1021|lluu精确解为精确解为可以观察到可以观察到采用采用Guass-Seidel迭代精确至迭代精确至当当x与与y方向步长减少到原来的方向步长减少到原来的1/2，误差减少到原来的，误差减少到原来的1/4，x方向与方向与y方向收敛阶均为方向收敛阶均为2 阶。阶。26复习复习:法线方向向量法线方向向量.1 方向向量方向向量(cos,cos),n (cos,cos,cos)n (,),ax y z 向量向量(,)|xyznaaa 那

15、么那么a的方向向量为的方向向量为:2 平面曲线平面曲线F(x,y)=0在点在点(x,y)处的法向量处的法向量:2222(,)yxxyxyFFFFFF 空间曲线空间曲线F(x,y,z)=0在点在点(x,y,z)处的法向量处的法向量:222222222(,)yzxxyzxyzxyzFFFFFFFFFFFF （三）边界条件的处理三）边界条件的处理(,)(,),(,),(,).ijijijijhu x yx yx yDux yD 1.矩形区域矩形区域(2)第三类边界条件第三类边界条件(1)第一类边界条件第一类边界条件第二类边界条件第二类边界条件(,),(,)0uux yx yDn hhhhDDDD在在

16、之之外外，即即在在2,1,1,1,1,1,0,0,0,2,1,1,1,1,1,0,0,0,01,2,2,01,22IjIjIjIjIjjjjjji JiJi Ji Ji JiiiiiuuujJhuuuhuuuiIkuuuk 四周增加一排节点四周增加一排节点2,1,1,2,Ijjii Juuuu 可用内点的差分格式在边界可用内点的差分格式在边界上成立得到的有关等式与边界离散相应的式子来消去上成立得到的有关等式与边界离散相应的式子来消去291coscosuuunxy 1 uy 1(0,1)n 即下边界时条件为即下边界时条件为1(,),(0)ux yy 例：例：,2 其中其中301 0(,)ixy在

17、下边界任取一边界点在下边界任取一边界点 1(,),ux yy (,0),ix0,1,iM(,0)(,0),iiuxxy 那么那么 1(,)ixy1(,)ixy 用中心差商代替一阶偏导数用中心差商代替一阶偏导数(,)(,)(,0),2iiiu xku xkxk 用用 ,iju代替代替iju用用 代替代替,1,1,02iiiuuk 0,1,iM 离散方法：离散方法：31在边界点在边界点 离散方程离散方程,例如五点格式例如五点格式,有有(,0),ix,1,1,02iiiuuk 0,1,iM 1,0,01,0,1,0,1,02222iiiiiiiuuuuuufhk 两式联立两式联立,则可以得到一个与则

18、可以得到一个与,1iu 无关的方程无关的方程1,0,01,0,1,02,0,0222222iiiiiiiuuuuuhfhk 1,iI(0,0),(1,0)I 处理稍有不同处理稍有不同.32,01,0,01,0,1,022142iiiiiiuuuufhh hk 特别特别1,0,01,0,1,02,0,022122222iiiiiiiuuuuuhfhh 0,01,00,11,10,11,1,hIIJIJuuuuuuu 这样这样,每一个边界点对应一个差分方程每一个边界点对应一个差分方程,将所有边界点和内点按照自然顺序排列将所有边界点和内点按照自然顺序排列,定义定义两边除以两边除以2(角点除以角点除以

19、4),则有则有,01,0,01,0,1,021 1112222iiiiiiuuuufhh 这些差分方程组可以写为一个代数方程组这些差分方程组可以写为一个代数方程组21hHugh 3321hHugh 分析一下系数矩阵分析一下系数矩阵BIIBIHIBIIB 21hHugh 对于第一类边界条件对于第一类边界条件,有有,1,2,iihi Juuuu I+2行行 H上边界上边界huI+2行行下边界下边界,1iu BBBBIIIIIIBI为为I+2维方阵维方阵,i Ju I+2I+2IIIII+2列列I+2列列,2iu 34例例:22221010316,(,)0,1 0,10,0.xxyyuux yDxy

20、uuuuuyxy 解解:布矩形网格布矩形网格(1)内点内点(2-5)四边界四边界(6)两个角点两个角点(7)分析其中三个块矩阵分析其中三个块矩阵35(1)内点内点 1,1,1,12123(2)16iji jiji ji ji juuuuuuh 1,2,3,1,2,3ij ,11,1,13831i jiji jiji juuuuu 将将h=1/4代入整理代入整理,有有36(2)左边界左边界00.xux 1,1,1,1,02jjjjuuuuh 由内点差分格式由内点差分格式增设虚点增设虚点,利用中心差商利用中心差商,得得(-1,1)(0,1)即即0,11,0,1,0,13831jjjjjuuuuu

21、(-1,2)(0,2)(-1,3)(0,3)(-1,0)(0,0)(-1,4)(0,4)0,10,1,0,138231jjjjuuuu 0,10,1,0,13314222jjjjuuuu 两边同除以两边同除以21,2,3i 1,0,01,0,111143222iiiiuuuu (3-4)上下边界同上下边界同(2)有有(5)右边界右边界10,xu 4,0,ju 0,1,2,3,4j ,31,4,41,41191132422iiiiuuuu 1,2,3i 10100,0 xxyyuuuuuyxy ,11,1,13831i jiji jiji juuuuu 内点：38(6)两个角点两个角点0,10,

22、01,00,13314222uuuu 对于对于(0,0)点点,按左边界离散得按左边界离散得按下边界离散得按下边界离散得1,00,01,00,111143222uuuu 按方程离散得按方程离散得0,11,00,01,00,13831uuuuu 三式联立三式联立,消去虚设点消去虚设点,得得0,01,00,18261uuu 两边同除以两边同除以4,得得0,01,00,11312224uuu 10100,0.xxyyuuuuuyxy 39对于对于(0,4)点点,同理可得同理可得0,30,41,4319112824uuu 令令0,03,04,00,13,14,10,43,44,4,huuuuuuuuuu

23、 那么那么20个方程按照自然顺序排列个方程按照自然顺序排列,则形成代数方程则形成代数方程:hHug 其中其中H 0 01 11 11 12 2E EK KK KE EK KK KE EK KK KE EK KK KE E01110gggggg E0K40(7)分析分析E0和和K(0,0):0,01,00,11312224uuu 1,0,01,0,111143222iiiiuuuu 1,2,3i (i,0):0,01,02,03,00,11,12,13,1,uuuuuuuu 132000002211400300221104003022100400032 E141分析分析E1和和K(0,1):1,

24、2,3i (i,1):0,01,02,03,00,11,12,13,10,21,22,23,2,uuuuuuuuuuuu 33000 4100000220300181003000030 018100300003 00180003 0,00,11,10,23314222uuuu ,01,1,11,1,23831iiiiiuuuuu 4201221142211422142E 14118118118E 2191821191242119124211924E 32333K 01111,4222Tg ，111 1,12Tg ，430,03,04,00,13,14,10,43,44,4,huuuuuuuuu

25、u 那么那么20个方程按照自然顺序排列个方程按照自然顺序排列,则形成代数方程则形成代数方程:hHug 其中其中H 0 01 11 11 12 2E EK KK KE EK KK KE EK KK KE EK KK KE E01110gggggg 2.一般区域一般区域1-1-11,0,1,0,1,(,)(,)(,)(,)|(,),(,),(,),(,),(,)(,)|(,)ijijhijijijijijijhijijhhxih iyjk jih jkxyi jDxyxyDxyxyxyxyDDxyxyDDDD：或或令令并并且且不不一一定定全全部部落落在在上上网格点网格点内部网格点内部网格点集合集合

26、边界节点集合边界节点集合1）直接转移法）直接转移法D(,)(,),(,)u x yx yx yD S在在D 上上(S)Su P不在不在,选与选与P 最靠近的网格线交点最靠近的网格线交点T()PuT D 上上第一边界第一边界2）线性插值）线性插值,TPk T,Q两点做线性插值两点做线性插值()()Pkuu Tu Qkk 11111()kkkkkkkkx xx xL xyyxxxx (,),(,)uux yx yDn uunx uuny（1）边界点）边界点P在在hD 上上外法线与坐标轴平行外法线与坐标轴平行()()()()()u Pu QP u PPh 第三边界第三边界外法线与坐标轴不平行外法线与

27、坐标轴不平行()cos(,)cos(,)u Pnuun xn yxy 12()()()()coscos()()()u Pu Qu Pu RP u PPhh (,),(,)uux yx yDn 12(coscos)uuxy （2）边界点）边界点 P不在不在hD 上上322311,1,12STh TRhPTh 1()()()PTuuQ u QQh 321231()()()()+RSPuuuQ u QQh（四）变系数方程（四）变系数方程(,)(,)(,)(,)(,),(,)0,(,)0,(,(,)|0 ,0 ,=(,)|=0,0)0(,)(,),(;,)=0,0uua x yb x yc x y u

28、f x yxxyyx yD a x yb x yc x yu x yxDx yx ay bDx yybxaxayx yDby 2211()+(),(,),(,)xijxijyijyijijijijijhijijijhaubuc ufx yDhkux yD 1.直接差分法直接差分法2.有限体积法有限体积法适合处理系数有间断，步长不等距问题，适合处理系数有间断，步长不等距问题，具有保持能量守恒等优点。具有保持能量守恒等优点。(,)ijP xy为内点为内点11213141(,),(,),(,),(,)ijijijijQ xyQxyQxyQxy为为P四邻点四邻点1234,NNNN,1,2,3,4iPQ

29、 i 为为中点中点(,)i j1i 1i 1j 1j 2N1N3N4N2Q1Q3Q4Q(,):,22 22ijiijjhhDx yxxxkkyyy (,)i j2N1N3N在阴影区域上积分在阴影区域上积分121111122222(,)(,)(,)(,)(,)ijjiiiijDyyua x ydxdyxxuua xyxya xyxy dyxx 由中矩形公式由中矩形公式4N11112222(,)(,)(,)(,)iiiijjjjuua xyxya xyxykxx (,)ijDub x ydxdyyy 11,11,2,221,11,12,221()()1()()ijijijijijiji jijij

30、i jijijiji ji jauuauuhbuubuuc ufk11112222(,)(,)(,)(,)iiiijjjjuub x yx yb x yxyhyy (,)(,)(,)(,)(,)ijijijijDijDc x y udxdyc xy u xy hkf x y dxdyf xy hk 同理同理（五）双调和方程（五）双调和方程44424224220,(,),(,),(,)(,),(,)(,)(,)|01,01uuuux yDxxyyuf x yx yDuug x yg x yx yDnnDx yxy 或或者者42244441,1210 xijxyijyijhkNuuuhhh +1,

31、1,+1,1+1,11,11,+11,1+2,2,+2,2208(+)2(+)+0ijijiji ji jijijijijijiji ji juuuuuuuuuuuuu1,22,20,21,31,12,32,10,30,113,21,21,41,00,10,10,20,20,30,31,01,01,21,20,2208(+)2(+)+,2uuuuuuuuuPuuuuufufufufuugh 1,22,21,31,12,32,13,21,40,20,30,11,00,2218()2()82()2uuuuuuuuffffhg（六）（六）特征值问题特征值问题(,)0(,)(,)|01,01uux y

32、Dux yDDx yxy 222sinsin,(,1,2,)=(+),(,1,2,)1pqpqupqp qpqp qhkN 21,+1,1,+1+(4)=0,(,),=0,1,1,=0,(,),ijiji ji jijijhijijhuuuuhuxyDi jNuxyD (,)222=sinsin,4=sin+sin22,=1,2,1p qijp iqjuNNpqhNNp qN 可用变量分离法可用变量分离法2,()p qp qO h五点差分五点差分 格式格式作业作业复习重点复习重点基本：基本：1、所有作业题。、所有作业题。2、三类方程求解的常用差分格式，其相容性以及适用范围。、三类方程求解的常用

33、差分格式，其相容性以及适用范围。3、分析经典差分格式稳定性（包括三层格式）；、分析经典差分格式稳定性（包括三层格式）；熟记所讲过问题差分格式的稳定条件，利用熟记所讲过问题差分格式的稳定条件，利用Lax等价定理等价定理 证明收敛性。证明收敛性。4、变分原理，有限元一维以及二维线性元的计算。、变分原理，有限元一维以及二维线性元的计算。在此基础上：在此基础上：5、掌握能量不等式判断稳定性的方法、掌握能量不等式判断稳定性的方法。6、构造有限差分格式的、构造有限差分格式的Taylor展开法，有限体积法。展开法，有限体积法。7、理解双曲型方程收敛的必要条件、理解双曲型方程收敛的必要条件C.F.L.条件。条件。8、了解如何补充初始条件以及用迎风格式处理边界条件。、了解如何补充初始条件以及用迎风格式处理边界条件。9、利用极值原理证明椭圆问题收敛性。、利用极值原理证明椭圆问题收敛性。