ch4 本征函数系与本征振动

实际上对于波动方程和热传导方程，还要考虑时间本征函数：
1. 对于波动方程，时间函数为cos(nπa t/l)或cos[(n+1/2)πa t/l) 2. 对于热传导方程，时间函数为exp(-λa2t) (文中该处p66有误)
考虑时间因子后的动画分别如程序ex302和ex303所示。 显然，对于波动方程，将形成驻波，固定端(第一类边界条件) 为波节，开口端(第二类边界条件)为波腹。 对于热传导方程，若两端温度相同(第一类边界条件)，对于任 何初始温度分布，最终都将趋于平衡态，
在波动问题中，还要考虑时间因子：T t sin kmnat 0
相应的动画演示程序为ex305
ex304
%ex305(p71) a0=1; a=2;b=1; x=0:0.05:a; y=0:0.025:b; [X,Y]=meshgrid(x,y); for m=1
for n=1:3 for i=1:50 k=sqrt((m*pi/a)^2+(n*pi/b)^2); t=i*0.02; Z=sin(k*a0*t)*sin(n*pi*X/a).*sin(m*pi*Y/b); T=['本征振动:','m=',int2str(m),',n=',int2str(n)]; mesh(X,Y,Z); axis([0,a,0,b,-1,1]); title('T'); p(:,i)=getframe; end; movie(p);

nm

kn2m

nπ
/
a2

mπ
/ b2;
其前四个模式如下页图所示。程序如下：
%ex304(p70) a=2;b=1; [m,n]=meshgrid(1:3); L=((n*pi/a).^2+(m*pi/b).^2) x=0:0.05:2; y=0:0.025:1; [X,Y]=meshgrid(x,y); w1=sin(pi*X/a).*sin(pi*Y/b); w2=sin(pi*X/a).*sin(2*pi*Y/b); w3=sin(2*pi*X/a).*sin(pi*Y/b); w4=sin(3*pi*X/a).*sin(pi*Y/b); figure(1); subplot(2,2,1);mesh(X,Y,w1); subplot(2,2,2);mesh(X,Y,w2); subplot(2,2,3);mesh(X,Y,w3); subplot(2,2,4);mesh(X,Y,w4);
4.2 二维本征值问题
4.2.1 矩形区域的本征模与本征振动
具有固定边界的矩形区域拉普拉斯方程的本征值问题：

2u


2u x2

2u y 2

ux,y;
u0, y ua, y ux,0 ux,b 0;
可通过分离变量求得本征函数：
令ux,y X xY y
%ex303(p67) function sfb x=0:0.01:1; t=10^-5*(0:10:2000); A=sf(1,t(1)); B=sf(2,t(1)); C=sf(3,t(1)); D=sf(4,t(1)); figure(1);subplot(2,1,1);
h1=plot(x,A,'r');hold on; h3=plot(x,B,'g');grid on; y1=max(abs(A));axis([0,1,-y1,y1]);

0;
X
l


0;


n 1/ 22 π2
/l2;
n 0,1,...
4. cos[(n+1/2)πx/l]
X X 0
X x cosn 1/ 2πx / l

X
0

0;
X
l


0;


n 1/ 22 π2
/l2;
n 0,1,...
r
r2
u r
1
sin


sin

u


1
sin 2
2u

2


k 2u;
可通过分离变量求得本征函数： 令ur,, RrY,

r2R 2rR

R

Y
1 sin



a2u
f
边值问题
热传导方程 (parabolic)：
u a2u f t
波动方程 (hyperbolic):
2u a2u f 2t
初边值问题
4.1 一维本征值问题
对于齐次问题，f=0，在一(二)维情况下，三种 方程分别变为：
二维拉普拉斯方程 uxx uyy 0 (0 x a;0 y b)
um,n , Jm kmn r m 其中kmn xnm / 0; xnm为m阶贝塞尔函数
的第n个零点
%ex106(p32-33) ; %贝塞耳(Bessel)函数\诺伊曼(Neumann)函数 clear; k=(0:3)'; m=20; x=.05:0.1:m; figure(1);subplot(2,1,1);y=besselj(k,x);
end; end;
4.2.2 园形区域的本征模与本征振动
Leabharlann Baidu
具有固定边界的圆形区域拉普拉斯方程的本征值问题：

2u


1




u


1
2
2u
2
k 2u;
u0 0;
可通过分离变量求得本征函数：令u, R


该方程与各种边值条件结合可以组成各种本征函数系
4.1.1 四种常见的正、余弦本征函数系
1.
sin(nπx/l)
X X 0
X x sin nπx / l

X
0

0; X l 0;
nπ / l2;
n
1,2,...
2.
cos(nπx/l)
X X 0
X x cosnπx / l

X
0

0;
X l
0;



nπ
/ l2;
n

0,1,...
3. sin[(n+1/2)πx/l]
X X 0
X x sin n 1/ 2πx / l

X
0
grid on; legend(‘k=1’,‘k=3’); title('驻波k=1、3') ; ym1=max(abs(w1));axis([0,1,-ym1,ym1]); subplot(2,1,2);h2=plot(x,w2,'r');hold on; h4=plot(x,w4,'g'); grid on; legend('k=2','k=4'); title('驻波k=2、4') ; ym4=max(abs(w4));axis([0,1,-ym4,ym4]) ; for n=2:length(t) w1=wf(1,t(n)); w2=wf(2,t(n)); w3=wf(3,t(n)); w4=wf(4,t(n)); set(h1,'ydata',w1); set(h3,'ydata',w3); drawnow; set(h2,'ydata',w2); set(h4,'ydata',w4); drawnow; end; function p=wf(k,t) x=0:0.001:1; a=1; p=sin(k*pi*x).*cos(k*pi*a*t);
二维拉普拉斯方程 X X 0;

Y


Y

0;
一维热传导方程 X X 0 T a2T 0 T ea2t
一维波动方程
X X 0

T


a2T

0

T

cos
at 0
这些方程的分离变量中都会出现方程：X X 0
plot(x,y);grid on; axis([0,m,-.5,1]); title('0-3阶Bessel函数'); subplot(2,1,2);y=bessely(k,x); plot(x,y);grid on; axis([0,m,-1,.6]); title('0-3阶Neumann函数');
ex306
图3-10 (m=0,n=0~3) ex307
图3-14 (m=1,n=0~3) ex308
图3-16 (m=2,n=0~3) ex309
4.3 三维本征值问题(球贝塞耳函数和勒让德函数)
4.3.1 球坐标系中的拉普拉斯方程
在球坐标系中拉普拉斯方程的本征值问题：
2u

1 r2

X X

0
k
2 x
X
X
a


0
X x sin nπx / a
X
X

Y Y




k
2 x

k
2 y
Y Y

0
k
2 y
y
Y b


0
Y y sin mπy / b
ux, y sin nπx / asin mπy / b;
4.1.2 本征函数系的图像及其运动
上节中的四个本征函数的图像如下页图。显然，端点为第一 类边界条件时，该端点为波节，端点为第二类边界条件时， 该端为波腹。程序如下：
%ex301(p65) clear;x=pi*(0:0.001:1); A=sin([1:4]'*x); B=cos([0:3]'*x); C=sin([1/2:7/2]'*x); D=cos([1/2:7/2]'*x); figure(1); 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);
%ex302(p65) function zb t=0:0.005:2.0; x=0:0.001:1; w1=wf (1,0); w2=wf(2,0); w3=wf (3,0); w4=wf(4,0); figure(1); subplot(2,1,1);h1=plot(x,w1,'r'); hold on; h3=plot(x,w3,'g');
一维热传导方程 ut a2uxx 0 (0 x l)
一维波动方程 utt a2uxx 0 (0 x l)
设函数 u(x,t) 具有变量分离形式，即它可表示为：
u(x, y) X (x)Y ( y) 或 u(x,t) X (x)T (t)
则上述三个方程可以写为：
第四章 本征函数系与本征振动
4.1 一维本征值问题；正余弦函数系 4.2 二维本征值问题；贝塞尔函数系 4.3 三维本征值问题；球谐函数系(勒让德函数系) 4.4 用偏微分方程工具箱(pdetool)求解本征值问题；
4.0 PDE 分类
Laplacian
算子：


2 x2

2 y 2
拉普拉斯方程 （elliptic）：
subplot(2,1,2); h2=plot(x,B);hold on; h4=plot(x,D,'g');grid on; y2=max(abs(B));axis([0,1,-y2,y2]); for n=2:length(t) A=sf(1,t(n)); set(h1,'ydata',A); C=sf(3,t(n)); set(h3,'ydata',C); drawnow; B=sf(2,t(n)); set(h2,'ydata',B); D=sf(4,t(n)); set(h4,'ydata',D); drawnow; end function wtx=sf(k,t) x=0:0.01:1; wtx=sin(k*pi*x)*exp(-(2*pi*k)^2*t);
sin

Y


1
sin
2Y

2R R
R





k
2
2R R k 2 m2 R 0
m2 0

Rk0 0
2

R
Jm
k
n m
r
m eim