C++矩阵计算效率

Question

我们正在尝试优化我们的 C++ 代码，我们有以下矩阵计算(使用 Eigen 库)

#include<Eigen/Dense>

int main(){

   MatrixXd P = MatrixXd::Random(30,30); // a  random double 30 x 30 matrix P
   MatrixXd M = MatrixXd::Random(30,30); // a  random double 30 x 30 matrix M
   Matrix<double, 30, 30> I; 
   I.setIdentity(); // I is an 30 x 30 identity matirx

   P = (I-M)*P

   return 0;

   }

其中都是nxn矩阵，I是单位矩阵。 我们发现重写上面的矩阵计算

   P= (I- M)*P

作为

   P = P-M*P

在使用 gcc 6.2 编译器的 Linux ubuntu 系统中获得约 4-8 倍的加速。我意识到编译器可能对单位矩阵和 I*P = P 这一事实一无所知，但我仍然无法理解是什么让效率提高了这么多。任何人都知道做出如此重大改进的可能原因吗？

Answer 1

首先，I.identity();不存在。您想要的是 I.setIdentity() 或 P = (MatrixXd::Identity(30,30)-M)*P。如果您使用第一个选项，Eigen 肯定需要对 I 和 M 进行完整的 30x30 减法(编译器很难看到与您的等效项第二个表达)。总的来说，这将导致两个临时文件(一个用于差异，一个用于产品)。

如果您实际使用了 I.Identity()，那么您将像调用成员函数一样调用静态函数，您的编译器至少应该警告您这一点。这实际上不会修改 I 并且您最终会在 I 中得到未初始化的值，这可能会包含一些 NaN 或非正规值，这两者都可能对浮点性能不利。当然你的结果是错误的。

总的来说，我认为编写方程式的最简单方法是

P -= M*P;

或

MatrixXd Pnew = P - M*P;