在不同机器上具有不同输出的并行代码

Question

当我在两台不同的机器上运行下面的代码时，我得到不同的输出，其中一台输出正确（sum = sum2），另一台则不正确。

我不知道为什么

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

int main(){

    const int NX=1000;    
    const int NY=1000;

    float x[NX+2];          
    float y[NX+2];          
    float u[NX+2][NY+2];    

    float x2;   // 
    float y2;
    float sum;
    float sum2;

    for (int i=0; i<NX+2; i++){
      for (int j=0; j<NY+2; j++){
        x2      = i;
        y2      = j;
        u[i][j] = x2+ y2;
        sum += u[i][j];
      }
    }
    for (int i=0; i<NX+2; i++){
      #pragma omp parallel for
      for (int j=0; j<NY+2; j++){
        x2      = i;
        y2      = j;
        u[i][j] = x2+ y2;
      }
    }

    for (int i=0; i<NX+2;i++){
      for (int j=0; j<NY+2; j++){
        sum2 += u[i][j];
      }
    }

    printf("%f \n", sum);
    printf("%f", sum2);
}

Answer 1

你需要初始化的值

float sum;
float sum2;

否则当操作：

sum += u[i][j];

和

sum2 += u[i][j];

转到undefined behaviour。这就是为什么您会看到两种不同的结果。

将两个变量都设置为零：

float sum = 0;
float sum2 = 0;

使用（至少）标志 -Wall 编译您的代码。如果您这样做了，您会看到以下警告：

main.c:17:7: warning: 'sum2' may be used uninitialized in this function [-Wmaybe-uninitialized]
   17 | float sum2;
      |       ^~~~
main.c:16:7: warning: 'sum' may be used uninitialized in this function [-Wmaybe-uninitialized]
   16 | float sum;
      |       ^~~

性能方面而不是并行化内部循环：

for (int i=0; i<NX+2; i++){
   #pragma omp parallel for
   for (int j=0; j<NY+2; j++){
      x2      = i;
      y2      = j;
      u[i][j] = x2+ y2;
    }
}

您应该分析使用 OpenMP collapse 选项并行化两个循环时发生的情况

#pragma omp parallel for collapse(2)
for (int i=0; i<NX+2; i++){       
   for (int j=0; j<NY+2; j++){
      u[i][j] = i + j;
    }
}

即使collapse 子句不是一个意见（例如， 它更慢），在性能方面，并行化外循环而不是内循环仍然会更好。首先，您避免了创建并行区域NX+2 次的开销。其次，由于外循环在列上迭代，而内循环在行上迭代，将第一个循环的迭代在线程之间划分可降低false-sharing的可能性。

此外，您还可以并行化其他两个循环。但是，您需要使用 OpenMP reduction clause 来避免在 sum 和 sum2 变量更新期间出现竞争条件。

最终代码如下所示：

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

int main(){

   const int NX=1000;    
   const int NY=1000;
   
   float u[NX+2][NY+2];    
   float sum = 0;
   float sum2 = 0;

   #pragma omp parallel for reduction(+:sum)
   for (int i=0; i<NX+2; i++){
     for (int j=0; j<NY+2; j++){
       sum += i+j;
     }
   }
   #pragma omp parallel for
   for (int i=0; i<NX+2; i++){
     for (int j=0; j<NY+2; j++){
        u[i][j] = i+j;
     }
   }

   #pragma omp parallel for reduction(+:sum2)
   for (int i=0; i<NX+2;i++){
     for (int j=0; j<NY+2; j++){
       sum2 += u[i][j];
     }
   }

   printf("%f \n", sum);
   printf("%f", sum2);
}