function hqr_demo(A)
% HQR Compute the QR factorization of A using
% Householder reflections.
% Input: m x n matrix A. There are no restrictions
% on the relative sizes of m and n.
% Output: m x n matrix Q and n x n upper triangular
% matrix R such that A = QR

[m, n] = size(A);

R = A;
Q = eye(m,m);
q = min(m-1, n);

disp('Initial matrix:');
A
spy(A);
pause;

for i=1:q
	[R(i:m,i:n), u] = hzero1(R(i:m,i:n));
	fprintf('Zero-out all elements below (%d,%d).\n',i,i);
	disp('The Householder reflection is:');
	Hu = eye(m);
	Hu(i:m,i:m) = eye(m-i+1) - (2*u*u')/(u'*u)
	disp('The current value of R is:');
	R
	disp('The current value of Q is:');
	Q(1:m,i:m) = hmulp(Q(1:m,i:m),u)

	% set almost-zero entries to true zero
	% display matrix and wait for a key
	R
% 	R(find(abs(R)<1.e-12))=0;
	R(abs(R)<1.e-12)=0;
	spy(R);
	pause;
	disp('_______________________________________________');
end
disp('Final values of Q and R');
Q
R
disp('Q*R =');
Q*R
disp('A');
A