对称正定矩阵的特征有效逆

Question

在 Eigen 中，如果我们有对称正定矩阵 A 那么我们可以通过

计算 A 的逆矩阵

A.inverse();

或

A.llt().solve(I);

其中 I 是与 A 大小相同的单位矩阵。但是有没有更有效的方法来计算对称正定矩阵的逆？

例如，如果我们将 A 的 Cholesky 分解写为 A = LL^{T}，那么 L^{-T} L^{-1} 是 A 的逆函数，因为 A L^{-T} L^{-1} = LL^{T} L^{-T} L^{-1} = I（并且其中L^{-T}表示L).

的转置的逆

所以我们可以得到 A 的 Cholesky 分解，计算它的倒数，然后得到那个倒数的叉积来找到 A 的倒数。但我的直觉是，计算这些显式步骤会比使用上面的 A.llt().solve(I) 慢。

在有人问之前，我确实需要一个显式逆 - 它是 Gibbs 采样器的一部分的计算。

Answer 1

对于 A.llt().solve(I)，您假设 A 是 SPD 矩阵并应用 Cholesky 分解来求解方程 Ax=I。求解方程的数学过程与您的显式方法完全相同。因此，如果您正确执行每一步，性能应该是相同的。

另一方面，对于 A.inverse()，您正在进行一般矩阵求逆，它对大矩阵使用 LU 分解。因此性能应该低于A.llt().solve(I);.

Answer 2

您应该针对您的特定问题分析代码以获得最佳答案。我正在对代码进行基准测试，同时尝试使用 googletest 库和 this repo:

评估这两种方法的可行性

#include <gtest/gtest.h>

#define private public
#define protected public

#include <kalman/Matrix.hpp>
#include <Eigen/Cholesky>
#include <chrono>
#include <iostream>

using namespace Kalman;
using namespace std::chrono;

typedef float T;
typedef high_resolution_clock Clock;

TEST(Cholesky, inverseTiming) {
    Matrix<T, Dynamic, Dynamic> L;
    Matrix<T, Dynamic, Dynamic> S;
    Matrix<T, Dynamic, Dynamic> Sinv_method1;
    Matrix<T, Dynamic, Dynamic> Sinv_method2;
    int Nmin = 2;
    int Nmax = 128;
    int N(Nmin);

    while (N <= Nmax) {

        L.resize(N, N);
        L.setRandom();
        S.resize(N, N);
        // create a random NxN SPD matrix
        S = L*L.transpose();
        std::cout << "\n";
        std::cout << "+++++++++++++++++++++++++ N = " << N << " +++++++++++++++++++++++++++++++++++++++" << std::endl;
        auto t1 = Clock::now();
        Sinv_method1.resize(N, N);
        Sinv_method1 = S.inverse();
        auto dt1 = Clock::now() - t1;
        std::cout << "Method 1 took " << duration_cast<microseconds>(dt1).count() << " usec" << std::endl;
        auto t2 = Clock::now();
        Sinv_method2.resize(N, N);
        Sinv_method2 = S.llt().solve(Matrix<T, Dynamic, Dynamic>::Identity(N, N));
        auto dt2 = Clock::now() - t2;
        std::cout << "Method 2 took " << duration_cast<microseconds>(dt2).count() << " usec" << std::endl;
        for(int i = 0; i < N; i++)
        {
            for(int j = 0; j < N; j++)
            {
                EXPECT_NEAR( Sinv_method1(i, j), Sinv_method2(i, j), 1e-3 );
            }
        }

        N *= 2;
        std::cout << "+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++" << std::endl;
        std::cout << "\n";
    }
}

上面的例子告诉我的是，对于我的大小问题，使用 method2 的加速可以忽略不计，而缺乏准确性（使用 .inverse() 调用作为基准）是显而易见的。

对称正定矩阵的特征有效逆

Eigen efficient inverse of symmetric positive definite matrix

c++

eigen

matrix-inverse