【Eigen】判断矩阵正定性及对称性的方法

finish time:2025-05-23 12:42

Cholesky分解可以用来判断一个矩阵是否是对称正定的，但并不能单独用来判断矩阵是否对称。具体原因如下：

Cholesky分解的前提条件
Cholesky分解要求矩阵必须是对称正定的（对于实矩阵）或Hermitian正定的（对于复矩阵）。
如果矩阵能够成功进行Cholesky分解，那么它一定是对称（或Hermitian）且正定的。
判断对称性的局限性
Cholesky分解本身并不能直接判断矩阵是否对称，因为：
- 如果矩阵不对称，Cholesky分解会直接失败（例如，在数值计算中抛出错误）。
- 但失败的原因可能是矩阵不对称，也可能是矩阵不正定。因此，仅凭Cholesky分解无法区分这两种情况。
要明确判断矩阵是否对称，需要直接检查矩阵是否满足 $ A = A^T $（实矩阵）或 $ A = A^H $（复矩阵）。

因此，最好先检查对称性，再验证正定性：在尝试Cholesky分解前，应先通过定义检查矩阵是否对称（例如计算 $ \|A - A^T\| $ 是否接近零）。如果矩阵对称，再通过Cholesky分解或其他方法（如特征值检查）验证正定性。

import numpy as np

def is_symmetric(A, tol=1e-8):
    return np.allclose(A, A.T, atol=tol)

def is_positive_definite(A):
    try:
        np.linalg.cholesky(A)
        return True
    except np.linalg.LinAlgError:
        return False

# 测试矩阵
A = np.array([[4, 1, 1], [1, 2, 1], [1, 1, 3]])  # 对称正定
B = np.array([[1, 2], [3, 4]])                   # 不对称

print("A是对称的:", is_symmetric(A))  # 输出 True
print("A是正定的:", is_positive_definite(A))  # 输出 True
print("B是对称的:", is_symmetric(B))  # 输出 False

+ 总结

Cholesky分解可以间接推断对称性，但仅适用于正定矩阵。对于非正定矩阵，分解失败无法确定是否对称。
直接检查 $ A = A^T $是判断对称性的可靠方法。

Python稀疏矩阵

在 Python 中，稀疏矩阵（如 scipy.sparse 格式）的对称性和正定性判断需要特殊处理，因为直接使用 A == A.T 或 np.linalg.cholesky 可能会遇到效率或兼容性问题。

1. 判断稀疏矩阵的对称性

由于稀疏矩阵通常以压缩格式存储（如 CSR、CSC），直接比较 A == A.T 可能不高效。推荐方法是计算 (A - A.T) 的范数，检查 A - A.T 的 Frobenius 范数是否接近 0

import numpy as np
from scipy.sparse import csr_matrix
from scipy.sparse.linalg import norm

def is_symmetric_sparse(A, tol=1e-8):
    if not (A.shape[0] == A.shape[1]):  # 非方阵一定不对称
        return False
    return norm(A - A.T, 'fro') < tol  # Frobenius 范数

# 示例
A = csr_matrix([[4, 1, 0], [1, 2, 1], [0, 1, 3]])  # 对称
B = csr_matrix([[1, 2], [3, 4]])                   # 不对称

print("A 对称:", is_symmetric_sparse(A))  # True
print("B 对称:", is_symmetric_sparse(B))  # False

对于小矩阵，可以逐元素比较

def is_symmetric_sparse_elementwise(A, tol=1e-8):
    if not (A.shape[0] == A.shape[1]):
        return False
    return (A != A.T).nnz == 0  # 非零元素位置是否相同

print("A 对称:", is_symmetric_sparse_elementwise(A))  # True

2. 判断稀疏矩阵的正定性

Cholesky 分解要求矩阵 对称正定，但 scipy.sparse 没有内置的稀疏 Cholesky 分解。

方法1是转换为稠密矩阵，仅适用于小矩阵；缺点是大稀疏矩阵会占用大量内存。

from scipy.linalg import cholesky

def is_positive_definite_sparse_dense(A):
    try:
        cholesky(A.toarray())  # 转为稠密矩阵
        return True
    except np.linalg.LinAlgError:
        return False

print("A 正定:", is_positive_definite_sparse_dense(A))  # True

更好的方法是 scipy.sparse.linalg 的迭代方法

检查矩阵是否对称（见上文）,然后计算最小特征值是否 > 0：

from scipy.sparse.linalg import eigsh

def is_positive_definite_sparse(A, tol=1e-8):
    if not is_symmetric_sparse(A, tol):
        return False
    # 计算最小特征值（使用 shift-invert 模式加速）
    min_eig = eigsh(A, k=1, which='SA', return_eigenvectors=False)[0]
    return min_eig > tol

print("A 正定:", is_positive_definite_sparse(A))  # True

+ 完整代码

import numpy as np
from scipy.sparse import csr_matrix
from scipy.sparse.linalg import norm, eigsh

def is_symmetric_sparse(A, tol=1e-8):
    if not (A.shape[0] == A.shape[1]):
        return False
    return norm(A - A.T, 'fro') < tol

def is_positive_definite_sparse(A, tol=1e-8):
    if not is_symmetric_sparse(A, tol):
        return False
    min_eig = eigsh(A, k=1, which='SA', return_eigenvectors=False)[0]
    return min_eig > tol

# 测试
A = csr_matrix([[4, 1, 0], [1, 2, 1], [0, 1, 3]])  # 对称正定
B = csr_matrix([[1, 2], [2, 1]])                   # 对称但不正定
C = csr_matrix([[1, 2], [3, 4]])                   # 不对称

print("A 对称:", is_symmetric_sparse(A))  # True
print("A 正定:", is_positive_definite_sparse(A))  # True
print("B 对称:", is_symmetric_sparse(B))  # True
print("B 正定:", is_positive_definite_sparse(B))  # False
print("C 对称:", is_symmetric_sparse(C))  # False

任务	推荐方法
判断对称性	`norm(A - A.T, 'fro') < tol` 或 `(A != A.T).nnz == 0`
判断正定性	先检查对称性，再用 `eigsh(A, which='SA')` 计算最小特征值是否 > 0
小矩阵	转为稠密矩阵后用 `np.linalg.cholesky`

对于 大规模稀疏矩阵，推荐使用 特征值方法（eigsh） 避免内存问题。

eigsh 适合大规模稀疏矩阵（如 1000x1000 以上）。
cholesky 适合小矩阵（如 100x100 以下）。

Eigen稀疏矩阵

在 C++ 中使用 Eigen 库判断 CSC（Compressed Sparse Column）稀疏矩阵的对称性和正定性，可以结合特征值计算和 Cholesky 分解。

1. 判断 CSC 稀疏矩阵的对称性

Eigen 的 SparseMatrix 默认是 列优先（CSC），我们可以通过比较 A 和 A.transpose() 来判断对称性。

方法：计算 `A - Aᵀ` 的范数

A.norm() 计算 Frobenius 范数（所有元素的平方和开根号）。
A.makeCompressed() 确保矩阵是 CSC 格式。

#include <iostream>
#include <Eigen/Sparse>
#include <Eigen/Dense>
using namespace Eigen;

bool isSymmetric(const SparseMatrix<double>& A, double tol = 1e-8) {
    if (A.rows() != A.cols()) return false;  // 非方阵不对称    
    SparseMatrix<double> AT = A.transpose();
    SparseMatrix<double> diff = A - AT;

    // 计算 Frobenius 范数
    double norm = diff.norm();
    return norm < tol;
}

int main() {
    // 示例：对称 CSC 矩阵
    SparseMatrix<double> A(3, 3);
    A.insert(0, 0) = 4; A.insert(0, 1) = 1;
    A.insert(1, 0) = 1; A.insert(1, 1) = 2; A.insert(1, 2) = 1;
    A.insert(2, 1) = 1; A.insert(2, 2) = 3;
    A.makeCompressed();  // 转换为 CSC 格式

    std::cout << "A 是否对称: " << isSymmetric(A) << std::endl;  // 输出 1 (true)

    // 示例：不对称矩阵
    SparseMatrix<double> B(2, 2);
    B.insert(0, 0) = 1; B.insert(0, 1) = 2;
    B.insert(1, 0) = 3; B.insert(1, 1) = 4;
    B.makeCompressed();

    std::cout << "B 是否对称: " << isSymmetric(B) << std::endl;  // 输出 0 (false)
    return 0;
}

2. 判断 CSC 稀疏矩阵的正定性

正定性要求：

矩阵对称（Hermitian）。
所有特征值 > 0。

方法 1：使用 Cholesky 分解（适用于小矩阵）

SimplicialLLT 是 Eigen 提供的 稀疏 Cholesky 分解，适用于 小规模稀疏矩阵。
如果分解失败（info() != Success），则矩阵不正定。

#include <Eigen/Cholesky>

bool isPositiveDefinite(const SparseMatrix<double>& A, double tol = 1e-8) {
    if (!isSymmetric(A, tol)) return false;  // 先检查对称性

    SimplicialLLT<SparseMatrix<double>> llt(A);
    if (llt.info() == Success) {
        return true;  // Cholesky 分解成功，说明正定
    } else {
        return false;
    }
}
int main() {
    SparseMatrix<double> A(3, 3);  // 正定矩阵
    A.insert(0, 0) = 4; A.insert(0, 1) = 1;
    A.insert(1, 0) = 1; A.insert(1, 1) = 2; A.insert(1, 2) = 1;
    A.insert(2, 1) = 1; A.insert(2, 2) = 3;
    A.makeCompressed();

    std::cout << "A 是否正定: " << isPositiveDefinite(A) << std::endl;  // 1 (true)

    SparseMatrix<double> B(2, 2);  // 不正定
    B.insert(0, 0) = 1; B.insert(0, 1) = 2;
    B.insert(1, 0) = 2; B.insert(1, 1) = 1;
    B.makeCompressed();

    std::cout << "B 是否正定: " << isPositiveDefinite(B) << std::endl;  // 0 (false)
    return 0;
}

方法 2：计算最小特征值（适用于大矩阵）

SelfAdjointEigenSolver 适用于 对称矩阵，计算所有特征值。
大矩阵优化：如果矩阵很大，可以用 ArpackSupport 或 Spectra 计算部分特征值（如最小特征值）。

#include <Eigen/Eigenvalues>
bool isPositiveDefiniteEigenvalues(const SparseMatrix<double>& A, double tol = 1e-8) {
    if (!isSymmetric(A, tol)) return false;

    // 计算最小特征值（使用 ARPACK 或 Spectra）
    SelfAdjointEigenSolver<SparseMatrix<double>> eigensolver(A);
    if (eigensolver.info() != Success) {
        std::cerr << "特征值计算失败！" << std::endl;
        return false;
    }
    double min_eig = eigensolver.eigenvalues().minCoeff();
    return min_eig > tol;
}
int main() {
    SparseMatrix<double> A(3, 3);  // 正定矩阵
    A.insert(0, 0) = 4; A.insert(0, 1) = 1;
    A.insert(1, 0) = 1; A.insert(1, 1) = 2; A.insert(1, 2) = 1;
    A.insert(2, 1) = 1; A.insert(2, 2) = 3;
    A.makeCompressed();

    std::cout << "A 是否正定 (特征值法): " << isPositiveDefiniteEigenvalues(A) << std::endl;  // 1 (true)

    return 0;
}

完整代码

#include <iostream>
#include <Eigen/Sparse>
#include <Eigen/Cholesky>
#include <Eigen/Eigenvalues>

using namespace Eigen;

bool isSymmetric(const SparseMatrix<double>& A, double tol = 1e-8) {
    if (A.rows() != A.cols()) return false;
    SparseMatrix<double> diff = A - A.transpose();
    return diff.norm() < tol;
}

bool isPositiveDefinite(const SparseMatrix<double>& A, double tol = 1e-8) {
    if (!isSymmetric(A, tol)) return false;

    // 方法 1: Cholesky 分解（小矩阵）
    SimplicialLLT<SparseMatrix<double>> llt(A);
    if (llt.info() == Success) return true;

    // 方法 2: 计算最小特征值（大矩阵）
    SelfAdjointEigenSolver<SparseMatrix<double>> eigensolver(A);
    if (eigensolver.info() != Success) return false;
    return eigensolver.eigenvalues().minCoeff() > tol;
}

int main() {
    SparseMatrix<double> A(3, 3);  // 对称正定
    A.insert(0, 0) = 4; A.insert(0, 1) = 1;
    A.insert(1, 0) = 1; A.insert(1, 1) = 2; A.insert(1, 2) = 1;
    A.insert(2, 1) = 1; A.insert(2, 2) = 3;
    A.makeCompressed();

    std::cout << "A 是否对称: " << isSymmetric(A) << std::endl;  // 1 (true)
    std::cout << "A 是否正定: " << isPositiveDefinite(A) << std::endl;  // 1 (true)

    return 0;
}

任务	推荐方法
判断对称性	计算 `A - Aᵀ` 的范数是否接近 0
判断正定性	1. 小矩阵：`SimplicialLLT`（Cholesky 分解） 2. 大矩阵：计算最小特征值
适用库	Eigen (`SparseMatrix`, `SimplicialLLT`, `SelfAdjointEigenSolver`)

这样，你可以高效地在 C++ 中使用 Eigen 判断 CSC 稀疏矩阵 的对称性和正定性！ 🚀