10

我是 OpenMp 编程的新手。我编写了一个简单的 c 程序来将矩阵与向量相乘。不幸的是，通过比较执行时间，我发现 OpenMP 比 Sequential 方式慢得多。

这是我的代码（这里的矩阵是 N*N int，向量是 N int，结果是 N long long）：

#pragma omp parallel for private(i,j) shared(matrix,vector,result,m_size)
for(i=0;i<m_size;i++)
{  
  for(j=0;j<m_size;j++)
  {  
    result[i]+=matrix[i][j]*vector[j];
  }
}

这是顺序方式的代码：

for (i=0;i<m_size;i++)
        for(j=0;j<m_size;j++)
            result[i] += matrix[i][j] * vector[j];

当我使用 999x999 矩阵和 999 向量尝试这两种实现时，执行时间为：

顺序：5439 毫秒并行：11120 毫秒

我真的不明白为什么 OpenMP 比顺序算法慢得多（慢 2 倍以上！）谁能解决我的问题？

4

3 回答 3

19

您的代码部分受到所谓的错误共享的影响，这对于所有缓存一致的系统都是典型的。简而言之，result[]数组的许多元素都适合同一个高速缓存行。当线程作为操作符的结果i写入时，保存该部分的缓存行变脏。然后，缓存一致性协议使其他内核中该缓存行的所有副本无效，并且它们必须从上层缓存或主存储器刷新其副本。作为的数组，一个高速缓存行（x86 上为 64 字节）包含 8 个元素，此外result[i]+=result[]resultlong longresult[i]同一缓存行中还有 7 个其他数组元素。因此，两个“相邻”线程可能会不断争夺缓存行的所有权（假设每个线程在单独的核心上运行）。

为了减轻您的情况下的错误共享，最简单的方法是确保每个线程都获得一个迭代块，其大小可被缓存行中的元素数整除。例如，您可以应用schedule(static,something*8)wheresomething应该足够大，这样迭代空间就不会被分割成太多块，但同时它应该足够小，以便每个线程都获得一个块。例如，对于m_size等于 999 和 4 个线程，您可以将该schedule(static,256)子句应用于parallel for构造。

代码运行速度变慢的另一个部分原因可能是，当启用 OpenMP 时，编译器可能会在分配共享变量时不愿意应用一些代码优化。OpenMP 提供了所谓的宽松内存模型，其中允许每个线程中共享变量的本地内存视图不同，并且flush提供构造以同步视图。volatile但是，如果编译器不能证明其他线程不需要访问去同步的共享变量，编译器通常会将共享变量视为隐含的。你的情况是其中之一，因为result[i]只分配给和价值result[i]从未被其他线程使用。在串行情况下，编译器很可能会创建一个临时变量来保存内部循环的结果，并且只会result[i]在内部循环完成后分配给。在并行情况下，它可能会决定这会result[i]在其他线程中创建一个临时的不同步视图，因此决定不应用优化。仅作记录，GCC 4.7.1 在-O3 -ftree-vectorize启用和不启用 OpenMP 的情况下执行临时变量技巧。

于 2013-05-05T02:45:57.367 回答
3

因为当 OpenMP 在线程之间分配工作时，会进行大量管理/同步，以确保共享矩阵和向量中的值不会以某种方式损坏。即使它们是只读的：人类很容易看到，但您的编译器可能看不到。

出于教学原因要尝试的事情：

0）如果matrix和vector不是会发生什么shared？

1）首先并行化内部“j-loop”，保持外部“i-loop”串行。走着瞧吧。

2）不要在 中收集总和result[i]，而是在变量中收集和，并且仅在内循环完成后将temp其内容分配result[i]给以避免重复的索引查找。不要忘记temp在内部循环开始之前初始化为 0。

于 2013-05-04T07:13:52.243 回答
0

我这样做是参考 Hristo 的评论。我尝试使用时间表（静态，256）。对我来说，这无助于更改默认的块大小。也许它甚至使情况变得更糟。我打印出线程号及其索引，无论是否设置时间表，很明显 OpenMP 已经选择了远离彼此的线程索引，因此错误共享似乎不是问题。对我来说，这段代码已经对 OpenMP 起到了很好的推动作用。

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

void loop_parallel(const int *matrix, const int ld, const int*vector, long long* result, const int m_size) {
    #pragma omp parallel for schedule(static, 250)
    //#pragma omp parallel for
    for (int i=0;i<m_size;i++) {
        //printf("%d %d\n", omp_get_thread_num(), i);
        long long sum = 0;
        for(int j=0;j<m_size;j++) {
            sum += matrix[i*ld +j] * vector[j];
        }
        result[i] = sum;
    }
}

void loop(const int *matrix, const int ld, const int*vector, long long* result, const int m_size) {
    for (int i=0;i<m_size;i++) {
        long long sum = 0;
        for(int j=0;j<m_size;j++) {
            sum += matrix[i*ld +j] * vector[j];
        }
        result[i] = sum;
    }
}

int main() {
    const int m_size = 1000;
    int *matrix = new int[m_size*m_size];
    int *vector = new int[m_size];
    long long*result = new long long[m_size];
    double dtime;

    dtime = omp_get_wtime();
    loop(matrix, m_size, vector, result, m_size);
    dtime = omp_get_wtime() - dtime;
    printf("time %f\n", dtime);

    dtime = omp_get_wtime();
    loop_parallel(matrix, m_size, vector, result, m_size);
    dtime = omp_get_wtime() - dtime;
    printf("time %f\n", dtime);

}
于 2013-05-05T18:58:21.680 回答