c++ - 快速单精度矩阵乘以 vector 积

Question

考虑一个大小为 48x16 的 float 矩阵 A 和一个大小为 1x48 的 float vector b。

请建议一种在常见桌面处理器 (i5/i7) 上尽可能快地计算 b×A 的方法。

背景。上述产品处于紧密循环中，因此其快速计算至关重要。目前我正在使用以下朴素算法:

inline void Critical(const float A[48][16], const float b[48], float x[16] ) const {
    for (int u = 0; u < 48; ++u) {
        for (int i = 0; i < 16; ++i) {
            x[i] += A[u][i] * b[u];
        }
    }
}

我曾尝试将乘法卸载到 MKL 的 SGEMV，然后再卸载到 SGEMM，但无济于事。简单的实现在 i7 4800MQ 上仍然运行得更快。

编辑 1。

采用静态分配的本征值与朴素算法的速度大致相同。

我试过打开优化的 GCC5、ICC 和 VC2015U3(/O3、快速数学、mtune=native 等)。 ICC 似乎在 Linux 和 Windows 上都能生成最快的代码。

EDIT2。

A的元素很小，max(|A_ui|) = 256。同理max(|b_u|) = 1.0。只要算法比原始算法快，合理的近似解也是可以接受的。

Answer 1

MKL 通常有很大的开销，因此对于小矩阵的性能很差。另一方面，Eigen 具有固定大小的矩阵优化，在小矩阵上表现良好。您还需要正确的编译选项才能使用 Eigen 获得最大性能。

#include <iostream>
#include <Eigen/Eigen>
#include <omp.h>

inline void Critical(const float A[48][16], const float b[48], float x[16]) {
  for (int i = 0; i < 16; ++i) {
    x[i] = 0;
  }
  for (int u = 0; u < 48; ++u) {
    for (int i = 0; i < 16; ++i) {
      x[i] += A[u][i] * b[u];
    }
  }
}

int main() {
  float a[48][16] = { 0 };
  float b[48] = { 0 };
  float x[16] = { 0 };
  Eigen::Matrix<float, 48, 16> ma;
  Eigen::Matrix<float, 1, 48> mb;
  Eigen::Matrix<float, 1, 16> mx;

  ma.setRandom();
  mb.setRandom();
  for (int i = 0; i < 48; ++i) {
    for (int j = 0; j < 16; ++j) {
      a[i][j] = ma(i, j);
    }
    b[i] = mb(i);
  }

  double t;
  int n = 10000000;

  t = omp_get_wtime();
  for (int i = 0; i < n; ++i) {
    Critical(a, b, x);
  }
  t = omp_get_wtime() - t;
  std::cout << "for-loop time: " << t << std::endl;

  t = omp_get_wtime();
  for (int i = 0; i < n; ++i) {
    mx = mb * ma;
  }
  t = omp_get_wtime() - t;
  std::cout << "eigen time: " << t << std::endl;

  Eigen::Map < Eigen::Matrix<float, 1, 16> > native_x(x);
  std::cout << "error: " << (mx - native_x).norm() << std::endl;

  return 0;
}

使用 g++ 5.2.1 编译时

$ g++ -fopenmp -O3 -DNDEBUG -I~/program/include/eigen3 -o test/gemv test/gemv.cpp && test/gemv
for-loop time: 2.53004
eigen time: 1.17458
error: 1.49636e-06

用icpc 16.0.2编译时

$ icpc -fopenmp -fast -DNDEBUG -I~/program/include/eigen3 -o test/gemv test/gemv.cpp && test/gemv
for-loop time: 1.03432
eigen time: 1.01054
error: 1.40769e-06

icpc 在 fop-loops 上使用自动矢量化，因此性能与 Eigen 相同。