我需要用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 01
0 2 0 1 0
3成功
1 0 0 0 1 00.333333 0
11 0 0 0 1.41421 01 0
1.414211 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";