c++ eigen A.inverse()*B 不等于 A.ldlt().solve(B)

Question

我想计算两个给定矩阵的乘积的轨迹，比如 A 和 B，Trace(AInv * B) 其中 * 是正则矩阵乘积，AInv 是逆矩阵A 的（对称且正定）和 B 是对称的。

解决方案 1：显式计算逆

注意到 Trace(AInv * B) 等同于取 AInv 和 B 的乘积之和：

     double sol1 = (A.inverse().cwiseProduct(B)).sum();

解决方案 2：使用 Eigen 库中的 ldlt 分解

     double sol2 = (A.selfadjointView<Lower>().ldlt().solve(B)).trace();

从理论上讲，这些解决方案应该是相同的，但在我的测试中，它们并不相同。听起来我好像错过了什么。由于 .ldlt().solve() 不是用来计算矩阵逆而是求解线性系统，我的问题是：.ldlt() 是否执行任何类型的归一化？如果不是，我做错了什么？

非常感谢！

Answer 1

计算 sol1 的语句是错误的：您需要转置其中一个操作数或使用矩阵-矩阵乘积：正确版本：

double sol1 = (A.inverse().cwiseProduct(B.transpose())).sum();
double sol1 = (A.inverse().lazyProduct(B)).diagonal().sum();
double sol1 = (A.inverse().lazyProduct(B)).trace();
double sol1 = (A.inverse() * B).diagonal().sum();
double sol1 = (A.inverse() * B).trace();

请注意，在 Eigen 中，当您编写 (A*B).diagonal() 时，只会计算 A*B 的对角线元素；而不是非对角线元素。

一般来说，不建议显式计算矩阵的逆，使用 A.lu().solve(B) 或 A.ldlt().solve(B) 会给你更准确的结果，而且速度也会更快，因为，除非 A很小(2,3,4)，A.inverse()相当于A.lu().solve(I)。将来，Eigen 很可能会重写这样的表达式：

A.inverse() * B

如：

A.lu().solve(B)

无论如何都是为了你。

c++ eigen A.inverse()*B 不等于 A.ldlt().solve(B)

c++ eigen A.inverse()*B not equal to A.ldlt().solve(B)

eigen

solver

matrix-inverse