有限元matlab_精品文档.ppt

1、基于基于MatlabMatlab的偏微分的偏微分方程数值解方程数值解求数值解方法求数值解方法差分方法差分方法有限元方法有限元方法MATLABMATLAB的的pedpepedpe函数函数MATLABMATLAB的的PDEtoolPDEtool工具箱工具箱 偏微分方程分类偏微分方程分类椭圆偏微分方程椭圆偏微分方程双曲线偏微分方程双曲线偏微分方程抛物线偏微分方程抛物线偏微分方程4 椭圆偏微分方程特例椭圆偏微分方程特例拉普拉斯方程拉普拉斯方程 拉普拉斯方程是最简单的椭圆偏微分方程，以下以拉普拉斯方程为例，讲述椭圆偏微分方程的的数值解法。拉普拉斯方程形式如下：5 椭圆偏微分方程边界条件椭圆偏微分方程边界

2、条件椭圆偏微分方程边界条件有以下三种提法：其中第一种提法最为普遍，下面以第一种边界条件，拉普拉斯方程为例介绍椭圆偏微分方程常用的五点差分格式和工字型差分格式的解法。6 五点差分格式五点差分格式五点差分格式最常用的格式，其形式如下：注意：这里的边界为矩形区域。注意：这里的边界为矩形区域。7 五点差分格式算法五点差分格式算法 注意注意：要保证：要保证x方向和方向和y方向上的网格步长相等才方向上的网格步长相等才能使用上面的公式。能使用上面的公式。8 五点差分格式在五点差分格式在MATLABMATLAB中实现中实现function u=peEllip5(nx,minx,maxx,ny,miny,max

3、y)%x方向的节点数：nx%求解区间x的左端：minx%求解区间x的右端：maxx%y方向的节点数：ny%求解区间y的左端：miny%求解区间y的右端：maxy%求解区间上的数值解：uformat long;hx=(maxx-minx)/(nx-1);hy=(maxy-miny)/(ny-1);u0=zeros(nx,ny);for j=1:ny u0(j,1)=EllIni2Uxl(minx,miny+(j-1)*hy);u0(j,nx)=EllIni2Uxr(maxx,miny+(j-1)*hy);endfor j=1:nx u0(1,j)=EllIni2Uyl(minx+(j-1)*hx

4、,miny);u0(ny,j)=EllIni2Uyr(minx+(j-1)*hx,maxy);end%边界条件的离散9 五点差分格式在五点差分格式在MATLABMATLAB中实现中实现A=-4*eye(nx-2)*(ny-2),(nx-2)*(ny-2);b=zeros(nx-2)*(ny-2),1);for i=1:(nx-2)*(ny-2);if mod(i,nx-2)=1 if i=1 A(1,2)=1;A(1,nx-1)=1;b(1)=-u0(1,2)-u0(2,1);else if i=(ny-3)*(nx-2)+1 A(i,i+1)=1;A(i,i-nx+2)=1;%注意边界节点的

5、离散方式 b(i)=-u0(ny-1,1)-u0(ny,2);else A(i,i+1)=1;A(i,i-nx+2)=1;A(i,i+nx-2)=1;b(i)=-u0(floor(i/(nx-2)+2,1);end end else if mod(i,nx-2)=0 if i=nx-210 五点差分格式在五点差分格式在MATLABMATLAB中实现中实现 A(i,i-1)=1;%注意边界节点的离散方式 A(i,i+nx-2)=1;b(i)=-u0(1,nx-1)-u0(2,nx);else if i=(ny-2)*(nx-2)A(i,i-1)=1;A(i,i-nx+2)=1;b(i)=-u0(

6、ny-1,nx)-u0(ny,nx-1);else A(i,i-1)=1;A(i,i-nx+2)=1;A(i,i+nx-2)=1;b(i)=-u0(floor(i/(nx-2)+1,nx);end end else if i1&i(ny-3)*(nx-2)&i0,a0格式不同是为了满足差分格式的稳定性，若第一个式子a0 for j=1:(n+M)u0(j)=IniU(minx+(j-M-1)*h);%向左延拓M个节点的函数值 endelse for j=1:(n+M)u0(j)=IniU(minx+(j-1)*h);%向左延拓M个节点的函数值 endendu1=u0;for k=1:M if

7、a0 for i=(k+1):n+M u1(i)=-dt*a*(u0(i)-u0(i-1)/h+u0(i);end else for i=1:n+M-k u1(i)=-dt*a*(u0(i+1)-u0(i)/h+u0(i);end 一维对流方程一维对流方程迎风格式迎风格式20 一维对流方程一维对流方程迎风格式算例迎风格式算例end u0=u1;endif a0 u=u1(M+1):M+n);else u=u1(1:n);endformat long;21 一维对流方程一维对流方程迎风格式算例迎风格式算例 然后在MATLAB窗口输入下列命令：u=peHypbYF(1,0.005,101,0,1,

8、100);22 一维对流方程一维对流方程迎风格式算例结果迎风格式算例结果t=0时，时，u,x分布图分布图t=0.5时，时，u,x分布图分布图23 一维对流方程一维对流方程拉克斯拉克斯-弗里德里希斯格式弗里德里希斯格式24 一维对流方程一维对流方程拉克斯拉克斯-弗里德里希斯格式弗里德里希斯格式25 拉克斯拉克斯-弗里德里希斯格式算例弗里德里希斯格式算例26 拉克斯拉克斯-弗里德里希斯格式算例弗里德里希斯格式算例27 拉克斯拉克斯-弗里德里希斯格式算例弗里德里希斯格式算例28 一维对流方程一维对流方程拉克斯拉克斯-温德洛夫格式温德洛夫格式29 一维对流方程一维对流方程拉克斯拉克斯-温德洛夫格式温德

9、洛夫格式30 一维对流方程一维对流方程拉克斯拉克斯-温德洛夫格式算例温德洛夫格式算例31 一维对流方程一维对流方程拉克斯拉克斯-温德洛夫格式算例温德洛夫格式算例32 一维对流方程一维对流方程拉克斯拉克斯-温德洛夫格式算例温德洛夫格式算例33 一维对流方程一维对流方程比姆比姆-沃明格式沃明格式34 一维对流方程一维对流方程比姆比姆-沃明格式沃明格式35 一维对流方程一维对流方程比姆比姆-沃明格式算例沃明格式算例36 一维对流方程一维对流方程多步格式多步格式 多步格式也有多种，这里只简单介绍其中几种格式。包括多步格式也有多种，这里只简单介绍其中几种格式。包括Richtmyer多多步格式、拉克斯步格

10、式、拉克斯-温德洛夫多步格式、温德洛夫多步格式、MacCormack多步格式。多步格式。37 一维对流方程一维对流方程多步格式多步格式 38 一维对流方程一维对流方程多步格式算例多步格式算例 Richtmyer多步格式算出的结果并不理多步格式算出的结果并不理想，不但左边有波动，而且光滑性也不好。拉想，不但左边有波动，而且光滑性也不好。拉克斯克斯-温德洛夫多步格式算出的结果比较不错，温德洛夫多步格式算出的结果比较不错，虽然左边有点小波动，但是初始函数的宽度和虽然左边有点小波动，但是初始函数的宽度和高度都保持的不错。高度都保持的不错。MacCormack多步格式多步格式求得的结果和拉克斯求得的结果

11、和拉克斯-温德洛夫多步格式算出温德洛夫多步格式算出的结果差不多。的结果差不多。39 双曲线偏微分方程双曲线偏微分方程二维对流方程二维对流方程40 二维对流方程二维对流方程拉克斯拉克斯-弗里德里希斯格式弗里德里希斯格式41 二维对流方程二维对流方程拉克斯拉克斯-弗里德里希斯格式弗里德里希斯格式42 二维对流方程二维对流方程拉克斯拉克斯-弗里德里希斯格式算例弗里德里希斯格式算例u=peHypb2LF(1,1,0.005,101,0,1,101,0,1,100);43 二维对流方程二维对流方程拉克斯拉克斯-弗里德里希斯格式算例弗里德里希斯格式算例 结果与初始值对比，可以看出，拉克斯结果与初始值对比，

12、可以看出，拉克斯-弗里德里希斯格式算出的弗里德里希斯格式算出的结果非常好。结果非常好。44 二维对流方程二维对流方程近似分裂格式近似分裂格式近似分裂格式也是一种不错的格式，其结果也非常接近理论值。近似分裂格式也是一种不错的格式，其结果也非常接近理论值。45 抛物线偏微分方程抛物线偏微分方程扩散方程扩散方程 在实际应用中遇到的抛物线偏微分方程主要是扩散方程。扩散方程有很强在实际应用中遇到的抛物线偏微分方程主要是扩散方程。扩散方程有很强的物理背景，例如不用物质之间的扩散过程、热传递过程、波传播等过程都可的物理背景，例如不用物质之间的扩散过程、热传递过程、波传播等过程都可以用扩散过程来描述。下面以扩

13、散方程为例介绍几种差分格式。以用扩散过程来描述。下面以扩散方程为例介绍几种差分格式。46 扩散方程扩散方程显式格式显式格式47 扩散方程扩散方程显式格式显式格式48 扩散方程扩散方程显式格式算例显式格式算例u=peParabExp(1,0.005,101,0,1,100)49 扩散方程扩散方程显式格式算例显式格式算例 我们知道显示格式虽然简单，但其精度很差，而且求得的解容易出现震荡。我们知道显示格式虽然简单，但其精度很差，而且求得的解容易出现震荡。次算例的结果如下：次算例的结果如下：数值结果震荡的非常厉害，说明显式格式在这种条件下不稳定。数值结果震荡的非常厉害，说明显式格式在这种条件下不稳定。

14、50 扩散方程扩散方程跳点格式跳点格式相同的算例数值结果同样震荡的非常厉害，说明跳点格式在这种条件下也不稳定。相同的算例数值结果同样震荡的非常厉害，说明跳点格式在这种条件下也不稳定。51 扩散方程扩散方程隐式格式隐式格式52 扩散方程扩散方程隐式格式隐式格式53 扩散方程扩散方程隐式格式隐式格式54 扩散方程扩散方程隐式格式算例隐式格式算例用隐式格式求解下面扩散方程的初值问题：u=peParabImp(1,0.005,101,0,1,0,1,100)55 扩散方程扩散方程隐式格式算例隐式格式算例结果是一条比较稳定的直线。结果是一条比较稳定的直线。56 扩散方程扩散方程克拉克克拉克-尼科尔森格式

15、尼科尔森格式57 扩散方程扩散方程克拉克克拉克-尼科尔森格式尼科尔森格式58 扩散方程扩散方程克拉克克拉克-尼科尔森格式尼科尔森格式59 扩散方程扩散方程克拉克克拉克-尼科尔森格式算例尼科尔森格式算例60 扩散方程扩散方程克拉克克拉克-尼科尔森格式算例尼科尔森格式算例61 扩散方程扩散方程加权隐式格式加权隐式格式 通过算例我们可以知道，通过取不同参数，用加权隐式格式算得的结果通过算例我们可以知道，通过取不同参数，用加权隐式格式算得的结果差别不大，只是达到稳态的时间不一样而已。差别不大，只是达到稳态的时间不一样而已。62 差分方法小结差分方法小结 以上我们介绍了差分方法在椭圆型、抛物型和双曲型偏

16、微分方程中的应用，用差分格式求解偏微分方程的基本步骤是一样的，首先把连续的问题离散化，建立差分格式，然后根据差分格式对求解区域进行网格剖分，最后求解方程。下面将简单介绍有限元方法。另外变分法、边界元法、混合有限元法和多重网格法等也是偏微分方程数值求解方法，有兴趣的同学可以参考相关书籍。63 有限元方法有限元方法介绍介绍 有限元方法是数值求解偏微分方程边值问题的一种方法，此方法首先于20世纪50年代初由工程师提出，并用于求解简单的结构问题。有限元方法是这一种系统的数值方法，并奠定其数学基础，是在60年代中期以冯康先生为代表的中国学者与西方学者独立并行完成的。有限元方法不同于差分方法，主要有以下三大特点：(1)从数学物理问题的变分原理出发，而不是从微分方程出发，因此事从问题的整体描述而不是从问题的局部描述出发。(2)对所考虑问题的区域(以二维情形为例)作三角形（或其他简单多边形）剖分，而不是仅作矩形剖分。(3)用剖分区域上的简单函数(如分片多项式)去逼近原问题的解，而不是只在剖分节点上的数值逼近。64 有限元方法有限元方法一维边值问题算例一维边值问题算例用有限元方法求解如下一维边值问题：解