disp(' ');
disp(' spring_mass_system.m  ver 1.0  February 24, 2012');
disp(' ');
disp(' by Tom Irvine ');
disp(' ');
disp(' This program assembles the mass and stiffness matrices');
disp(' for a spring-mass system. It also solves the corresponding');
disp(' generalized eigenvalue problem. ');
%
clear mass;
clear stiff;
clear length;
clear k;
clear T;
clear fn;
clear omegan;
clear ModeShapes;
clear MST;
%
disp(' ');
disp(' Each dof should have exactly one mass. ');
disp(' ');
%
disp(' Select units:  1=English  2=metric ');
iu=input(' ');
%
if(iu==1)
    disp(' ');
    disp(' Select input mass unit: 1=lbm  2=lbf sec^2/in ');
    imass=input(' ');
end
disp(' ');
%
if(iu==1)
    if(imass==1)
        disp(' The mass file must have two columns:  dof number & mass(lbm) ');
    else
        disp(' The mass file must have two columns:  dof number & mass(lbf sec^2/in) ');
    end
else
    disp(' The mass file must have two columns:  dof number & mass(kg) ');
end
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
disp(' Select mass file input method ');
disp('   1=external ASCII file ');
disp('   2=file preloaded into Matlab ');
disp('   3=Excel file ');
file_choice = input('');
%
clear THM;
if(file_choice==1)
        [filename, pathname] = uigetfile('*.*');
        filename = fullfile(pathname, filename);
        fid = fopen(filename,'r');
        THM = fscanf(fid,'%g %g',[2 inf]);
        THM=THM';
end
if(file_choice==2)
        FS = input(' Enter the matrix name:  ','s');
        THM=evalin('caller',FS);
end
if(file_choice==3)
        [filename, pathname] = uigetfile('*.*');
        xfile = fullfile(pathname, filename);
%        
        THM = xlsread(xfile);
%         
end
%
ndof=length(THM(:,1));
%
mass=zeros(ndof,ndof);
stiff=zeros(ndof,ndof);
%
start_dof=min(THM(:,1));
%
for i=1:ndof
    j=THM(i,1);
%    
    if(start_dof==0)
        j=j+1;
    end
%
    mass(j,j)=THM(i,2);
end
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
disp('  ');
disp(' Each spring is connect to two dofs. ');
disp(' ');
%
disp(' Select stiffness file input method ');
disp('   1=external ASCII file ');
disp('   2=file preloaded into Matlab ');
disp('   3=Excel file ');
file_choice = input('');
%
%
if(iu==1)
    disp(' The stiffness file must have three columns:  dof numbers 1 & 2, stiffness(lbf/in) ');
else
    disp(' The stiffness file must have three columns:  dof numbers 1 & 2, stiffness(N/m) ');
end
%
clear THM;
if(file_choice==1)
        [filename, pathname] = uigetfile('*.*');
        filename = fullfile(pathname, filename);
        fid = fopen(filename,'r');
        THM = fscanf(fid,'%g %g  %g',[3 inf]);
        THM=THM';
end
if(file_choice==2)
        FS = input(' Enter the matrix name:  ','s');
        THM=evalin('caller',FS);
end
if(file_choice==3)
        [filename, pathname] = uigetfile('*.*');
        xfile = fullfile(pathname, filename);
%        
        THM = xlsread(xfile);
%         
end
%
%
n=length(THM(:,1));
%
for i=1:n
%
    n1=THM(i,1);
    n2=THM(i,2);
%
    if(start_dof==0)
        n1=n1+1;
        n2=n2+1;
    end
%
    k=THM(i,3);
%
    T=zeros(2,ndof);
    T(1,n1)=1;
    T(2,n2)=1;
%
    KL=[k -k;-k k];
%
    stiff=stiff+(T'*KL*T);
%
end
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
if(iu==1 && imass==1)
  mass=mass/386;
end
%
mass
%
stiff
%
[fn,omegan,ModeShapes,MST]=Generalized_Eigen(stiff,mass,1);
%