111xxx.m

三自由度无阻尼系统振动的MATLAB分析

%Calculate the natural frequencies,the normalized mode shapes and the FRFS
clear
close
a=0.001;b=0.001;
M=[100,0,0;0,100,0;0,0,100];
K1=[2,-1,0;-1,2,-1;0,-1,2];
K=1.0e+5*K1;
C=a*M+b*K;           %比例阻尼
[V,D]=eig(K,M);
for j=1:1:3
    omega(j)=sqrt(D(j,j));
    for k=1:1:3
        Y(k,j)=V(k,j)/V(3,j);
    end
end
omega
Y
w=0:0.1:70;                   %频率取值范围和步长
H=zeros(701,3,3);
for j=1:1:701
    Z=K-w(j)^2*M+w(j)*C*i;    %动刚度矩阵
    HH=inv(Z);                %频响函数矩阵        
    H(j,:,:)=HH;
end;
Hr=real(H);Hi=imag(H);      
Habs=abs(H);Hang=(180/pi)*angle(H);
figure(1)
subplot(2,1,1),plot(w,Hr(:,1,1)),grid on
title('实频曲线'),xlabel('\omega'),ylabel('Hr')
subplot(2,1,2),plot(w,Hi(:,1,1)),grid on
title('虚频曲线'),xlabel('\omega'),ylabel('Hi')
figure(2)
subplot(2,1,1),plot(w,Habs(:,1,1)),grid on
title('幅频曲线'),xlabel('\omega'),ylabel('|H|')
subplot(2,1,2),plot(w,Hang(:,1,1)),grid on
title('相频曲线'),xlabel('\omega'),ylabel('\psi')
figure(3)
plot(Hr(:,1,1),Hi(:,1,1)),grid on
title('乃奎斯特圆'),xlabel('Re(H)'),ylabel('Im(H)')
%画跨点 的频响
figure(4)
subplot(2,1,1),plot(w,Hr(:,2,1)),grid on
title('实频曲线'),xlabel('\omega'),ylabel('Hr')
subplot(2,1,2),plot(w,Hi(:,2,1)),grid on
title('虚频曲线'),xlabel('\omega'),ylabel('Hi')
figure(5)
subplot(2,1,1),plot(w,Habs(:,2,1)),grid on
title('幅频曲线'),xlabel('\omega'),ylabel('|H|')
subplot(2,1,2),plot(w,Hang(:,2,1)),grid on
title('相频曲线'),xlabel('\omega'),ylabel('\psi')
figure(6)
plot(Hr(:,2,1),Hi(:,2,1)),grid on
title('乃奎斯特圆'),xlabel('Re(H)'),ylabel('Im(H)')
%画跨点 的频响
figure(7)
subplot(2,1,1),plot(w,Hr(:,3,1)),grid on
title('实频曲线'),xlabel('\omega'),ylabel('Hr')
subplot(2,1,2),plot(w,Hi(:,3,1)),grid on
title('虚频曲线'),xlabel('\omega'),ylabel('Hi')
figure(8)
subplot(2,1,1),plot(w,Habs(:,3,1)),grid on
title('幅频曲线'),xlabel('\omega'),ylabel('|H|')
subplot(2,1,2),plot(w,Hang(:,3,1)),grid on
title('相频曲线'),xlabel('\omega'),ylabel('\psi')
figure(9)
plot(Hr(:,3,1),Hi(:,3,1)),grid on
title('乃奎斯特圆'),xlabel('Re(H)'),ylabel('Im(H)')