Eigen框架下如何直接基于稀疏矩阵求解Ax=b并获取核空间？

我需要求解稀疏矩阵方程 A*x = b 的解 x，同时获取矩阵 A 的核空间（零空间）。目前我通过将稀疏矩阵转换为稠密矩阵实现了需求，示例代码如下：

#include <iostream>
#include <Eigen/Dense>
#include <Eigen/SparseQR>
int main() {
    // 这是一个示例问题，实际使用的矩阵更大更稀疏。
    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()函数的类。

当然可以直接基于稀疏矩阵完成这两项操作！Eigen提供了针对稀疏矩阵的分解工具，完全不用把大稀疏矩阵转成稠密格式浪费内存，下面给你具体的实现方案：

用SparseQR一步完成求解+核空间计算 #

你已经在使用SparseQR求解Ax=b了，其实这个类本身就可以用来获取核空间，不用额外转稠密矩阵。SparseQR会对矩阵做QR分解，其中正交矩阵Q的后几列正好对应A的核空间。

具体步骤是：

1. 用rank()方法获取矩阵的秩，核空间的维度就是列数 - 秩
2. 通过matrixQ().block()提取Q中对应的列，得到核空间的基

调整后的代码示例：

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

int main() {
    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;

    // 初始化SparseQR求解器，用COLAMD排序减少填充
    Eigen::SparseQR<Eigen::SparseMatrix<double>, Eigen::COLAMDOrdering<int>> solver;
    solver.compute(A);

    // 求解Ax=b
    Eigen::VectorXd x = solver.solve(b);
    std::cout << "Solution:\n" << x << std::endl;

    // 计算并输出核空间
    int rank = solver.rank();
    int kernel_dim = A.cols() - rank;
    // 提取Q中从rank列开始的kernel_dim列，作为核空间的基
    Eigen::MatrixXd kernel = solver.matrixQ().block(0, rank, A.rows(), kernel_dim);
    std::cout << "Kernel (sparse-based):\n" << kernel << std::endl;

    return 0;
}

这里需要注意：matrixQ()返回的是稀疏格式的正交矩阵，但block()提取的核空间基是稠密矩阵——不过这没关系，因为核空间的维度通常远小于矩阵本身的规模，内存开销完全可控。

其他可选的稀疏方案 #

如果SparseQR不符合你的场景需求，还可以考虑这些选项：

• SparseLU：适合方阵的LU分解，能高效求解方程，但不能直接获取核空间，需要额外做零空间投影的处理，步骤稍复杂。
• 迭代求解器（BiCGSTAB/GMRES）：针对超大规模稀疏矩阵，适合求解方程，但核空间计算需要结合投影方法，实现起来更繁琐。
• BDCSVD：虽然是稠密SVD，但支持分块处理，适合内存装不下的超大矩阵，但本质还是依赖稠密操作，仅当核空间维度很小时实用。

小提示 #

• 稀疏矩阵的核空间计算依赖分解的数值稳定性，SparseQR用COLAMD排序减少填充，在大多数工程场景下足够稳定。
• 如果你的矩阵是正定对称矩阵，可以用SimplicialLDLT或SimplicialLLT快速求解方程，但这些分解类不支持直接获取核空间，需要配合其他方法。

内容的提问来源于stack exchange，提问作者Dominik Mokriš