OpenMP 和 C++：私有变量

Question

我对 OpenMP 和 C++ 很陌生，也许正因为如此，我遇到了一些非常基本的问题。

我正在尝试做一个所有变量都是私有的静态调度（以防万一，为了验证获得的结果与非并行的结果相同）。

当我看到诸如 bodies 之类的变量时，我不知道它们来自哪里，因为它们以前没有定义。

是否可以将所有出现的变量，如bodies，都定义为私有的？怎么做到的

  std::vector<phys_vector> forces(bodies.size());

  size_t i, j; double dist, f, alpha;


  #pragma omp parallel for schedule(static) private(i, j, dist, f, alpha)
  for (i=0; i<bodies.size(); ++i) {
    for (j = i+1; j<bodies.size(); ++j) {
      dist = distance(bodies[i], bodies[j]);
      if (dist > param.min_distance()) {
        f = attraction(bodies[i], bodies[j], param.gravity(), dist);
        alpha = angle(bodies[i],bodies[j]);
        phys_vector deltaf{ f * cos(alpha) , f * sin(alpha) };
        forces[i] += deltaf;
        forces[j] -= deltaf;
      }
    }
  }
  return forces;
}

PS：以当前代码，执行结果与非并行执行不同

Answer 1

需要重申的是，您的 bodies 变量不是随机出现的，而是随机出现的。您应该确切地找出它声明的位置以及定义的内容。但是，因为您只访问 bodies 的元素并且从不更改它们，所以这个变量应该是 shared 无论如何，所以不是您的问题。

您的实际问题来自 forces 变量。您必须确保不同的线程不会更改同一个 j 的变量 forces[j]。如果您遵循循环的逻辑，则可以确保 forces[i] 仅由不同的线程访问，因此它们之间没有争用。但是对于同一个 j 的 forces[j] 可以很容易地通过并行 i 循环的不同迭代进行修改。您需要做的是 reduce on your array 按照 Whosebug link.

的答案之一

Answer 2

NoseKnowsAll 已正确识别您的问题。

我想详细解释一下为什么会出现这个问题。你可以用这样的方形循环来完成这个：

#pragma omp parallel for
for(int i=0; i<n; i++) {
    if(i==j) continue;
    phys_vector sum = 0;
    for(int j=0; j<n; j++) {
        //calculate deltaf
        sum += deltaf;
    }
    forces[i] = sum;
}

它使用 n*(n-1) 次迭代并且易于并行化。

但是因为 force(i,j) = -force(j,i) 我们可以在一半的迭代中做到这一点，n*(n-1)/2，使用三角循环（这就是你所做的）：

for(int i=0; i<n; i++) {
    phys_vector sum = 0;
    for(int j=i+1; j<n; j++) {
        //calculate deltaf
        sum += deltaf;
        forces[j] -= deltaf;
    }
    forces[i] = sum;
}

问题是，当您进行此优化时，它会使并行化外循环变得更加困难。有两个问题：写入 forces[j] 并且迭代不再分布均匀，即第一个线程运行比最后一个线程进行了更多的迭代。

简单的解决方案是将内循环并行化

#pragma omp parallel
for(int i=0; i<n; i++) {
    phys_vector sum = 0;
    #pragma omp for
    for(int j=i+1; j<n; j++) {
        //calculate deltaf
        sum += deltaf;
        forces[j] -= deltaf;
    }
    #pragma omp critical
    forces[i] += sum;
}

这使用了总共 n*(n-1)/2 次迭代中的 n*nthreads 次关键操作。因此，随着 n 变大，关键操作的成本会变小。您可以为每个线程使用私有 forces 向量并将它们合并到关键部分，但我认为这不是必需的，因为关键操作在外循环而不是内循环。

这是一个解决方案，它 fuses the triangular loop 允许每个线程运行经过相同的迭代次数。

unsigned n = bodies.size();
unsigned r = n*(n-1)/2;
#pragma omp parallel
{
    std::vector<phys_vector> forces_local(bodies.size());
    #pragma omp for schedule(static)
    for(unsigned k=0; k<r; k++) {
        unsigned i  = (1 + sqrt(1.0+8.0*k))/2;
        unsigned j = i - k*(k-1)/2;
        //calculate deltaf
        forces_local[i] += deltaf;
        forces_local[j] -= deltaf;
    }
    #pragma omp critical
    for(unsigned i=0; i<n; i++) forces[i] += forcs_local[i];
}

我对以前融合三角形的方法不满意（因为它需要使用浮点数和 sqrt 函数）所以我想出了一个基于 this answer.[=34= 的更简单的解决方案]

这会将三角形映射到矩形，反之亦然。首先，我转换为宽度为 n 但宽度为 n*(n-1)/2 的矩形（与三角形相同）。然后我计算矩形的（行，列）值，然后映射到三角形（跳过对角线）我使用以下公式

//i is the row, j is the column of the rectangle
if(j<=i) {
    i = n - i - 2;
    j = n - j - 1;
}

让我们选择一个例子。考虑以下 n=5 个三角形环对

(0,1), (0,2), (0,3), (0,4)
       (1,2), (1,3), (1,4)
              (2,3), (2,4)
                     (3,4)

将其映射到矩形变为

(3,4), (0,1), (0,2), (0,3), (0,4)
(2,4), (2,3), (1,2), (1,3), (1,4)

具有偶数值的三角形循环以相同的方式工作，尽管它可能不那么明显。例如 n = 4.

(0,1), (0,2), (0,3)
       (1,2), (1,3)
              (2,3)

这变成了

(2,3), (0,1), (0,2), (0,3)
(1,2), (1,3)

这不完全是一个矩形，但映射效果相同。我本可以将其映射为

 (0,1), (0,2), (0,3)
 (2,3), (1,2), (1,3)

这是一个矩形，但我需要两个公式来计算奇数和偶数三角形的大小。

这是使用矩形到三角形映射的新代码。

unsigned n = bodies.size();
#pragma omp parallel
{
    std::vector<phys_vector> forces_local(bodies.size());
    #pragma omp for schedule(static)
    for(unsigned k=0; k<n*(n-1)/2; k++) {
        unsigned i = k/n;
        unsigned j = k%n;
        if(j<=i) {
            i = n - i - 2;
            j = n - j - 1;
        }
        //calculate deltaf
        forces_local[i] += deltaf;
        forces_local[j] -= deltaf;
    }
    #pragma omp critical
    for(unsigned i=0; i<n; i++) forces[i] += forcs_local[i];
}

OpenMP 和 C++：私有变量

OpenMP and C++: private variables

c++

multithreading

openmp