OpenMP 循环运行代码比串行循环慢

Question

我正在运行这个简洁的小重力模拟，在串行执行中它需要 4 分钟多一点，当我在 a 中并行化一个循环时，它会增加到大约 7 分钟，如果我尝试并行化更多循环，它会增加到超过20分钟。我发布了一个稍微缩短的版本，没有进行一些初始化，但我认为它们并不重要。我发布了 7 分钟版本，但是有一些我想在循环中添加并行化的 cmets。感谢您帮助我处理我乱七八糟的代码。

#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <string.h>
#include <omp.h>

#define numb 1000
int main(){
  double pos[numb][3],a[numb][3],a_local[3],v[numb][3];
  memset(v, 0.0, numb*3*sizeof(double));
  double richtung[3];
  double t,deltat=0.0,r12 = 0.0,endt=10.;
  unsigned seed;
  int tcount=0;
  #pragma omp parallel private(seed) shared(pos)
  {
    seed = 25235 + 16*omp_get_thread_num();
    #pragma omp for 
    for(int i=0;i<numb;i++){
      for(int j=0;j<3;j++){
        pos[i][j] = (double) (rand_r(&seed) % 100000 - 50000);
      }
    }
  }
  for(t=0.;t<endt;t+=deltat){
    printf("\r%le", t);
    tcount++;
    #pragma omp parallel for shared(pos,v)
    for(int id=0; id<numb; id++){
      for(int l=0;l<3;l++){
        pos[id][l] = pos[id][l]+(0.5*deltat*v[id][l]);
        v[id][l] = v[id][l]+a[id][l]*(deltat);
      }
    }
    memset(a, 0.0, numb*3*sizeof(double));
    memset(a_local, 0.0, 3*sizeof(double));
    #pragma omp parallel for private(r12,richtung) shared(a,pos)
    for(int id=0; id <numb; ++id){
      for(int id2=0; id2<id; id2++){
        for(int k=0;k<3;k++){
          r12 += sqrt((pos[id][k]-pos[id2][k])*(pos[id][k]-pos[id2][k]));
        }
        for(int k=0; k<3;k++){
          richtung[k] = (-1.e10)*(pos[id][k]-pos[id2][k])/r12;
          a[id][k] += richtung[k]/(((r12)*(r12)));
          a_local[k] += (-1.0)*richtung[k]/(((r12)*(r12)));
          #pragma omp critical
          {
            a[id2][k] += a_local[k];
          }
        }
        r12=0.0;
      }
    }
    #pragma omp parallel for shared(pos)
    for(int id =0; id<numb; id++){
      for(int k=0;k<3;k++){
        pos[id][k] = pos[id][k]+(0.5*deltat*v[id][k]);
      }
    }
    deltat= 0.01;
  }
  return 0;
}

我正在使用 g++ -fopenmp -o test_grav test_grav.c 编译代码，我只是在 shell 中测量时间 time ./test_grav。 当我使用 get_numb_threads() 获取它显示的线程数 4. top 还显示超过 300%（有时 ~380%）的 cpu 使用率。有趣的小事实，如果我在时间循环之前启动并行区域（意味着最外层的 for 循环）并且没有任何实际的 #pragma omp for 它相当于为每个主要创建一个并行区域（最外层循环的三秒）环形。所以我认为这是一个优化的事情，但我不知道如何解决它。谁能帮帮我？

编辑：我使示例可验证并降低了numb 之类的数字，以使其更好地可测试，但问题仍然存在。即使我按照 TheQuantumPhysicist 的建议移除了临界区域，也没有那么严重。

Answer 1

我认为关键部分是问题的原因。考虑将所有关键部分置于并行化循环之外，并在并行化结束后运行它们。

试试这个：

#pragma omp parallel shared(a,pos)
{
#pragma omp for private(id2,k,r12,richtung,a_local) 
for(id=0; id <numb; ++id){
    for(id2=0; id2<id; id2++){
        for(k=0;k<3;k++){
            r12 += sqrt((pos[id][k]-pos[id2][k])*(pos[id][k]-pos[id2][k]));
        }
        for(k =0; k<3;k++){
            richtung[k] = (-1.e10)*(pos[id][k]-pos[id2][k])/r12;
            a[id][k] += richtung[k]/(((r12)*(r12))+epsilon);
            a_local[k]+= richtung[k]/(((r12)*(r12))+epsilon)*(-1.0);
        }
    }
}
}
for(id=0; id <numb; ++id){
    for(id2=0; id2<id; id2++){
        for(k=0;k<3;k++){
            a[id2][k] += a_local[k];
        }
    }
}

关键部分会导致锁定和阻塞。如果您可以使这些部分保持线性，您将在性能上大获全胜。

请注意，我说的是一种句法解决方案，我不知道它是否适用于您的情况。但要明确：如果您的系列中的每个点都依赖于下一个点，那么并行化不是您的解决方案；至少使用 OpenMP 进行简单的并行化。