c++ 基于OpenMP的COO矩阵向量积并行化

bvjveswy 于 2023-08-09 发布在其他

关注(0)|答案(1)|浏览(125)

我正在尝试使用C++中的OpenMP来提高稀疏矩阵向量积的性能。我以COO格式存储稀疏矩阵，也就是说，我有一个3个数组的结构体，对应于稀疏矩阵的每个非零项。对于结构体的每个索引，我可以找到非零项的行索引、列索引和值。此外，我可以通过以下方式指定函数要使用的线程数

export OMP_NUM_THREADS=n

字符串
其中n是我想使用的线程数。目前，我的代码如下

void ompMatvec(const Vector& x, Vector& y) const {

    std::string envName = "OMP_NUM_THREADS";
    std::string thread_count = getenv(envName); 
    int thread_count2 = atoi(thread_count.c_str()); 
    Vector y2(y.numRows()); 
    size_type k; 
    #pragma omp parallel num_threads(thread_count2) 
    #pragma omp for 
      for (k = 0; k < arrayData.size(); ++k) {
        y(rowIndices[k]) += arrayData[k] * x(colIndices[k]);
      }
  }

型
然而，当我测量性能时，我发现我的速度不是太高。我将上面的并行化函数与：

void matvec(const Vector& x, Vector& y) const {
    for (size_type k = 0; k < arrayData.size(); ++k) {
      y(rowIndices[k]) += arrayData[k] * x(colIndices[k]);
    }
  }

型
我想提一下，我已经创建了一个Vector类，其中包含私有成员函数.numRows（），它本质上提供了向量的长度。我也在4核上运行代码。是否存在可提高OpenMP API性能的实施变更？或者它受到我的程序运行的内核数量的限制？
任何和所有的建议都非常感谢。谢谢你，谢谢
更新：尝试避免上述竞态条件：

void ompMatvec(const Vector& x, Vector& y) const {

    std::string envName = "OMP_NUM_THREADS";
    std::string thread_count = getenv(envName); 
    int thread_count2 = atoi(thread_count.c_str());  
    size_type k; 
    #pragma omp parallel num_threads(thread_count2) \
    default(none) private(k) 
    Vector y2(y.numRows());
    #pragma omp for 
      for (k = 0; k < arrayData.size(); ++k) {
        y2(rowIndices[k]) += arrayData[k] * x(colIndices[k]);
      }
    #pragma omp critical
      for(k = 0; k < y.numRows(); ++k){
        y(k) += y2(k); 
      }
  }

型

c++

来源：https://stackoverflow.com/questions/43974180/coo-matrix-vector-product-parallelization-with-openmp

1条答案

按热度按时间

esyap4oy1#

正如@Zulan在他的评论中所解释的那样，ompMatvec的第一个实现将产生竞争条件并简单地计算垃圾。
据我所知，有两种方法可以用OpenMP并行化COO matvec。
第一种是使用原子操作来强制同步。然后，COO matvec的可能实现由下式给出：

void dcoomv_atomic(SparseMatrixCOO *A, double *x, double *y) {
  for (int i=0; i<A->m; i++) y[i] = 0.;
  #pragma omp parallel for
  for (int i=0; i<A->nnz; i++) {
    #pragma omp atomic
    y[A->row_idx[i]] += A->val[i] * x[A->col_idx[i]];
  }
}

字符串
与第一个实现不同的是，这将起作用。然而，同步的成本非常大。在我的例子中，当我用16个线程运行这段代码时，我没有得到任何加速，而是比串行实现慢了7倍。因此，这种方法在实践中是无用的。
第二，就像@祖蓝在评论中提到的，可以做减分。OpenMP 4.5+支持C/C++的数组缩减。具有数组缩减的COO matvec的实现给出如下：

void dcoomv_builtin_array_reduction(SparseMatrixCOO *A, double *x, double *y) {
  for (int i=0; i<A->m; i++) y[i] = 0.;
  #pragma omp parallel for reduction(+:y[:A->m])
  for (int i=0; i<A->nnz; i++) {
    y[A->row_idx[i]] += A->val[i] * x[A->col_idx[i]];
  }
}

型
这段代码编译得很好。但是，在执行过程中出现了一个分段错误。很明显，我做错了什么。如果有人知道我应该如何实际使用OpenMP的内置数组缩减功能，请告诉我。请注意，我使用英特尔C编译器icc 2021.6.0使用-fopenmp标志进行编译，该编译器支持OpenMP 4.5和OpenMP 5.0的子集。
由于我很难使OpenMP的内置数组减少工作，因此我可以从头开始实现数组减少，如下所示：

void dcoomv_array_reduction_from_scratch(SparseMatrixCOO *A, double *x, double *y) {
  for (int i=0; i<A->m; i++) y[i] = 0.;
  double *YP;
  #pragma omp parallel 
  {
    int P = omp_get_num_threads();
    int p = omp_get_thread_num();
    #pragma omp single
    {
      YP = (double*)mkl_malloc(A->m * P * sizeof(double), sizeof(double));
      for (int i=0; i<A->m*P; i++) YP[i] = 0.;
    }
    #pragma omp for
    for (int i=0; i<A->nnz; i++) {
      YP[p * A->m + A->row_idx[i]] += A->val[i] * x[A->col_idx[i]];
    }
    #pragma omp for
    for (int i=0; i<A->m; i++) {
      for (int p=0; p<P; p++) {
        y[i] +=  YP[A->m * p + i];
      }
    }
  }
  mkl_free(YP);
}

型
这个代码可以工作。然而，它也是低效的，因为当使用16个线程时，与串行实现相比，我得到的不是加速，而是大约3倍的减速。
和你一样，我也觉得这个问题很有趣，因为根据我自己的基准测试经验，MKL函数mkl_sparse_d_mv似乎对于CSR和BSR存储格式是多线程的，但对于COO和CSC则不是。
然而，在玩了足够多的游戏之后，我明白了对COO（和CSC）matvec进行多线程处理并不简单，或者至少，它不容易转换为加速，即使以额外的内存分配为代价。也许这就是为什么MKL开发人员决定不为COO和CSC稀疏矩阵格式并行化mkl_sparse_d_mv的原因。

赞(0）回复(0）举报 2023-08-09

我来回答

c++ 基于OpenMP的COO矩阵向量积并行化

1条答案

相关问题

热门标签

最新问答