c++ - 使用openmp进行LU分解

我有一个问题:LU分解算法的并行版本与序列同时运行:

void lup_od_omp(double* a, int n){

int i,j,k;

for(k = 0; k < n - 1; ++k)
{
    #pragma omp parallel for shared(a,n,k) private(i,j)
    for(i = k + 1; i < n; i++)
    {
        a[i*n + k] /= a[k*n + k];
        for(j = k + 1; j < n; j++)
        {
            a[i*n + j] -= a[i*n + k]*a[k*n + j];
        }
    }
}}

也许我做错了什么？

最佳答案

由于您仅在两个内核上工作，因此并行化实际上可能会妨碍矢量化器。 SSE2上的矢量化将为您提供每个操作2倍的数据带宽，在AVX上为4倍。

双线程有很多同步开销，这可能会减慢您的速度，尤其是在向量化宽松的情况下。同样由于某种原因，除非调用#pragma omp使其实际使用线程，否则omp_set_num_threads不会启动任何线程。

与向量化也有关的另一件事是，并非所有编译器都知道a[i*n + j]旨在处理二维数组，因此最好一开始就声明它。

这是对代码的略微优化，可以在Xeon上很好地运行:

void lup_od_omp(int n, double (*a)[n]){
    int i,k;

    for(k = 0; k < n - 1; ++k) {
        // for the vectoriser
        for(i = k + 1; i < n; i++) {
            a[i][k] /= a[k][k];
        }

        #pragma omp parallel for shared(a,n,k) private(i) schedule(static, 64)
        for(i = k + 1; i < n; i++) {
            int j;
            const double aik = a[i][k]; // some compilers will do this automatically
            for(j = k + 1; j < n; j++) {
                a[i][j] -= aik * a[k][j];
            }
        }
    }
}

3000x3000 icc -O2数组的运行时:

Your code sequential:  0:24.61 99%  CPU
Your code 8 threads :  0:05.21 753% CPU
My   code sequential:  0:18.53 99%  CPU
My   code 8 threads :  0:05.42 766% CPU

在另一台机器上，我在AVX上进行了测试(256位 vector ，每个运算符4倍):

My code on AVX sequential :  0:09.45 99%  CPU
My code on AVX 8 threads  :  0:03.92 766% CPU

如您所见，我在向量化器方面做了一些改进，但对于并行部分并没有做太多事情。