Poisson方程九点差分格式_米瑞琪
Laplace九点差分格式

1.5609 1.8067 2.1132 2.5190 3.0488
1.6172 1.8649 2.1590 2.5661 3.1403
1.6393 1.8897 2.0944 2.4093 3.1553
精确解:
1.3610 1.5029 1.6031 1.6589 1.6688
1.6078 1.7754 1.8939 1.9598 1.9714
共 10 页 第3页
Laplace 方程九点差分格式
的光滑解.
3.2 矩阵形式 Au b; 其中
A1 A2
A2 A1 A2
A
,
ቤተ መጻሕፍቲ ባይዱ
A2
A1
A
2
A2 A1 ( N 1)*(N 1)
20 4
4
20 4
A1
4 20
4
,
4
20( N 1)*(N 1)
4 1
1 4 1
A2
1
Laplace 方程九点差分格式
开始
k=2
输入区间,节点数N 给出迭代次数n
对向量b赋值, 对向量u赋初值u0; 对矩阵A1,A2赋值
m=1
m<n+1 N
输出数值解u1
k<N-1
Y
N
A1*u1(k)=A2*u1(k-1)A2*u0(k+1);
K=k+1
A1*u1(1)=A2*u2(2)+b(1);
中南林业科技大学
本科课程设计说明书
学 院:
理学院
专业年级: 2008 级信息与计算科学二班
课 程:
科学计算课程设计
论文题目: Laplace 方程九点差分格式
二维泊松方程的差分格式有限差分法

有限差分法(Finite Differential Method)是基于差分原理的一种
数值计算法。其基本思想:将场域离散为许多小网格,应用差分原理,将
求解连续函数的泊松方程的问题转换为求解网格节点上 的差分方程组的
问题。
1. 二维泊松方程的差分格式
二维静电场边值问题:
2
x 2
2
y 2
F
(1)
f (s)
(2)
L
通常将场域分成足够小的正方形网格, 网格线之间的距离为h,节点0,1,2,3,4上
的电位分别用0 ,1,和2 ,表3 示。4
设函数 在x0处可微 , 则沿x方向在 x0处的泰勒公式展开为
x
n (K )
Kn )
0
1 4
(1
2
3
4)
若场域离散为矩形网格, 差分格式为:
1•
2
1 h12
(1
2)
1 h2 2
( 2
4
)
(
1 h12
1 h2 2
)20
F
2.边界条件的离散化处理 ⑴第一类边界条件 给边界离散节点直接赋已知电位值。
⑵对称边界条件 合理减小计算场域, 差分格式为
•
0
1 4
(21
2
4
h2F)
⑶第二类边界条件 边界线与网格线相重合的差分格式:
(3)
将 x 和x1 分x别3 代入式(3),得
1
0
h(
x
)0
1 2!
h
2
(
2
x 2
)0
1 3!
h
3
(
3
x3
第2章 地球物理中常用数值解法基本原理-有限差分法

等代入
U f x, y ,
2
x, y
Ui , j 1 2Uij Ui , j 1 h
2 2
hUij
Ui 1, j 2Uij Ui 1, j h
2 1
fij
2 截断误差为 O h
第二节 椭圆型偏微分方程的有限差分解法
2.1 差分格式 ——九点差分格式
2 1 U xi , y j 2 U xi 1 , y j U xi , y j h1 h1 2 x 2! x 3 4 5 U x , y U x , y 1 1 1 U xi , y j 5 i j i j 3 4 6 h h h O h 1 1 1 1 3 4 5 3! x 4! x 5! x 2 1 U xi , y j 2 U xi 1 , y j U xi , y j h1 h1 2 x 2! x 3 4 5 1 U xi , y j 3 1 U xi , y j 4 1 U xi , y j 5 6 h h h O h 1 1 1 1 3! x3 4! x 4 5! x 5
如果两个节点满足
i1 i2 j1 j2 1 ,称其为相邻节点。
非正则内点
正则内点——邻点都在区域内;
第二节 椭圆型偏微分方程的有限差分解法
2.1 差分格式 ——九点差分格式 对正则内点,
U xx , U yy
U xi , y j x U xi , y j x
2
特征方程
2 dy a12 a12 a11a22 dx a11
2 a12 a11a22 0
泊松方程隐格式和显格式

泊松方程隐格式和显格式
泊松方程是一种偏微分方程,通常用于描述物理现象,例如电荷分布、热传导等。
在数值分析中,我们常常需要将泊松方程离散化,以便在计算机上求解。
隐格式和显格式是两种常见的离散化方法。
隐格式离散化方法:
隐格式方法将偏微分方程转化为一个非线性方程组,然后使用迭代法(如牛顿法)求解。
具体来说,对于泊松方程:
Δφ = f
隐格式离散化方法将该方程转化为:
(Δ_i,j φ_i,j) = f_i,j
其中Δ_i,j 是离散化的拉普拉斯算子,φ_i,j 是未知的电势,f_i,j 是已知的源项。
通过迭代求解该非线性方程组,可以得到电势的近似解。
显格式离散化方法:
显格式方法将偏微分方程转化为一个线性方程组,然后直接求解。
具体来说,对于泊松方程:
Δφ = f
显格式离散化方法将该方程转化为:
(Δ_i,j - Δ_i+1,j) φ_i,j = f_i,j - f_i+1,j
其中Δ_i,j 是离散化的拉普拉斯算子,φ_i,j 是未知的电势,f_i,j 是已知的源项。
通过求解该线性方程组,可以得到电势的近似解。
总结:
隐格式离散化方法适用于非线性问题,而显格式离散化方法适用于线性问题。
在数值分析中,根据问题的性质选择合适的离散化方法非常重要。
椭圆型方程的差分解法

椭圆型方程的差分解法1.引言考虑问题①二维Poisson 方程2222(,)u u f x y x y ⎛⎫∂∂-+= ⎪∂∂⎝⎭, (,)x y ∈Ω 其中Ω为2R 中的一个有界区域,其边界Γ为分段光滑曲线。
在Γ上u 满足下列边界条件之一:⑴(,)u x y αΓ=(第一边值条件), ⑵(,)ux y n βΓ∂=∂(第二边值条件), ⑶(,)uku x y n γΓ∂+=∂(第三边值条件), (,),(,),(,),(,),(,)f x y x y x y x y k x y αβγ都是连续函数,0k ≥.2.差分格式将区间[,]a b 作m 等分,记为11()/,,0i h b a m x a ih i m =-=+≤≤;将区间[,]c d 作n 等分,记为22()/,,0i h d c n y c jh j n =-=+≤≤.称1h 为x 方向的步长,2h 为y 方向的步长。
2.1 Poisson 方程五点差分格式参考单如图所示:以(,)i j x y 为中心沿y 方向Taylor 展开:41)(),j u y o h +①41)(),j u y o h +②41(),u h21(),o h ③22(),o h ④(,),i j ij f x y R -=+(,),i j f x y -=○6 j+1考虑到边值条件(,)(,)u x y x y αΓ=,构成差分格式:11112212(,)2(,)(,)(,)2(,)(,)(,),(,)(,),i j i j i j i j i j i j i j u x y u x y u x y u x y u x y u x y f x y h h u x y x y α+-+-Γ⎧-+-+⎛⎫-+=⎪ ⎪⎨⎝⎭⎪=⎩○72.2 Poisson 方程九点差分格式由上式 ③ + ④ 得:11112212442221244222222122222(,)2(,)(,)(,)2(,)(,)(,)1(,)()12(,)(,)1(,)12i j i j i j i j i j i j h i j i j iji j i j i j u x y u x y u x y u x y u x y u x y u x y h h u u u x y h h o h x y u x y u x y u x y h h x y x y +-+--+-+=+⎡⎤∂∂=∆+++⎢⎥∂∂⎣⎦⎛⎫∂∂⎛⎫∂∂=∆+++- ⎪ ⎪ ⎪∂∂∂∂⎝⎭⎝⎭422212222242222212122222(,)()12(,)(,)(,)1(,)()1212i j i j i j i j i j u x y h h o h x y f x y f x y u x y h h f x y h h o h x y x y ∂++∂∂⎛⎫∂∂∂+=--+-+ ⎪ ⎪∂∂∂∂⎝⎭○8 又()41122222211111112212311111(,)(,)2(,)(,)()1[(,)2(,)(,)2(,)2(,)(,)(,)2(,)(,)]()i j xx i j xx i j xx i j i j i j i j i j i j i j i j i j i j u x y u x y u x y u x y o h x y h u x y u x y u x y u x y u x y u x y h h u x y u x y u x y o h +-+++-++-+----∂-+=+∂∂=-+--++-++ 则得到:222222121121112112222221211212122222221112111211()(,)(210)(,)()(,)(210)(,)20()(,)(210)(,)(210)(,)()(,)()(,)i j i j i j i j i j i j i j i j i j h h u x y h h u x y h h u x y h h u x y h h u x y h h u x y h h u x y h h u x y h h u x y ---+--++-+++-++--++-+++-+--+-+2212222241222,12(,)(,)1(,)()12i j i j i j h hf x y f x y f x y h h o h x y ⎛⎫∂∂=--++ ⎪ ⎪∂∂⎝⎭○9 舍去截断误差得到逼近Poisson 方程的九点差分方程○10:()()2212,11,,11,1,11,11,11,122122212(,)[42]121(,)(,),12i j i j i j i j i j i j i j i j i j i j ij xx i j yy i j h h u x y u u u u u u u u u h h f h f x y h f x y -++-+---++-++-∆--+++++++''''=++考虑到边值条件(,)(,)u x y x y αΓ=,构成差分格式○11:()()2212,11,,11,1,11,11,11,122122212(,)[42]121(,)(,),12(,)(,),i j i j i j i j i j i j i j i j i j i j ijxx i j yy i j h h u x y u u u u u u u u u h h f h f x y h f x y u x y x y α-++-+---++-+Γ⎧+-∆--+++++++⎪⎪⎪''''=++⎨⎪⎪=⎪⎩3.格式求解3.1 Poisson 方程五点差分格式记122,1,j j j m j m j u u u u u --⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎣⎦,0.j n ≤≤ 矩阵格式改写为:11,11j j j j Du Cu Du f j m -+++=≤≤-,其中2221212222112122221121222112(1)111211112111121112m h h h h h h h C h h h h h h h -⎡⎤⎛⎫+-⎢⎥ ⎪⎝⎭⎢⎥⎢⎥⎛⎫⎢⎥-+- ⎪⎢⎥⎝⎭⎢⎥=⎢⎥⎢⎥⎛⎫⎢⎥-+- ⎪⎢⎥⎝⎭⎢⎥⎛⎫⎢⎥-+ ⎪⎢⎥⎝⎭⎣⎦,22222222(1)1111m h h D h h -⎡⎤-⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥-⎢⎥⎣⎦,10212212111(,)(,)(,)(,)1(,)(,)j j j j m j m j m j m f x y x y h f x y f f x y f x y x y h ---⎡⎤+Φ⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥+Φ⎢⎥⎣⎦, 可进一步写为:110222211(1)*(1).n n n n n n m u f Du C D u f D C D u f DC D u f Du D C -------⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-⎣⎦⎣⎦⎣⎦3.2 Poisson 方程九点差分格式记122,1,j j j m j m j u u u u u --⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎣⎦,0.j n ≤≤ 矩阵格式改写为:11,11j j j j Du Cu Du f j m -+++=≤≤-,其中2222121222222212121222222212121222221212(1)20()(210)(210)20()(210)(210)20()(210)(210)20()m h h h h h h h h h h C h h h h h h h h h h -⎡⎤+-⎢⎥-+-⎢⎥⎢⎥=⎢⎥-+-⎢⎥⎢⎥-+⎣⎦, 2222211222222212211222222212211222221221(1)(210)()()(210)()()(210)()()(210)m h h h h h h h h h h D h h h h h h h h h h -⎡⎤--+⎢⎥-+--+⎢⎥⎢⎥=⎢⎥-+--+⎢⎥⎢⎥-+-⎣⎦,22121022221211(,)(210)(,)(,)(,)(,)(210)(,)j j j j m j m j m j m f x y h h x y f x y f f x y f x y h h x y ---⎡⎤--Φ⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥+-Φ⎣⎦, 可进一步写为:110222211(1)*(1).n n n n n n m u f Du C Du f D C D u f DC D u f Du D C -------⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-⎣⎦⎣⎦⎣⎦4.数值例子4.1 Poisson 方程五点差分格式计算如下问题:22220,01,01,(0,)sin cos ,(2,)(sin cos ),01,(,0),(,1)(sin1cos1),0 1.x x u u x y x y u y y y u y e y y y u x e u x e x ⎛⎫∂∂-+=<<<< ⎪∂∂⎝⎭=+=+≤≤==+<<其精确解为:(,)(sin cos ).x u x y e y y =+,11,1,,1,222222122112112()(,),i j i j i j i j i j i j u u u u u f x y h h h h h h -+-++=++++ 考虑到本例中h1=h2,则有2,11,1,,1,(,),4i j i j i j i j i j i j u u u u h f x y u -+-+++++=利用Gauss-Seidel 迭代方法对k=0,1,2,……,计算112,11,1,,11(,),41,2,....,1;1,2,...., 1.k k k k i j i j i j i j i j k ij u u u u h f x y u i m j n ++--+++++++==-=-表1 部分结点处的精确解和取不同步长时所得的数值解表2 取不同步长时部分结点处数值解的误差绝对值图1 取h=1/4时所得的数值解曲线图2 取h=1/4时所得的误差曲线图3 取h=1/16时所得的数值解曲线图4 取h=1/16时所得的误差曲线图5 取h=1/64时所得的数值解曲线图6 精确解曲线图7 取h=1/64时所得的误差曲线4.2 Poisson 方程九点差分格式计算如下问题:22220,01,01,(0,)sin cos ,(2,)(sin cos ),01,(,0),(,1)(sin1cos1),0 1.x x u u x y x y u y y y u y e y y y u x e u x e x ⎛⎫∂∂-+=<<<< ⎪∂∂⎝⎭=+=+≤≤==+<<其精确解为(,)(sin cos ).x u x y e y y =+222222221212121112122222222121112111211211222211120()(,)12(,)()(,)(102)(,)()(,)()(,)()(,)(102)(,)(102)(,)(10i j i j i j i j i j i j i j i j i j h h u x y h h f x y h h u x y h h u x y h h u x y h h u x y h h u x y h h u x y h h u x y h ----++++--++=+++-+++++++-+-+2212)(,)i j h u x y +-考虑到本例中h1=h2,则有,11,1,,11,11,11,11,1,4(),20i j i j i j i j i j i j i j i j i j u u u u u u u u u -+-+--++-++-+++++++=利用Gauss-Seidel 迭代方法对k=0,1,2,……,计算1111,11,1,,11,11,11,11,11,4(),201,2,....,1;1,2,...., 1.k k k k k k k k i j i j i j i j i j i j i j i j k i j u u u u u u u u u i m j n ++++-+-+--++-++-++++++++==-=-表1 部分结点处的精确解和取不同步长时所得的数值解表2 取不同步长时部分结点处数值解的误差绝对值表3 取不同步长时部分结点处数值解的最大误差图1 取h=1/4时所得的数值解曲线图2 取h=1/16时所得的数值解曲线图3 取h=1/64时所得的数值解曲线图4 取h=1/4时所得的误差曲线图5 取h=1/16时所得的误差曲线图6 取h=1/64时所得的误差曲线5.结论观察Poisson方程五点格式,方程以较快速度迭代收缩。
九点差分格式五点差分格式matlab

九点差分格式五点差分格式matlab九点差分格式概述九点差分格式是一种数值计算方法,用于求解偏微分方程。
它是有限差分法的一种,通常用于求解二维泊松方程或热传导方程。
原理九点差分格式是通过在目标点周围的八个邻居点上进行数值计算来估计目标点的值。
其基本形式为:f(x,y) = (1/6h^2) * [4f(x+h,y)+4f(x-h,y)+4f(x,y+h)+4f(x,y-h)-20f(x,y)] + (1/3h^2) * [f(x+h,y+h)+f(x-h,y-h)+f(x+h,y-h)+f(x-h,y+h)]其中,h为网格尺寸。
优缺点九点差分格式的优点是精度较高,适用于求解高精度问题。
但由于需要使用八个邻居点进行计算,因此其计算量较大。
应用九点差分格式通常应用于求解二维泊松方程或热传导方程等偏微分方程问题。
五点差分格式概述五点差分格式是一种数值计算方法,也是有限差分法的一种。
它通常用于求解二维泊松方程或热传导方程等偏微分方程问题。
原理五点差分格式是通过在目标点周围的四个邻居点上进行数值计算来估计目标点的值。
其基本形式为:f(x,y) = (1/h^2) * [4f(x+h,y)+4f(x-h,y)+4f(x,y+h)+4f(x,y-h)-20f(x,y)]其中,h为网格尺寸。
优缺点五点差分格式的优点是计算量较小,适用于求解大规模问题。
但由于只使用了四个邻居点进行计算,因此其精度较低。
应用五点差分格式通常应用于求解二维泊松方程或热传导方程等偏微分方程问题。
由于其计算量较小,因此也适用于求解大规模问题。
Matlab中的实现Matlab是一种常用的数值计算软件,可以很方便地实现九点差分格式和五点差分格式。
以下是两种方法的Matlab代码:九点差分格式:function [u] = nine_point(f, h)[m, n] = size(f);u = zeros(m, n);for i = 2:m-1for j = 2:n-1u(i,j) = (1/6*h^2)*(4*f(i+1,j)+4*f(i-1,j)+4*f(i,j+1)+4*f(i,j-1)-20*f(i,j)) + (1/3*h^2)*(f(i+1,j+1)+f(i-1,j-1)+f(i+1,j-1)+f(i-1,j+1));endendend五点差分格式:function [u] = five_point(f, h)[m, n] = size(f);u = zeros(m, n);for i = 2:m-1for j = 2:n-1u(i,j) = (1/h^2)*(4*f(i+1,j)+4*f(i-1,j)+4*f(i,j+1)+4*f(i,j-1)-20*f(i,j));endendend以上是九点差分格式和五点差分格式的基本原理、优缺点和应用,以及在Matlab中的实现方法。
二维扩散方程的9点格式有限近似解法

二维扩散方程的9点格式有限近似解法扩散方程,又称拉普拉斯方程,是一个基本的偏微分方程,用于描述物质、信息向外扩散的过程,应用范围广泛,比如热传导、流体力学、化学反应等。
二维扩散方程是扩散方程研究的重要类型。
9点格式有限近似解法是求解二维扩散方程的一种方法,拥有精确性高、计算量小的优点,被广泛应用于对二维扩散方程的数值模拟。
本文的目的是研究9点格式有限近似解法求解二维扩散方程的方法。
文章将从以下几个方面进行阐述:1)9点格式有限近似解法的数学模型;2)9点格式有限近似解法的实现流程;3)9点格式有限近似解法的优缺点。
二、二维扩散方程及9点格式有限近似解法的数学模型二维扩散方程的数学模型可以用下面的方程表示:$$frac{partial P}{partialt}=Dleft(frac{partial^2P}{partialx^2}+frac{partial^2P}{partial y^2}right)$$其中,$P(x,y,t)$ 为扩散的物质的浓度,$D$ 为扩散系数,$x,y$别为方程的两个空间变量,$t$为时间变量。
经过空间离散,可以将上式离散化为下面的六阶线性方程组:$$frac{P_{i,j}^{n+1}-P_{i,j}^n}{Delta t}=Dleft(frac{P_{i-1,j}^{n+1}+P_{i+1,j}^{n+1}-2P_{i,j}^{n+1}}{D eltax^2}+frac{P_{i,j-1}^{n+1}+P_{i,j+1}^{n+1}-2P_{i,j}^{n+1}}{Delta y^2}right)$$其中$P_{i,j}^n$表示在空间变量$x,y$分别取$x_i,y_j$,在时间变量$t$取$t_n$时的物质浓度,$Delta t,Delta x,Delta y$分别表示时间步长,x轴、y轴空间步长。
上述方程可以采用9点格式近似解法处理:$$P_{i,j}^{n+1}=frac{DDelta t}{Delta x^2+Delta y^2}left[(Delta x^2-Deltay^2)(P_{i-1,j}^{n+1}+P_{i+1,j}^{n+1})+(Delta y^2-Deltax^2)(P_{i,j-1}^{n+1}+P_{i,j+1}^{n+1})+4Delta x^2Deltay^2P_{i,j}^nright]$$三、 9点格式有限近似解法的实现流程(一)按照X轴和Y轴的方向分别对边界值进行赋值;(二)计算第一次迭代时候的内部点值;(三)多次迭代求解,直到收敛;(四)检查计算结果是否稳定,若不稳定则重复多次迭代;(五)最终得到稳定的计算结果。
Poisson方程九点差分格式_米瑞琪

Poisson方程九点差分格式_米瑞琪数值实验报告I实验名称Poisson方程九点差分格式实验时间2016年 4 月 15 日姓名米瑞琪班级信息1303学号04成绩一、实验目的,内容1、理解Poisson方程九点差分格式的构造原理;2、理解因为网格点的不同排序方式造成的系数矩阵格式的差异;3、学会利用matlab的spdiags(),kron()函数生成系数矩阵;二、算法描述针对一个Poisson方程问题:在Poisson方程五点差分格式的基础上,采用Taylor展开分析五点差分算子的截断误差,可以得到:为了提高算子截断误差的精度,在(1)式中配凑出了差分算子的形式,将原Poisson方程代入(1)式有:考虑,有:将(3)代回(2)可得得到Poisson方程的九点差分格式:在计算机上实现(4)式,需要在五点差分格式的基础上在等式两端分别增加一部分,将等式左侧新增的部分写成紧凑格式,有:对于该矩阵,可以看成是两个矩阵的组合:以及则生成这两个矩阵可以采用Kroncker生成,方法类似于五点差分格式。
对于右端添加的关于f(x,y)的二阶导数,可以采用中心差分格式进行近似代替,即:写成相应的紧凑格式有:该式中的矩阵又可以分解为两个矩阵的和:%计算误差u_real=@(x,y)exp(pi*(x+y))*sin(pi*x).*sin(pi*y);for i=1:N1-1u_m((i-1)*(N2-1)+1:i*(N2-1))=u_real(x(i),y);endu_v=u_m';err_d=max(abs(u_d-u_v));sol=reshape(u_d,N2-1,N1-1);mesh(X,Y,sol)四. 数值结果针对课本P93给出的问题,分别采用步长,将计算出的误差列表如下:步长五点差分格式误差九点差分格式误差可见采用九点差分格式可以进一步缩小误差,达到更高阶的精度。
五. 计算中出现的问题,解决方法及体会在生成九点差分格式的时候,等号右端涉及到了对f的二阶偏导,我最初利用符号函数定义了f,随后求出其二阶偏导(仍然是符号函数)之后带入网格点,求f二阶偏导的精确解,但是代入过程相当繁琐,运行速度非常慢,最终我改变策略,选用f关于x,y的二阶中心差分格式替代精确值,最终得到了相对满意的结果。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数值实验报告I
实验名称Poisson方程九点差分格式实验时间2016年 4 月 15 日姓名米瑞琪班级信息1303学号04成绩
一、实验目的,内容
1、理解Poisson方程九点差分格式的构造原理;
2、理解因为网格点的不同排序方式造成的系数矩阵格式的差异;
3、学会利用matlab的spdiags(),kron()函数生成系数矩阵;
二、算法描述
针对一个Poisson方程问题:
在Poisson方程五点差分格式的基础上,采用Taylor展开分析五点差分算子的截断误差,可以得到:
为了提高算子截断误差的精度,在(1)式中配凑出了差分算子的形式,将原Poisson方程代入(1)式有:
考虑,有:
将(3)代回(2)可得
得到Poisson方程的九点差分格式:
在计算机上实现(4)式,需要在五点差分格式
的基础上在等式两端分别增加一部分,将等式左侧新增的部分写成紧凑格式,有:
对于该矩阵,可以看成是两个矩阵的组合:
以及
则生成这两个矩阵可以采用Kroncker生成,方法类似于五点差分格式。
对于右端添加的关于f(x,y)的二阶导数,可以采用中心差分格式进行近似代替,即:
写成相应的紧凑格式有:
该式中的矩阵又可以分解为两个矩阵的和:
%计算误差
u_real=@(x,y)exp(pi*(x+y))*sin(pi*x).*sin(pi*y);
for i=1:N1-1
u_m((i-1)*(N2-1)+1:i*(N2-1))=u_real(x(i),y);
end
u_v=u_m';
err_d=max(abs(u_d-u_v));
sol=reshape(u_d,N2-1,N1-1);
mesh(X,Y,sol)
四. 数值结果
针对课本P93给出的问题,分别采用步长,将计算出的误差列表如下:
步长五点差分格式误差九点差分格式误差
可见采用九点差分格式可以进一步缩小误差,达到更高阶的精度。
五. 计算中出现的问题,解决方法及体会
在生成九点差分格式的时候,等号右端涉及到了对f的二阶偏导,我最初利用符号函数定义了f,随后求出其二阶偏导(仍然是符号函数)之后带入网格点,求f二阶偏导的精确解,但是代入过程相当繁琐,运行速度非常慢,最终我改变策略,选用f关于x,y的二阶中心差分格式替代精确值,最终得到了相对满意的结果。
教
师
评
语
指导教师:年月日。