%% BASIC / INEFFICIENT HOUSEHOLDER
% This could be made more efficient by 
%   * Not explicitly forming each H_k
%   * Not explicitly forming H'
% It could be made more accurate by 
%   * choosing the appropriate sign of norm(a,2)

%% Setup (example problem)

A = [4 -5 ; 0 6 ; -4 2 ; 2 5]
b = [1 2 3 9]';

m = 4;
n = 2;

I = eye(m);

%% Step 1 (1st column)
a = A(:,1)
v = a - norm(a,2)*I(:,1) % Note: Inaccurate!
H1 = I - 2*(v*v')/(v'*v) % Note: Inefficient!
A1 = H1*A                % Note: Inefficient!

%% Step 2 (2nd column)

a = [0 ; A1(2:end,2)]
v = a - norm(a,2)*I(:,2) % Note: Inaccurate!
H2 = I - 2*(v*v')/(v'*v) % Note: Inefficient!
A2 = H2*H1*A             % Note: Inefficient!

%% Check
R = A2(1:n,1:n)
H = H2*H1;
Q = H'    

[Qa, Ra] = qr(A)

%% Compute LS solution
c = Q'*b; % Note: Inefficient!
c1 = c(1:n);
x = R\c1

test = A\b


