Reputation: 4114
Simple matrix multiplication in RcppEigen
gives a matrix of correct dimensions, but the first element is repeated to all the elements.
Rcpp
source code using RcppEigen
for matrix multiply:
#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);
}
Calling from 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
Upvotes: 2
Views: 138
Reputation: 368509
You are going too fast. Use intermediate results -- your product is now what you think it is. And even with the rest of what you have, you are getting confused about when Eigen::Map
makes sense and when not. You also do not need as<>
and wrap()
-- RcppEigen
takes care of that.
#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>
Upvotes: 4