首页 > 解决方案 > C++ Cholesky 分解

问题描述

我需要使用 C++ 重写一些 MatLab 代码。

在 Matlab 代码中,我们调用该函数chol来计算上三角矩阵。

对于 C++ 部分,我开始关注 Eigen。但是,我很难获得与 Matlab 的功能等效的chol功能。

我尝试使用LDLTEigen 类:

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

int main() {

  MatrixXd matA(2, 2);
  matA << 1, 2, 3, 4;

  MatrixXd matB(4, 4);
  matB << matA, matA/10, matA/10, matA;
  matB = matB*matB.transpose();

  Eigen::LDLT<MatrixXd> tmp(matB);
  MatrixXd U = tmp.matrixU();
  cout << U << endl;
}

但结果与 Matlab 代码不同:

matB = [  1   2 0.1 0.2
          3   4 0.3 0.4
        0.1 0.2   1   2
        0.3 0.4   3   4];
matB = matB*matB';
D = chol(matB);

标签: c++matlabeigen

解决方案


使用您的代码示例和Matlab 文档,我在使用 LLT 而不是 LDLT(在线)时得到相同的结果:

#include <iostream>
#include <Eigen/Dense>
using namespace Eigen;
using std::cout;

int main()
{
  MatrixXd matA(3,3);
  matA << 1, 0, 1, 0, 2, 0, 1, 0, 3;
  cout << matA << "\n\n";
  Eigen::LDLT<MatrixXd> tmp(matA);
  cout << ((tmp.info() == Success) ? "succeeded" : "failed") << "\n\n";
  MatrixXd U = tmp.matrixL();
  cout << U << "\n\n";
  // Using LLT instead
  cout << MatrixXd(matA.llt().matrixL()) << "\n\n";
  cout << MatrixXd(matA.llt().matrixU()) << "\n\n";
}

输出:

1 0 1
0 2 0
1 0 3

成功了

1 0 0
0 1 0
0.333333 0 1

1 0 0
0 1.41421 0
1 0 1.41421

1 0 1
0 1.41421 0
0 0 1.41421


推荐阅读