问题描述
我刚刚开始使用我的矢量化代码。我的矩阵向量乘法代码没有被gcc自动向量化，我想知道为什么。 This pastebin contains the output from -fopt-info-vec-missed。
我无法理解输出所告诉我的内容，无法查看其与我编写的代码的匹配程度。
例如，我看到许多行说not enough data-refs in basic block，但是我无法通过Google搜索在线找到很多有关此的详细信息。我还看到存在与内存对齐有关的问题，例如UnkNown misalignment,naturally aligned和vector alignment may not be reachable。我所有的内存分配都是使用double为malloc类型分配的，我相信可以保证为该类型对齐。
环境：在WSL2上使用gcc进行编译
gcc -v: gcc version 7.5.0 (Ubuntu 7.5.0-3ubuntu1~18.04)
#include <time.h>
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>
#define N 4000 // Matrix size will be N x N
#define T 1
//gcc -fopenmp -g vectorisation.c -o main -O3 -march=native -fopt-info-vec-missed=missed.txt
void doParallelcomputation(double *A,double *V,double *results,unsigned long matrixSize,int numThreads)
{
omp_set_num_threads(numThreads);
unsigned long i,j;
#pragma omp parallel for simd private(j)
for (i = 0; i < matrixSize; i++)
{
// double *AHead = &A[i * matrixSize];
// double tmp = 0;
for (j = 0; j < matrixSize; j++)
{
results[i] += A[i * matrixSize + j] * V[j];
// also tried tmp += A[i * matrixSize + j] * V[j];
}
// results[i] = tmp;
}
}
void genRandVector(double *S,unsigned long size)
{
srand(time(0));
unsigned long i;
for (i = 0; i < size; i++)
{
double n = rand() % 5;
S[i] = n;
}
}
void genRandMatrix(double *A,unsigned long size)
{
srand(time(0));
unsigned long i,j;
for (i = 0; i < size; i++)
{
for (j = 0; j < size; j++)
{
double n = rand() % 5;
A[i*size + j] = n;
}
}
}
int main(int argc,char *argv[])
{
double *V = (double *)malloc(N * sizeof(double)); // v in our A*v = parV computation
double *parV = (double *)malloc(N * sizeof(double)); // Parallel computed vector
double *A = (double *)malloc(N * N * sizeof(double)); // NxN Matrix to multiply by V
genRandVector(V,N);
doParallelcomputation(A,V,parV,N,T);
free(parV);
free(A);
free(V);
return 0;
}
解决方法
添加double *restrict results来保证输入/输出不重叠，没有 OpenMP但带有-ffast-math。 https://godbolt.org/z/qaPh1v
您需要具体告知OpenMP减少量，以使其放松FP-数学关联性。 （-ffast-math对OpenMP矢量化器没有帮助）。有了这些，我们就能得到您想要的：
#pragma omp simd reduction(+:tmp)
仅使用restrict而没有-ffast-math或-fopenmp时，您会得到全部垃圾：它执行SIMD FP乘法，但是将4x vaddsd的内容解压缩为标量累加器，根本无法掩盖FP延迟。
使用restrict和-fopenmp（无快速运算符），它仅执行标量FMA。
使用restrict和-ffast-math（未使用-fopenmp 或#pragma进行注释），可以很好地自动矢量化：vfmadd231pd ymm循环，随机播放/在外部添加水平和。 （但不并行化）。 https://godbolt.org/z/f36oG3
使用restrict和-ffast-math（使用 -fopenmp）仍然不会自动矢量化。 OpenMP矢量化器有所不同，也许没有利用快速运算法则，而是需要您告诉它缩减的方法吗？
还要注意，在数据布局中，要并行化的循环（外部）与要通过SIMD矢量化的循环（内部）不同。内部的点积循环位于连续的内存中，因此读取它们最有意义，而不是尝试将4个不同列中的数据SIMD混洗到一个向量中以累积4个result[i+0..3]结果为1个向量。
但是，将外部循环展开4来将每个V[j+0..3]与来自4个不同列的数据一起使用会提高计算强度（每个FMA接近于1个负载，而不是2个）
（只要V[] 和一行矩阵适合L1d缓存，就可以了。如果不是，那实际上就很糟糕，应该被缓存阻塞。或者实际上如果展开外部循环，则显示矩阵的4行。）
还请注意，double tmp = 0;是个好主意：您当前的版本已添加到result[i]中，并在编写前先进行阅读。在将其用作纯输出之前，需要先进行零初始化。
自动vec自动par版本：
我认为这是正确的； asm看起来像是自动并行化以及自动矢量化了内部循环。
void doParallelComputation(double *restrict A,double *restrict V,double *restrict results,unsigned long matrixSize,int numThreads)
{
omp_set_num_threads(numThreads);
unsigned long i,j;
#pragma omp parallel for private(j)
for (i = 0; i < matrixSize; i++)
{
// double *AHead = &A[i * matrixSize];
double tmp = 0;
// TODO: unroll outer loop and cache-block it.
#pragma omp simd reduction(+:tmp)
for (j = 0; j < matrixSize; j++)
{
//results[i] += A[i * matrixSize + j] * V[j];
tmp += A[i * matrixSize + j] * V[j]; //
}
results[i] = tmp; // write-only to results,not adding to old value.
}
}
使用OpenMPified帮助函数doParallelComputation._omp_fn.0:中的矢量化内部循环编译（Godbolt）
# gcc7.5 -xc -O3 -fopenmp -march=skylake
.L6:
add rdx,1 # loop counter; newer GCC just compares the end-pointer
vmovupd ymm2,YMMWORD PTR [rcx+rax] # 32-byte load
vfmadd231pd ymm0,ymm2,YMMWORD PTR [rsi+rax] # 32-byte memory-source FMA
add rax,32 # pointer increment
cmp rdi,rdx
ja .L6
然后，循环后水平效率的总和；不幸的是，OpenMP矢量化器不如“正常” -ftree-vectorize矢量化器那样聪明，但是这需要-ffast-math才能在这里做任何事情。