//программа решения задачи НДС шарнирно-соединенной фермы
//разработал Губин М.В.

//начальные условия
E=2*10^6; s=1;

//матрица координат узлов
xy=[1 1;
    8.072 8.072;
    15.144 1];
    
//определяем число строк в матрице xy (числов узлов фермы)
n=size(xy,1);

//определяем длину элемента
for i=1:(n-1)
    l(i,1)=sqrt((xy(i+1,1)-xy(i,1))^2+(xy(i+1,2)-xy(i,2))^2);
end

//определяем матрицу углов для элементов
for i=1:(n-1)
    g(i,1)=atan((xy(i+1,2)-xy(i,2))/(xy(i+1,1)-xy(i,1)))*360/(2 * %pi);
end

//определяем матрицу A
for i=1:(n-1)
    A(i,1)=E*s/l(i,1);
end

//матрица сил
F=  [
    ];

//матрица кинематических связей
a= [1 0 0;
    0 1 0;
    0 1 0;
    0 0 1];

//вывод значения х на экран
function [x] = print_x(x)
    disp (x);
endfunction

//определение общей матрицы жесткости
function [K] = Mk_opt(A,g)
k=zeros(2*n,2*n);
 for p=0:n-2
    Mkg=[cosd(g(p+1)) -sind(g(p+1));-sind(g(p+1)) cosd(g(p+1))];
    Mk=A(p+1)*[Mkg -Mkg;-Mkg Mkg];
    for i=1:4
        for j=1:4
            k(i+2*p,j+2*p)=k(i+2*p,j+2*p)+Mk(i,j);
        end
    end
 end   
K=k;
endfunction

//start programm
print_x("координаты узлов x,y =");
print_x(xy);

print_x("количество узлов n =");
print_x(n);

print_x("длина сегментов l(i) =");
print_x(l);

print_x("углы для сегментов g(i) =");
print_x(g);

//вычисляем татрицу жесткости
K=Mk_opt(A,g);
print_x("общая матрица жесткости K =");
print_x(K);