RcppArmadillo на нескольких ядрах процессора

Question

У меня есть следующая функция RccpArmadillo, которая отлично работает, если я выполняю ее на одном ядре процессора. Но если я использую несколько ядер, то R сработает. Все остальные функции Rcpp, которые я создал до сих пор, отлично работают на нескольких ядрах (с foreach), только RccpArmadillo кажется проблематичным. Любые идеи, как это исправить?

cppFunction('double augmentedDickeyFullerCpp(NumericVector a, NumericVector b, double gamma, double mu, int lags) { 

    if (gamma < 0) { 
     return 0; 
    } 

    int n = a.size()-1; 
    int lags2 = lags + 1; 
    // first rows, then columns 
    NumericMatrix x(n-lags2,lags2); 
    NumericMatrix zdifflag(n-lags2+1,lags2); 
    NumericVector diff(n); 
    NumericVector zdiff(n-lags2+1); 
    NumericVector residuals(n+1); 
    residuals[0] = a[0] - gamma * b[0] - mu; 


    // residuals a is y and b is x 
    for(int i = 1; i < n+1; i++) { 
     residuals[i] = a[i] - gamma * b[i] - mu; 
     diff[i-1] = residuals[i] - residuals[i-1]; 
    } 
    for(int i = 0; i < n-lags2+1; i++) { 
     zdifflag[0,i] = residuals[i+lags2-1]; 
    } 

    for(int j = 0; j < n-lags2+1; j++) { 
     for(int i = 0; i < lags2; i++) { 
      x(j,i) = diff[j+lags2-1-i]; 
      if (i > 0) { 
       zdifflag(j,i) = x(j,i); 
      } 
     } 
     zdiff[j] = x(j,0); 
    } 

    int length = zdifflag.nrow(), k = zdifflag.ncol(); 

    arma::mat X(zdifflag.begin(), length, k, false); // reuses memory and avoids extra copy 
    arma::colvec y(zdiff.begin(), zdiff.size(), false); 
    arma::colvec coef = arma::solve(X, y); // fit model y ~ X 
    arma::colvec res = y - X*coef;   // residuals 

    // std.errors of coefficients 
    //arma::colvec res = y - X*coef[0]; 
    // sqrt(sum(residuals^2)/(length - k)) 
    double s2 = std::inner_product(res.begin(), res.end(), res.begin(), 0.0)/(length - k);               
    arma::colvec std_err = arma::sqrt(s2 * arma::diagvec(arma::pinv(arma::trans(X)*X))); 

    return coef[0]/std_err[0];  

}',depends = "RcppArmadillo", includes="#include <RcppArmadillo.h>") 

Answer 1

Обычно я рекомендую помещать код в небольшой пакет и каждый параллельный рабочий загружает пакет. Это, как известно, работает как в последовательном, так и в параллельном режимах, тогда как для функции ad-hoc может быть слишком хрупкой для параллельного выполнения, полагаясь на cppFunction().