稀疏矩阵的多线程程序

Multithreaded Program for Sparse Matrices

本文关键字：多线程 程序      更新时间：2023-10-16

我是多线程的新手。我正在尝试设计一个解决稀疏矩阵的程序。在我的代码中，我将向量向量点积和 Matix 向量积多次称为子例程，以达到最终解决方案。我正在尝试使用开放 MP（尤其是上述两个子例程）并行化代码。我也有顺序代码，我不打算在它们之间并行化。

我的问题是如何处理调用子例程时创建的线程。我应该在每个子例行调用结束时设置障碍吗？

另外，我应该在哪里设置线程数？ 

Mat_Vec_Mult(MAT,x0,rm);
#pragma omp parallel for schedule(static)
for(int i=0;i<numcols;i++)  
    rm[i] = b[i] - rm[i];
#pragma omp barrier

#pragma omp parallel for schedule(static)
for(int i=0;i<numcols;i++)  
    xm[i] = x0[i];
#pragma omp barrier

double* pm = (double*) malloc(numcols*sizeof(double));

    #pragma omp parallel for schedule(static)
for(int i=0;i<numcols;i++)  
    pm[i] = rm[i];
#pragma omp barrier
scalarProd(rm,rm,numcols);

谢谢

编辑：

对于标量点积，我使用以下代码段：

double scalarProd(double* vec1, double* vec2, int n){
double prod = 0.0;
int chunk = 10; 
int i;
//double* c = (double*) malloc(n*sizeof(double));
omp_set_num_threads(4);
// #pragma omp parallel shared(vec1,vec2,c,prod) private(i)
#pragma omp parallel
{
    double pprod = 0.0;
    #pragma omp for
    for(i=0;i<n;i++) {
        pprod += vec1[i]*vec2[i];
    }
    //#pragma omp for reduction (+:prod)
    #pragma omp critical
    for(i=0;i<n;i++) {
        prod += pprod;
    }
}

return prod;
}

我现在在我的共轭梯度函数中添加了时间计算代码，如下所示：

start_dotprod = omp_get_wtime();
rm_rm_old = scalarProd(rm,rm,MAT->ncols);
    run_dotprod = omp_get_wtime() - start_dotprod;
fprintf(timing,"Time taken by rm_rm dot product : %lf n",run_dotprod);

观察到的结果 ： 点积所花费的时间 顺序版本 ： 0.000007s 并行版本 ： 0.002110

我正在我的英特尔 I7 笔记本电脑上的 Linux 操作系统上使用 gcc -fopenmp 命令进行简单的编译。

我目前正在使用大小为 n = 5000 的矩阵。

由于同一个点积被多次调用，直到实现收敛（大约 80k 次），我的总体速度大大降低了。

请提出一些改进建议。任何帮助都非常感谢！

老实说，我建议在更高层次上并行化。我的意思是尽量减少您使用的#pragma omp parallel数量。每次尝试在线程之间拆分工作时，都会产生 OpenMP 开销。尽可能避免这种情况。

所以在你的情况下，至少我会尝试：

Mat_Vec_Mult(MAT,x0,rm);
double* pm = (double*) malloc(numcols*sizeof(double)); // must be performed once outside of parallel region
// all threads forked and created once here
#pragma omp parallel for schedule(static)
for(int i = 0; i < numcols; i++) {
    rm[i] = b[i] - rm[i]; // (1)
    xm[i] = x0[i];        // (2) does not require (1)
    pm[i] = rm[i];        // (3) requires (1) at this i, not (2)
}  
// implicit barrier at the end of omp for
// implicit join of all threads at the end of omp parallel
scalarProd(rm,rm,numcols);

请注意，我是如何证明无论如何，您的循环之间实际上不需要任何障碍。

如果你的大部分时间都花在这个计算阶段，你肯定会看到相当大的改进。但是，我有理由相信，您的大部分时间都花在Mat_Vec_Mult()上，也许也花在scalarProd()上，因此您将节省的时间可能很少。

**

编辑 **

根据您的编辑，我看到了一些问题。（1） 在测试算法性能时，始终使用 -O3 进行编译。（2）您将无法改善需要.000007秒才能完成的操作的运行时间;这几乎是瞬间的。这又回到了我之前所说的：尝试在更高的级别上并行化。CG方法本质上是一种顺序算法，但肯定有研究论文详细描述了并行CG。（3） 标量产品的实现不是最佳的。事实上，我怀疑你对矩阵向量积的实现也不是。我个人会做以下几点：

double scalarProd(double* vec1, double* vec2, int n) {
    double prod = 0.0;
    int i;
    // omp_set_num_threads(4); this should be done once during initialization somewhere previously in your program
    #pragma omp parallel for private(i) reduction(+:prod)
    for (i = 0; i < n; ++i) {
        prod += vec1[i]*vec2[i];
    }
    return prod;
}
（

4）有整个库（LAPACK，BLAS等）具有高度优化的矩阵向量，向量向量等操作。任何线性代数库都必须建立在它们之上。因此，我建议在开始在此处重新创建轮子并尝试实现自己的轮子之前，先考虑使用这些库之一来执行两个操作。

相关文章：