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

Posted 2023-02-16

【中文标题】OpenMP 循环运行代码比串行循环慢

【英文标题】：OpenMP loop runs code slower than serial loop

【发布时间】：2017-03-06 08:55:58

【问题描述】：

我正在运行这个简洁的小重力模拟，在串行执行中它需要 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 的建议移除了临界区域，也没有那么严重。

【问题讨论】：

临界区看起来很邪恶。难道你不能只重新运行循环并将关键部分保留在外部而不进行并行化吗？

临界区通过我的加速a[id2][0,1,2] 解决了一个竞争条件，就像减少到一个数组上一样。而且我确实需要 id2 循环，所以我正在写正确的力量 int a。

richtung 在哪里定义？如果它是一个数组或一个指针，它会产生很大的不同。如果它是一个数组，那么 OpenMP 将为每个线程创建私有数组（如您所愿）。如果它是一个指针，那么您只会为每个线程获得一个私有指针。您是否检查过并行版本是否得到相同的答案？您在循环结束时执行 `r12=0.0;`，因此每个线程的 r12 的初始值是未定义的。

numb 有多大。您需要做足够的工作来克服 OpenMP 开销。

@Haemiltoen 我了解，但您可能不了解我的替代解决方案。只需创建另一个循环在并行循环之外，不要在那里使用 OpenMP。这对你有用吗？如果是这样，它肯定会比你那里的任何东西都要好得多。关键部分意味着您正在使用锁定和互斥锁，根据定义，这会减慢您的代码速度。

【参考方案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 进行简单的并行化。

【讨论】：

是的，这将解决临界区问题，而不是 a_local[k] 我将不得不使用 a[id][k]*(-1.0) 这意味着我运行两个循环使用 (1000^2 -1000)/2 迭代而不是一个.

@Haemiltoen 那么这肯定更快。去做吧，看看会发生什么。请记住，迭代次数并不重要。重要的是你在循环中做了什么。

@TheQuantimPhysicist 好吧，现在它的运行速度比临界区域快，但 5 分钟仍然比串行慢 1 分钟。但是我不能并行化第二个循环，否则我会得到一个竞争条件，我必须通过某种临界区域或减少到一个应该相同的数组来解决它。但是感谢您的帮助，我想我学到了一些关于 OpenMP 的知识。

@Haemiltoen 不客气。如果此问题已关闭，请通过复选标记选择上面的答案。

@Haemiltoen 你为什么去掉复选标记？怎么回事？