R – C/C++ programming Katia Oleinik koleinik@bu.edu Scientific Computing and Visualization Boston University

1. R – C/C++ programmingKatia Oleinikkoleinik@bu.eduScientific Computing and VisualizationBoston University http://www.bu.edu/tech/research/training/tutorials/list/

2. R-C/C++ programming Goal – performance enhancement. Benefits – use of existing C/C++ libraries and memory management Base R package provides 3 types of interfaces between R and C/C++ .C() .Call() .External() – used to create R packages There are other R packages that provide interface between R and C/C++ (and other languages such as FORTRAN and Python): Rcpp

3. R-C/C++ programming .C() interface exC1.c /* exC1.c – example C function to be called from R */ void exampleC1(int *iVec){ iVec[0] = 7; return; } Important: Function returns no values – it is VOID All the values that need to be changed/returned by a function must be passed through its arguments.

4. R-C/C++ programming .C() interface katana:~ % R CMD SHLIB exC1.c gcc -std=gnu99 -I/usr/local/IT/R-2.13.2/lib64/R/include -I/usr/local/include -fpic -g -O2 -c exC1.c -o exC1.o gcc -std=gnu99 -shared -L/usr/local/lib64 -o exC1.so exC1.o katana:~ % Important: In linux (and R) environment commands are case sensitive!

5. R-C/C++ programming .C() interface ># load C function to R workspace > dyn.load("exC1.so") > Note: In windows after the function is compiled it will be named exC1.dll

6. R-C/C++ programming .C() interface ># load C function to R workspace > dyn.load("exC1.so") > # create a vector > iv <- 1:3

7. R-C/C++ programming .C() interface ># load C function to R workspace > dyn.load("exC1.so") > # create a vector > iv <- 1:3 > # call c-function > out <- .C("exampleC1", newVec = as.integer(iv))

8. R-C/C++ programming .C() interface ># load C function to R workspace > dyn.load("exC1.so") > # create a vector > iv <- 1:3 > # call c-function > out <- .C("exampleC1", newVec = as.integer(iv)) > out $newVec [1] 7 2 3

9. R-C/C++ programming .C() interface • Note: • Allocate memory to the vectors passed to .C in R by creating vectors of the right length • The first argument to .C is a character string of the C function name • The rest of the arguments are R objects to be passed to the C function. • All arguments should be coerced to the correct R storage mode to prevent mismatching of types that can lead to errors • C returns a list object • The second argument in this example is given a name newVec. This name is used to access the component in the returned list object.

10. R-C/C++ programming .C() interface Note: R has to allocate memory for the arrays passed to and from C. R has to pass objects of correct type R copies its arguments prior to passing them to C and then creates a copy of the values passed back from C.

11. R-C/C++ programming exC2.c /* exC2.c – example C function to be called from R */ /* normalize the vector */ include<math.h> include <string.h> void exampleC2(char **c, double *A, double *B, int *ierr){ double len = 0; /*local variable – vector length */ inti; for (i=0; i<3; i++) len += pow( A[i]), 2); /* check if the vector is degenerate */ if ( len < 0.000001){ ierr[0] = -1; /*error – null vector */ stncpy(c,“Error”, 5); return; } /* calculate output vector len = pow(len, 0.5); for (i=0; i<3; i++)B[i] = A[i] / len ; ierr[0] = 0; strncpy(c,“OK”, 2); return; }

12. R-C/C++ programming .C() interface katana:~ % R CMD SHLIB exC2.c gcc -std=gnu99 -I/usr/local/IT/R-2.13.2/lib64/R/include -I/usr/local/include -fpic -g -O2 -c exC2.c -o exC2.o gcc -std=gnu99 -shared -L/usr/local/lib64 -o exC2.so exC2.o

13. R-C/C++ programming .C() interface ># load C function to R workspace > dyn.load("exC2.so") > # create error vector > ierr_in <- 0 > # create input vector > A_in<- c(2, 3, 6) > # create output vector > B_in<- c(0, 0, 0) > # create message vector (make sure it is long enough!) > C_in<- c("")

14. R-C/C++ programming .C() interface > # execute C function > out <- .C("exampleC2", + C_out= as.character(C_in), + A_out= as.numeric(A_in), + B_out= as.numeric(B_in), + ierr_out= as.integer(ierr_in))

15. R-C/C++ programming .C() interface > out $C_out [1] "OK " $A_out [1] 2 3 6 $B_out [1] 0.2857143 0.4285714 0.8571429 $ierr_out [1] 0

16. R-C/C++ programming .C() interface > # create input vector > A_in <- c(0, 0, 0) > # execute C function > out <- .C("exampleC2", "exampleC2", + C_out= as.character(C_in), + A_out= as.numeric(A_in), + B_out= as.numeric(B_in), + ierr_out = as.integer(ierr_in))

17. R-C/C++ programming .C() interface > out $C_out [1] "error " $A_out [1] 0 0 0 $B_out [1] 0 0 0 $ierr_out [1] -1

18. R-C/C++ programming .C() interface Note: To compile more than one C file: R CMD SHLIB file1.c file2.c file3.c The resulting file will be named file1.so

19. R-C/C++ programming • .Call() interface • does not copy arguments before and after calling c-function • it is possible to find the length of the input vector inside c-function • an easier access to wide-range of R – objects • NA (missing values) handling • Access to vectors’ attributes

20. R-C/C++ programming .Call() interface – passing a value exC3.c /* exC3.c – example C function to be called from R with .Call interface*/ /* access R object (scalar value) inside c-function */ include<R.h> /* 2 standard includes for .Call interface) */ include <Rdefines.h> SEXP exampleC3 ( SEXPiValue ){ return (R_NilValue); /* “void” function must return “NULL” value */ } • Note: • All objects passed between R and C/C++ are of type SEXP – SimpleEXPression. • 2 standard includes needed for .Call interface • If function is void it should return R_NilValue object.

21. R-C/C++ programming .Call() interface – passing a value exC3.c /* exC3.c – example C function to be called from R with .Call interface*/ /* access R object (scalar value) inside c-function */ include<R.h> include <Rdefines.h> SEXP exampleC3 ( SEXPiValue ){ intlocal_iValue; /* convert R object to c-accessible variable */ local_iValue= INTEGER_VALUE(iValue); return (R_NilValue); }

22. R-C/C++ programming .Call() interface – passing a value exC3.c /* exC3.c – example C function to be called from R with .Call interface*/ /* access R object (scalar value) inside c-function */ include<R.h> include <Rdefines.h> SEXP exampleC3 ( SEXPiValue ){ intlocal_iValue; /* convert R object to c-accessible variable */ local_iValue= INTEGER_VALUE(iValue); /* print value of the local variable*/ printf(" In exampleC3 iValue = %d\n", local_iValue); return (R_NilValue); }

23. R-C/C++ programming .Call() interface – passing a value ># load C function to R workspace – same as before > dyn.load("exC3.so") ># call C function > out<-.Call("exampleC3", 7) In exampleC3 iValue = 7 ># explore output > out NULL

24. R-C/C++ programming .Call() interface – passing a vector exC4.c /* exC4.c - example C function to be called from R */ /* normalize the vector and return its length */ include<R.h> include <Rdefines.h> include <Rmath.h> SEXP exampleC4 ( SEXP Vector ){ SEXP rLen; return (rLen); /* return a value */ } • Note: • Rmath.h include provides access to many R-functions include rnorm(),rgamma(), etc. • Function should return SEXP object. • .Call() interface allows for changing the function arguments – be careful!

25. R-C/C++ programming exC4.c SEXP exampleC4 ( SEXP Vector ){ SEXP rLen; /* output value – length of a vector */ double * pVector; /* local variable - pointer to the input vector */ double vLen = 0; /* local variable to calculate intermediate values */ intlen; /* local variable – size of the input vector */ inti; /* local variable – loop index */ return (rLen); /* return a value */ }

26. R-C/C++ programming exC4.c SEXP exampleC4 ( SEXP Vector ){ SEXP rLen; double * pVector; double vLen = 0; intlen; inti; /* get the pointer to the vector */ pVector=NUMERIC_POINTER(Vector); return (rLen); /* return a value */ } Note: Use INTEGER_POINTER()and CHARACTER_POINTER() to get pointer to integer and character arrays respectfully

27. R-C/C++ programming exC4.c SEXP exampleC4 ( SEXP Vector ){ SEXP rLen; double * pVector; double vLen = 0; intlen; inti; /* get the pointer to the vector */ pVector=NUMERIC_POINTER(Vector); /* number of elements in the array */ len = length(Vector); return (rLen); /* return a value */ } Note: We can get the size of the input R-vector !

28. R-C/C++ programming exC4.c SEXP exampleC4 ( SEXP Vector ){ SEXP rLen; double * pVector; double vLen = 0; intlen; inti; pVector=NUMERIC_POINTER(Vector); len= length(Vector); /* allocate storage for integer variable (array works also!) */ PROTECT(rLen= NEW_NUMERIC(1)); UNPROTECT(1); return (rLen); /* return a value */ } • Note: • To allocate integer and character arrays use NEW_INTEGER(len)and NEW_CHARACTER(len) functions respectfully • PROTECT() and UNPROTECT() command must be balanced!

29. R-C/C++ programming exC4.c SEXP exampleC4 ( SEXP Vector ){ SEXP rLen; double * pVector; double vLen = 0; intlen; inti; pVector=NUMERIC_POINTER(Vector); len= length(Vector); PROTECT(rLen= NEW_NUMERIC(1)); /* calculate the length */ for( i=0; i < len; i++)vLen += pow(pVector[i], 2); if ( vLen > 0.000001){ vLen = pow( vLen,0.5 ); /* Here we are working with a pointer - it WILL change R vector */ for( i=0; i < len; i++ )pVector[i] /= vLen; } /* copy the value of local variable into R-object */ REAL(rLen)[0] = vLen; UNPROTECT(1); return (rLen); /* return a value */ }

30. R-C/C++ programming .Call() interface – passing an array ># load C function to R workspace – same as before > dyn.load("exC4.so") ># define and input array > A_in<- c( 2, 3, 6) ># call C function > out<-.Call("exampleC4", A_in) ># input array changed !!! > A_in [1] 0.2857143 0.4285714 0.8571429 > out [1] 7

31. This tutorial has been made possible by Scientific Computing and Visualization group at Boston University. Katia Oleinikkoleinik@bu.edu