多自由度阻尼系统固有振型的MATLAB求解程序

function [v,w,zeta]=vbr_sf(m,d,k)
%vbr_sf  vbr_sf(m,d,k)
%        [v,w,zeta]=vbr4(m,d,k)
%        function vbr_sf finds the mode shapes and natural frequencies of
%        a linear second order matrix equation.  
%        [v,w]=vbr_sf(m,k) finds the mode shapes and natural frequencies
%        for the undamped case.
if nargin==2
  k=d;
  [v,w]=eig(m\k);
  w=sqrt(w);
end
if nargin==3
  if norm(d/m*k-k/m*d) < 1e-8*norm(k/m*d)
    disp('Damping is proportional, eigenvectors are real.')
    [v,w]=eig(m\k);
    w=sqrt(w);
    zeta=(v'*m*v)\(v'*d*v)/2/w;
   else
    disp('Damping is non-proportional, eigenvectors are complex.')
    a=[0*k eye(length(k));-m\k -m\d];
    [v,w1]=eig(a);
    w=abs(w1);
    zeta=-real(w1)/w;
  end
end
w=diag(w);zeta=diag(zeta);

author:yanyongju