逆迭代法(1)

%逆迭代算法求第一振型

m=400;

k=10000;

M=m.*[1 0 0;0 1 0;0 0 1];%质量矩阵

K=k*[2 -1 0;-1 2 -1;0 -1 2];%刚度矩阵

F=inv(K)

D=F*M;

a=0;

X=[1 1 1]';%迭代初始向量

for j=0:100000000

j %第j次迭代

X

Y=D*X

p=a;

a=(Y'*K*Y)/(Y'*M*Y)

if abs((a-p)/a)<1e-7

b1=X %b1为第一振型

a1=a %a1为第一阶频率

break

else

X=Y/sqrt([Y'*M*Y])

end

pause

end

%第二振型

M=400.*[1 0 0;0 1 0;0 0 1]; %质量矩阵

K=1e4.*[2 -1 0;-1 2 -1;0 -1 2]; %刚度矩阵

b1=[0.0250 0.0353 0.0250]'; %第一振型

a1=14.645 %第一振型的频率

I=eye(3);

m1=(b1)'*M*b1;

S=I-(1/m1).*b1*(b1)'*M;

X=[1 1 1]';

X=S*X; %迭代初始向量

a=0;

D=F*M;

D=D-(1/a1.*m1).*b1*(b1)'*M;

for j=0:100000000

j %第j次迭代

X

Y=D*X

p=a;

a=(Y'*K*Y)/(Y'*M*Y)

if abs((a-p)/a)<1e-7

b2=X %b1为第一振型

a2=a %a1为第一阶频率

break

else

X=Y/sqrt([Y'*M*Y])

end

pause

end

%第三振型

M=400.*[1 0 0;0 1 0;0 0 1]; %质量矩阵

K=1e4.*[2 -1 0;-1 2 -1;0 -1 2]; %刚度矩阵

b1=[0.0250 0.0353 0.0250]'; %第一振型

b2=[ 0.0249 -0.0355 0.0249]' %第二振型

a1=14.645 %第一振型的频率

a2=86.3538 %第二振型的频率

I=eye(3);

m1=(b1)'*M*b1;

m2=(b2)'*M*b2;

S=I-(1/m1).*b1*(b1)'*M-(1/m2).*b2*(b2)'*M;

X=[1 1 1]';

X=S*X; %迭代初始向量

a=0;

D=F*M;

D=D-(1/a1.*m1).*b1*(b1)'*M-(1/a2.*m2).*b2*(b2)'*M;

for j=0:100000000

j %第j次迭代

X

Y=D*X

p=a;

a=(Y'*K*Y)/(Y'*M*Y)

if abs((a-p)/a)<1e-7

b1=X %b1为第一振型

a1=a %a1为第一阶频率

break

else

X=Y/sqrt([Y'*M*Y])

end

pause

end

本文来源：https://www.wddqw.com/doc/f38b6f2dba68a98271fe910ef12d2af90242a8ce.html