上机练习专题5：数学物理方程的计算机求解和可视化

1、数学物理建模与计算机辅助设计上机练习专题5：数学物理方程的计算机求解和可视化Page 2对代数表达式进行可视化二维本征值问题矩形区域的本征模与本征振动边长为b和c的的四周固定的矩形膜振动的本征值问题为采用分离变量法可以得到本征模和本征值为Page 3对代数表达式进行可视化二维本征值问题绘制前4个本征函数的图形%P70_1.ma=2; b=1;m,n=meshgrid(1:3);L=(m*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*

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); subplot(2,2,4); mesh(X,Y,w22);Page 4对代数表达式进行可视化二维本征值问题动态可视化：设具有时间因子%P71_1.mb=2; c=1; x=0:0.02:b;y=0:0.02:c; X,Y=

3、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*pi/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.

4、avi)Page 5PDE工具箱求解数学物理方程PDE：Partial Differential EquationPDE工具箱功能设定二维定解区域、边界条件以及方程的形式和系数用有限元法(FEM)生成网格、方程离散化并求出数值解对所求解的可视化，图像和动画PDE工具箱可以求解的4种偏微分方程椭圆型方程抛物型方程双曲形方程特征值方程Page 6PDE工具箱求解数学物理方程PDE工具箱给定的两种边界条件Dilichlet条件Neumann条件用PDE工具箱求解二维本征值问题Page 7PDE工具箱求解数学物理方程求解二维本征值问题PDE工具箱求解数学物理方程1.求解定解问题Page 8取PDE工具

5、箱求解数学物理方程2.求矩形薄膜的横向振动。定解问题为Page 9取Page 10差分法求解数学物理方程导数的差分公式在x附近将函数f(x)展开成泰勒公式前差公式后差公式二阶导数的差分公式(中心差分公式)Page 11差分法求解数学物理方程用数值方法求解微分方程的时候步骤：(1)变量空间网格化(2)微分方程变差分方程(3)定解条件设置(4)整理得到递推公式Page 12差分法求解数学物理方程一维波动问题两端固定的弦的振动问题设定解条件为解的代数表达式为Page 13差分法求解数学物理方程两端固定的弦的振动问题的求解步骤(1)变量空间网格化(2)微分方程变差分方程Page 14差分法求解数学物理

6、方程两端固定的弦的振动问题的求解步骤(3)定解条件设置思考：1. 求时间分量的递推公式需要几个初始值？2. 对坐标分量呢?3. 两个初始值时，另一个初始值如何处理？Page 15差分法求解数学物理方程两端固定的弦的振动问题的求解步骤(4)整理得到递推公式Page 16差分法求解数学物理方程%P164_1.mN=4001; dx=0.0024;dt=0.0005; c=dt*dt/dx/dx;x=linspace(0,1,420);u(1:420,1)=0;u(181:240,1)=0.05*sin(pi*x(181:240)*7);u(2:419,2)=u(2:419,1)+c/2*(u(3:

7、420,1)-2*u(2:419,1)+u(1:418,1);h=plot(x,u(:,1),linewidth,3);axis(0,1,-0.05,0.05);set(h,erasemode,xor,markersize,18);for k=2:N set(h,XData,x,YData,u(:,2); drawnow; u(2:419,3)=2*u(2:419,2)-u(2:419,1)+c*(u(3:420,2)-2*u(2:419,2)+u(1:418,2); u(2:419,1)=u(2:419,2); u(2:419,2)=u(2:419,3);end两端固定的弦的振动问题程序的编制关键在于：(1) 正确的递推公式(2) 正确设定定解条件Page 17本专题小结直接对代数表达式求解并进行