openmp中的并行for循环答案

【问题标题】：Parallel for loop in openmpopenmp中的并行for循环
【发布时间】：2012-07-31 04:01:34
【问题描述】：

我正在尝试并行化一个非常简单的 for 循环，但这是我很长时间以来第一次尝试使用 openMP。我对运行时间感到困惑。这是我的代码：

#include <vector>
#include <algorithm>

using namespace std;

int main () 
{
    int n=400000,  m=1000;  
    double x=0,y=0;
    double s=0;
    vector< double > shifts(n,0);


    #pragma omp parallel for 
    for (int j=0; j<n; j++) {

        double r=0.0;
        for (int i=0; i < m; i++){

            double rand_g1 = cos(i/double(m));
            double rand_g2 = sin(i/double(m));     

            x += rand_g1;
            y += rand_g2;
            r += sqrt(rand_g1*rand_g1 + rand_g2*rand_g2);
        }
        shifts[j] = r / m;
    }

    cout << *std::max_element( shifts.begin(), shifts.end() ) << endl;
}

我用

编译

g++ -O3 testMP.cc -o testMP  -I /opt/boost_1_48_0/include

也就是说，没有“-fopenmp”，我得到了这些时间：

real    0m18.417s
user    0m18.357s
sys     0m0.004s

当我使用“-fopenmp”时，

g++ -O3 -fopenmp testMP.cc -o testMP  -I /opt/boost_1_48_0/include

我得到了这些数字：

real    0m6.853s
user    0m52.007s
sys     0m0.008s

这对我来说没有意义。如何使用八核只能导致 3 倍性能提升？我的循环编码是否正确？

【问题讨论】：

标签： c++ multithreading performance parallel-processing openmp

【解决方案1】：

您应该为x 和y 使用OpenMP reduction 子句：

#pragma omp parallel for reduction(+:x,y)
for (int j=0; j<n; j++) {

    double r=0.0;
    for (int i=0; i < m; i++){

        double rand_g1 = cos(i/double(m));
        double rand_g2 = sin(i/double(m));     

        x += rand_g1;
        y += rand_g2;
        r += sqrt(rand_g1*rand_g1 + rand_g2*rand_g2);
    }
    shifts[j] = r / m;
}

使用reduction，每个线程在x 和y 中累积自己的部分和，最后将所有部分值相加以获得最终值。

Serial version:
25.05s user 0.01s system 99% cpu 25.059 total
OpenMP version w/ OMP_NUM_THREADS=16:
24.76s user 0.02s system 1590% cpu 1.559 total

见 - 超线性加速:)

【讨论】：

【解决方案2】：

让我们尝试了解如何使用 OpenMP 并行化简单的 for 循环

#pragma omp parallel
#pragma omp for
    for(i = 1; i < 13; i++)
    {
       c[i] = a[i] + b[i];
    }

假设我们有3 可用线程，这就是将会发生的事情

首先

线程被分配一组独立的迭代

最后

线程必须在工作共享结构结束时等待

【讨论】：

【解决方案3】：

因为这个问题的关注度很高，所以我决定添加一点 OpenMP 背景来帮助那些访问它的人

#pragma omp parallel 与threads 团队创建一个并行区域，其中每个线程执行parallel region 所包含的整个代码块。从OpenMP 5.1 可以阅读更正式的描述：

当一个线程遇到一个并行构造时，一组线程被创建 来执行并行区域 (..)。这遇到并行构造的线程成为主线程新团队的线程，持续时间的线程数为零新的平行区域。 新团队中的所有主题，包括主线程，执行区域。 创建团队后，团队中的线程数在持续时间内保持不变那个平行区域。

#pragma omp parallel for 创建了一个parallel region（如前所述），并且将使用default chunk size 和default schedule 为该区域的threads 分配它所包含的循环的迭代这是通常 static。但是请记住，default schedule 在 OpenMP 标准的不同具体实现之间可能会有所不同。

您可以从OpenMP 5.1 阅读更正式的描述：

worksharing-loop 结构指定一个或一个的迭代更多相关的循环将由线程中的线程并行执行团队在他们的隐含任务的背景下。 迭代是分布在团队中已经存在的线程中执行工作共享循环区域所在的并行区域绑定。

Moreover,

并行循环结构是指定并行循环的快捷方式包含带有一个或多个相关联的循环构造的构造循环，没有其他语句。

或者非正式地，#pragma omp parallel for 是构造函数#pragma omp parallel 和#pragma omp for 的组合。在您的情况下，这意味着：

#pragma omp parallel for 
for (int j=0; j<n; j++) {

    double r=0.0;
    for (int i=0; i < m; i++){

        double rand_g1 = cos(i/double(m));
        double rand_g2 = sin(i/double(m));     

        x += rand_g1;
        y += rand_g2;
        r += sqrt(rand_g1*rand_g1 + rand_g2*rand_g2);
    }
    shifts[j] = r / m;
}

将创建一组线程，并为每个线程分配最外层循环的迭代块。

为了更清楚地说明，4 线程与 chunk_size=1 和 static schedule 的 #pragma omp parallel for 将导致类似：

在代码方面，循环将被转换为逻辑上类似于：

for(int i=omp_get_thread_num(); i < n; i+=omp_get_num_threads())
{  
    c[i]=a[i]+b[i];
}

在哪里omp_get_thread_num()

omp_get_thread_num 例程返回线程号，在调用线程的当前团队。

和omp_get_num_threads()

返回当前团队中的线程数。在一个顺序程序 omp_get_num_threads 的部分返回 1。

或者换句话说，for(int i = THREAD_ID; i < n; i += TOTAL_THREADS)。其中THREAD_ID 范围从0 到TOTAL_THREADS - 1，TOTAL_THREADS 表示在并行区域上创建的团队的线程总数。

掌握了这些知识，并查看您的代码，可以看到您对变量“x”和“y”的更新有一个竞态条件。这些变量在线程之间共享并在并行区域内更新，即：

     x += rand_g1;
     y += rand_g2;

要解决这种竞争条件，您可以使用 OpenMP' reduction 子句：

指定每个线程私有的一个或多个变量是并行结束时归约操作的主题地区。

简而言之，reduction 子句将为每个线程创建变量“x”和“y”的私有副本，并在并行区域结束时对所有“x”和“y”变量进行求和进入初始线程中的原始“x”和“y”变量。

#pragma omp parallel for reduction(+:x,y)
for (int j=0; j<n; j++) {

    double r=0.0;
    for (int i=0; i < m; i++){

        double rand_g1 = cos(i/double(m));
        double rand_g2 = sin(i/double(m));     

        x += rand_g1;
        y += rand_g2;
        r += sqrt(rand_g1*rand_g1 + rand_g2*rand_g2);
    }
    shifts[j] = r / m;
}

【讨论】：

【解决方案4】：

您最多可以实现（！）是线性加速。现在我不记得哪个是来自 linux 的时间，但我建议您使用 time.h 或（在 c++ 11 中）“chrono”并直接从程序测量运行时。最好将整个代码打包成一个循环，运行 10 次，然后平均得到 prog 的大约运行时间。

此外，您还遇到了 x,y 的问题 - 它不符合并行编程中数据局部性的范式。

【讨论】：

“你最多能实现的（！）是线性加速。” - 错误的！使用适当的本地数据访问模式，对于像这样的令人尴尬的并行问题，通常会观察到超线性加速，因为更多的数据适合组合的 CPU 缓存。