Eigen中的矩阵乘法非常缓慢

Question

我实现了一个Gauss-Newton优化过程，其中涉及通过解决线性系统来计算增量，即Hx = b。 H矩阵的计算方式为H = J.transpose() * W * J，而b则是由b = J.transpose() * (W * e)计算得出的，其中e为误差向量。这里的Jacobian是一个n x 6的矩阵，其中n在数千个，并且在迭代过程中保持不变，而W是一个n x n的对角重量矩阵，其将在迭代过程中改变（某些对角元素将被设置为零）。然而，我遇到了速度问题。

当我不添加权重矩阵W时，即H = J.transpose()*J和b = J.transpose()*e，我的高斯牛顿过程可以在0.02秒内快速运行30次迭代。但是当我添加定义在迭代循环之外的W矩阵时，它变得非常缓慢（30次迭代需要0.3〜0.7秒），而我不明白这是我的编码问题还是通常需要这么长时间。

这里的所有内容都是Eigen矩阵和向量。

我使用Eigen库中的.asDiagonal()函数从逆方差的向量中定义了我的W矩阵，然后在计算H和b时使用它。然后速度就变得非常慢了。我希望能获取一些有关此巨大减速的潜在原因的提示。

编辑：

这里只有两个矩阵。Jacobian肯定是密集的。权重矩阵是通过函数vec.asDiagonal()从向量生成的，该函数来自密集库，因此我认为它也是密集的。

代码非常简单，唯一导致时间变化的区别是添加了权重矩阵。以下是代码片段：

for (int iter=0; iter<max_iter; ++iter) {
    // obtain error vector
    error = ...  
    // calculate H and b - the fast one
    Eigen::MatrixXf H = J.transpose() * J;
    Eigen::VectorXf b = J.transpose() * error;
    // calculate H and b - the slow one
    Eigen::MatrixXf H = J.transpose() * weight_ * J;
    Eigen::VectorXf b = J.transpose() * (weight_ * error);
    // obtain delta and update state
    del = H.ldlt().solve(b);
    T <- T(del)   // this is pseudo code, meaning update T with del
}

这是一个类中的函数，为了调试目的，权重矩阵现在被定义为一个类变量，可以被函数访问，在函数被调用之前就已经定义好了。

- CathIAS

---

T <- T(del)？？ - erip

使用del命令更新T。 - CathIAS

可能是一个重载的运算符，但这看起来很可疑...非常像 T < -T(del) 或者 T 小于 -T(del)。 - erip

你有这个问题的 MCVE 吗？ - harold

@AviGinsburg 这是一个好想法，但我现在基本上是在乘以一个身份矩阵W。H只有6×6，所以解决起来不应该是问题。 - CathIAS

Answer 1

我猜测weight_被声明为密集的MatrixXf？如果是这样，那么请在您使用weight_的任何地方都用w.asDiagonal()替换它，或者将后者作为asDiagonal表达式的别名：

auto weight = w.asDiagonal();

这样，Eigen就会知道weight是一个对角矩阵，并且计算将会按预期进行优化。

- ggael

Answer 2

因为矩阵乘法只是对角线，所以您可以将其更改为使用系数乘法，如下所示：

MatrixXd m;
VectorXd w;
w.setLinSpaced(5, 2, 6);

m.setOnes(5,5);

std::cout << (m.array().rowwise() * w.array().transpose()).matrix() << "\n";

同样地，矩阵向量乘积可以写成：

(w.array() * error.array()).matrix()

这可以避免矩阵中的零元素。没有一个我可以依据的最小可重现示例，你的结果可能会有所不同...

- Avi Ginsburg