椭圆型方程差分法
- 格式:pdf
- 大小:244.03 KB
- 文档页数:34
关于解二维椭圆型高精度差分方程的交替方向迭代法交替方向迭代法是一种常用于解二维椭圆型差分方程的方法。
在许多实际应用中,需要对这种方程进行高精度的求解。
交替方向迭代法可以有效地解决这个问题,并且在计算效率和内存利用方面都具有优势。
交替方向迭代法的基本思想是将二维椭圆型差分方程化为两个一维的差分方程,然后分别对其进行迭代求解。
具体地,将二维椭圆型差分方程表示为:$-\frac{\partial^2 u}{\partial 某^2}-\frac{\partial^2u}{\partial y^2}=f(某,y)$。
其中,$u(某,y)$是要求解的未知函数,$f(某,y)$是已知的函数。
接下来,我们可以将上式分别用某轴和y轴方向进行差分,得到两个一维的差分方程:$-\frac{u_{i-1,j}-2u_{i,j}+u_{i+1,j}}{\Delta 某^2}-\frac{u_{i,j-1}-2u_{i,j}+u_{i,j+1}}{\Delta y^2}=f_{i,j}$。
$-\frac{u_{i-1,j}-2u_{i,j}+u_{i+1,j}}{\Delta 某^2}-\frac{u_{i,j-1}-2u_{i,j}+u_{i,j+1}}{\Delta y^2}=f_{i,j}$。
其中,$u_{i,j}$表示在坐标$(i,j)$处的未知函数值,$\Delta 某$和$\Delta y$分别表示在某轴和y轴方向的差分间隔。
接下来就是交替方向迭代法的核心部分。
首先,我们需要对第一个方程中的某方向进行求解,可以通过迭代求解来不断逼近最终的解。
假设我们用$u^{k}_{i,j}$表示在第k次迭代时,在坐标$(i,j)$处的近似解。
那么对于每个坐标$(i,j)$,我们可以通过以下公式来更新近似解:$u^{k+1}_{i,j}=\frac{1}{2}(\frac{u^{k}_{i-1,j}+u^{k}_{i+1,j}}{\Delta 某^2}+\frac{u^{k}_{i,j-1}+u^{k+1}_{i,j+1}}{\Delta y^2}-f_{i,j})$。
椭圆型方程的差分方法差分方法是一种数值计算方法,使用近似的差商来表示微分方程。
椭圆型方程是一类常见的偏微分方程,具有重要的数学和物理应用。
在本文中,我们将介绍椭圆型方程的差分方法,并讨论其优点和缺点。
一、椭圆型方程的差分近似L[u]=-∂(p∂u/∂x)/∂x-∂(q∂u/∂y)/∂y+r(x,y)u=f(x,y)其中,L[u]是一个偏微分算子,u(x,y)是未知函数,p(x,y),q(x,y),r(x,y),f(x,y)是已知函数。
椭圆型方程的解通常在一个区域Ω上求解。
差分方法的主要思想是用网格来离散化区域Ω,将连续的偏微分方程转化为离散的代数方程。
对于椭圆型方程,我们可以选择矩形网格,其中Ω可以被划分为N*M个小矩形,并且网格的步长为Δx和Δy。
假设我们要在网格点(xi, yj)处求解未知函数的值uij,其中i和j分别表示网格的行索引和列索引。
我们可以使用中心差分法来近似x和y方向的偏导数,从而得到离散形式的椭圆型方程:L[u] ≈ -(p(xi+1/2, yj)(ui+1,j - ui-1,j)/Δx^2 + p(xi,yj+1/2)(ui,j+1 - ui,j-1)/Δy^2) + q(xi,yj)uij = f(xi,yj)其中,p(xi+1/2, yj)和p(xi, yj+1/2)分别表示在(xi+1/2, yj)和(xi, yj+1/2)处的系数。
可以通过有限差分方式计算出这些系数。
将上述公式在每个网格点(xi, yj)处形成一个方程,从而得到一个线性方程组。
通过求解这个线性方程组,我们可以得到网格点上的未知函数值。
二、椭圆型方程差分方法的优点和缺点差分方法是一种简单有效的数值计算方法,具有以下优点:1.可以处理任意形状的区域Ω:差分方法可以适应不规则网格和复杂区域,因此适用于各种几何形状的椭圆型方程求解。
2.数值稳定性:差分方法可以确保数值解的稳定性,避免数值上的不稳定问题。
3.线性时间复杂度:差分方法的计算复杂度通常是线性的,即解方程的时间随着网格点数的增加而线性增加。
《微分方程数值解》大作业(一)——椭圆型方程编程计算:采用五点差分格式求如下椭圆型方程2222uu x y f (x,y),(x,y);∂∂∂∂--=∈Ω其中f (x,y)、Ω及边条件为:1. f (x,y)0,= (1,2)(0,1)Ω=⨯, 且边条件如下:222u(x,0)2ln x,u(x,1)ln(x 1)1x 2;u(1,y)ln(1y ),u(2,y)ln(4y ),0y 1.⎧==+<<⎪⎨=+=+<<⎪⎩, 问题存在精确解为: 22(,)ln()u x y x y =+2.f (x,y)4,=- (0,1)(0,2)Ω=⨯,且边条件如下:2222u(x,0)x ,u(x,2)(x 2)0x 1;u(0,y)y ,u(1,y)(y 1),0y 2.⎧==-<<⎪⎨==-<<⎪⎩, 问题存在精确解为: 2(,)()u x y x y =-3.f (x,y)cos(x y)cos(x y),=++- (0,)(0,)2πΩ=π⨯,且边条件如下:u(x,0)cos x,u(x,)00x ;2u(0,y)cos y,u(,y)cos y,0y .2π⎧==<<π⎪⎪⎨π⎪=π=-<<⎪⎩, 问题存在精确解为: (,)cos cos u x y x y =.代码:主函数1,差分解function g=fivepoints(x1,x2,y1,y2,M,N)%变步长法h=(x2-x1)/M; %横轴步长k=(y2-y21/N; %纵轴步长m=M-1;n=N-1;h1=h^2;r=h1/k^2; %五点中的上下两个点的系数t=2+2*r; %五点中的中心点的系数x=x1+(x2-x1)*(0:M)/M; %x,y向量表示横纵坐标y=y1+(y2-y1)*(0:N)/N;a=zeros(m*n,m*n);b=zeros(m*n,1);%初始化a,b矩阵,a为系数矩阵%内部的(m-2)*(n-2)个点for i=2:m-1for j=2:n-1a(i+(j-1)*m,:)=[zeros(1,i-1+(j-2)*m) -r zeros(1,m-2) -1 t -1 zeros(1,m-2) -r zeros(1,(n-j)*m-i)];b(i+(j-1)*m)=h1*f(x(i+1),y(j+1));endend%下边缘j=1;for i=2:m-1a(i+(j-1)*m,:)=[zeros(1,i-2) -1 t -1 zeros(1,m-2) -r zeros(1,(n-j)*m-i)];b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+r*bottom(x(i+1));end;%右边缘i=m;for j=2:n-1a(i+(j-1)*m,:)=[zeros(1,(j-1)*m-1) -r zeros(1,m-2) -1 t zeros(1,m-1) -r zeros(1,(n-j)*m-i)];b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+right(y(j+1));end%上边缘j=n;for i=2:m-1a(i+(j-1)*m,:)=[zeros(1,i-1+(j-2)*m) -r zeros(1,m-2) -1 t -1 zeros(1,m-i-1)];b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+r*top(x(i+1));end%左边缘i=1;for j=2:n-1a(i+(j-1)*m,:)=[zeros(1,i-1+(j-2)*m) -r zeros(1,m-1) t -1 zeros(1,m-2) -rzeros(1,(n-j)*m-i)];b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+left(y(j+1));end;%左下角的那个点i=1;j=1;a(1,:)=[t -1 zeros(1,m-2) -r zeros(1,(n-1)*m-1)];b(1)=h1*f(x(2),y(2))+r*bottom(x(2))+left(y(2));%右下角的那个点i=m;j=1;a(i+(j-1)*m,:)=[zeros(1,m-2) -1 t zeros(1,m-1) -r zeros(1,(n-2)*m)]; b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+r*bottom(x(i+1))+right(y(j+1)); %左上角的那个点i=1;j=n;a(i+(j-1)*m,:)=[zeros(1,(n-2)*m) -r zeros(1,m-1) t -1 zeros(1,m-2)]; b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+r*top(x(i+1))+left(y(j+1));%右上角的那个点i=m;j=n;a(i+(j-1)*m,:)=[zeros(1,(n-1)*m-1) -r zeros(1,m-2) -1 t];b(i+(j-1)*m)=h1*f(x(i+1),y(j+1))+r*top(x(i+1))+right(y(j+1));u=a\bab2,精确解:function g=ni(x1,x2,y1,y2,M,N)m=M-1;n=N-1;x=x1+(x2-x1)*(0:M)/M;y=y1+(y2-y1)*(0:N)/N;for i=1:mfor j=1:nu1(i+(j-1)*m)=f1(x(i+1),y(j+1))endend(1)辅助函数function g=f(x,y)g=0;function g=bottom(x)g=2*log(x);function g=right(y)g=log(4+y^2);function g=top(x)g=log(x^2+1);function g=left(y)g=log(1+y^2);function g=f1(x,y)g=log(x^2+y^2);运行fivepoints(1,2,0,1,4,4)u =数值解0.4847467147016780.8376456266975491.1390195099193150.5944295076643080.9158860659528741.1974022894530100.7539416986884711.0340668399966291.287784599003526a =4 -1 0 -1 0 0 0 0 0 -1 4 -1 0 -1 0 0 0 0 0 -1 4 0 0 -1 0 0 0 -1 0 0 4 -1 0 -1 0 0 0 -1 0 -1 4 -1 0 -1 0 0 0 -1 0 -1 4 0 0 -1 0 0 0 -1 0 0 4 -1 0 0 0 0 0 -1 0 -1 4 -1 0 0 0 0 0 -1 0 -1 4b =0.5069117244448540.8109302162163292.5210301235267010.2231435513142101.4469189829363251.3872704470929461.1786549963416462.919669266564466运行ni(1,2,0,1,4,4)u1 =精确解Columns 1 through 30.485507815781701 0.838329190404443 1.139434283188365 Columns 4 through 60.594707107746693 0.916290731874155 1.197703191312341 Columns 7 through 90.753771802376380 1.034073767530539 1.287854288306638 误差很小(2)辅助函数function g=f(x,y)g=-4;function g=bottom(x)g=x^2;function g=right(y)g=(y-1)^2;function g=top(x)g=(x-2)^2;function g=left(y)g=y^2;function g=f1(x,y)g=(x-y)^2;fivepoints(1,2,0,1,4,4)fivepoints(0,1,0,2,4,4)u =0.062500000000000-0.0000000000000000.0625000000000000.5625000000000000.2500000000000000.0625000000000001.5625000000000001.0000000000000000.562500000000000a =Columns 1 through 32.500000000000000 -1.000000000000000 0 -1.000000000000000 2.500000000000000-1.0000000000000000 -1.000000000000000 2.500000000000000 -0.250000000000000 0 00 -0.250000000000000 00 0 -0.2500000000000000 0 00 0 00 0 0Columns 4 through 6-0.250000000000000 0 00 -0.250000000000000 00 0 -0.2500000000000002.500000000000000 -1.000000000000000 0 -1.000000000000000 2.500000000000000-1.0000000000000000 -1.000000000000000 2.500000000000000 -0.250000000000000 0 00 -0.250000000000000 00 0 -0.250000000000000Columns 7 through 90 0 00 0 00 0 0-0.250000000000000 0 00 -0.250000000000000 00 0 -0.2500000000000002.500000000000000 -1.000000000000000 0 -1.000000000000000 2.500000000000000-1.0000000000000000 -1.000000000000000 2.500000000000000b =0.015625000000000-0.1875000000000000.1406250000000000.750000000000000-0.250000000000000-0.2500000000000002.7656250000000000.3125000000000000.390625000000000精确解ni(0,1,0,2,4,4)u1 =u1 =Columns 1 through 30.062500000000000 0 0.062500000000000 Columns 4 through 60.562500000000000 0.2500000000000000.062500000000000Columns 7 through 91.562500000000000 1.0000000000000000.562500000000000误差很小(3)辅助函数function g=f(x,y)g=cosd(x+y)+cosd(x-y);function g=bottom(x)g=cosd(x);function g=right(y)g=-cosd(y);function g=top(x)g=0;function g=left(y)g=cosd(y);function g=f1(x,y)g=cosd(x)*cosd(y);数值解Pi=3.1415926fivepoints(0,pi,0,pi/2,4,4)u =0.6578183624886530.000000024999241-0.6578183271343870.5049807980892560.000000019229497-0.5049807708946410.2736443626241530.000000010432161-0.273644347870850a =10 -1 0 -4 0 0 0 0 0 -1 10 -1 0 -4 0 0 0 0 0 -1 10 0 0 -4 0 0 0 -4 0 0 10 -1 0 -4 0 0 0 -4 0 -1 10 -1 0 -4 0 0 0 -4 0 -1 10 0 0 -4 0 0 0 -4 0 0 10 -1 0 0 0 0 0 -4 0 -1 10 -1 0 0 0 0 0 -4 0 -1 10b =4.5582604075302670.000000137720159-4.5582602127645491.323957*********0.000000023374742-1.3239570281549570.7165204234523470.000000012650320-0.716520405562093精确解ni(0,pi,0,pi/2,4,4)u1 =Columns 1 through 30.653281493003155 0.000000024755257-0.653281457993935Columns 4 through 60.500000013397448 0.000000018946853-0.499999986602551Columns 7 through 90.270598066826879 0.000000010253963-0.270598052325585误差很小注:(1)需要对数值解与精确解作比较,以及不同步长选取下的误差比较。
关于解二维椭圆型高精度差分方程的交替方向迭代法交替方向迭代法(alternating direction iterative method,简称ADI法)是一种常用于解二维椭圆型高精度差分方程的数值方法。
它的基本思想是将二维问题拆分为两个一维问题,然后依次迭代求解,从而达到解二维问题的目的。
一、AD方法的原理二维椭圆型方程可以表示为:>(A1u(i+1,j)+B1u(i,j+1)-(2A1+B1+C1)u(i,j)+B1u(i,j-1)+A1u(i-1,j))/Δx^2+(D1Δy^2)u(i,j)=f(i,j),式中i,j分别为空间坐标,f(i,j)为已知的外力场,D1为系数。
对于前一半时间步,我们考虑j方向上的一维问题,可将上述方程改写为:>(A1+B1)u(i+1,j)+(C1-2A1-B1)u(i,j)+A1u(i-1,j)=f1(i,j),式中f1(i,j)为已知。
对于后一半时间步,我们考虑i方向上的一维问题,可将上述方程改写为:>(A2+B2)u(i,j+1)+(C2-2A2-B2)u(i,j)+A2u(i,j-1)=f2(i,j),式中f2(i,j)为已知。
通过交替迭代这两个方程,逐步逼近真实解。
二、AD方法的步骤1.选取初始解u0(i,j),以及参数Δx,Δy和迭代终止准则ε。
2.开始迭代,设定初始误差e0=ε+1,n=0。
3. 使用前一半时间步的方程求解每个i所对应的j方向的一维问题,得到一个临时解ut(i,j)。
4.使用后一半时间步的方程求解每个j所对应的i方向的一维问题,得到一个新的解u(n+1)(i,j)。
5. 计算误差en=,u(n+1)-u(n),2,n=n+16. 如果en ≤ ε,迭代结束,得到近似解u(n+1);否则返回第3步。
7.算法结束。
三、AD方法的优势1.AD方法是一种分离的迭代方法,当网格较大时其计算量相对较小。
2.AD方法是一种隐式迭代方法,对稳定性和收敛性有较好的保证。
课程大纲课程编号(理学院)课程名称随机规划学时40基本预备知识 1. 概率统计2. 最优化理论与算法3. 随机过程授课方式讲授、研讨基本要求掌握随机规划模型的类型。
(3TKH 主要类型),了解分布问题中参数LP 及其最优值得表达式,了解Z(3 )的可测性及其概率分布,掌握简单分布问题的计算方法,了解逼近方法和最优值的数学期望的估计,掌握有补偿的二阶段问题和二阶段问题的数值解法,了解概率约束规划和随机拟次梯度法,了解上图收敛性。
教材及参考书《随机规划》,王全德编著,南京大学出版社,1990 年。
《随机线性规划》,Kall 著,王金德译,南京大学出版社。
讲授的主要内容:(每章后附学时数)1.随机规划的模型(6 学时)1.1分布问题,二阶段有补偿问题,概率约束问题;1.2多阶段有补偿问题和多阶段概率约束计划;1.3各类问题的统一形式与相互关系。
2.分布问题:(6 学时)2.1参数LP;2.2Z(3)的可测性;2.3最优化Z(3 )的概率分布;2.4简单分布问题的计算方法;2.5逼近方法与最优值的数学期望的估计。
3.有补偿二阶段问题(8 学时)3.1一般有补偿二阶段的问题;3.2具有固定补偿矩阵的情形;3.3具有完备和简单补偿矩阵的二阶段问题。
4.二阶段问题的数值解法(8 学时)4.1具有离散随机变量的二阶段问题的解法;4.2简单补偿问题的解法。
5.概率约束规划(6 学时)可行解集合的特性,约束函数的分析性质,数值解法,逼近方法。
6.随机拟次梯度法(* )(2 学时)7. 应用举例(2 学时)8. 上图收敛性(2 学时)注:(*)只做了解课程名称学时基本预备知识值代数601. 数学分析2. 线性代数3. 矩阵论4. 计算方法授课方式讲授基本要求1. 知道矩阵计算的基本工具,熟悉Vandermonde、Toeplitz 等方程组的解法及某些迭代法的收敛性,了解多项式加速技巧。
2.掌握不完全分解预先共轭梯度法,广义共轭剩余法,Lanczos 方法,求解特征值问题的同伦方法和分而治之法以及求解Jacobi 矩阵特征值反问题的正交约化法。
第4章 椭圆型方程的有限差分法§ 2 一维差分格式1、用积分插值法导出逼近微分方程的差分格式。
, d / du 、 du 上 , Lu=- (p )+r +qu=f,a<x<b, dx dx dxu(a)= a ,u(b)=如果系数p,q,r 以及右端f 光滑,则可用中矩形公式计算得解:考虑在[a,b ]内任一小区间[x ⑴,x (2)],将上式在此区间上积分得 X ⑵d dU曲 d ;(p(x)d ;)dxW(x ⑴)W(x ⑵)x(2)dU(1) r ——dxxdx x ⑵dU 」 xr dxx (2) x (1) qudx X (2)X (1) qudxX (2)x (1)fdXX ⑵X (1) fdX其中,W(x) p吧特别地,取[X ⑴,X ⑵]为对偶单元[x i 1/2 , x i 1/2 ],则W (x i 1/2)W(x i 1/2)x1/2dUr ——dxxi 1/2 qudxxi 1/ 2xi 1/ 2 fdx 。
Xi 1/ 2-J..将(1.2)改写成—dxW(x)再沿[x i 1/2, x i 1/2]积分,得U i U iX iW(x) dx , Xi 1p(x)矩形公式,得W 1/2a iU i U ih i a iNdx x1 p(x)]xi 1/ 2又qUdxXi 1/ 2h ih i2:d i U i , d ixi 1/ 2 x q(x)dxX1/2x1/2dU , r ——dx.U i 1 U i 1bi "T-2 h i h i 1x 1/2 xr(x)dxxi1/22ih i h i 1x1/ 2 f (x)dxxi 1/2(1.1) (1.2)利用中(1.3)(1.4)(1.5)(1.6)将(1.3) ~ (1.5)代入1.1),即得微分方程的差分格式a iU i U i 1h*h i h i 1)d i U b^Ui21-(h i h i 1) i则逼近阶为O(h 2)。
数学中的椭圆型方程数学中,椭圆型方程是一类非常重要且广泛应用的方程类型。
它们在许多领域中起着重要作用,包括物理学、工程学、生物学和经济学等。
本文将介绍椭圆型方程的基本概念、性质和一些常见的应用。
一、椭圆型方程的定义和性质椭圆型方程是指二阶偏微分方程的一种形式,通常表示为:\[a\frac{{\partial^2 u}}{{\partial x^2}} + b\frac{{\partial^2u}}{{\partial x \partial y}} + c\frac{{\partial^2 u}}{{\partial y^2}} = f(x,y)\]其中,\(a\)、\(b\)、\(c\)是与\(u\)相关的系数,\(f(x, y)\)是已知的函数。
椭圆型方程中的二阶导数对\(u\)的贡献是正的。
椭圆型方程具有以下性质:1. 线性性质:椭圆型方程是线性的,这意味着如果\(u_1\)和\(u_2\)是该方程的解,那么\(c_1u_1 + c_2u_2\)也是该方程的解,其中\(c_1\)和\(c_2\)是常数。
2. 正定性质:椭圆型方程中的系数满足\(b^2 - 4ac < 0\)时,方程被称为正定的。
正定性质保证了方程解的唯一性和稳定性。
3. 边界条件:对于椭圆型方程,需要指定边界条件才能得到唯一解。
常见的边界条件包括Dirichlet边界条件(给定边界上的函数值)、Neumann边界条件(给定边界上的法向导数值)和Robin边界条件(给定边界上的线性组合)。
二、椭圆型方程的应用1. 热传导方程:热传导方程是一种椭圆型方程,用于描述物体中的热传导过程。
它在工程学和物理学中具有广泛应用,例如分析热交换器、传热管和材料热扩散等问题。
2. 电势方程:电势方程是一种椭圆型方程,用于描述电场中的电势分布。
它在电磁学和电子学中起着重要作用,用于分析电场和电势的分布以及导体和介质之间的电荷传输。
3. 流体力学方程:流体力学方程也可以表达为椭圆型方程的形式。
椭圆型方程的差分解法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方程五点格式,方程以较快速度迭代收缩。
精选全文完整版(可编辑修改)用差分法解椭圆型偏微分方程-(Uxx+Uyy)=(pi*pi-1)e^xsin(pi*y)0<x<2; 0<y<1U(0,y)=sin(pi*y),U(2,y)=e^2sin(pi*y); 0=<y<=1U(x,0)=0, U(x,1)=0;0=<x<=2先自己去看一下关于五点差分法的理论书籍Matlab程序:unction[p e u x y k]=wudianchafenfa(h,m,n,kmax,ep)% g-s迭代法解五点差分法问题%kmax为最大迭代次数%m,n为x,y方向的网格数,例如(2-0)/0.01=200;%e为误差,p为精确解symstemp;u=zeros(n+1,m+1);x=0+(0:m)*h;y=0+(0:n)*h;for(i=1:n+1)u(i,1)=sin(pi*y(i));u(i,m+1)=exp(1)*exp(1)*sin(pi*y(i));endfor(i=1:n)for(j=1:m)f(i,j)=(pi*pi-1)*exp(x(j))*sin(pi*y(i));endendt=zeros(n-1,m-1);for(k=1:kmax)for(i=2:n)for(j=2:m)temp=h*h*f(i,j)/4+(u(i,j+1)+u(i,j-1)+u(i+1,j)+u (i-1,j))/4;t(i,j)=(temp-u(i,j))*(temp-u(i,j));u(i,j)=temp;endendt(i,j)=sqrt(t(i,j));if(k>kmax)break;endif(max(max(t))<ep)break;endendfor(i=1:n+1)for(j=1:m+1)p(i,j)=exp(x(j))*sin(pi*y(i));e(i,j)=abs(u(i,j)-exp(x(j))*sin(pi*y(i)));endEnd在命令窗口中输入:[p e u xyk]=wudianchafenfa(0.1,20,10,10000,1e-6) k=147surf(x,y,u) ;xlabel(‘x’);ylabel(‘y’);zlabel(‘u’);Title(‘五点差分法解椭圆型偏微分方程例1’)就可以得到下图surf(x,y,p)surf(x,y,e)[p e uxy k]=wudianchafenfa(0.05,40,20,10000,1e-6)[pe u x y k]=wudianchafenfa(0.025,80,40,10000,1e-6)为什么分得越小,误差会变大呢?我们试试运行:[pe u x yk]=wudianchafenfa(0.025,80,40,10000,1e-8)K=2164surf(x,y,e)误差变小了吧还可以试试[p e ux y k]=wudianchafenfa(0.025,80,40,10000,1e-10) K=3355误差又大了一点再试试[peu x y k]=wudianchafenfa(0.025,80,40,10000,1e-11)k=3952 误差趋于稳定总结:最终的误差曲面与网格数有关,也与设定的迭代前后两次差值(ep,看程序)有关;固定网格数,随着设定的迭代前后两次差值变小,误差由大比变小,中间有一个最小值,随着又增大一点,最后趋于稳定。
椭圆型方程的高阶隐式紧致差分方法研究的开题报告题目:椭圆型方程的高阶隐式紧致差分方法研究一、研究背景椭圆型方程在现实问题中具有广泛的应用,如热传导方程、泊松方程等。
为了求解这些方程,需要采用数值方法进行近似求解。
目前常用的数值方法有有限差分法、有限元法、谱方法等。
其中,有限差分法因其简单易行、计算速度快等优点,得到了广泛的应用。
在有限差分法中,采用隐式求解方法可以保证数值解的稳定性,避免数值解中的震荡现象。
同时,使用高阶紧致差分格式可以提高数值解的精度。
因此,将隐式求解方法和高阶紧致差分格式相结合,得到高精度、稳定的数值解,具有重要的研究价值和实际应用意义。
二、研究目的和内容本课题的研究目的是发展一种高阶隐式紧致差分方法,用于求解椭圆型方程。
具体内容如下:1.分析椭圆型方程的特点,确定所采用的数值方法;2.设计高阶紧致差分格式,并应用隐式求解方法,推导数值解的公式;3.编写相应的数值程序,并验证数值方法的稳定性和精度;4.对求解结果进行分析和验证,探究数值方法的优化方法。
三、研究方法本课题的研究方法主要包括理论分析和数值计算。
具体方法如下:1.分析椭圆型方程的特点,研究其数学性质;2.根据数学性质和所需求解的精度,确定采用的高阶紧致差分格式和隐式求解方法;3.推导数值解的公式,计算数值解;4.在计算过程中,对数值方法的稳定性和精度进行验证;5.对数值计算结果进行分析和验证,探究数值方法的优化方法。
四、研究意义本课题的研究结果可以具有以下意义:1.提高椭圆型方程的数值求解精度,满足不同实际问题对求解结果的要求;2.研究开发新的求解方法,丰富数值方法的研究成果;3.发展高性能计算技术,提高计算效率和求解速度;4.为实际问题的数值求解提供有效的数值方法和数值程序。