OpenMP 并行化 GEMM 的错误结果

Wrong result with OpenMP to parallelize GEMM

我知道 OpenMP 共享所有工作程序之间在外部范围内声明的所有变量。这就是我的问题的答案。但是我真的很困惑为什么函数 omp3 提供正确的结果而函数 omp2 提供错误的结果。

void omp2(double *A, double *B, double *C, int m, int k, int n) {
    for (int i = 0; i < m; ++i) {
#pragma omp parallel for
        for (int ki = 0; ki < k; ++ki) {
            for (int j = 0; j < n; ++j) {
                C[i * n + j] += A[i * k + ki] * B[ki * n + j];
            }
        }
    }
}
void omp3(double *A, double *B, double *C, int m, int k, int n) {
    for (int i = 0; i < m; ++i) {
        for (int ki = 0; ki < k; ++ki) {
#pragma omp parallel for
            for (int j = 0; j < n; ++j) {
                C[i * n + j] += A[i * k + ki] * B[ki * n + j];
            }
        }
    }
}

问题是这一行存在竞争条件:

C[i * n + j] += ...

不同的线程可以同时读取和写入相同的内存位置(C[i * n + j]),这会导致数据竞争。在 omp2 中可能会发生此数据竞争,但在 omp3 中不会发生。

解决方案(如@Victor Eijkhout 所建议)是重新排序循环,使用局部变量来计算最内层循环的总和。在这种情况下 C[i * n + j] 只更新一次,所以你摆脱了数据竞争并且最外层的循环可以并行化(这提供了最好的性能):

   #pragma omp parallel for
    for (int i = 0; i < m; ++i) {
        for (int j = 0; j < n; ++j) {
            double sum=0;            
            for (int ki = 0; ki < k; ++ki) {            
                sum += A[i * k + ki] * B[ki * n + j];
            }
            C[i * n + j] +=sum;
        }        
    }

请注意,您可以使用 collapse(2) 子句,这可能会提高性能。