c++ - Eigen:获取稀疏矩阵的核

Question

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

一种可能是convert A 到密集表示。

#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() .

Answer 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) = complement(Range(Aᵀ)) = complement(Range (Q₁)) = Null(Q₁ᵀ)).

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

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

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