%% EK 301: Preliminary Design
% Kara Walp, Nikita Vinay Kishan
%% 2.2: Algorithm
%% 2.2.1: Defining the Truss Parameters
% First we define a connection matrix C. This matrix has j rows and m columns, where the row represents the joint number and the column represents the member number (thus, the first column corresponds to member 1 and the first row corresponds to joint 1). We indicate the connection of member 1 to joints 1 and 2 by placing a ‘1‘ in column 1 of rows 1 and 2. If a given joint is not connected to a given member, we place a ‘0’ in the matrix location that corresponds to the joint and member numbers. This can be summarized as:
% C_j,m = {1, if member m connected to joint j, 0, elsewhere
% The end result should work out such that the sum of each column is 2, since each member is only connected to two joints. The sum of each row should equal the number of members attached to the corresponding joint. For the arrangement shown in Fig. 2, we would have the following connection matrix:

C = [1 1 0 0 0 0 0 0 0 0 0 0 0
     1 0 1 1 0 0 0 0 0 0 0 0 0
     0 0 0 1 1 1 0 0 0 0 0 0 0
     0 1 1 0 1 0 1 1 0 0 0 0 0
     0 0 0 0 0 1 1 0 1 1 0 0 0
     0 0 0 0 0 0 0 1 1 0 1 1 0
     0 0 0 0 0 0 0 0 0 1 1 0 1
     0 0 0 0 0 0 0 0 0 0 0 1 1]

% Next we construct a connection matrix for the support forces along each axis, where Sx and Sy are matrices with j rows and 3 columns. Note that for our statically-determinate truss, supported by one pin and one roller joint, we will have a total of three unknown reactions. In each matrix, for each unknown reaction force, put a ‘1’ in the column that corresponds to the joint j (there should be only a single entry of ‘1’; note that this is true in both matrices even though we know that for the loading conditions in this project, there are no support forces in the x direction) (and, another parenthetical note: be sure you understand and agree with the previous parenthetical remark!). For the example in Fig. 2 (where a pin and roller support is modeled), we get:

Sx = [1 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;]

Sy = [0 1 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 1;]

% The truss is also defined by the location of the joints; to capture these we construct two location vectors X and Y. Each has j elements corresponding to the relevant location of the jth joint; that is:
% X = [x1, x2, x3, x4,...], Y = [y1, y2, y3, y4,...];

C = [1 1 0 0 0 0 0 0 0 0 0 0 0
     1 0 1 1 0 0 0 0 0 0 0 0 0
     0 0 0 1 1 1 0 0 0 0 0 0 0
     0 1 1 0 1 0 1 1 0 0 0 0 0
     0 0 0 0 0 1 1 0 1 1 0 0 0
     0 0 0 0 0 0 0 1 1 0 1 1 0
     0 0 0 0 0 0 0 0 0 1 1 0 1
     0 0 0 0 0 0 0 0 0 0 0 1 1]

% Next we construct a connection matrix for the support forces along each axis, where Sx and Sy are matrices with j rows and 3 columns. Note that for our statically-determinate truss, supported by one pin and one roller joint, we will have a total of three unknown reactions. In each matrix, for each unknown reaction force, put a ‘1’ in the column that corresponds to the joint j (there should be only a single entry of ‘1’; note that this is true in both matrices even though we know that for the loading conditions in this project, there are no support forces in the x direction) (and, another parenthetical note: be sure you understand and agree with the previous parenthetical remark!). For the example in Fig. 2 (where a pin and roller support is modeled), we get:

Sx = [1 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;]

Sy = [0 1 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 0;
      0 0 1;]

% The truss is also defined by the location of the joints; to capture these we construct two location vectors X and Y. Each has j elements corresponding to the relevant location of the jth joint; that is:
% X = [x1, x2, x3, x4,...], Y = [y1, y2, y3, y4,...];

X= [0 0 4 4 8 8 12 12];
Y= [0 4 8 4 8 4 4 0];

% Finally, we define the load vector L which represents the known forces that act on each joint. This vector has 2j elements; the first j elements correspond to loads in the x direction for each of the j joints, and the last j elements represent loads along the y direction. In our example in Fig. 2, the vertical force W = −mgȷˆ for the load is placed on joint 2. In a standard force equilibrium equation, the force W = −mgȷˆ would be combined with the other vertical forces to equal zero. Once this known force is moved into the L matrix (see Eqn. 5 below), its sign is flipped and becomes positive. Therefore:

mg = 25
L = [0 0 0 0 0 0 0 0 0 0 0 -mg 0 0  0 0 ]'

%

[jointno, membno]= size(C);
xvecs=2*jointno;
yvecs=2*jointno;
[reacvecrow, reacveccol]= size(Sx)

Atot=zeros(xvecs,reacveccol);
sum=0;
for j = 1:membno
for i = 1:jointno
if (C(i,j) == 1)
for k = (i + 1):jointno
if C(k,j) == 1
Atot(i,j) = (X(k)-X(i))/(sqrt( ((X(k)-X(i))^2)+((Y(k)-Y(i))^2) ));
Atot(k,j) = (X(i)-X(k))/(sqrt( ((X(k)-X(i))^2)+((Y(k)-Y(i))^2) ));
Atot(i + jointno,j) = (Y(k)-Y(i))/(sqrt( ((X(k)-X(i))^2)+((Y(k)-Y(i))^2) ));
Atot(k + jointno,j) = (Y(i)-Y(k))/(sqrt( ((X(k)-X(i))^2)+((Y(k)-Y(i))^2) ));
sum=sum+(sqrt( ((X(k)-X(i))^2)+((Y(k)-Y(i))^2) ));
if (isnan(Atot(i,j)))
Atot(i,j) = 0;
end
if (isnan(Atot((i + jointno),j)))
Atot(i+jointno,j) = 0;
end
if (isnan(Atot(k,j)))
Atot(k,j) = 0;
end
if (isnan(Atot((k + jointno),j)))
Atot(k+jointno,j) = 0;
end
end
end

break;
end
end
end
clear i
clear j
clear k
%RXN Forces
for i = 1:jointno
for j = 1:reacveccol
if (Sx(i,j) == 1)
Atot(i,j + membno) = 1;
end
if (Sy(i,j) == 1)
Atot((i + jointno),(j + membno)) = 1;
end
end
end
Atot
clear i
clear j
clear k
for i = 1:length(L)
if (L(i) > 0)
loads = L(i);
end
end
% fprintf functions and calculating T
T = inv(Atot) * L; 
fprintf('EK301, Section A2, overloaded, Nikita Vinay Kishan and Kara Walp, 11/11/2022.\n');
fprintf('Load: %d oz\n',mg);
fprintf('Member forces in oz\n');
%Locating Tension (T) members or Compression (C) members
for i = 1:membno
if T(i) > 1  
memtype = 'T';
fprintf('m%d: %0.3f (%c)\n',i,T(i),memtype);
elseif ((T(i) < 1) && (T(i) ~= 0))
memtype = 'C';
%Output Magnitude
T(i)=T(i) * -1; 
fprintf('m%d: %0.3f (%c)\n',i,T(i),memtype);
else
fprintf('m%d: %1.0f\n',i,T(i));
end
end
%Print RXN Forces
fprintf('Reaction forces in oz: \n');   
if (reacveccol == 3)
fprintf('Sx1: %1.2f\n',T(membno + 1));
fprintf('Sy1: %1.2f\n',T(membno + 2));
fprintf('Sy2: %1.2f\n',T(membno + 3));
elseif (reacveccol == 4)
fprintf('Sx1: %1.2f\n',T(membno + 1));
fprintf('Sx2: %1.2f\n',T(membno + 2));
fprintf('Sy3: %1.2f\n',T(membno + 3));
fprintf('Sy4: %1.2f\n',T(membno + 4));
end
%Cost Calculation
cost=sum+(jointno * 10); 
fprintf('Cost of truss: $%0.2f\n',cost);
clear i
clear j
%Load/Cost Ratio Calculation
for i=1:(size(T) - reacveccol)
for j=1:(size(T) - reacveccol)
if (T(i) > T(j))
max=T(i);
end
end
end
ratio = max/cost;
fprintf('Theoretical max load/cost ratio in oz/$: %0.4f\n',ratio);



%% 2.2.2: Construting the Equilibrium Equations
%



%
% Finally, we define the load vector L which represents the known forces that act on each joint. This vector has 2j elements; the first j elements correspond to loads in the x direction for each of the j joints, and the last j elements represent loads along the y direction. In our example in Fig. 2, the vertical force W = −mgȷˆ for the load is placed on joint 2. In a standard force equilibrium equation, the force W = −mgȷˆ would be combined with the other vertical forces to equal zero. Once this known force is moved into the L matrix (see Eqn. 5 below), its sign is flipped and becomes positive. Therefore:

