矩阵sum_数值计算（十五）Householder变换矩阵QR分解

function  [ Q, R ] =householder()
%{
程序功能：
1、Householder变换又被称为初等反射矩阵或镜像变换
2、这一变换将一个向量变换为由一个超平面反射的镜像，是一种线性变换。
3、有某一个向量 x， 想通过一个标准正交矩阵 P 将 x 转换为 y，有什么方法可以求出矩阵 P？
    这里，我们有一个更加快捷的公式，即为 Householder Transformations.
4、Hx=y, x,y为两个不相等的n维向量
5、Householder变换将一个矩阵A分解为正交矩阵Q和上三角矩阵R.
6、u为n维实向量，且u'u=1。矩阵H=I-2uu' 称为houseHolder变换。
7、利用schmidt正交分解
%}
%     
    clear,clc
%     A=[1,2,2;1,0,2; 0,1,1];
%     A=[12,10,4;10,8,-5;4,-5,3];
    A=[0,3,1; 0,4,-2; 2,1,2];
    N=length(A);
    x=zeros(N,N);
    for k=1:N
        x(:,k)=A(:,k);
    end
    %正交化
    y=zeros(N,N);
    e=zeros(N,N);
    y(:,1)=x(:,1);
    for k=2:N
        s=0;
        for j=1:k-1
            s=s+vecMul(x(:,k),y(:,j))/vecMul(y(:,j),y(:,j))*y(:,j) ;%正交化
        end
        y(:,k)=x(:,k)-s
        
    end
    for k=1:N
        e(:,k)=standard(y(:,k)) ;  %归一化
    end
    Q=zeros(N,N);
    R=zeros(N,N);
    for i=1:N
        Q(:,i)=e(:,i);
%         R(:,i)=max( y(:,i)./x(:,i) );
    end
    %---验证----
    R=QA;
    
end

%向量归一化
function y=standard(x)
    s=sum(x.*x);
    y=x/sqrt(s);
end
function sys=vecMul(x,y)
        
%     sys=sqrt( sum(x.*y) );
    sys= sum(x.*y) ;
end