嵌套循环和openmp的麻烦

Question

我无法将openmp应用于这样的嵌套循环:

        #pragma omp parallel shared(S2,nthreads,chunk) private(a,b,tid)
    {
        tid = omp_get_thread_num();
        if (tid == 0)
        {
            nthreads = omp_get_num_threads();
            printf("\nNumber of threads = %d\n", nthreads);
        }
        #pragma omp for schedule(dynamic,chunk)
        for(a=0;a<NREC;a++){
            for(b=0;b<NLIG;b++){
                S2=S2+cos(1+sin(atan(sin(sqrt(a*2+b*5)+cos(a)+sqrt(b)))));
            }
        } // end for a
    } /* end of parallel section */

当我将串口与openmp版本进行比较时,最后一个给出了奇怪的结果.即使我删除了#pragma omp,openmp的结果也不正确,你知道为什么或者可以指出一个关于双循环和openmp的明确的教程吗？

Answer 1

这是竞争条件的典型例子.您的每个openmp线程都在同时访问和更新共享值,并且没有任何保证,某些更新不会丢失(充其量)或者最终的答案不会是胡言乱语(最坏的情况下).

有竞争条件的是他们敏感地依赖于时间; 在较小的情况下(例如,使用较小的NREC和NLIG),您有时可能会错过这个,但在更大的情况下,它最终会出现.

你得到错误答案的原因#pragma omp for是,只要你进入并行区域,你的所有openmp线程都会启动; 除非你使用像omp for(一个所谓的工作共享构造)这样的东西来分割工作,每个线程将在并行部分中执行所有操作 - 因此所有线程将执行相同的整个总和,所有线程都同时更新S2.

您必须小心OpenMP线程更新共享变量.OpenMP具有atomic允许您安全地修改共享变量的操作.下面是一个例子(不幸的是,你的例子对求和顺序非常敏感,很难看出发生了什么,所以我有点改变你的总和:).在mysumallatomic,每个线程都S2像以前一样更新,但这次安全地完成了:

#include <omp.h>
#include <math.h>
#include <stdio.h>

double mysumorig() {

    double S2 = 0;
    int a, b;
    for(a=0;a<128;a++){
        for(b=0;b<128;b++){
            S2=S2+a*b;
        }
    }

    return S2;
}


double mysumallatomic() {

    double S2 = 0.;
#pragma omp parallel for shared(S2)
    for(int a=0; a<128; a++){
        for(int b=0; b<128;b++){
            double myterm = (double)a*b;
            #pragma omp atomic
            S2 += myterm;
        }
    }

    return S2;
}


double mysumonceatomic() {

    double S2 = 0.;
#pragma omp parallel shared(S2)
    {
        double mysum = 0.;
        #pragma omp for
        for(int a=0; a<128; a++){
            for(int b=0; b<128;b++){
                mysum += (double)a*b;
            }
        }
        #pragma omp atomic
        S2 += mysum;
    }
    return S2;
}

int main() {
    printf("(Serial)      S2 = %f\n", mysumorig());
    printf("(All Atomic)  S2 = %f\n", mysumallatomic());
    printf("(Atomic Once) S2 = %f\n", mysumonceatomic());
    return 0;
}

但是,原子操作确实会损害并行性能(毕竟,重点是防止变量的并行操作S2!)所以更好的方法是进行求和,只做两次求和后的原子操作而不是做128*128次; 这是mysumonceatomic()例程,每个线程只会产生一次同步开销而不是每个线程16k次.

但这是一种常见的操作,不需要自己实施.可以使用OpenMP内置功能进行缩减操作(缩减是指计算列表总和,查找列表的最小值或最大值等操作,只能通过查看一次一个元素来完成到目前为止的结果和@ejd建议的下一个元素.OpenMP将工作并且更快(它的优化实现比您自己使用其他OpenMP操作所做的快得多).

如您所见,两种方法都有效:

$ ./foo
(Serial)      S2 = 66064384.000000
(All Atomic)  S2 = 66064384.000000
(Atomic Once) S2 = 66064384.00000

Answer 2

问题不在于双循环,而在于变量S2.尝试在for指令上添加一个减少条款:

#pragma omp用于调度(动态,块)减少(+:S2)