c++ - OpenMP 缩减为 C++ 模板指定大小的数组会导致未定义的行为

Question

我是 OpenMP 的新手，但我正在尝试使用它来加速对具有大量行和少量列的二维数组条目的某些操作。同时，我正在使用缩减来计算每列中所有数组值的总和。代码看起来像这样(我稍后会解释奇怪的部分):

template <unsigned int NColumns>
void Function(int n_rows, double** X, double* Y)
{
    #pragma omp parallel for reduction(+:Y[:NColumns])
    for (int r = 0; r < n_rows; ++r)
    {
        for (int c = 0; c < NColumns; ++c)
        {
            X[r][c] = some_complicated_stuff(X[r], X[r][c]);
            Y[c] += X[r][c];
        }
    }
}

澄清一下，X 是分配在堆上的 n_rows x NColumns 大小的二维数组，Y 是 NColumns 大小的一维数组大批。 some_complicated_stuff 实际上并没有作为一个单独的函数实现，但是我在该行中对 X[r][c] 所做的操作仅取决于 X[r] [c] 和一维数组 X[r] 中的其他值。

NColumns 作为模板参数而不是作为常规参数(如 n_rows)传入的原因是当 NColumns 是在编译时已知，编译器可以更积极地优化上述函数中的内部循环。我知道当程序运行时 NColumns 将成为少数值之一，所以稍后我有类似这样的代码:

cin >> n_cols;
double** X;
double Y[n_cols];

// initialise X and Y, etc. . .

for (int i = 0; i < n_iterations; ++i)
{
    switch (n_cols)
    {
        case  2: Function< 2>(X, Y); break;
        case 10: Function<10>(X, Y); break;
        case 60: Function<60>(X, Y); break;
        default: throw "Unsupported n_cols."; break;
    }
    // . . .
    Report(i, Y); // see GDB output below
}

通过测试，我发现将此 NColumns“参数”作为模板参数而不是普通函数参数实际上可以显着提高性能。然而，我还发现，一旦出现极难发生的情况(比如，大约每 10^7 次调用 Function)，程序就会挂起——更糟糕的是，它的行为有时会因一次运行而改变下一个程序。这种情况很少发生，以至于我在隔离错误时遇到了很多麻烦，但我现在想知道这是否是因为我在我的 OpenMP 缩减中使用了这个 NColumns 模板参数。

我注意到 similar StackOverflow question询问在缩减中使用模板类型，这显然会导致未指定的行为 - OpenMP 3.0 spec说

If a variable referenced in a data-sharing attribute clause has a type derived from a template, and there are no other references to that variable in the program, then any behavior related to that variable is unspecified.

在这种情况下，使用的不是模板类型本身，但我的情况大致相同。我是不是在这里搞砸了，还是错误更有可能出现在代码的其他部分？

我正在使用 GCC 6.3.0。

如果它更有帮助，这里是 Function 中的真实代码。 X 实际上是一个扁平化的二维数组； ww 和 min_x 在别处定义:

#pragma omp parallel for reduction(+:Y[:NColumns])
for (int i = 0; i < NColumns * n_rows; i += NColumns)
{
    double total = 0;
    for (int c = 0; c < NColumns; ++c)
        if (X[i + c] > 0)
            total += X[i + c] *= ww[c];

    if (total > 0)
        for (int c = 0; c < NColumns; ++c)
            if (X[i + c] > 0)
                Y[c] += X[i + c] = (X[i + c] < min_x * total ? 0 : X[i + c] / total);
}

只是为了让情节更复杂一点，我将 gdb 附加到挂起的程序的运行进程，这是回溯显示的内容:

#0  0x00007fff8f62a136 in __psynch_cvwait () from /usr/lib/system/libsystem_kernel.dylib
#1  0x00007fff8e65b560 in _pthread_cond_wait () from /usr/lib/system/libsystem_pthread.dylib
#2  0x000000010a4caafb in omp_get_num_procs () from /opt/local/lib/libgcc/libgomp.1.dylib
#3  0x000000010a4cad05 in omp_get_num_procs () from /opt/local/lib/libgcc/libgomp.1.dylib
#4  0x000000010a4ca2a7 in omp_in_final () from /opt/local/lib/libgcc/libgomp.1.dylib
#5  0x000000010a31b4e9 in Report(int, double*) ()
#6  0x3030303030323100 in ?? ()
[snipped traces 7-129, which are all ?? ()]
#130 0x0000000000000000 in ?? ()

Report() 是一个在程序主循环内调用但不在 Function() 内调用的函数(我已将其添加到上面的中间代码片段中) ，并且 Report() 不包含任何 OpenMP 编译指示。这是否说明了正在发生的事情？

请注意，可执行文件在进程开始运行和我将 GDB 附加到它之间发生了变化，这需要引用新的(已更改的)可执行文件。所以这可能意味着符号表被搞乱了。

Answer 1

我已经设法部分解决了这个问题。

其中一个问题是程序的行为具有不确定性。这是因为 (1) OpenMP 执行线程完成顺序的缩减，这是不确定的，并且 (2) 浮点加法是非关联的。我假设减少将按线程编号顺序执行，但事实并非如此。因此，任何使用浮点运算减少的 OpenMP for 构造都可能是不确定的，即使从一次运行到下一次运行的线程数相同，只要线程数更大比 2。关于此事的一些相关 StackOverflow 问题是 here和 here .

另一个问题是程序偶尔会挂起。我一直无法解决这个问题。在挂起的进程上运行 gdb 总是会在堆栈跟踪的顶部产生 __psynch_cvwait ()。每执行 10^8 次并行化 for 循环，它就会挂起。

希望对您有所帮助。