C语言 使用pragma omp parallel for计算pi的问题

Question

我写了下面的代码来计算圆周率的值，它可以工作:

#include <omp.h>
#include <stdio.h>
static long num_steps = 1000000;
double step;
#define NUM_THREADS 16
int main()
{
int i, nthreads;
double tdata, pi, sum[NUM_THREADS];
omp_set_num_threads(NUM_THREADS);
step = 1.0 / (double)num_steps;
tdata = omp_get_wtime();
#pragma omp parallel
{
int i, id, nthrds;
double x;
id = omp_get_thread_num();
nthrds = omp_get_num_threads();
if (id == 0)
nthreads = nthrds;
for (i = id, sum[id] = 0.0; i < num_steps; i = i + nthrds)
{
x = (i + 0.5) * step;
sum[id] = sum[id] + 4.0 / (1.0 + x * x);
}
}
tdata = omp_get_wtime() - tdata;
for (i = 0, pi = 0.0; i < nthreads; i++)
{
pi = pi + sum[i] * step;
}
printf("pi=%f and it took %f seconds", pi, tdata);
}

然后我了解到我可以使用#pragma omp parallel for，然后我不必手动将计算分解到不同的线程。所以我这样写:

#include <omp.h>
#include <stdio.h>
static long num_steps = 1000000;
double step;
#define NUM_THREADS 16
int main()
{
int i;
double tdata, pi, x, sum = 0.0;
omp_set_num_threads(NUM_THREADS);
step = 1.0 / (double)num_steps;
tdata = omp_get_wtime();
#pragma omp parallel for
{
for (i = 0; i < num_steps; i++)
{
x = (i + 0.5) * step;
sum = sum + 4.0 / (1.0 + x * x);
}
}
tdata = omp_get_wtime() - tdata;
pi = sum * step;
printf("pi = %f and compute time = %f seconds", pi, tdata);
}

然而，这不起作用并且输出错误的pi值。我做错了什么?谢谢你。

Answer 1

您的代码中有两个竞态条件会导致不正确的结果:

• 一个已经在评论中指出了，即。，线程间共享的变量sum()的更新)，可通过减量条款解决;
• 另一个是变量x的更新，它也是在线程之间共享的。这个竞争条件可以通过简单地将变量设置为私有来解决。

两种可能的解决方案:

1. 使用OpenMP的私有构造函数

   #pragma omp parallel for reduction(+:sum) private(x)
   for (i = 0; i < num_steps; i++)
   {
   x = (i + 0.5) * step;
   sum = sum + 4.0 / (1.0 + x * x);
   }
   

2. 在

   的并行范围内声明变量'x'

   #pragma omp parallel for reduction(+:sum)
   for (i = 0; i < num_steps; i++)
   {
   double x = (i + 0.5) * step;
   sum = sum + 4.0 / (1.0 + x * x);
   }
   

最终代码如下所示:

#include <omp.h>
#include <stdio.h>
static long num_steps = 1000000;
#define NUM_THREADS 16
int main()
{
double sum = 0.0;
omp_set_num_threads(NUM_THREADS);
double step = 1.0 / (double)num_steps;

double tdata = omp_get_wtime();

#pragma omp parallel for reduction(+:sum)
for (int i = 0; i < num_steps; i++)
{
double x = (i + 0.5) * step;
sum = sum + 4.0 / (1.0 + x * x);
}

tdata = omp_get_wtime() - tdata;
double pi = sum * step;
printf("pi = %f and compute time = %f seconds", pi, tdata);
}