GCC 暗示向量化

Question

我希望 GCC 向量化以下代码。-fopt-info 告诉我 GCC 当前没有。我相信问题是 W 的跨步访问或 k 的可能向后递增。请注意，height 和 width 是常量，index_type 当前设置为 unsigned long。

我删除了一些评论

114  for (index_type k=height-1;k+1>0;k--) {
116    for (index_type i=0;i<width;i++) {
117      Yp[k*width + i] = 0.0;                                                            
119      for (index_type j=0;j<width;j++) {                                                                            
121        Yp[k*width + i] += W[k*width*width + j*width + i]*Yp[(k+1)*width + j];
122      }
123      Yp[k*width + i] *= DF(Ap[k*width + i]);
124    }
125  }

我正在编译 gcc -O3 -ffast-math -fopt-info -std=c11 ./neural.c -o neural -lm

有什么好的方法可以把这个矢量化吗？你能告诉我更多信息吗？

我的索引方法是不是一个坏主意（即 k*width*width + ...）？我需要动态分配，并且我认为在内存中保持接近会更好地实现优化。

编辑：这可能有用

这些行的 -fopt-info-missed 输出

./neural.c:114:3: note: not vectorized: multiple nested loops.
./neural.c:114:3: note: bad loop form.
./neural.c:116:5: note: not vectorized: control flow in loop.
./neural.c:116:5: note: bad loop form.
./neural.c:119:7: note: step unknown.
./neural.c:119:7: note: reduction used in loop.
./neural.c:119:7: note: Unknown def-use cycle pattern.
./neural.c:119:7: note: not vectorized: complicated access pattern.
./neural.c:119:7: note: bad data access.
./neural.c:110:21: note: not vectorized: not enough data-refs in basic block.
./neural.c:110:58: note: not vectorized: not enough data-refs in basic block.
./neural.c:110:62: note: not vectorized: not enough data-refs in basic block.
./neural.c:117:18: note: not vectorized: not enough data-refs in basic block.
./neural.c:114:37: note: not vectorized: not enough data-refs in basic block.

编辑：

最小的例子是 HERE

我正在尝试使用 BLAS。在最小的例子中，它运行得更快，但在整个代码中它更慢......不知道为什么

编辑：

编译器正在优化代码。固定的。 BLAS 现在更快了。修复是针对整个代码，而不是最小的示例。

编辑：

与上次编辑link中的代码相同

#include <math.h>
#include <cblas.h>
#include <stdlib.h>
#include <stdio.h>

typedef float value_type;
typedef unsigned long index_type;

static value_type F(value_type v) {
  return 1.0/(1.0 + exp(-v));
}

static value_type DF(value_type v) {
  const value_type Ev = exp(-v);
  return Ev/((1.0 + Ev)*(1.0 + Ev));
}

#ifndef WITH_BLAS

static void get_Yp(const value_type * __restrict__ Ap, const value_type * __restrict__ W,
           value_type * __restrict__ Yp, const value_type * __restrict__ Dp,
           const index_type height, const index_type width) {
  for (index_type i=0;i<width;i++) {
    Yp[height*width + i] = 2*DF(Ap[height*width + i])*(Dp[i] - F(Ap[height*width + i]));
  }

  for (index_type k=height-1;k+1>0;k--) {
    for (index_type i=0;i<width;i++) {
      Yp[k*width + i] = 0.0;
      for (index_type j=0;j<width;j++) {
    Yp[k*width + i] += W[k*width*width + j*width + i]*Yp[(k+1)*width + j];
      }
      Yp[k*width + i] *= DF(Ap[k*width + i]);
    }
  }
}

#else

static void get_Yp(const value_type * __restrict__ Ap, const value_type * __restrict__ W,
           value_type * __restrict__ Yp, const value_type * __restrict__ Dp,
           const index_type height, const index_type width) {
  for (index_type i=0;i<width;i++) {
    Yp[height*width + i] = 2*DF(Ap[height*width + i])*(Dp[i] - F(Ap[height*width + i]));
  }

  for (index_type k=height-1;k+1>0;k--) {
    cblas_sgemv(CblasRowMajor, CblasTrans, width, width, 1,
        W+k*width*width, width, Yp+(k+1)*width, 1, 0, Yp+k*width, 1);
    for (index_type i=0;i<width;i++)
      Yp[k*width + i] *= DF(Ap[k*width + i]);
  }
}

#endif

int main() {
  const index_type height=10, width=10000;

  value_type *Ap=malloc((height+1)*width*sizeof(value_type)),
    *W=malloc(height*width*width*sizeof(value_type)),
    *Yp=malloc((height+1)*width*sizeof(value_type)),
    *Dp=malloc(width*sizeof(value_type));

  get_Yp(Ap, W, Yp, Dp, height, width);
  printf("Done %f\n", Yp[3]);

  return 0;
}

Answer 1

根据我的经验，要求 GCC 正确矢量化是一件很麻烦的事情。特别是如果您希望充分利用现代硬件（例如 AVX2）。我在我的代码中处理了很多矢量化问题。我的解决方案：根本不要尝试使用 GCC 进行矢量化。而是根据 BLAS（基本线性代数子程序）调用和 link 使用为现代硬件实现 BLAS 的 OpenBLAS 库来制定您想要矢量化的所有操作。 OpenBLAS 还在共享内存节点（使用 pthreads 或 OpenMP）上提供并行化，这取决于您的应用程序 - 对于进一步提高执行速度非常重要：通过这种方式，您可以将矢量化与（共享内存）并行化结合起来。

此外，我建议您调整内存以充分利用 AVX and/or AVX2 等。 IE。不要使用 malloc 或 new 分配内存，使用 memalign 或 aligned_alloc（取决于您的系统支持）。例如，如果您打算使用 AVX2，您应该对齐您的分配，以便地址是 64 字节的倍数（8 * 8 双精度）。

Answer 2

j-loop 是 可向量化 SIMD 缩减循环，具有 "width" 元素的恒定步幅。您可以使用现代编译器对其进行矢量化。此代码可使用英特尔编译器进行矢量化处理，并且在某些情况下应由 GCC 进行矢量化处理。
首先，Reduction是"vectorizable"trueloop-carried依赖的特例。所以你不能安全地向量化它，除非 "reduction" 模式是 (a) auto-recognized 由编译器（不是那么容易，严格来说不是 valid/expected 行为）或 (b) 由开发人员传达给明确使用 OpenMP 或类似标准的编译器。

为了"communicate"为了编译器有一个减少 - 你需要使用 #pragma omp simd reduction (+ : variable_name) 在 j-loop.

之前

仅从 OpenMP4.0 开始支持。所以你必须使用支持 OpenMP4.x 的 GCC 版本。引自 https://gcc.gnu.org/wiki/openmp："GCC 4.9 supports OpenMP 4.0 for C/C++, GCC 4.9.1 also for Fortran"

我还会临时使用局部变量来累积减少量（OpenMP4.0 需要这样使用减少量变量）：

 tmpSUM = 0; 
 #pragma omp simd reduction (+: tmpSUM) 
 for (index_type j=0;j<width;j++) {                                                                            
        tmpSUM += W[k*width*width + j*width + i]*Yp[(k+1)*width + j];
      }
 Yp[k*width + i] = tmpSUM

我还建议使用 signed int 而不是 unsigned，因为 unsinged 归纳变量对于所有现代向量化器来说都非常糟糕，至少会引入额外的开销。如果使用 unsigned 是主要原因之一，我不会感到惊讶，"confused" GCC。
现在，您可以 not-satisfied 接受我的回复，因为它说明了它应该如何工作（以及它在 ICC/ICPC 中的工作方式）。它没有考虑 GCC 的特定细微差别（对于减少来说似乎是反对的），正如在 GCC 优化报告中看到的那样。

所以，如果你仍然局限于 GCC，我建议：

确保它足够新鲜 GCC（4.9 或更高版本）
使用带符号的归纳变量并仍然尝试对临时本地 tmp SUM 进行 omp simd 缩减（因为它应该启用更高级的向量化技术）
如果 none 以上有帮助，请查看 "weird" 此处描述的内容（与您的情况非常相似）： What do gcc's auto-vectorization messages mean? 或者考虑使用其他 compiler/compiler 版本。
1. 最后评论：您的代码中的访问模式以及更普遍的 const-stride 是否如此糟糕？回答：对于某些 platforms/compilers const-stride 来说 而不是 会破坏你的表现。但是，理想情况下你需要更多 cache-friendly 算法。检查 Not sure how to explain some of the performance results of my parallelized matrix multiplication code 。如果您的代码确实 memory-bound 并且您没有时间自己处理内存优化，则另一种选择是考虑 MKL 或 BLAS（如其他回复所建议的那样）。

Answer 3

我在 GCC 5.3.1、Clang 3.7.1、ICC 13.0.1 和 MSVC 2015 中测试了以下代码

void foo(float *x) 
    unsigned i;
    for (i=0; i<1024; i++) x[0] += x[1024 + i];
}

我对 GCC、Clang 和 ICC 使用 -Ofast，对 MSVC 使用 /O2 /fp:fast。查看程序集表明只有 ICC 设法对循环进行矢量化。

但是，使用相同的编译选项，所有编译器都对以下代码进行了向量化

void foo2(float *x) {
    float sum = 0.0f;
    unsigned i;
    for (i=0; i<1024; i++) sum += x[1024 + i];
    x[0] = sum;
}

我不确定为什么只有 ICC 矢量化 foo。我似乎很清楚没有依赖性。

GCC 不展开循环（和 -funroll-loops does not help）。 MSVC 展开两次，Clang 展开四次，ICC 展开八次。 Intel 处理器，因为至少 Core2 至少有 3 个周期的添加延迟，所以展开四次或更多次将比两次或 none 更好地工作。

无论如何使用

value_type sum = 0.0;
for (index_type j=0;j<width;j++) {
    sum += W[k*width*width + j*width + i]*Yp[(k+1)*width + j];
}
Yp[k*width + i] = sum*DF(Ap[k*width + i]);

矢量化你的循环。

GCC 暗示向量化

GCC Hinting at Vectorization

gcc

vectorization

auto-vectorization