1

我正在尝试编写一个获取 mxn 矩阵并 QR 分解它的程序。

它还没有完成，但我遇到了一个问题。我尝试使用 wikipedia http://en.wikipedia.org/wiki/QR_decomposition中显示的示例运行我的程序 

A=[12,-51,4;6,167,-68;-4,24,-41]

他们所谓的 Q1、Q2 ......我称之为 Qtemp。每次我计算 Qtemp 时，我都会打印它以查看我得到与维基百科相同的结果。对于 Q1 我会，但对于 Q2 我不会。

他们的 Q2 和我的 Q2 值相同，但符号不同。他们有一个+我有一个-，他们有一个-，我有一个+。

这是我的代码：

    Q=eye(m);
R=A;
for i=1:min(m-1,n)
    ei=zeros(n,1);
    ei(i,1)=1;
    x=A(:,i);
    for j=1:i-1
        x(j,1)=0;
    end
    u=x-norm(x)*ei;
    v=u/norm(u);
    Qtemp=eye(m)-2*(v*v');
    A=Qtemp*A;
    disp(Qtemp);
end

我实际上只是复制了他们的算法并将其转换为代码，但第二个 Qtemp 的输出仍然很糟糕。

4

1 回答 1

1

看起来您并没有在每次迭代中减少块的大小。一切似乎都是相同的函数m和n（您没有在代码中定义）。请参阅维基百科页面上定义 A' 并使用它来构建 Q 2的行（仅较低的三分之二）。下面是我的一些代码，适用于执行 3×3 矩阵的 QR 分解，可能会有所帮助。请特别注意，第二个块仅适用于A(:,2)and q(2:3,:)：

function [q,r]=qr3(A)

u = A(:,1);
u(1) = u(1)-(1-2*(u(1)<0))*norm(u); % Flip < to > to match sign convention of qr
u = u/norm(u);
u(~isfinite(u)) = sqrt(3)/3;
q = -2*(u*u');
q([1 5 9]) = q([1 5 9])+1;

u = q(2:3,:)*A(:,2);
u(1) = u(1)-(1-2*(u(1)<0))*norm(u); % Flip < to > to match sign convention of qr
u = u/norm(u);
u(~isfinite(u)) = sqrt(2)/2;
q(:,2:3) = q(:,2:3)*[1-2*u(1)^2 -2*u(1)*u(2);
                     -2*u(1)*u(2) 1-2*u(2)^2];
r = triu(q'*A);

qr上面的代码和维基百科上详述的方法使用了与 Matlab函数不同的符号约定。有关如何翻转标志，请参阅我在代码中的注释。

于 2013-08-04T22:08:07.947 回答