 % resavanje sistema Ax=b
A=[3.1 1.5 1; 1.5 2.5 0.5; 1 0.5 4.2];
b=[10.83 9.2 17.1]';

% prvi korak: svodimo na 0 sve potrebne elemente vektora prve kolone matrice A
x=A(:,1);
sigma = norm(x);
beta = inv(sigma*(sigma + abs(x(1))));
k = -sigma*sign(x(1));
u=x - k*[1 0 0]';
H1=eye(3)-beta*u*u';
T1=H1;
A1=T1*A;

% drugi korak: svodimo na 0 sve potrebne elemente vektora druge kolone matrice A
x=A1([2:end],2);
sigma = norm(x);
beta = inv(sigma*(sigma + abs(x(1))));
k = -sigma*sign(x(1));
u=x -k*[1 0]';
H2=eye(2)-beta*u*u';
T2=[eye(1) zeros(1,2); zeros(2,1) H2];
A2=T2*A1;

% Matrica Q iz QR dekompozicije matrice A
Q=T1*T2;
% provera: Q*A2 = A

%resavanje sistema sa gornje trougaonom matricom (ako je dozvoljen operator \):
x=A2\(Q'*b)
