The key problem is that when you are multiplying the matrix by a scalar in C++, you're using Rcpp
's syntactic sugar for *
, which is vectorized. For whatever reason, it doesn't understand how to return a matrix (I haven't looked at the documentation extensively).
If we instead multiply each element of each matrix by the scalar, you get the expected results:
FooList.R
// [[Rcpp::depends(RcppArmadillo)]]
#include <RcppArmadillo.h>
using namespace Rcpp;
using namespace arma;
// I got this template from here: http://stackoverflow.com/a/18014655/1315767
template <typename WHAT>
class ListOf : public List {
public:
template <typename T>
ListOf( const T& x) : List(x){}
WHAT operator[](int i){ return as<WHAT>( ( (List*)this)->operator[]( i) ) ; }
} ;
// [[Rcpp::export]]
List FooList(NumericVector fi1, ListOf<NumericMatrix> Ct){
List TempList(Ct.size());
NumericMatrix ct(2,2);
for(int i=0; i<Ct.size(); i++){
ct = Ct[i] ;
for (int j=0; j < ct.nrow(); j++) {
for (int k=0; k < ct.ncol(); k++) {
ct(j, k) *= fi1[i]; // Multiply each element of the matrix by the scalar in fi1
}
}
TempList[i] = ct;
}
return TempList;
}
Interactive Session:
> sourceCpp("FooList.cpp")
> A <- replicate(3,matrix(1:4, 2), simplify=FALSE) # a list of matrices
> vec <- 0.5 * c(1:3) # a vector
> FooList(vec, A)
[[1]]
[,1] [,2]
[1,] 0.5 1.5
[2,] 1.0 2.0
[[2]]
[,1] [,2]
[1,] 1 3
[2,] 2 4
[[3]]
[,1] [,2]
[1,] 1.5 4.5
[2,] 3.0 6.0