%% Chapter 4 Part 1

%% Constants
a = 1/sqrt(2);
b = a^-1;
%% Symbolics
syms f1 f2 f3 f4 f5 f6 f7 f8 f9 f10 f11 f12;
syms fh ful fur;

%% Homogeneous System
% A is a coefficent matrix to the non-homogeneous linear system Ax=B.
% Each row vector is from one of the joint equations.
% The row vector order is [f1 f2 f3 f4 f5 f6 f7 f8 f9 f10 f11 f12 f13 fh ful fur]
A = [a  1  0  0  0  0  0  0  0  0  0  0  0 -1  0  0;
     a  0  0  0  0  0  0  0  0  0  0  0  0  0 -1  0;
     0  1  0  0  0 -1  0  0  0  0  0  0  0  0  0  0;
     0  0  1  0  0  0  0  0  0  0  0  0  0  0  0  0;
     a  0  0 -1 -a  0  0  0  0  0  0  0  0  0  0  0;
     a  0  1  0  a  0  0  0  0  0  0  0  0  0  0  0;
     0  0  0  1  0  0  0 -1  0  0  0  0  0  0  0  0;
     0  0  0  0  0  0  1  0  0  0  0  0  0  0  0  0;
     0  0  0  0  a  1  0  0 -a -1  0  0  0  0  0  0;
     0  0  0  0  a  0  1  0  a  0  0  0  0  0  0  0;
     0  0  0  0  0  0  0  0  0  1  0  0 -1  0  0  0;
     0  0  0  0  0  0  0  0  0  0  1  0  0  0  0  0;
     0  0  0  0  0  0  0  1  a  0  0 -a  0  0  0  0;
     0  0  0  0  0  0  0  0  a  0  1  a  0  0  0  0;
     0  0  0  0  0  0  0  0  0  0  0  a  1  0  0  0;
     0  0  0  0  0  0  0  0  0  0  0  a  0  0  0 -1;]
% B is the solution vector.
B = [0; 0; 0; 12; 0; 0; 0; 0; 0; 16; 0; 20; 0; 0; 0; 0;]

%% Solution
% The non-homogeneous linear system Ax=B describes the truss freebody.
x = linsolve(A, B);
%% Printing answers
n = length(x);
% Forces 1..13
fprintf('Force solutions:\n');
for k = 1:n-3
    fprintf('f%d\t=\t%.4f\n', k, x(k));
end
% Forces up and horizontal
fprintf('fh\t=\t%.4f\n', x(14));
fprintf('ful\t=\t%.4f\n', x(15));
fprintf('fur\t=\t%.4f\n', x(16));
%% Force analysis
v = x(1:13);                        % Only consider member forces
i = find(v == max(v));              % Max tension
j = find(v == min(v));              % Max compression
m = find(abs(v) == max(abs(v)));    % Max stress
k = find(abs(v) == min(abs(v)));    % Min stress
fprintf('Member(s) under greatest tension:\n'); disp(i)
fprintf('Member(s) under greatest compression:\n'); disp(j)
fprintf('Member(s) under greatest stress:\n'); disp(m)
fprintf('Member(s) under least stess:\n'); disp(k)