% u''''=f ---- u, u' is 0 a peremen

% A merevségi mátrix és terhelési vektor konstrukciója
n=30;
h=1/(n+1);
Q=[12/h^3,6/h^2,-12/h^3,6/h^2;
    6/h^2,4/h,-6/h^2,2/h;
    -12/h^3,-6/h^2,12/h^3,-6/h^2;
    6/h^2,2/h,-6/h^2,4/h];

A=zeros(2*n,2*n);
b=zeros(2*n,1);
A(1:2,1:2)=Q(3:4,3:4);
A(2*n-1:2*n,2*n-1:2*n)=Q(1:2,1:2);
for k=1:n-1
    b(2*k-1)=-2*h;
    A(2*k-1:2*k+2,2*k-1:2*k+2)=A(2*k-1:2*k+2,2*k-1:2*k+2)+Q;
end
b(2*n-1)=-2*h;

% Az egyenletrendszer megoldása 
u=A\b;

% A megoldás ábrázolása (csak lineáris interpolációval)
plot([0:h:1],[0;u(1:2:end);0])