kondal
kondal

Reputation: 161

How to write Rcpp function for simple matrix multiplication in R

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

Answers (2)

Dirk is no longer here
Dirk is no longer here

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

Andrie
Andrie

Reputation: 179558

You have to keep in mind that 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

Related Questions