function main
format long
E = 2e11;
L_el = 100e-3;
b = 10e-3;
h = 20e-3;
I = (b * h^3) / 12;

F = 10;
M = F * L_el;
f = [0; 0; 2 * F; 0; 4 * F; 0; 0; -2 * M; 0; 0];

u_flag = [0; 0; 1; 1; 1; 0; 1; 1; 0; 0];

N_el = 4;
N_node = 5;

N_dof = N_node * 2;

node_cnt = 4;

K_el = zeros(node_cnt);
K_glob = zeros(N_dof);

M_index = [1 2 3 4; 3 4 5 6; 5 6 7 8; 7 8 9 10];

L = [L_el, 1.5 * L_el, L_el, L_el];

for i = 1:N_el
    K_el = K_calc(E, L(i), I);
    for j = 1:node_cnt
        for k = 1:node_cnt
           K_glob(M_index(i, j), M_index(i, k)) = K_glob(M_index(i, j), M_index(i, k))+ K_el(j, k);
        end
    end
end

for i = 1:length(u_flag)
    if ~u_flag(i, 1)
        K_glob(:, i) = 0;
        K_glob(i ,:) = 0;
        K_glob(i, i) = 1;
    end
end

final_K_glob = K_glob;

u = pinv(final_K_glob) * f

end

function K_el_loc = K_calc(E, L, I)


K_el_loc = [12*E*I/L^3, 6*E*I/L^2, -12*E*I/L^3, 6*E*I/L^2;
    6*E*I/L^2, 4*E*I/L, -6*E*I/L^2, 2*E*I/L;
    -12*E*I/L^3,-6*E*I/L^2, 12*E*I/L^3, -6*E*I/L^2;
    6*E*I/L^2, 2*E*I/L, -6*E*I/L^2, 4*E*I/L];

end