QR分解中的正交化输出稍微不准确的正交化矩阵答案

【问题标题】：Orthogonalization in QR Factorization outputting slightly innaccurate orthogonalized matrixQR分解中的正交化输出稍微不准确的正交化矩阵
【发布时间】：2017-02-20 01:20:14
【问题描述】：

我正在为 QR 分解编写代码，但由于某种原因，我的正交方法无法按预期工作。基本上，我的 proj() 方法是输出随机投影。这是代码：

apmatrix<double> proj(apmatrix<double> v, apmatrix<double> u)   
//Projection of u onto v
{
//proj(v,u) = [(u dot v)/(v dot v)]*v
    double a = mult(transpose(u,u),v)[0][0], b = mult(transpose(v,v),v)[0][0], c = (a/b);
    apmatrix<double>k;
    k.resize(v.numrows(),v.numcols());
    for(int i = 0; i<v.numrows(); i++)
    {
        for(int j = 0; j<v.numcols(); j++)
        {
            k[i][j]=v[i][j]*c;
        }
    }
    return k;
}

我使用手动矩阵输入单独测试了该方法，它似乎工作正常。这是我的正交方法：

apmatrix<double> orthogonal(apmatrix<double> A)   //Orthogonal
{
    /*
    n = (number of columns of A)-1
    x = columns of A
    v0 = x0
    v1 = x1 - proj(v0,x1)
    vn = xn - proj(v0,xn) - proj(v1,xn) - ... - proj(v(n-1),xn)
    V = {v1, v2, ..., vn} or [v0 v1 ... vn]
    */
    apmatrix<double> V, x, v;
    int n = A.numcols();
    V.resize(A.numrows(),n);
    x.resize(A.numrows(), 1);
    v.resize(A.numrows(),1);
    for(int i = 0; i<A.numrows(); i++)
    {
        x[i][0]=A[i][1];
        v[i][0]=A[i][0];
        V[i][0]=A[i][0];
    }
    for (int c = 1; c<n; c++)   //Iterates through each col of A as if each was its own matrix
    {
        apmatrix<double>vn,vc; //vn = Orthogonalized v (avoiding matrix overwriting of v); vc = previously orthogonalized v
            vn=x;
        vc.resize(v.numrows(), 1);
        for(int i=0; i<c; i++)   //Vn = an-(sigma(t=1, n-1, proj(vt, xn))
        {
            for(int k = 0; k<V.numrows(); k++)
                vc[k][0] = V[k][i]; //Sets vc to designated v matrix
            apmatrix<double>temp = proj(vc, x);
            for(int j = 0; j<A.numrows(); j++)
            {
                vn[j][0]-=temp[j][0]; //orthogonalize matrix
            }
        }
        for(int k = 0; k<V.numrows(); k++)
        {
            V[k][c]=vn[k][0]; //Subtracts orthogonalized col to V
            v[k][0]=V[k][c]; //v is redundant. more of a placeholder
        }
        if((c+1)<A.numcols())  //Matrix Out of Bounds Checker
        {
            for(int k = 0; k<A.numrows(); k++)
            {
                vn[k][0]=0;
                vc[k][0]=0;
                x[k][0]=A[k][c+1]; //Moves x onto next v
            }
        }
    }
    system("PAUSE");
    return V;
}

出于测试目的，我一直在使用二维数组：[[1,1,4],[1,4,2],[1,4,2],[1,1,0]]。每列都是它自己的 4x1 矩阵。矩阵应分别输出为：[1,1,1,1]T、[-1.5,1.5,1.5,-1.5]T 和 [2,0,0,-2]T。现在发生的情况是第一列输出正确（它是同一个矩阵），但第二列和第三列输出的结果可能相似但不等于它们的预期值。

同样，每次我调用正交方法时，它都会输出不同的东西。我认为这是由于 proj() 方法中输入的数字，但我不完全确定。

apmatrix 来自 AP 大学理事会，当时他们教 cpp。它类似于 Java 中的向量或 ArrayList。

Here 是指向 apmatrix.cpp 和文档或条件（可能更有用）apmatrix.h 的链接。 Here 是完整代码的链接（我添加了视觉标记以查看计算机在做什么）。

假设所有自定义方法都按预期工作是公平的（也许矩阵回归除外，但这无关紧要）。并确保在尝试分解之前使用 enter 方法输入矩阵。代码可能效率低下，部分原因是我不久前自学了 cpp，并且我一直在尝试不同的方法来修复我的代码。感谢您的帮助！

【问题讨论】：

价值观与应有的价值观有何不同？
@AhmedFasih 每个元素始终在预期值的 +-1 范围内。一次迭代返回：[-1, 2, 2, -1]T 用于第二个矩阵和 [2.8, -1, -1, 1.2]T。有时第二个矩阵实际上是正确的值。另外，请注意，第二个矩阵中所有元素的总和等于预期矩阵的总和。第三个矩阵也是如此，除了其中两个元素应该是 0。
好吧，如果不是数值舍入错误，那么我怀疑是缓冲区溢出或内存访问类型问题，因为该算法是完全确定的，但显示的是随机输出。让它通过 valgrind 或 Clang 的一些消毒剂？
这是 Scipy 为您的输入矩阵返回的 Q 矩阵（请注意，您可以将每列除以标量并且仍然具有正交基）：array([[-0.5, 0.5, 0.7, -0. ], [-0.5, -0.5, -0. , 0.7], [-0.5, -0.5, 0. , -0.7], [-0.5, 0.5, -0.7, 0. ]])。每个子列表都是 Numpy 约定中的一行，因此查看每个子列表的第一个元素以查看第一列。前三个匹配您的三个预期向量...
@AhmedFasih 今天做了更多测试后，我发现这实际上是一些内存问题。我发现由于某种原因，如果在循环中声明变量或 apmatrix 对象，初始化，然后重复该循环，内存不会完全擦除存储在该变量或对象中的值。这在我的代码中的两个地方都有说明。无论出于何种原因，我必须在 proj 方法中将双精度数 a、b 和 c 设置为 0，在 mult 方法中将 apmatrixdh 设置为 0，否则它们会在下一次迭代中存储一些值。非常感谢您的帮助！

标签： c++ linear-algebra qr-code projection orthogonal

【解决方案1】：

如cmets中所说：

@AhmedFasih 今天做了更多测试后，我发现这实际上是一些>内存问题。我发现由于某种原因，如果一个变量或 apmatrix 对象 > 在循环中声明、初始化，然后重复该循环，则 >memory 不会完全擦除存储在该变量或对象中的值。这 > 在我的代码中的两个地方注明。无论出于何种原因，我必须在 proj 方法中将 >doubles a、b 和 c 设置为 0，在 >mult 方法中将 apmatrixdh 设置为 0，否则它们会在下一次迭代中存储一些值。非常感谢您的帮助！

【讨论】：