Блог

Простое умножение матрицы в RcppEigen дает матрицу правильных размеров, но первый элемент повторяется для всех элементов.

Rcpp исходный код, использующий RcppEigen для умножения матрицы:

 #include <Rcpp.h>
#include <RcppEigen.h>


using namespace Rcpp;

// [[Rcpp::depends(RcppEigen)]]
// [[Rcpp::export]]
NumericMatrix myfun(const NumericVector a ,  const NumericVector b) {


        const Eigen::Map<Eigen::MatrixXd> a_eig(Rcpp::as<Eigen::Map<Eigen::MatrixXd> >(a));
        const Eigen::Map<Eigen::MatrixXd> b_eig(Rcpp::as<Eigen::Map<Eigen::MatrixXd> >(b));

        return wrap(a_eig * b_eig);
}
  

Вызов из R :

 a=matrix(data=1,nrow=10,ncol=1)
b=sample(100,10)
a * b


      [,1]
 [1,]   67
 [2,]   59
 [3,]   19
 [4,]   68
 [5,]   83
 [6,]    4
 [7,]   28
 [8,]   88
 [9,]   97
[10,]   43

myfun(a,b)

      [,1]
 [1,]   67
 [2,]   67
 [3,]   67
 [4,]   67
 [5,]   67
 [6,]   67
 [7,]   67
 [8,]   67
 [9,]   67
[10,]   67
  

Вы идете слишком быстро. Используйте промежуточные результаты — ваш продукт теперь такой, каким вы его считаете. И даже с остальным, что у вас есть, вы путаетесь в том, когда Eigen::Map имеет смысл, а когда нет. Вам также не нужно as<> и wrap() — RcppEigen позаботится об этом.

Код

 #include <Rcpp.h>
#include <RcppEigen.h>

// [[Rcpp::depends(RcppEigen)]]

// [[Rcpp::export]]
Eigen::MatrixXd myfun(Eigen::Map<Eigen::MatrixXd> a,
                      Eigen::Map<Eigen::MatrixXd> b) {
  Eigen::MatrixXd res = a * b;
  return res;
}
  

Вывод

  R> Rcpp::sourceCpp("~/git/stackoverflow/55797031/answer.cpp")
 R> myfun(matrix(2,2,2), matrix(3,2,2))
      [,1] [,2]
 [1,]   12   12
 [2,]   12   12
 R>
 R> myfun(matrix(as.numeric(1:4),2,2), matrix(as.numeric(4:1),2,2))
      [,1] [,2]
 [1,]   13    5
 [2,]   20    8
 R>