function [Q,R]=cgs(A)
[m,n]=size(A);
Q = A;
R = zeros(n);
for i=1:n
for k = 1:i-1
 R(k,i) = Q(:,k)'*A(:,i);
 Q(:,i) = Q(:,i)-R(k,i)*Q(:,k);
end
 R(i,i) = norm(Q(:,i));
 Q(:,i) = Q(:,i)/R(i,i);
end
