波动方程课件.pptVIP

波动方程全文共56页，当前为第32页。 7.2 二维波动问题 7.2.1 矩形膜的振动 矩形膜的四周边界固定，膜的初位移已知。定解问题是： 显式差分格式：(取 ) 波动方程全文共56页，当前为第33页。 %ex621; (p183) 二维矩域波动方程的差分解 clear; N=80; M=40; K=500; a=2; b=1; c=1; A=1; dx=a/N; dy=b/M; dt=0.01; t=dt*(0:K); C=(c*dt/dx)^2; x=dx*(0:N); y=dx*(0:M); [Y,X]=meshgrid(y,x); u1=zeros(N+1,M+1); u2=u1; u3=u1; I=2:N; J=2:M; u1=A*X.*(X-a).*sin(pi*Y/b); %初位移； u2(I,J)=u1(I+1,J)+u1(I,J+1)-4*u1(I,J)+u1(I-1,J)+u1(I,J-1); u2(I,J)=u1(I,J)+C/2*u2(I,J); figure(1);mesh(X,Y,u1); title(二维波动);axis([0 2 0 1 -1 1]);p(1)=getframe; for k=2:K+1; u3(I,J)=u2(I+1,J)+u2(I,J+1)-4*u2(I,J)+u2(I-1,J)+u2(I,J-1); u3(I,J)=2*u2(I,J)-u1(I,J)+C*u3(I,J); u1=u2; u2=u3; if mod(k,5)==0; mesh(X,Y,u3);axis([0 2 0 1 -1 1]); p(k/5+1)=getframe; end; end; movie(p); 波动方程全文共56页，当前为第34页。 波动方程全文共56页，当前为第35页。 7.2.2 圆膜的振动 半径为ρ0的圆膜，边界固定，膜的初始形状为旋转抛物面，而初始速度为0，求膜的振动。定解问题是： 显式差分格式：(取 ) 波动方程全文共56页，当前为第36页。 %ex622; (p186) 二维圆域波动方程的差分解 clear; N=90; K=500; a=1; r0=1; u0=0.4; dr=r0/N; dt=0.01; t=dt*(0:K); C=(a*dt/dr)^2; r=dr*(0:N); fai=2*pi*(0:N)/N; [R,Fai]=meshgrid(r,fai);[X,Y]=pol2cart(Fai,R); u1=zeros(1,N+1); u2=u1; u3=u1; I=2:N; u1=u0*(1-(r/r0).^2); %初位移； u2(I)=u1(I+1)-2*u1(I)+u1(I-1)+(u1(I+1)-u1(I-1))/(2*(I-1)); u2(I)=u1(I)+C/2*u2(I); u2(1)=u2(2); figure(1);mesh(X,Y,ones(N+1,1)*u1); title(二维波动);axis([-r0 r0 -r0 r0 -.5 .5]);p(1)=getframe; for k=2:K+1; u3(I)=u2(I+1)-2*u2(I)+u2(I-1)+(u2(I+1)-u2(I-1))/(2*(I-1)); u3(I)=2*u2(I)-u1(I)+C*u3(I); u3(1)=u3(2); u1=u2; u2=u3; if mod(k,5)==0; mesh(X,Y,ones(N+1,1)*u3); axis([-r0 r0 -r0 r0 -.5 .5]);p(k/5+1)=getframe; end;end; movie(p); 波动方程全文共56页，当前为第37页。 波动方程全文共56页，当前为第38页。 半径为ρ0的匀质圆柱，高为L，上、下底面固定，侧面自由，初始位移为零, 初始速度为u0ρ2z, 求柱体内的振动情况。 定解问题为： 其解析解为： 7.3 三维波动问题 7.3.1 柱体内的振动 波动方程全文共56页，当前为第39页。 显式差分格式：(取 ) 波动方程全文共56页，当前为第40页。 %ex615; (p187) 柱体内振