容器元素上的 OpenMP 缩减

Question

我有一个嵌套循环，外部迭代很少，内部迭代很多。在内部循环中，我需要计算一个总和，所以我想使用 OpenMP 缩减。外循环在容器上，因此减少应该发生在该容器的元素上。这是一个最小的人为示例：

#include <omp.h>
#include <vector>
#include <iostream>

int main(){
    constexpr int n { 128 };

    std::vector<int> vec (4, 0);
    for (unsigned int i {0}; i<vec.size(); ++i){

        /* this does not work */
        //#pragma omp parallel for reduction (+:vec[i])
        //for (int j=0; j<n; ++j)
        //  vec[i] +=j;

        /* this works */
        int* val { &vec[0] };
        #pragma omp parallel for reduction (+:val[i])
        for (int j=0; j<n; ++j)
            val[i] +=j;

        /* this is allowed, but looks very wrong. Produces wrong results
         * for std::vector, but on an Eigen type, it worked. */
        #pragma omp parallel for reduction (+:val[i])
        for (int j=0; j<n; ++j)
            vec[i] +=j;
    }
    for (unsigned int i=0; i<vec.size(); ++i) std::cout << vec[i] << " ";
    std::cout << "\n";

    return 0;
}

问题是，如果我将缩减子句写成 (+:vec[i])，我会得到错误 ‘vec’ does not have pointer or array type，它的描述性足以找到解决方法。但是，这意味着我必须引入一个新变量并稍微更改代码逻辑，而且我发现看代码应该做什么不太明显。

我的主要问题是，是否有 better/cleaner/more 标准方法来编写容器元素的缩减。

我还想知道上述代码中显示的第三种方法有点为什么以及如何工作。我实际上正在使用 Eigen 库，该变体在其容器上似乎工作得很好（尽管尚未对其进行广泛测试），但在 std::vector 上，它产生的结果介于零和实际结果 (8128)。我认为它应该有效，因为 vec[i] 和 val[i] 应该都评估为取消引用同一地址。但是很遗憾，显然不是。

我正在使用 OpenMP 4.5 和 gcc 9.3.0。

Answer 1

对于第二个示例，与其存储指针然后始终访问同一元素，不如使用局部变量：

    int val = vec[i];
    #pragma omp parallel for reduction (+:val)
    for (int j=0; j<n; ++j)
        val +=j;
    vec[i] = val;

对于第三个循环，我怀疑问题是因为 reduction 子句命名了一个变量，但是您从不在循环中用该名称更新该变量，因此编译器看不到任何减少.使用 Eigen 可能会使代码分析起来更复杂一些，从而导致循环工作。

Answer 2

我会分三部分回答你的问题：

1.在上面的示例中使用 std::vec 执行 OpenMP 缩减的最佳方法是什么？

i) 使用你的方法，即创建一个指针 int* val { &vec[0] };

ii) 声明一个新的共享变量，如@1201ProgramAlarm 已回答。

iii) 声明一个 user defined reduction（这在您的简单情况下并不适用，但请参阅下面的 3. 以获得更有效的模式）。

2。为什么第三个循环不起作用，为什么它与 Eigen 一起起作用？

就像前面的答案一样，您告诉 OpenMP 对内存地址 X 执行归约和，但您正在对内存地址 Y 执行加法，这意味着归约声明被忽略并且您的加法受到通常的线程竞争条件。

您并没有真正提供关于您的 Eigen 企业的太多细节，但这里有一些可能的解释：

i) 你并没有真正使用多线程（检查 n = Eigen::nbThreads( )）

ii) 您没有禁用 Eigen 自己的并行性，这可能会破坏您自己对 OpenMP 的使用，例如EIGEN_DONT_PARALLELIZE 编译器指令。

iii) 存在竞争条件，但您没有看到它，因为 Eigen 操作需要更长的时间，您使用的线程数量较少并且只写入少量值 => 线程干扰的发生率较低相互产生错误的结果。

3。我应该如何使用 OpenMP 并行化此场景（技术上不是您明确提出的问题）？

您不应仅并行化内部循环，而应同时并行化两者。您拥有的序列代码越少越好。在这种情况下，每个线程都有自己的 vec 向量的私有副本，在所有元素都被各自的线程求和后，它会减少。此解决方案最适合您提供的示例，但如果您使用非常大的向量和非常多的线程（或 RAM 非常有限），则可能运行出现 RAM 问题。

#pragma omp parallel for collapse(2) reduction(vsum : vec)
for (unsigned int i {0}; i<vec.size(); ++i){
    for (int j = 0; j < n; ++j) {
        vec[i] += j;
    }
}

其中 vsum 是用户定义的缩减，即

#pragma omp declare reduction(vsum : std::vector<int> : std::transform(omp_out.begin(), omp_out.end(), omp_in.begin(), omp_out.begin(), std::plus<int>())) initializer(omp_priv = decltype(omp_orig)(omp_orig.size()))

在你使用它的函数之前声明 reduction，你就可以开始了

容器元素上的 OpenMP 缩减

OpenMP reduction on container elements

c++

openmp

reduction