c++ - 如何以比自动矢量化更好的性能进行手动代码矢量化以进行边缘检测

Question

我一直在关注 this coursera 类(class)，在某些时候给出了下面的代码，讲师声称矢量化是通过包含 #pragma omp simd 来完成的。内外之间for循环，因为引导矢量化很难。我如何自己矢量化类(class)中使用的代码，有没有办法比我简单地添加 #pragma omp simd 获得更好的性能？继续前进？

template<typename P>
void ApplyStencil(ImageClass<P> & img_in, ImageClass<P> & img_out) {

  const int width  = img_in.width;
  const int height = img_in.height;

  P * in  = img_in.pixel;
  P * out = img_out.pixel;

  for (int i = 1; i < height-1; i++)
    for (int j = 1; j < width-1; j++) {
      P val = -in[(i-1)*width + j-1] -   in[(i-1)*width + j] - in[(i-1)*width + j+1] 
    -in[(i  )*width + j-1] + 8*in[(i  )*width + j] - in[(i  )*width + j+1] 
    -in[(i+1)*width + j-1] -   in[(i+1)*width + j] - in[(i+1)*width + j+1];

      val = (val < 0   ? 0   : val);
      val = (val > 255 ? 255 : val);

      out[i*width + j] = val;
    }

}

template void ApplyStencil<float>(ImageClass<float> & img_in, ImageClass<float> & img_out);

我正在使用 进行编译gcc 与 -march=native -fopenmp 的标志AVX512 支持 Skylake 处理器。

❯ gcc -march=native -Q --help=target|grep march
  -march=                           skylake

❯ gcc -march=knl -dM -E - < /dev/null | egrep "SSE|AVX" | sort
#define __AVX__ 1
#define __AVX2__ 1
#define __AVX512CD__ 1
#define __AVX512ER__ 1
#define __AVX512F__ 1
#define __AVX512PF__ 1
#define __SSE__ 1
#define __SSE2__ 1
#define __SSE2_MATH__ 1
#define __SSE3__ 1
#define __SSE4_1__ 1
#define __SSE4_2__ 1
#define __SSE_MATH__ 1
#define __SSSE3__ 1

Answer 1

这是一些未经测试的概念验证实现，它使用 4 个添加、1 个 fmsub 和每个数据包 3 个负载(而不是 9 个负载、7 个添加、1 个 fmsub，用于直接实现)。我省略了钳位(对于 float 图像至少看起来不寻常，而对于 uint8 它什么也不做，除非您将 P val = ... 更改为 auto val = ... ，正如彼得在评论中注意到的那样)--但您可以轻松自己添加。

此实现的想法是将左右像素 ( x0_2 ) 以及所有 3 个 ( x012 ) 相加，并将它们从 3 个连续行 ( a012 + b0_2 + c012 ) 中相加，然后从中间像素乘以8.
在每个循环结束时删除 a012 的内容并移动 bX至 aX和 cX至 bX为下一次迭代。
applyStencil函数只是为每列 16 个像素应用第一个函数(从 col = 1 开始，最后只对最后 16 列执行可能重叠的计算)。如果您的输入图像少于 18 列，您需要以不同方式处理(可能通过屏蔽加载/存储)。

#include <immintrin.h>

void applyStencilColumn(float const *in, float *out, size_t width, size_t height)
{
  if(height < 3) return; // sanity check
  float const* last_in = in + height*width;
  __m512 a012, b012, b0_2, b1;
  __m512 const eight = _mm512_set1_ps(8.0);
  {
    // initialize first rows:
    __m512 a0 = _mm512_loadu_ps(in-1);
    __m512 a1 = _mm512_loadu_ps(in+0);
    __m512 a2 = _mm512_loadu_ps(in+1);
    a012 = _mm512_add_ps(_mm512_add_ps(a0,a2),a1);
    in += width;
    __m512 b0 = _mm512_loadu_ps(in-1);
    b1 = _mm512_loadu_ps(in+0);
    __m512 b2 = _mm512_loadu_ps(in+1);
    b0_2 = _mm512_add_ps(b0,b2);
    b012 = _mm512_add_ps(b0_2,b1);
    in += width;
  }
  // skip first row for output:
  out += width;

  for(; in<last_in; in+=width, out+=width)
  {
    // precalculate sums for next row:
    __m512 c0 = _mm512_loadu_ps(in-1);
    __m512 c1 = _mm512_loadu_ps(in+0);
    __m512 c2 = _mm512_loadu_ps(in+1);
    __m512 c0_2 = _mm512_add_ps(c0,c2);
    __m512 c012 = _mm512_add_ps(c0_2, c1);

    __m512 outer = _mm512_add_ps(_mm512_add_ps(a012,b0_2), c012);
    __m512 result = _mm512_fmsub_ps(eight, b1, outer);

    _mm512_storeu_ps(out, result);
    // shift/rename registers (with some unrolling this can be avoided entirely)
    a012 = b012;
    b0_2 = c0_2; b012 = c012; b1 = c1;
  }
}

void applyStencil(float const *in, float *out, size_t width, size_t height)
{
  if(width < 18) return; // assert("special case of narrow image not implemented");

  for(size_t col = 1; col < width - 18; col += 16)
  {
    applyStencilColumn(in + col, out + col, width, height);
  }
  applyStencilColumn(in + width - 18, out + width - 18, width, height);
}

可能的改进(留作练习):

applyStencilColumn可以作用于 32、48、64、... 像素的列以获得更好的缓存位置(只要您有足够的寄存器)。当然，这使得实现这两个功能稍微复杂一些。

如果您展开 for(; in<last_in; in+=width) 的 3(或 6、9、...)次迭代循环，就不需要实际移动寄存器(加上展开的一般好处)。

如果您的宽度是 16 的倍数，您可以确保至少商店大部分对齐(第一列和最后一列除外)。

您可以同时迭代少量行(通过向主函数添加另一个外部循环并以固定高度调用 applyStencilColumn。确保行集之间有适当的重叠量。(理想的行数可能取决于图像的大小。

您也可以始终添加 3 个连续像素，但将中心像素乘以 9( 9*b1-outer )。然后(通过一些簿记工作)您可以添加 row0+(row1+row2)和 (row1+row2)+row3获取 row1和 row2中间结果(有 3 个而不是 4 个加法)。不过，水平地做同样的事情看起来更复杂。

当然，您应该始终对任何自定义 SIMD 实现与编译器从通用实现生成的内容进行测试和基准测试。