For the qr factorization using classic Gram-Schmidt algorithm, I found the 2 different implementations below. The first one uses the
for loop to compute the upper triangular matrix
R, but the second one uses the matrix-vector multiplication. Since they are mathematically equivalent, but they produce different results. I cannot find why it is so, can anyone give me some indications? Furthermore, I often find that whenever I understand an algorithm theoretically, if I implemented it, I will have some trouble because of different numerical results. This is why I always ask questions about how, actually, to implement the algorithm with specific code. Because even if I know the algorithm, I still cannot generate the correct results. Thanks.
% test the classical Gram-Schmidt algorithm clc;clear;format compact; A = hilb(7); %% method 1 R(1,1) = norm(A(:,1)); Q(:,1) = A(:,1)/R(1,1); m = size(A,2); for j=2:m for i=1:j-1 R(i,j) = Q(:,i)'*A(:,j);% the difference between 2 methods from below end q_hat = A(:,j)-Q(:,1:j-1)*R(1:j-1,j); R(j,j) = norm(q_hat); Q(:,j) = q_hat/R(j,j); end %% method 2 R1(1,1) = norm(A(:,1)); Q1(:,1)=A(:,1)/R1(1,1); for j=2:m R1(1:j-1,j) = Q1(:,1:j-1)'*A(:,j); % the difference between 2 methods from above q_hat = A(:,j)-Q1(:,1:j-1)*R1(1:j-1,j); R1(j,j) = norm(q_hat); Q1(:,j) = q_hat/R(j,j); end %% compare norm(A-Q*R) norm(Q'*Q-eye(m)) norm(A-Q1*R1) norm(Q1'*Q1-eye(m))
My results are as follows:
ans = 0 ans = 0.3369 ans = 8.9238e-10 ans = 0.1890