%
% Seminar No. 3 - truss 2D
%

clc

% array of coordinates
xz=[0, 6;
    16, 0;
    4, 6;
    8, 6;
    12, 6;
    16, 6;
    4, 3;
    8, 3;
    12, 3;
    16, 3;];

% connectivity:
con = [1, 3;
       3, 4;
       4, 5;
       5, 6;
       7, 8;
       8, 9;
       9, 10;
       7, 1;
       7, 3;
       8, 3;
       7, 4;
       8, 4;
       9, 4;
       8, 5;
       9, 5;
       10, 5;
       9, 6;
       10, 6;
       2, 9;
       2, 10];

% array of code numbers
global lm;
lm=[1,2,5,6;
    5,6,7,8;
    7,8,9,10;
    9,10,11,12;
    13,14,15,16;
    15,16,17,18;
    17,18,19,20;
    1,2,13,14;
    5,6,13,14;
    5,6,15,16;
    13,14,7,8;
    7,8,15,16;
    7,8,17,18;
    15,16,9,10;
    9,10,17,18;
    9,10,19,20;
    17,18,11,12;
    11,12,19,20;
    17,18,3,4;
    19,20,3,4];
% number of elements
nelem=20;
EA =210*0.04*1.e9;
% declaration and zeoring:
f=zeros(20,1);
k=zeros(20);
r=zeros(20,1);

% localization
for i=1:nelem
    ke = truss2d([xz(con(i,1),1) xz(con(i,1),2) xz(con(i,2),1) xz(con(i,2),2)], EA);
    k=assembly(k,ke,i);
end

% prescribed displacement and forces
r(3)=0.0005;
f(7)=-20000;
f(8)=-50000;
f(9)=-20000;
f(10)=-50000;

% solution of dispalcements
kuu=k(5:20,5:20);
fu=f(5:20);
u=kuu\(fu-k(5:20,1:5)*r(1:5));

% resulting vector of displacements
r=[0;0;r(3);0;u];

% printig displacement matrix
fprintf('Displacement vector:\n');
fprintf('--------------------\n');
for i=1:10
    fprintf(1,'%2d\t\t{%12e, %12e}\n',i,r(2*i-1),r(2*i));
end


