% *****************  EIGENVALUE PROBLEM ******************
% Program to solve for natural frequencies (eigenvalues) and
%   resulting mode shapes (eigenvectors) in simple structures;
% Uses eig() function;
% Author: Gautam Wagle and N.J. Salamon, Apr 2002


m=1;
k=10;
x=[3 2 1 0]; % Reverse x because MATLAB's origin is at bottom;
             %   Then the mode shapes plot as we picture them;

stiff=[2*k -k 0;-k 2*k -k;0 -k k]; % stifness matrix is reduced 
                                   % by eliminating first row/col;
mass = m*eye(3); % eye is the identity matrix;
% To solve using MATLABs function eig(A,B);

[Vector,Value] = eig(stiff,mass) % Vector= eigenvector
                                 % Value = eigenvalue
plotvector= [0 0 0;Vector]; % expand Vector to include node 1 with
                            % zero displacements
for i=1:3;
    figure; 
    plot(plotvector(:,i),x);
end;

% To solve using polynomial expansion

syms lambda; % specifying a symbolic lambda to get poly exprn.

characteristic_poly = det(stiff - lambda*mass)

% NOTE: The above poly will give a polynomial expression in lambda
%       The next step will be to find the roots of the polynomial
%       equation using roots(A) where A is a row vector of 
%       coefficients of lambda in the equation
%
% eg. If characteristic_poly = 1000-600*lambda+50*lambda^2-lambda^3
%       then A = [-1 50 -600 1000] and  Value = roots(A)