From: Alireza Mahani <alireza.s.mahani_at_gmail.com>

Date: Thu, 14 Jul 2011 08:21:07 -0700 (PDT)

# mvMultiply.r

# comparing performance of matrix multiplication in R (using '%*%' operator)

vs. calling compiled code (using .Call function)

# y [m x 1] = A [m x n] %*% x [n x 1]

}

*# mvMultiply.cc
*

#include <Rinternals.h>

#include <R.h>

}

Date: Thu, 14 Jul 2011 08:21:07 -0700 (PDT)

Below are the results that I get, along with discussion (tR and tCall are in sec):

INCLUDE_DATAPREP,ROWMAJOR,tR,tCall

F,F,13.536,13.875 F,T,13.824,14.299 T,F,13.688,18.167 T,T,13.982,30.730

Interpretation: The execution time for the .Call line is nearly identical to the call to R operator '%*%'. Two data preparation lines for the .Call method create the overhead:

A <- t(A) (~12sec, or 12usec per call)

dim(A) <- dim(A)[1] * dim(A)[2] (~4sec, or 4usec per call)

While the first line can be avoided by providing options in c++ function (as is done in the BLAS API), I wonder if the second line can be avoided, aside from the obvious option of rewriting the R scripts to use vectors instead of matrices. But this defies one of the primary advantages of using R, which is succinct, high-level coding. In particular, if one has several matrices as input into a .Call function, then the overhead from matrix-to-vector transformations can add up. To summarize, my questions are:

1- Do the above results seem reasonable to you? Is there a similar penalty
in R's '%*%' operator for transforming matrices to vectors before calling
BLAS functions?

2- Are there techniques for reducing the above overhead for developers
looking to augment their R code with compiled code?

Regards,

Alireza

# mvMultiply.r

# comparing performance of matrix multiplication in R (using '%*%' operator)

vs. calling compiled code (using .Call function)

# y [m x 1] = A [m x n] %*% x [n x 1]

rm(list = ls())

system("R CMD SHLIB mvMultiply.cc")

dyn.load("mvMultiply.so")

INCLUDE_DATAPREP <- F

ROWMAJOR <- F #indicates whether the c++ function treats A in a row-major or
column-major fashion

m <- 100 n <- 10 N <- 1000000

diffVec <- array(0, dim = N)

tR <- 0.0

tCall <- 0.0

for (i in 1:N) {

A <- runif(m*n); dim(A) <- c(m,n) x <- runif(n) t1 <- proc.time()[3] y1 <- A %*% x tR <- tR + proc.time()[3] - t1 if (INCLUDE_DATAPREP) { t1 <- proc.time()[3] } if (ROWMAJOR) { #since R will convert matrix to vector in a column-majorfashion, if the c++ function expects a row-major format, we need to transpose A before converting it to a vector

A <- t(A) } dim(A) <- dim(A)[1] * dim(A)[2] if (!INCLUDE_DATAPREP) { t1 <- proc.time()[3] } y2 <- .Call("matvecMultiply", as.double(A), as.double(x), as.logical(c(ROWMAJOR))) tCall <- tCall + proc.time()[3] - t1 diffVec[i] <- max(abs(y2 - y1))

}

cat("Data prep time for '.Call' included: ", INCLUDE_DATAPREP, "\n") cat("C++ function expects row-major matrix: ", ROWMAJOR, "\n") cat("Time - Using '%*%' operator in R: ", tR, "sec\n") cat("Time - Using '.Call' function: ", tCall, "sec\n") cat("Maximum difference between methods: ", max(diffVec), "\n")

dyn.unload("mvMultiply.so")

#include <Rinternals.h>

#include <R.h>

extern "C" {

SEXP matvecMultiply(SEXP A, SEXP x, SEXP rowmajor) {

double *rA = REAL(A), *rx = REAL(x), *ry; int *rrm = LOGICAL(rowmajor); int n = length(x); int m = length(A) / n; SEXP y; PROTECT(y = allocVector(REALSXP, m)); ry = REAL(y); for (int i = 0; i < m; i++) { ry[i] = 0.0; for (int j = 0; j < n; j++) { if (rrm[0] == 1) { ry[i] += rA[i * n + j] * rx[j]; } else { ry[i] += rA[j * m + i] * rx[j]; } } } UNPROTECT(1); return(y);

}

}

-- View this message in context: http://r.789695.n4.nabble.com/Performance-of-C-and-Call-functions-vs-native-R-code-tp3665017p3667896.html Sent from the R devel mailing list archive at Nabble.com. ______________________________________________ R-devel_at_r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-develReceived on Thu 14 Jul 2011 - 15:25:18 GMT

This quarter's messages: by month, or sorted: [ by date ] [ by thread ] [ by subject ] [ by author ]

*
Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.
Archive generated by hypermail 2.2.0, at Tue 19 Jul 2011 - 15:50:10 GMT.
*

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-devel.
Please read the posting
guide before posting to the list.
*