使用RcppArmadillo与OpenMP并行时R会话崩溃问题求助

崩溃核心原因排查 #

• 第一个最直接的问题：并行区域变量作用域错误，出现多线程数据竞争
  你定义的temp1、temp2、temp2a、t都在for循环外部，OpenMP默认会将这些变量判定为所有线程共享变量，多个线程同时修改同一个临时矩阵/变量，会直接导致内存越界、非法访问，触发R会话崩溃。
• 第二个可优化的问题：循环内存在大量重复冗余计算，既拖慢效率也增加并行风险
  你在循环内每次都计算temp2 = dat.rows(find(dat.col(3) == j && dat.col(2) == h))，这个判断条件和循环变量i完全无关，完全可以提前到循环外计算一次即可，不需要每次循环反复执行。
• 额外注意：R的所有API都是单线程不安全的，并行区域内绝对不能调用任何R API、创建Rcpp对象（你当前代码的返回List是在并行区外创建的，暂时没有这个问题，但后续修改要注意）

修复后的代码 #

// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::plugins(openmp)]]

#include <RcppArmadillo.h>
#include <iostream>
#include<algorithm>

using namespace arma;

// [[Rcpp::export]]
Rcpp::List A_hj_Comp(const arma::mat& X, int h, int j, double v, double u, int ncores=1){
    arma::mat dat = X.rows(find(X.col(2)==h));
    arma::mat ts = dat.rows( find(dat.col(3) == j && (dat.col(3) != dat.col(2))) );
    arma::vec ts1 = ts.col(5);
    arma::vec ts_1 = sort(unique(ts1));
    arma::vec ts2 = ts_1.elem( find(ts_1 > v && ts_1 <= u) );
    
    int nn=ts2.size();
    arma::vec res1(ts2.size()); 
    arma::vec res2(ts2.size());
    arma::vec A = arma::zeros<arma::vec>(nn+1); 
    // 提前计算循环内不变的temp2，避免重复计算
    arma::mat temp2 = dat.rows(find(dat.col(3) == j && dat.col(2) == h));
    
    #if defined(_OPENMP)
        #pragma omp parallel for num_threads(ncores)
    #endif
    for(arma::uword i = 0; i < ts2.size(); i++) {
        // 所有临时变量定义在循环内部，自动成为线程私有变量，避免数据竞争
        double t = ts2(i);
        arma::mat temp1 = dat.rows(find(dat.col(4) < t && dat.col(5) >= t));
        res1(i) = sum(1/temp1.col(6));
        arma::mat temp2a = temp2.rows(find(temp2.col(5) == t));
        res2(i) = sum(1/temp2a.col(6));
    }
    A.subvec(1,nn) = cumsum(res2/res1); 
    
    Rcpp::List res; 
    res["ts"] = ts2;
    res["A"] = A;
    return(res);
}

Windows平台额外适配说明 #

你使用的是Windows平台的mingw版本R，编译时保留原来的环境变量设置即可：

Sys.setenv("PKG_CXXFLAGS" = "-fopenmp")
Sys.setenv("PKG_LIBS" = "-fopenmp")
sourceCpp("case1.cpp")

如果依然出现偶发崩溃，可以尝试将OpenMP的线程数限制在8核以内，部分旧版本mingw的OpenMP实现对高核心数的支持存在已知bug。

内容的提问来源于stack exchange，提问作者Sam