在对数时间内平行减少

Question

给定n部分和,可以将log2并行步骤中的所有部分和相加.例如,假设有八个线程与八个部分和:s0, s1, s2, s3, s4, s5, s6, s7.这可以在这样的log2(8) = 3连续步骤中减少;

thread0     thread1    thread2    thread4
s0 += s1    s2 += s3   s4 += s5   s6 +=s7
s0 += s2    s4 += s6
s0 += s4

我想用OpenMP做这个,但我不想使用OpenMP的reduction子句.我想出了一个解决方案,但我认为可以使用OpenMP的task子句找到更好的解决方案.

这比标量加法更通用.让我选择一个更有用的情况:一个数组减少(见这里,这里,并在这里为更多关于阵列减少).

假设我想在阵列上进行数组缩减a.下面是一些代码,它们为每个线程并行填充私有数组.

int bins = 20;
int a[bins];
int **at;  // array of pointers to arrays
for(int i = 0; i<bins; i++) a[i] = 0;
#pragma omp parallel
{
    #pragma omp single   
    at = (int**)malloc(sizeof *at * omp_get_num_threads());        
    at[omp_get_thread_num()] = (int*)malloc(sizeof **at * bins);
    int a_private[bins];
    //arbitrary function to fill the arrays for each thread
    for(int i = 0; i<bins; i++) at[omp_get_thread_num()][i] = i + omp_get_thread_num();
}

在这一点上,我有一个指向每个线程数组的指针数组.现在我想将所有这些数组加在一起并将最终总和写入a.这是我提出的解决方案.

#pragma omp parallel
{
    int n = omp_get_num_threads();
    for(int m=1; n>1; m*=2) {
        int c = n%2;
        n/=2;
        #pragma omp for
        for(int i = 0; i<n; i++) {
            int *p1 = at[2*i*m], *p2 = at[2*i*m+m];
            for(int j = 0; j<bins; j++) p1[j] += p2[j];
        }
        n+=c;
    }
    #pragma omp single
    memcpy(a, at[0], sizeof *a*bins);
    free(at[omp_get_thread_num()]);
    #pragma omp single
    free(at);
}

让我试着解释这段代码的作用.我们假设有八个线程.让我们定义+=运算符来表示对数组求和.例如s0 += s1是

for(int i=0; i<bins; i++) s0[i] += s1[i]

那么这段代码就行了

n   thread0     thread1    thread2    thread4
4   s0 += s1    s2 += s3   s4 += s5   s6 +=s7
2   s0 += s2    s4 += s6
1   s0 += s4

但是这个代码并不像我想的那样理想.

一个问题是存在一些需要所有线程同步的隐式障碍.这些障碍不是必要的.第一个障碍是在填充阵列和进行缩减之间.第二个障碍是在#pragma omp for减少的声明中.但是我不能使用nowait这个方法来删除障碍.

另一个问题是有几个线程不需要使用.例如,有八个线程.减少的第一步只需要四个线程,第二步需要两个线程,最后一步只需要一个线程.但是,这种方法将涉及减少中的所有八个线程.虽然,其他线程无论如何都没有做多少,应该直接进入障碍并等待,这可能不是什么大问题.

我的直觉是使用omp task子句可以找到更好的方法.不幸的是,我对这个task条款几乎没有经验,到目前为止我所做的所有努力都比我现在失败的更好.

有人可以提出一个更好的解决方案来使用例如OpenMP的task条款来减少对数时间吗？

---

我找到了一种解决障碍问题的方法.这会异步减少.唯一剩下的问题是它仍然将不参与减少的线程放入繁忙的循环中.这个方法使用类似堆栈的东西将关键部分中的指针推送到堆栈(但从不弹出它们)(这是关键部分没有隐式障碍的关键之一.堆栈是串行操作但是并行减少.

这是一个有效的例子.

#include <stdio.h>
#include <omp.h>
#include <stdlib.h>
#include <string.h>

void foo6() {
    int nthreads = 13;
    omp_set_num_threads(nthreads);
    int bins= 21;
    int a[bins];
    int **at;
    int m = 0;
    int nsums = 0;
    for(int i = 0; i<bins; i++) a[i] = 0;
    #pragma omp parallel
    {
        int n = omp_get_num_threads();
        int ithread = omp_get_thread_num();
        #pragma omp single
        at = (int**)malloc(sizeof *at * n * 2);
        int* a_private = (int*)malloc(sizeof *a_private * bins);

        //arbitrary fill function
        for(int i = 0; i<bins; i++) a_private[i] = i + omp_get_thread_num();

        #pragma omp critical (stack_section)
        at[nsums++] = a_private;

        while(nsums<2*n-2) {
            int *p1, *p2;
            char pop = 0;
            #pragma omp critical (stack_section)
            if((nsums-m)>1) p1 = at[m], p2 = at[m+1], m +=2, pop = 1;
            if(pop) {
                for(int i = 0; i<bins; i++) p1[i] += p2[i];
                #pragma omp critical (stack_section)
                at[nsums++] = p1;
            }
        }

        #pragma omp barrier
        #pragma omp single
        memcpy(a, at[2*n-2], sizeof **at *bins);
        free(a_private);
        #pragma omp single
        free(at);
    }
    for(int i = 0; i<bins; i++) printf("%d ", a[i]); puts("");
    for(int i = 0; i<bins; i++) printf("%d ", (nthreads-1)*nthreads/2 +nthreads*i); puts("");
}

int main(void) {
    foo6();
}

我觉得可以找到一个更好的方法,使用的任务不会在繁忙的循环中使用线程.

Answer 1

实际上,使用递归的分而治之的方法实现干净利落的任务非常简单.这几乎是教科书代码.

void operation(int* p1, int* p2, size_t bins)
{
    for (int i = 0; i < bins; i++)
        p1[i] += p2[i];
}

void reduce(int** arrs, size_t bins, int begin, int end)
{
    assert(begin < end);
    if (end - begin == 1) {
        return;
    }
    int pivot = (begin + end) / 2;
    /* Moving the termination condition here will avoid very short tasks,
     * but make the code less nice. */
#pragma omp task
    reduce(arrs, bins, begin, pivot);
#pragma omp task
    reduce(arrs, bins, pivot, end);
#pragma omp taskwait
    /* now begin and pivot contain the partial sums. */
    operation(arrs[begin], arrs[pivot], bins);
}

/* call this within a parallel region */
#pragma omp single
reduce(at, bins, 0, n);

据我所知,没有不必要的同步,并且关键部分没有奇怪的轮询.它也可以自然地使用不同于您的排名数据的数据大小.我发现它非常干净,易于理解.所以我确实认为这比你的两个解决方案都要好.

但让我们来看看它在实践中的表现*.为此,我们可以使用Score-p和Vampir:

*_{bins=10000所以减少实际上需要一点时间.在没有涡轮增压的24核Haswell系统上执行.gcc 4.8.4 , -O3. 我在实际执行周围添加了一些缓冲区来隐藏初始化/后处理}

该图显示了应用程序中任何线程在水平时间轴上发生的情况.从上到下的树实现:

1. omp for 环
2. omp critical 任务.
3. omp task

这很好地展示了具体实现如何实际执行.现在看来for循环实际上是最快的,尽管有不必要的同步.但是这种性能分析仍然存在许多缺陷.例如,我没有固定线程.在实践中,NUMA(非统一内存访问)非常重要:核心是否在它自己的套接字的自己的缓存/内存中有这些数据？这是任务解决方案变得不确定的地方.在简单的比较中不考虑重复之间非常显着的差异.

如果还原操作在运行时变为可变,那么任务解决方案将变得比同步的for循环更好.

该critical解决方案有一些有趣的方面,被动线程不连续的等待,所以他们会更容易消耗CPU资源.这可能对性能有害,例如在turbo模式的情况下.

请记住,task通过避免立即返回的产生任务,该解决方案具有更多优化潜力.这些解决方案的执行方式也高度依赖于特定的OpenMP运行时.英特尔的运行时似乎对任务更糟糕.

我的建议是:

• 以最佳算法复杂度实施最易维护的解决方案
• 测量代码的哪些部分对运行时实际重要
• 根据实际测量结果分析瓶颈是什么.根据我的经验,它更多的是关于NUMA和调度,而不是一些不必要的障碍.
• 根据您的实际测量值执行微观优化

线性解决方案

下面是线性的时间表proccess_data_v1,从这个问题.

OpenMP 4减少

所以我想到了OpenMP减少.棘手的部分似乎是在at没有副本的情况下从循环内的数组中获取数据.我用初始化worker数组,NULL并且只是第一次移动指针:

void meta_op(int** pp1, int* p2, size_t bins)
{
    if (*pp1 == NULL) {
        *pp1 = p2;
        return;
    }
    operation(*pp1, p2, bins);
}

// ...

// declare before parallel region as global
int* awork = NULL;

#pragma omp declare reduction(merge : int* : meta_op(&omp_out, omp_in, 100000)) initializer (omp_priv=NULL)

#pragma omp for reduction(merge : awork)
        for (int t = 0; t < n; t++) {
            meta_op(&awork, at[t], bins);
        }

令人惊讶的是,这看起来并不太好:

_{顶部是icc 16.0.2,底部是gcc 5.3.0,两者都有-O3.}

两者似乎都实现了序列化的缩减.我试图调查gcc/ libgomp,但我不会立即明白发生了什么.从中间代码/反汇编,它们似乎将最终合并包装在GOMP_atomic_start/ end- 中,这似乎是一个全局互斥.同样icc将调用包装到operationa中kmpc_critical.我认为进行昂贵的自定义减少操作没有太多优化.传统的减少可以通过硬件支持的原子操作来完成.

注意每个operation都是更快的,因为输入是在本地缓存的,但由于序列化,它总体上较慢.由于差异很大,这不是一个完美的比较,早期的截图是不同的gcc版本.但趋势很明显,我也有关于缓存效果的数据.