c++ - OpenMP 还原同步错误

Question

我正在尝试将一个循环与相互依赖的循环并行化，我已经尝试过缩减并且代码有效，但结果是错误的，我认为缩减适用于总和但不适用于数组的更新在正确的循环中，是否有一种方法可以使循环并行化以获得正确的结果？

#pragma omp parallel for reduction(+: sum)
for (int i = 0; i < DATA_MAG; i++)
{
    sum += H[i];

    LUT[i] = sum * scale_factor;
}

Answer 1

reduction 子句为团队中的每个线程创建了 sum 的私有(private)拷贝，就像私有(private)子句已用于 sum 一样。在 for 循环之后，每个私有(private)拷贝的结果与 sum 的原始共享值相结合。由于共享 sum 仅在 for 循环之后更新，因此您不能在 for 循环内依赖它。

在这种情况下你需要做一个前缀和。不幸的是 parallel prefix-sum using threads是大 DATA_MAG 的内存带宽限制，它主要由小 DATA_MAG 的 OpenMP 开销决定。但是，在小型和大型之间可能存在一些最佳点，您可以在其中看到使用线程的一些好处。

但在您的情况下，DATA_MAG 仅为 256，这非常小，无论如何都不会从 OpenMP 中受益。您可以做的是使用 SIMD。然而，据我所知，OpenMP 的 simd 构造对于前缀和来说还不够强大。但是，您可以通过像这样手动起诉内在函数来做到这一点

#include <stdio.h>
#include <x86intrin.h>

#define N 256

inline __m128 scan_SSE(__m128 x) {
  x = _mm_add_ps(x, _mm_castsi128_ps(_mm_slli_si128(_mm_castps_si128(x), 4)));
  x = _mm_add_ps(x, _mm_castsi128_ps(_mm_slli_si128(_mm_castps_si128(x), 8)));
  return x;
}

void prefix_sum_SSE(float *a, float *s, int n, float scale_factor) {
  __m128 offset = _mm_setzero_ps();
  __m128 f = _mm_set1_ps(scale_factor);
  for (int i = 0; i < n; i+=4) {
    __m128 x = _mm_loadu_ps(&a[i]);
    __m128 out = scan_SSE(x);
    out = _mm_add_ps(out, offset);
    offset = _mm_shuffle_ps(out, out, _MM_SHUFFLE(3, 3, 3, 3));
    out = _mm_mul_ps(out, f); 
    _mm_storeu_ps(&s[i], out);
  }
}

int main(void) {
  float H[N], LUT[N];
  for(int i=0; i<N; i++) H[i] = i;
  prefix_sum_SSE(H, LUT, N, 3.14159f);
  for(int i=0; i<N; i++) printf("%.1f ", LUT[i]); puts("");
  for(int i=0; i<N; i++) printf("%.1f ", 3.14159f*i*(i+1)/2); puts("");

}

参见 here有关使用 SSE 和 AVX 的 SIMD 前缀和的更多详细信息。