《椭圆型方程的差分解法(共20页).docx》由会员分享，可在线阅读，更多相关《椭圆型方程的差分解法(共20页).docx（20页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、精选优质文档-倾情为你奉上椭圆型方程的差分解法1.引言考虑问题二维Poisson方程， 其中为中的一个有界区域，其边界为分段光滑曲线。在上满足下列边界条件之一：（第一边值条件），（第二边值条件），（第三边值条件），都是连续函数，.2.差分格式将区间作m等分，记为；将区间作n等分，记为.称为方向的步长，为方向的步长。2.1 Poisson方程五点差分格式参考单如图所示：j+1jj-1i+1i-1i以为中心沿方向Taylor展开：由 + 得：整理得到：同理，以为中心沿方向Taylor展开：代入原方程得：得到截断误差：其中是原方程光滑解，舍去截断误差得到逼近Poisson方程的五点差分方程：考虑到边
2、值条件，构成差分格式：2.2 Poisson方程九点差分格式由上式 + 得：又则得到：舍去截断误差得到逼近Poisson方程的九点差分方程：考虑到边值条件，构成差分格式：3. 格式求解3.1 Poisson方程五点差分格式记，矩阵格式改写为：，其中，可进一步写为：3.2 Poisson方程九点差分格式记，矩阵格式改写为：，其中，可进一步写为：4.数值例子4.1 Poisson方程五点差分格式计算如下问题：其精确解为：考虑到本例中h1=h2，则有利用Gauss-Seidel迭代方法对k=0,1,2,计算表1 部分结点处的精确解和取不同步长时所得的数值解(h1,h2)(x,y)(1/4,1/4)(
3、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.2.2.1.2.2.(1/8,1/8)1.2.2.1.2.2.(1/16,1/16)1.2.2.1.2.2.(1/32,1/32)1.2.2.1.2.2.(1/64,1/64)1.2.2.1.2.2.精确解1.2.2.1.2.2.表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-04
4、0.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-06表3 取不同
5、步长时部分结点处数值解的最大误差(h1,h2)(1/4,1/4)0.0016*(1/8,1/8)4.2732e-043.(1/16,1/16)1.0885e-043.(1/32,1/32)2.6148e-054.(1/64,1/64)1.1928e-052.图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方程九点差分格式计算如下问题：其精确解为考虑到本例中h1=h2，则有利用Gauss-S
6、eidel迭代方法对k=0,1,2,计算表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.2.2.1.2.2.(1/8,1/8)1.2.2.1.2.2.(1/16,1/16)1.2.2.1.2.2.(1/32,1/32)1.2.2.1.2.2.(1/64,1/64)1.2.2.1.2.2.精确解1.2.2.1.2.2.表2 取不同步长时部分结点处数值解的误差绝对值(h1,h2)(x,y)(1/4,1/4)(1/2,1/4)(3/4,1/4)(
7、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)4.0638e-085.3998e-083.5875e-083.3749e-084.4845e-082.9794e-08(1/32,1/32)1.3096e-051.9001e-051.7247e-051.4215e-052.0586e-051.8853e-05(1/64,1/64
8、)6.1867e-068.3780e-064.7959e-065.5548e-067.4932e-064.0549e-06表3 取不同步长时部分结点处数值解的最大误差(h1,h2)(1/4,1/4)5.1568e-09*(1/8,1/8)5.1568e-091.(1/16,1/16)6.9859e-080.(1/32,1/32)1.0268e-060.(1/64,1/64)1.5435e-050.图1 取h=1/4时所得的数值解曲线图2 取h=1/16时所得的数值解曲线图3 取h=1/64时所得的数值解曲线图4 取h=1/4时所得的误差曲线图5 取h=1/16时所得的误差曲线图6 取h=1/6
9、4时所得的误差曲线5.结论观察Poisson方程五点格式，方程以较快速度迭代收缩。对于Poisson方程九点格式，一开始有较快速度和较好结果，但是与五点格式比起来，迭代速度因步长增加而扩散，结果很不理想，且运行速度与五点格式相比慢了许多，与预计结果有失偏颇。参考文献：1 孙志忠. 偏微分方程数值解法M. 北京：科学出版社， 20052 李荣华. 偏微分方程数值解法M. 北京：高等教育出版社， 20043 易昆南. 基于数学建模的数学实验M. 长沙：中国铁道出版社， 2014源代码：%结点取值，结点误差及最大误差输出函数functionp,e,u,x,y,k=five(h,m,n,kmax,ep
10、)%p代表数值解，h代表步长,m,n,分别代表x、y轴结点数,kmax代表最大迭代次数，ep代表最大误差,e用于存储精确解syms temp;u=zeros(n+1,m+1);x=0+(0:m)*h;y=0+(0:n)*h;for(j=1:n+1) %边值条件 u(j,1)=sin(y(j)+cos(y(j); u(j,m+1)=exp(1)*(sin(y(j)+cos(y(j);endfor(i=1:m+1) u(1,i)=exp(x(i); u(n+1,i)=exp(x(i)*(sin(1)+cos(1);endt=zeros(n-1,m-1);for(k=1:kmax) %存储迭代结果f
11、or(j=2:n)for(i=2:m) temp=(u(j,i+1)+u(j,i-1)+u(j+1,i)+u(j-1,i)/4; t(j,i)=(temp-u(j,i)*(temp-u(j,i); u(j,i)=temp;endend t(j,i)=sqrt(t(j,i);if(kkmax) %结果满足kmax次停止迭代break;endif(max(max(t)e); %输出最大误差 e=abs(u(j,i)-p(j,i);endif(x(i)=1/4&y(j)=1/4)|(x(i)=1/2&y(j)=1/4)|(x(i)=3/4&y(j)=1/4)|(x(i)=1/4&y(j)=3/4)|(x(i)=1/2&y(j)=3/4)|(x(i)=3/4&y(j)=3/4) %输出结点处函数值和误差 disp(u(j,i); disp(abs(u(j,i)-p(j,i);endendenddisp(e);surf(x,y,u); %近似解曲线surf(x,y,p); %精确解曲线surf(x,y,e); %误差曲线专心-专注-专业