MATLAB计算四类数学物理方程的举例求解ppt课件.ppt
《MATLAB计算四类数学物理方程的举例求解ppt课件.ppt》由会员分享,可在线阅读,更多相关《MATLAB计算四类数学物理方程的举例求解ppt课件.ppt(172页珍藏版)》请在三一办公上搜索。
1、数学物理建模与计算机辅助设计,第5章 四类数学物理方程的求解举例,本章内容,5.1 求解本征值型数学物理方程 5.2 求解稳定型数学物理方程5.3 求解热传导型数学物理方程5.4 求解波动型数学物理方程,本征值问题简介,用分离变量法解数学物理方程的时候,最后都会归结到求解本征值问题上去。在利用本征函数系展开法求解数学物理方程的时候,需要对所用的本征函数系有较好的理解在本部分中讲介绍一维和二维本征值问题和相应的本征函数系,5.1 求解本征值型的数学物理方程,一维本征值问题,一维波动方程 令则分离变量可得,5.1 求解本征值型的数学物理方程,一维本征值问题,一维热传导方程令分离变量可得,5.1 求
2、解本征值型的数学物理方程,一维本征值问题,二维拉普拉斯方程令分离变量可得,5.1 求解本征值型的数学物理方程,一维本征值问题,特点:分离变量后的方程都是一元一阶或二阶常微分方程一阶常微分方程表明变量是衰减的,具有阻尼,即解中具有e的负指数因子二阶常微分方程表明变量是振荡的,将具有一次或多次倍频的正弦或余弦的形式,即为本征函数。这儿的频率既可以指空间频率,也可指时间频率。下面讲解二阶常微分方程获得的本征函数系,5.1 求解本征值型的数学物理方程,一维本征值问题的本征函数系,本征函数系求解本征值问题注意到这里的两个边界上都是第一类齐次边界条件。其本征值和本征函数分别是,5.1 求解本征值型的数学物
3、理方程,一维本征值问题的本征函数系,本征函数系求解本征值问题注意到这里的两个边界上都是第二类齐次边界条件。其本征值和本征函数分别是,5.1 求解本征值型的数学物理方程,一维本征值问题的本征函数系,本征函数系求解本征值问题注意到这里的两个边界中左边界是第一类边界条件,右边界上是第二类齐次边界条件。其本征值和本征函数分别是,5.1 求解本征值型的数学物理方程,一维本征值问题的本征函数系,本征函数系求解本征值问题注意到这里的两个边界中左边界是第二类边界条件,右边界上是第一类齐次边界条件。其本征值和本征函数分别是,5.1 求解本征值型的数学物理方程,本征函数系图像,x=0:0.001:1; (p65)
4、A=sin(pi*1:4*x);B=cos(pi*(0:3)*x);C=sin(pi*(1/2:7/2)*x);D=cos(pi*(1/2:7/2)*x);subplot(4,1,1); plot(x,A);subplot(4,1,2); plot(x,B);subplot(4,1,3); plot(x,C);subplot(4,1,4); plot(x,D);,采用了矢量化技术,5.1 求解本征值型的数学物理方程,本征函数系运动图像(振动问题),function zb %存储于文件zb.m中 (p65)t=0:0.005:2.0; x=0:0.001:1;ww1=wfun(1,0); ww2
5、=wfun(2,0);ww3=wfun(3,0); ww4=wfun(4,0);ymax1=max(abs(ww1);figuresubplot(2,1,1); h1=plot(x,ww1,r,linewidth,5); hold on;h3=plot(x,ww3,g,linewidth,5); axis(0 1 ymax1 ymax1);ymax4=max(abs(ww4);subplot(2,1,2); h2=plot(x,ww2,r,linewidth,5); hold on;h4=plot(x,ww4,g,linewidth,5); axis(0 1 ymax4 ymax4);,5.1
6、 求解本征值型的数学物理方程,本征函数系运动图像(振动问题),for n=2:length(t) ww1=wfun(1,t(n); set(h1,ydata,ww1); ww3=wfun(3,t(n); set(h3,ydata,ww3); drawnow ww2=wfun(2,t(n); set(h2,ydata,ww2); ww4=wfun(4,t(n); set(h4,ydata,ww4); drawnowendfunction wtx=wfun(k,t) %调用函数存储于wfun.m中x=0:0.001:1; a=1;wtx=cos(k*pi*a*t).*sin(k*pi*x);,5.
7、1 求解本征值型的数学物理方程,本征函数系运动图像(热传导问题),function sfb %存放于sfb.m中 (p67)x=0:0.01:1; t=1e-5:0.0001:0.005;ww1=sfbf(1,t(1); ww2=sfbf(2,t(1);ww3=sfbf(3,t(1); ww4=sfbf(4,t(1);ymax3=max(abs(ww3);figuresubplot(2,1,1); h1=plot(x,ww1,r,linewidth,3); hold on;h3=plot(x,ww3,g,linewidth,3); axis(0 1 ymax3 ymax3);ymax4=max
8、(abs(ww4);subplot(2,1,2); h2=plot(x,ww2,r,linewidth,3); hold on;h4=plot(x,ww4,g,linewidth,3); axis(0 1 ymax4 ymax4);,5.1 求解本征值型的数学物理方程,本征函数系运动图像(热传导问题),for n=2:length(t) ww1=sfbf(1,t(n); set(h1,ydata,ww1); ww3=sfbf (3,t(n); set(h3,ydata,ww3); drawnow ww2=sfbf (2,t(n); set(h2,ydata,ww2); ww4=sfbf (4,
9、t(n); set(h4,ydata,ww4); drawnowendfunction cht=sfbf (k,t) %调用函数存储于wfun.m中x=0:0.01:1;cht=sin(k*pi*x).*exp(-k2*pi2*22*t);,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,矩形区域的本征模与本征振动边长为b和c的的四周固定的矩形本征模的本征值问题为采用分离变量法可以得到本征模和本征值为,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,下面绘制前4个本征函数的图形 (p70_1.m)a=2; b=1;m,n=meshgrid(1:3);L=(m
10、*pi./b).2+(n*pi./b).2);x=0:0.01:a; y=0:0.01:b;X,Y=meshgrid(x,y);w11=sin(pi*Y./b).*sin(pi*X./a); w12=sin(2*pi*Y./b).*sin(pi*X./a); w21=sin(pi*Y./b).*sin(2*pi*X./a);w22=sin(pi*Y./b).*sin(3*pi*X./a);figuresubplot(2,2,1); mesh(X,Y,w11); subplot(2,2,2); mesh(X,Y,w12);subplot(2,2,3); mesh(X,Y,w21); subplo
11、t(2,2,4); mesh(X,Y,w22);,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,下面绘制前4个本征函数的运动图形:波动图形波动问题中的时间因子的形式为这里的不妨取其中的正弦项加以展示,function p71_1 %p71_1.mb=2; c=1;x=0:0.02:b;y=0:0.02:c;X,Y=meshgrid(x,y);Z=zeros(51,51);p=moviein(2*3*60);for m=1:2 for n=1:3 for i=1:60 a=sqrt(m*pi/c).2+(n*p
12、i/b).2); Z=sin(a*i*.02*pi)*sin(m*pi*Y./c).*sin(n*pi*X./b); mesh(X,Y,Z) t=本征振动:,m=,int2str(m), n=,int2str(n); title(t); axis(0 b 0 c -1 1); p(:,(m-1)*3+(n-1)*60+i)=getframe; end endendMOVIE2AVI(p,D:A.avi),5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域
13、情况,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,5.1 求解本征值型的数学物理方程,二维本征值问题:矩形区域情况,5.1 求解本征值型的数学物理方程,二维本征值问题:圆形区域情况,半径为 周边固定的圆形膜的本征问题是本征模是其中下面讨论第一类其次边界条件下的本征值情形,5.1 求解本征值型的数学物理方程,二维本征值问题:圆形区域情况,m=0,n=0,1,2,3的贝塞尔函数的前4个本征函数图形r=0:0.01:1; %(p71_1.m)y=inline(besselj(0,x),x);x(1)=fzero(
14、y,1);plot(r,y(x(1)*r)hold onfor k=1:3 x(k+1)=fzero(y,x(k)+3); plot(r,y(x(k+1)*r);end,5.1 求解本征值型的数学物理方程,二维本征值问题:圆形区域情况,圆域内m=0,n=0,1,2,34个本征模图形 (p75_1.m)XX,YY=meshgrid(-1:0.01:1);Q,R=cart2pol(XX,YY);R(find(R1)=NaN;y=inline(besselj(0,x),x);x(1)=fzero(y,3);meshc(XX,YY,y(x(1)*R)for k=1:3 x(k+1)=fzero(y,x
15、(k)+2); figure meshc(XX,YY,y(x(k+1)*R);end,5.1 求解本征值型的数学物理方程,二维本征值问题:圆形区域情况,圆域内m=1,n=0,1,2,3的4个本征模图形 (p76_1.m)XX,YY=meshgrid(-1:0.01:1);Q,R=cart2pol(XX,YY);R(find(R1)=NaN;y=inline(besselj(1,x),x);x(1)=fzero(y,3);meshc(XX,YY,y(x(1)*R ).*sin(Q)for k=1:3 x(k+1)=fzero(y,x(k)+2); figure meshc(XX,YY,y(x(k
16、+1)*R ).*sin(Q);end,5.1 求解本征值型的数学物理方程,二维本征值问题:圆形区域情况,圆域内m=2,n=0,1,2,3的4个本征模图形XX,YY=meshgrid(-1:0.01:1);Q,R=cart2pol(XX,YY);R(find(R1)=NaN;y=inline(besselj(2,x),x);x(1)=fzero(y,4);meshc(XX,YY,y(x(1)*R ).*sin(2*Q)for k=1:3 x(k+1)=fzero(y,x(k)+2); figure meshc(XX,YY,y(x(k+1)*R ).*sin(2*Q);end,5.1 求解本征值
17、型的数学物理方程,二维本征值问题:圆形区域情况,用PDE工具箱求解步骤命令窗口中键入pdetool进入PDE工具箱界面以原点为中心绘一个半径为1的圆选择方程类型为Eigenmodes,参数c=1,a=0,d=1边界条件设置为固定(默认)细分网格设置求解参数中本征值搜索范围为0100打开作图对话框选择需要作图的类型和本征值,并作图,二维本征值问题:球函数的图形,实数形式的球函数 其中 是连带勒让德函数。复数形式的球函数,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,球函数的图形球函数是在球面上的二元函数,它是从球函数方程的本征问题中得到的本征函数系。球函数的图形是空间图形,为
18、了画出其图形,必须指定球的半径对复数形式的球函数,必须对其实部和虚部分别作图可以用角变量作为平面上的x、y轴的变量画出球函数图,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,l=3; m=2; R=4; A=3; delta=pi/40; %(p81_1.m)theta0=0:delta:pi; phi=0:2*delta:2*pi;phi,theta=meshgrid(phi, theta0);Ymn=legendre(l,cos(theta0); Ymn=Ymn(m+1,:);L=size(theta,1); yy=repmat(Ymn,1,L);Reyy=yy.*co
19、s(m*phi); Imyy=yy.*sin(m*phi);ReM=max(max(abs(Reyy);Rerho=R+A*Reyy/ReM;Rer=Rerho.*sin(theta); Rex=Rer.*cos(phi);Rex=Rer.*sin(phi); Rez=Rerho.*cos(theta);,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,subplot(1,2,1); surf(Rex,Rey,Rez); light; lighting phong;axis(square); axis(-5 5 -5 5 -5 5); axis (off); view(40,
20、30)title(实球谐函数);ImM=max(max(abs(Imyy);Imrho=R+A*Imyy/(ImM+eps*(ImM=0);Imr=Imrho.*sin(theta); Imx=Imr.*cos(phi);Imx=Imr.*sin(phi); Imz=Imrho.*cos(theta);subplot(1,2,2); surf(Imx,Imy,Imz); light; lighting phong;axis(square); axis(-5 5 -5 5 -5 5); axis (off); view(40,30)title(虚球谐函数);,5.1 求解本征值型的数学物理方程,
21、二维本征值问题:球函数的图形,figuresubplot(3,1,1); surf(theta,phi,Reyy);xlabel(theta); ylabel(phi);subplot(3,1,2); surf(theta,phi,Imyy);xlabel(theta); ylabel(phi);subplot(3,1,3); surf(theta,phi,(Reyy.2+Imyy.2);xlabel(theta); ylabel(phi);,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,5.1 求解本征值型
22、的数学物理方程,二维本征值问题:球函数的图形,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,5.1 求解本征值型的数学物理方程,二维本征值问题:球函数的图形,5.1 求解本征值型的数学物理方程,二维本征值问题:拉普拉斯方程,矩形区域的拉普拉斯方程求解方法:利用PDE工具箱,5.1 求解本征值型的数学物理方程,二维本征值问题:拉普拉斯方程,采用Heat transfer类型的时候,此类方程属于椭圆型方程采用Generic scalar类型的时候,此类方程为,5.1 求解本征值型的数学物理方程,二维本征值问题:拉普拉斯方程,5.1 求解本征值型的数学物理方程,二维拉普拉斯方程,
23、问题1:阳光照射的圆柱的温度分布 半径为a,表面熏黑的均匀长圆柱,在温度为零的空气中受阳光的照射,阳光垂直于柱轴,热流强度为q,试求柱内稳定的温度分布定解问题是其中, 5.2 求解稳定型的数学物理方程,二维拉普拉斯方程,问题的解析解为a=1; H=1.5; k=0.3; h=0.2; q=5; N=20; %(p91_1.m)p=0:0.05:1; theta=0:pi/30:2*pi;Th,P=meshgrid(theta,p); X,Y=pol2cart(Th,P);z=0;u0=q/H/pi+1/(k+H*a)*q/2.*P.*sin(Th);for n=1:N zz=2*q/pi/a(
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- MATLAB 计算 数学 物理 方程 举例 求解 ppt 课件
链接地址:https://www.31ppt.com/p-2002398.html