Using MATLAB, the implementation process of QR was written based on householder change
1, Householder change algorithm
function [H, V, beta] = Householder (x)
% x:inout param. x is a vector which the size is n*1
% v and beta:is PA Ram which construct h matrix
% h is hoseholder matrix. H = I-beta*v*v '
%derive parameter V and beta
v = zeros (size (x));
Beta = zeros (size (x));
%houserholder algorithm begin
X_len = length (x);
X_max = max (x);
x = X./x_max;
Zgama = x (2:end) ' *x (2:end);
V (1) = 1;
V (2:end) = x (2:end);
If Zgama = = 0
beta = 0;
else
alpha = sqrt (x (1) ^2 + Zgama);
if X (1) <=0
V (1) = x (1)-alpha;
Else
V (1) =-zgama./(x (1) + alpha);
End
beta = 2*v (1) ^2./(Zgama + V (1) ^2);
v = v./v (1);
End
%beta = 2./(v ' *v);
H = Eye (x_len,x_len)-beta*v*v ';
End
2, QR decomposition
A = rand (300,20);
Tic
A_hang = size (a,1);
A_lie = size (a,2);
H = Cell (A_lie, 1);
Hs = Eye (A_hang, A_hang);
R2 = A;
Q2 = Eye (A_hang,a_hang);
For i = 1:a_lie
[Hi, VI, Betai] = Householder (R2 (i:end,i));
%h{i} = Blkdiag (eye (i-1), Hi);
R2 (i:end,i:end) = HI*R2 (i:end,i:end);
Q2 = Q2*blkdiag (Eye (i-1), Hi);
End
Q2 (ABS (Q2) < 1e-10) = 0;
R2 (ABS (R2) < 1e-10) = 0;
TOC
[Q1, R1] = QR (A);%matlab internal algorithm (very very fast)
q1*r1-a
q2*r2-a erroq
= Q1 + Q2
ErroR = R1 + R2
Mainly write the implementation of the QR process, the algorithm is easy to understand. But for yourself matrix data forms are also converted.
For example: sparse matrices do not compute the transformation of 0 elements.