почему Rcpp медленнее, чем R при использовании умножения матриц?

Question

Чтобы ускорить мой пакет, который включает много матричных вычислений, я использую Rcpp, чтобы переписать весь код.Однако некоторые функции работают даже медленнее, чем раньше.Я использую microbenchmark для анализа и нахожу, что умножение матриц в Rcpp происходит медленнее.Почему это произойдет?И как ускорить мой пакет?Большое спасибо.Код Rcpp выглядит следующим образом:

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
NumericMatrix mmult(const NumericMatrix& a, const NumericMatrix& b){
if (a.ncol() != b.nrow()) stop ("Incompatible matrix dimensions");
NumericMatrix out(a.nrow(),b.ncol());
NumericVector rm1, cm2;
for (int i = 0; i < a.nrow(); ++i) {
  rm1 = a(i,_);
  for (int j = 0; j < b.ncol(); ++j) {
    cm2 = b(_,j);
    out(i,j) = std::inner_product(rm1.begin(), rm1.end(), cm2.begin(), 0.);
  }
}
return out;}

Код R выглядит следующим образом:

X = matrix(rnorm(10*10,1),10,10)
Y = matrix(rnorm(10*10,1),10,10)


microbenchmark(
  mmult(X,Y),
  X%*%Y)

Результат:

Unit: microseconds
    expr    min      lq      mean median     uq      max neval
 mmult(X, Y) 45.720 48.9860 126.79228 50.385 51.785 6368.512   100
 X %*% Y  5.599  8.8645  12.85787  9.798 10.730  153.486   100

Answer 1

Это противоположное, но ожидаемое результат из того, что видели для умножения матрицы на вектор .Здесь R использует BLAS для выполнения всей тяжелой работы, которая может даже работать параллельно.Вы отбрасываете все оптимизированное управление памятью, выполняемое в библиотеке BLAS, используя умножение наивной матрицы.

Вместо того, чтобы пытаться заново изобретать низкоуровневые вещи, такие как умножение матриц, вы можете попытаться реализовать более крупные части вашегокод, использующий что-то вроде RcppArmadillo, который использует ту же библиотеку BLAS, что и R, но также (и не только!) предлагает удобный синтаксис.