没有线程本地副本的 OpenMP reduce 数组

Question

我希望 OpenMP 将 large 数组缩减为较小的动态数组。例如，哪里

large[] = {1, 1, 1, 2, 2, 2, 3, 3, 3};

// OpenMP reduce length-3 sublists of large, to produce:

small[] = {3, 6, 9};

我的要求与 this question 类似，但我有重要的附加限制：

必须支持 OpenMP 3.1（所以我不能像 this answer 那样使用 OpenMP 4.5 的数组缩减）
不能拥有 small 的线程私有副本（如 this answer），因为 small 可以与 large 一样大，因此线程私有副本可能会导致堆栈溢出
large 数组必须按顺序迭代，以便进行良好的缓存（因为它可能 非常大）

其他一些细节：

large的求和元素不一定相邻
这种减少发生在一个函数内。 small（输出）数组由用户预先分配，并通过引用传递。理想情况下，该函数不应分配 small
small 和 large 的长度都是 2 的幂（例如 2、4、8、16 ...)。 small 的每个元素都是从相同数量的 large 元素（2 的另一个幂）中减去的。

为清楚起见，这里是示例串行伪代码：

void myfunc(double* small, double* large, int lLen, // other args) {

    for (int i=0; i<lLen; i++) {
    
        int sInd = // determined by other args
        
        small[sInd] += large[i];
    }
}

这是一个示例实现，由于它使用了 OpenMP 4.5 的数组缩减，因此不合格。此外，它还不受欢迎地使用了 small 的本地副本（需要使用 reduction）。

void myfunc(double* small, int sLen, double* large, int lLen, // other args) {

    double smallLocal[sLen];
    int i, sInd;

    #pragma omp parallel private (i) for
    for (i=0; i<sLen; s++)
        smallLocal[i] = 0;
    
    #pragma omp parallel private (i, sInd) reduction (+:smallLocal) for
    for (i=0; i<largeLen; i++) {

        sInd = // determined by other args
        
        smallLocal[sInd] += large[i];
    }

    #pragma omp parallel private (i) for
    for (i=0; i<sLen; s++)
        small[i] = smallLocal[i];
}

我希望我能在 OpenMP 3.1 中实现同样的事情，甚至可能通过自己管理数组元素缩减来避免 smallLocal 的本地副本。我该怎么做？

Answer 1

这里有一个不太令人满意的使用锁的解决方案，它需要一个与 small 一样大的本地 locks 数组（可能需要 malloc'd）。

void myfunc(double* small, int sLen, double* large, int lLen, // other args) {

    int i, sInd;

    // clear small
    #pragma omp parallel private (i) for
    for (i=0; i<sLen; s++)
        small[i] = 0;

    // initialise locks
#ifdef _OPENMP
    omp_lock_t locks[sLen];
    #pragma omp parallel private (i) for
    for (i=0; i<sLen; s++)
        omp_init_lock(&locks[i]);
#endif 
    
    // main loop
    #pragma omp parallel private (i, sInd) for
    for (i=0; i<largeLen; i++) {

        sInd = // determined by other args
        
        // update small with locks
#ifdef _OPENMP
        omp_set_lock(&locks[sInd]);
        small[sInd] += large[i];
        omp_unset_lock(&locks[sInd]);
#else 
        small[sInd] += large[i];
#endif
    }

    // destroy locks
#ifdef _OPENMP
    #pragma omp parallel private (i) for
    for (i=0; i<sLen; s++)
        omp_destroy_lock(&locks[i]);
#endif 
}

正如预期的那样，这比使用 OpenMP 4.5 的数组缩减要慢：

对于 sLen = 8、lLen = 1048576，使用锁 26x 较慢
对于 sLen = 32768、lLen = 1048576，使用锁 1.5x 较慢

需要说明的是，OpenMP 4.5 的数组缩减只兼容本地栈数组，不兼容动态内存。这将缩减数组的大小限制为 ~1 MiB，或 2^17 双精度数 (lLen = 131072)。相反，这个锁解决方案可以处理任何大小，因为 locks 可以 malloc'd。

Answer 2

你可以做这种在科学代码中很常见的事情：

#pragma omp parallel for
for(int i = 0; i < N/3; i++)
{
    for (int j = 0; j < 3; j++)
    {
        small[i] += large[3*i + j];
    }
}

Answer 3

您可以使用 OpenMPs atomic update 结构，它已经存在于 OpenMP 3.1 中：

void myfunc(double* small, double* large, int lLen, // other args) {

    #pragma omp parallel for
    for (int i=0; i<lLen; i++) {
    
        int sInd = // determined by other args
        #pragma omp atomic update
        small[sInd] += large[i];
    }
}

这应该比使用锁更快。

没有线程本地副本的 OpenMP reduce 数组

OpenMP reduce array without thread-local copies

c

openmp

reduction