Eigen：获取稀疏矩阵的核答案

【问题标题】：Eigen: Obtain the kernel of a sparse matrixEigen：获取稀疏矩阵的核
【发布时间】：2019-02-19 12:29:17
【问题描述】：

给定一个稀疏矩阵A和一个向量b，我想得到方程A * x = b的解x以及A的核。

一种可能性是convertA 到密集表示。

#include <iostream>
#include <Eigen/Dense>
#include <Eigen/SparseQR>

int main()
{
    // This is a toy problem. My actual matrix
    // is of course bigger and sparser.
    Eigen::SparseMatrix<double> A(2,2);
    A.insert(0,0) = 1;
    A.insert(0,1) = 2;
    A.insert(1,0) = 4;
    A.insert(1,1) = 8;
    A.makeCompressed();

    Eigen::Vector2d b;
    b << 3, 12;

    Eigen::SparseQR<Eigen::SparseMatrix<double>,
                    Eigen::COLAMDOrdering<int> > solver;
    solver.compute(A);
    std::cout << "Solution:\n" << solver.solve(b) << std::endl;

    Eigen::Matrix2d A_dense(A);
    std::cout << "Kernel:\n" << A_dense.fullPivLu().kernel() << std::endl;
    return 0;
}

是否可以直接在稀疏表示中做同样的事情？除了FullPivLu，我在任何地方都找不到函数kernel()。

【问题讨论】：

懒得自己实现它，但是应该计算A.transpose()的SparseQR并通过将qr.matrixQ()乘以Q的最后一列来获得qr.rank()的最后一列一个稀疏（或密集）的单位矩阵。
@chtz 你不想把你的评论变成答案吗？如果对它的工作原理有一点解释，我很乐意接受。

标签： c++ sparse-matrix linear-algebra eigen

【解决方案1】：

我认为@chtz 的答案几乎是正确的，除了我们需要采用最后的 A.cols() - qr.rank() 列。这是一个数学推导。

假设我们对矩阵 Aᵀ 进行 QR 分解

Aᵀ * P = [Q₁ Q₂] * [R; 0] = Q₁ * R

其中P是置换矩阵，因此

Aᵀ = Q₁ * R * P⁻¹。

我们可以看到 Range(Aᵀ) = Range(Q₁ * R * P⁻¹) = Range(Q₁)（因为 P 和 R 都是可逆的）。

由于 Aᵀ 和 Q₁ 具有相同的范围空间，这意味着 A 和 Q₁ᵀ 也将具有相同的零空间，即 Null(A) = Null(Q₁ᵀ)。（这里我们使用 Range(M) 和 Null(Mᵀ) 对任何矩阵 M 互为补码的性质，因此 Null(A) = 补码（Range(Aᵀ)) = 补码（Range(Q₁)) =空（Q₁ᵀ））。

另一方面，由于矩阵 [Q₁ Q₂] 是正交矩阵，Null(Q₁ᵀ) = Range(Q₂)，因此 Null(A) = Range(Q₂)，即内核(A) = Q₂。

由于 Q₂ 是正确的 A.cols() - qr.rank() 列，您可以调用 rightCols(A.cols() - qr.rank()) 来检索 A 的内核。

有关内核空间的更多信息，您可以参考https://en.wikipedia.org/wiki/Kernel_(linear_algebra)

【讨论】：