椭圆型方程的差分解法

1椭圆型方程的差分解法1.引言考虑问题①二维Poisson方程2222(,)uufxyxy，(,)xy其中为2R中的一个有界区域，其边界为分段光滑曲线。在上u满足下列边界条件之一：⑴(,)uxy（第一边值条件），⑵(,)uxyn（第二边值条件），⑶(,)ukuxyn（第三边值条件），(,),(,),(,),(,),(,)fxyxyxyxykxy都是连续函数，0k.2.差分格式将区间[,]ab作m等分，记为11()/,,0ihbamxaihim；将区间[,]cd作n等分，记为22()/,,0ihdcnycjhjn.称1h为x方向的步长，2h为y方向的步长。2.1Poisson方程五点差分格式参考单如图所示：2j+1j3以(,)ijxy为中心沿y方向Taylor展开：2344111111(,)(,)(,)(,)(,)(,)(),2624ijijxijxxijxxxijxxxxijhhhuxyuxyhuxyuxyuxyuxyoh①2344111111(,)(,)(,)(,)(,)(,)(),2624ijijxijxxijxxxijxxxxijhhhuxyuxyhuxyuxyuxyuxyoh②由①+②得：42411111(,)(,)2(,)(,)(,)(),12ijijijxxijxxxxijhuxyuxyuxyhuxyuxyoh整理得到：21121121(,)2(,)(,)(,)(,)(),12ijijijxxijxxxxijuxyuxyuxyhuxyuxyohh③同理，以(,)ijxy为中心沿y方向Taylor展开：21122222(,)2(,)(,)(,)(,)(),12ijijijyyijyyyyijuxyuxyuxyhuxyuxyohh④代入原方程得○5：11112212(,)2(,)(,)(,)2(,)(,)(,),ijijijijijijijijuxyuxyuxyuxyuxyuxyfxyRhh得到截断误差：44222212441(,)(,)()(),12ijijhijijuuRuxyuxyhhohohxy其中u是原方程光滑解，舍去截断误差得到逼近Poisson方程的五点差分方程：11112212(,)2(,)(,)(,)2(,)(,)(,),ijijijijijijijuxyuxyuxyuxyuxyuxyfxyhh○6考虑到边值条件(,)(,)uxyxy，构成差分格式：11112212(,)2(,)(,)(,)2(,)(,)(,),(,)(,),ijijijijijijijuxyuxyuxyuxyuxyuxyfxyhhuxyxy○742.2Poisson方程九点差分格式由上式③+④得：11112212442221244222222122222(,)2(,)(,)(,)2(,)(,)(,)1(,)()12(,)(,)1(,)12ijijijijijijhijijijijijijuxyuxyuxyuxyuxyuxyuxyhhuuuxyhhohxyuxyuxyuxyhhxyxy422212222242222212122222(,)()12(,)(,)(,)1(,)()1212ijijijijijuxyhhohxyfxyfxyuxyhhfxyhhohxyxy○8又41122222211111112212311111(,)(,)2(,)(,)()1[(,)2(,)(,)2(,)2(,)(,)(,)2(,)(,)]()ijxxijxxijxxijijijijijijijijijijuxyuxyuxyuxyohxyhuxyuxyuxyuxyuxyuxyhhuxyuxyuxyoh则得到：222222121121112112222221211212122222221112111211()(,)(210)(,)()(,)(210)(,)20()(,)(210)(,)(210)(,)()(,)()(,)ijijijijijijijijijhhuxyhhuxyhhuxyhhuxyhhuxyhhuxyhhuxyhhuxyhhuxy2212222241222,12(,)(,)1(,)()12ijijijhhfxyfxyfxyhhohxy○9舍去截断误差得到逼近Poisson方程的九点差分方程○10：2212,11,,11,1,11,11,11,122122212(,)[42]121(,)(,),12ijijijijijijijijijijijxxijyyijhhuxyuuuuuuuuuhhfhfxyhfxy考虑到边值条件(,)(,)uxyxy，构成差分格式○11：52212,11,,11,1,11,11,11,122122212(,)[42]121(,)(,),12(,)(,),ijijijijijijijijijijijxxijyyijhhuxyuuuuuuuuuhhfhfxyhfxyuxyxy3.格式求解3.1Poisson方程五点差分格式记122,1,jjjmjmjuuuuu，0.jn矩阵格式改写为：11,11jjjjDuCuDufjm，其中2221212222112122221121222112(1)111211112111121112mhhhhhhhChhhhhhh，622222222(1)1111mhhDhh，10212212111(,)(,)(,)(,)1(,)(,)jjjjmjmjmjmfxyxyhfxyffxyfxyxyh，可进一步写为：110222211(1)*(1).nnnnnnmufDuCDufDCDufDCDufDuDC3.2Poisson方程九点差分格式记122,1,jjjmjmjuuuuu，0.jn矩阵格式改写为：11,11jjjjDuCuDufjm，其中2222121222222212121222222212121222221212(1)20()(210)(210)20()(210)(210)20()(210)(210)20()mhhhhhhhhhhChhhhhhhhhh，2222211222222212211222222212211222221221(1)(210)()()(210)()()(210)()()(210)mhhhhhhhhhhDhhhhhhhhhh，722121022221211(,)(210)(,)(,)(,)(,)(210)(,)jjjjmjmjmjmfxyhhxyfxyffxyfxyhhxy，可进一步写为：110222211(1)*(1).nnnnnnmufDuCDufDCDufDCDufDuDC4.数值例子4.1Poisson方程五点差分格式计算如下问题：22220,01,01,(0,)sincos,(2,)(sincos),01,(,0),(,1)(sin1cos1),01.xxuuxyxyuyyyuyeyyyuxeuxex其精确解为：(,)(sincos).xuxyeyy,11,1,,1,222222122112112()(,),ijijijijijijuuuuufxyhhhhhh考虑到本例中h1=h2，则有2,11,1,,1,(,),4ijijijijijijuuuuhfxyu利用Gauss-Seidel迭代方法对k=0,1,2,……,计算112,11,1,,11(,),41,2,....,1;1,2,....,1.kkkkijijijijijkijuuuuhfxyuimjn8表1部分结点处的精确解和取不同步长时所得的数值解(h1,h2)(x,y)(1/4,1/4)(1/2,1/4)(3/4,1/4)(1/4,3/4)(1/2,3/4)(3/4,3/4)(1/4,1/4)1.5626222.0065892.5760251.8157532.3314922.993189(1/8,1/8)1.5620992.0056822.5751831.8150312.3305222.992321(1/16,1/16)1.5618362.0054462.5750131.8148082.3302692.992092(1/32,1/32)1.5617942.0053382.5794591.8146762.3302202.992003(1/64,1/64)1.5616742.0053582.5749371.8147432.3302762.992011精确解1.5617812.0054372.5746941.8147492.3301832.992014表2取不同步长时部分结点处数值解的误差绝对值(h1,h2)(x,y)(1/4,1/4)(1/2,1/4)(3/4,1/4)(1/4,3/4)(1/2,3/4)(3/4,3/4)(1/4,1/4)8.4094e-040.00120.00119.0500e-040.00130.0012(1/8,1/8)2.1769e-043.1498e-042.8202e-042.3451e-043.3891e-043.0635e-04(1/16,1/16)5.4874e-057.9388e-057.1220e-055.9144e-058.5448e-057.7433e-05(1/32,1/32)1.3096e-051.9001e-051.7247e-051.4215e-052.0586e-051.8853e-05(1/64,1/64)6.1867e-068.3780e-064.7959e-065.5548e-067.4932e-064.0549e-069表3取不同步长时部分结点处数值解的最大误差(h1,h2)()Eh(2)/()EhEh(1/4,1/4)0.0016*(1/8,1/8)4.2732e-043.744267(1/16,1/16)1.0885e-043.925769(1/32,1/32)2.6148e-054.162842(1/64,1/64)1.1928e-052.192153图1取h=1/4时所得的数值解曲线10图2取h=1/4时所得的误差曲