C++ 乔列斯基因式分解



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

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

对于C++部分,我开始关注本征。 但是,我很难获得与Matlab的chol函数相当的功能。

我尝试使用本征的LDLT类:

#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);

使用您的代码示例和 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 << "nn";
Eigen::LDLT<MatrixXd> tmp(matA);
cout << ((tmp.info() == Success) ? "succeeded" : "failed") << "nn";
MatrixXd U = tmp.matrixL();
cout << U << "nn";
// Using LLT instead
cout << MatrixXd(matA.llt().matrixL()) << "nn";
cout << MatrixXd(matA.llt().matrixU()) << "nn";
}

输出:

1 0

1
0 2 0 1 0
3

成功

1 0 0 0 1 00.333333 0

1

1 0 0 0 1.41421 01 0

1.41421

1 0 1 0 1.41421 0 00
1.41421

根据 Eigen::LDTL 的文档,第二个模板参数_UpLo默认为Lower,但您省略了该参数并希望计算上三角矩阵。

所以你的类实例化应该看起来像这样(不知道这里是否正确的特征定义是否Upper(:

Eigen::LDLT<MatrixXd, Upper> tmp(matB);

编辑@chtz是对的 - 使用Upper不会给你预期的结果,因为LDLT类是用于具有枢轴的鲁棒cholesky分解。 所以在 除了@Avi的正确答案之外,您还可以使用正确的类进行标准 cholesky 分解:

Eigen::LLT<MatrixXd> tmp(matA);
cout <<  MatrixXd(tmp.matrixU()) << "nn";

相关内容

  • 没有找到相关文章

最新更新