删除单个构造会导致执行不正确

Removing single construct results in incorrect execution

本文关键字:执行 不正确 单个构 删除      更新时间:2023-10-16

此代码按预期工作:

#include <iostream>
#include <cmath>
#include <omp.h>
//https://stackoverflow.com/questions/37970024/jacobi-relaxation-in-mpi
#define max(a, b) (a)>(b)?(a):(b)
const int m = 2001;
const int n = 1500;
const int p = 4;
double v[m + 2][m + 2];
double x[m + 2];
double y[m + 2];
double _new[m + 2][m + 2];
double maxdiffA[p + 1];
int icol, jrow;
int main() {
omp_set_num_threads(p);
double h = 1.0 / (n + 1);
double start = omp_get_wtime();
#pragma omp parallel for private(icol) shared(x, y, v, _new)
for (icol = 0; icol <= n + 1; ++icol) {
x[icol] = y[icol] = icol * h;
_new[icol][0] = v[icol][0] = 6 - 2 * x[icol];
_new[n + 1][icol] = v[n + 1][icol] = 4 - 2 * y[icol];
_new[icol][n + 1] = v[icol][n + 1] = 3 - x[icol];
_new[0][icol] = v[0][icol] = 6 - 3 * y[icol];
}
const double eps = 0.01;

#pragma omp parallel private(icol, jrow) shared(_new, v, maxdiffA)
{
while (true) { //for [iters=1 to maxiters by 2]
#pragma omp single
for (int i = 0; i < p; i++) maxdiffA[i] = 0;
#pragma omp for
for (icol = 1; icol <= n; icol++)
for (jrow = 1; jrow <= n; jrow++)
_new[icol][jrow] =
(v[icol - 1][jrow] + v[icol + 1][jrow] + v[icol][jrow - 1] + v[icol][jrow + 1]) / 4;
#pragma omp for
for (icol = 1; icol <= n; icol++)
for (jrow = 1; jrow <= n; jrow++)
v[icol][jrow] = (_new[icol - 1][jrow] + _new[icol + 1][jrow] + _new[icol][jrow - 1] +
_new[icol][jrow + 1]) / 4;
#pragma omp for
for (icol = 1; icol <= n; icol++)
for (jrow = 1; jrow <= n; jrow++)
maxdiffA[omp_get_thread_num()] = max(maxdiffA[omp_get_thread_num()],
fabs(_new[icol][jrow] - v[icol][jrow]));
#pragma omp barrier
double maxdiff = 0.0;
for (int k = 0; k < p; ++k) {
maxdiff = max(maxdiff, maxdiffA[k]);
}

if (maxdiff < eps)
break;
#pragma omp single
std::cout << maxdiff << std::endl;
}
}
double end = omp_get_wtime();
printf("start = %.16lfnend = %.16lfndiff = %.16lfn", start, end, end - start);
return 0;
}

它输出

1.12454
<106 iterations here>
0.0100436
start = 1527366091.3069999217987061
end = 1527366110.8169999122619629
diff = 19.5099999904632568

但是如果我删除

#pragma omp single
std::cout << maxdiff << std::endl;

该程序似乎运行无限长,或者我得到

start = 1527368219.8810000419616699
end = 1527368220.5710000991821289
diff = 0.6900000572204590

为什么会这样?

你在 while 循环的开头覆盖maxdiffA- 这必须与最后的读取maxdiffA隔离以检查条件。否则,一个线程可能已经在另一个线程有机会读取它们之前重置了这些值。循环末尾的omp single构造由于构造末尾的隐式屏障omp single充当隔离。但是,在omp single构造开始时没有障碍。此外,"一大堆代码">也不是一个安全屏障。因此,如果没有有效的隐式屏障,则必须使用#pragma omp barrier保护重置代码的入口。

也就是说,我强烈建议重新构建代码以具有共享退出条件,该条件也在single构造中计算。这使得所有线程进程同时退出while循环更加清晰。否则,代码定义不明确。