Reputation: 161
I have wrote a Rcpp code to compute element wise matrix multiplication in R. But when try to run this code R stops working and its exiting. How to correct this function?
Thanks in advance.
library(Rcpp)
func <- 'NumericMatrix mmult( NumericMatrix m , NumericMatrix v, bool byrow=true )
{
if( ! m.nrow() == v.nrow() ) stop("Non-conformable arrays") ;
if( ! m.ncol() == v.ncol() ) stop("Non-conformable arrays") ;
NumericMatrix out(m) ;
for (int i = 1; i <= m.nrow(); i++)
{
for (int j = 1; j <= m.ncol(); j++)
{
out(i,j)=m(i,j) * v(i,j) ;
}
}
return out ;
}'
cppFunction( func )
m1<-matrix(1:4,2,2)
m2<-m1
r1<-mmult(m1,m2)
r2<-m1*m2
Upvotes: 1
Views: 999
Reputation: 368519
The (at least to me) obvious choice is to use RcppArmadillo:
R> cppFunction("arma::mat matmult(arma::mat A, arma::mat B) { return A % B; }",
+ depends="RcppArmadillo")
R> m1 <- m2 <- matrix(1:4,2,2)
R> matmult(m1,m2)
[,1] [,2]
[1,] 1 9
[2,] 4 16
R>
as Armadillo is strongly typed, and has an element-by-element multiplication operator (%
) which we use in the one-liner it takes.
Upvotes: 4
Reputation: 179558
You have to keep in mind that c++ uses 0 indexed arrays. (See Why does the indexing start with zero in 'C'? and Why are zero-based arrays the norm? .)
So you need to define your loop to run from 0
to m.nrow() - 1
Try this:
func <- '
NumericMatrix mmult( NumericMatrix m , NumericMatrix v, bool byrow=true )
{
if( ! m.nrow() == v.nrow() ) stop("Non-conformable arrays") ;
if( ! m.ncol() == v.ncol() ) stop("Non-conformable arrays") ;
NumericMatrix out(m) ;
for (int i = 0; i < m.nrow(); i++)
{
for (int j = 0; j < m.ncol(); j++)
{
out(i,j)=m(i,j) * v(i,j) ;
}
}
return out ;
}
'
Then I get:
> mmult(m1,m2)
[,1] [,2]
[1,] 1 9
[2,] 4 16
> m1*m2
[,1] [,2]
[1,] 1 9
[2,] 4 16
Upvotes: 4