
clear all
clc

NodeNum=xlsread('Truss2D-Input.xlsx','General','C2');
ElemNum=xlsread('Truss2D-Input.xlsx','General','C4');
PrescDOFNum=xlsread('Truss2D-Input.xlsx','General','C6');

NodeCoord=xlsread('Truss2D-Input.xlsx','NodeCoord');
ElemConec=xlsread('Truss2D-Input.xlsx','ElemConec');
ElemProp=xlsread('Truss2D-Input.xlsx','ElemProp');
PrescDOF=xlsread('Truss2D-Input.xlsx','BoundaryConditions');
Force=xlsread('Truss2D-Input.xlsx','ExternalForces');


% Computing Ttal Stiffness Matrix, Using the 'StiffnessMat' Function
Kmat=StiffnessMat(NodeNum,ElemNum,ElemProp,NodeCoord,ElemConec)

% Computing Modified Stiffness Matrix
KmatMod=Kmat;
for i=1:PrescDOFNum
    KmatMod(PrescDOF(i,1),:)=0;
    KmatMod(:,PrescDOF(i,1))=0;
    KmatMod(PrescDOF(i,1),PrescDOF(i,1))=1;
end

% Computing Modified Forces Vector
for i=1:NodeNum
    ForceVec(2*i-1)=Force(i,1);
    ForceVec(2*i)=Force(i,2);
end
ForceMod=ForceVec';
for i=1:PrescDOFNum
    ForceMod=ForceMod-Kmat(:,PrescDOF(i,1))*PrescDOF(i,2);
end
for i=1:PrescDOFNum
    ForceMod(PrescDOF(i,1))=Kmat(PrescDOF(i,1),PrescDOF(i,1))*PrescDOF(i,2);
end

Diplacement=inv(KmatMod)*ForceMod

Reaction=Kmat*Diplacement-ForceVec'