矩陣分解QR， Householder變換__矩陣分解

使用matlab，基於householder變化寫了QR的實現過程

1、Householder變化演算法

function [ H, v, beta ] = householder( x )% x : inout param. x is a vector which size is n*1% v and beta : is param which construct H matrix% H is hoseholder Matrix. H = I - beta*v*v' %derive parameter v and betav = zeros(size(x));beta = zeros(size(x));%Houserholder Algorithm beginx_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分解

A = rand(300,20);ticA_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);endQ2(find(abs(Q2) < 1e-10)) = 0;R2(find(abs(R2) < 1e-10)) = 0;toc[Q1, R1] = qr(A); %matlab內部演算法（非常非常的快）Q1*R1 - AQ2*R2 - AerroQ = Q1 + Q2erroR = R1 + R2

主要寫了QR的實現過程，演算法很容易明白。但是針對你自己矩陣資料形式還要轉換。

例如：疏鬆陣列就不要計算0元素的變換啦。