通过多线程为稀疏矩阵赋值时的段错误

Segfault when assign values to sparse matrix by multi-threads

我打算用从一系列步骤派生的值填充一个稀疏矩阵,以提高效率,OpenMP 用于加速这些进程,我发现它在使用 1 个线程时工作正常,但遇到了段错误对于多线程,我准备了一个简单的demo代码来重现错误,真心希望有人能帮我一个忙

#include <RcppArmadillo.h>
#include <omp.h>

// [[Rcpp::plugins(cpp11)]]
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::depends(bigmemory, BH)]]

using namespace std;
using namespace Rcpp;
using namespace arma;

// [[Rcpp::export]]
arma::sp_mat test(arma::vec x, int n, int threads = 1){
    omp_set_num_threads(threads);
    arma::sp_mat m(n, n);
    #pragma omp parallel for schedule(dynamic) 
    for(int i = 0; i < n; i++){
        for(int j = 0; j < n; j++){
            m(i, j) = x[i * n + j];
        }
    }
    return m;
}

# run 
a<-test(sample(c(0,1,2),100*100,rep=T), n=100, threads=1)
a<-test(sample(c(0,1,2),100*100,rep=T), n=100, threads=10)

tl;dr: 你不能。

长话短说:(Rcpp)Armadillo 的主要优势之一是在底层操作之上非常好且一致 API,这些操作介于令人困惑和难以使用之间。缺点之一是我们很容易忽视底层数据结构。

Dense 矩阵(基本上总是)固定的内存块。本质上,一个大小为行 x 列的向量。这就是使我们能够在 R 和 (Rcpp)Armadillo 之间进行高效 'zero copy' 传输的原因。它还允许我们在非重叠块上同时工作。这很重要,例如 RcppParallel 充分利用了它。 OpenMP 在这里工作。

稀疏 矩阵是(我在这里简化)具有相互依赖性的动态 list/vector 类型。所以并发工作根本行不通。伤心。但事实就是如此。一旦您更仔细地查看稀疏矩阵的通用数据结构(如 e.g. R 的 Matrix 包所做的那样),它就会变得清晰。例如 this wikipedia piece 是一个非常体面和详尽的介绍。