From: Thomas Lumley <tlumley_at_u.washington.edu>

Date: Mon 25 Dec 2006 - 03:52:49 GMT

https://stat.ethz.ch/mailman/listinfo/r-devel Received on Mon Dec 25 14:56:44 2006

Date: Mon 25 Dec 2006 - 03:52:49 GMT

*>
*

> I have data that is stored in R lists. I would like to pass the list to

*> C code using the .Call function, perform mathematical operations on the
**> elements of the lists (which might be other lists), and return results
**> in a new list.
**>
**> As a simple example, suppose I have a list of n matrices (dimension
**> k-by-k) that I want to invert, and return the result in another list.
**> I wrote the C code below that reads the list as a SEXP. Then, for each
**> of the i=1..n matrices, I convert the R vector to a gsl matrix, perform
**> the linear algebra operations and populate the new list with the
**> inverted matrix.
**>
**> The problem is that I can't figure out the best way to get the real
**> elements into the list SEXP to return it to R. When I run the code
**> below, the compiler (R CMD SHLIB sexp1.c) gives me:
**>
**> -bash-3.00$ R CMD SHLIB sexp1.c
**> gcc -I/usr/lib64/R/include -I/usr/lib64/R/include -I/usr/include -I/usr/local/include -fpic -O2 -g -c sexp1.c -o sexp1.o
**> sexp1.c: In function `invMatList':
**> sexp1.c:37: warning: passing arg 3 of `SET_VECTOR_ELT' from incompatible pointer type gcc -shared -L/usr/local/lib64 -o sexp1.so sexp1.o -L/usr/lib64 -lgsl -lgslcblas -L/usr/lib -L/usr/lib64/R/lib -lR
**>
*

<snip>

> SEXP r = PROTECT(allocVector(VECSXP,n));

*> double rm[k][k];
*

This allocates memory for one k x k C array wherever C allocates temporary storage. You need to allocate n vectors of length k*k on the R heap to create a list of n matrices in R. So you need to use allocVector(k*k, REALSXP) to allocate space for each element of the list, inside the loop, and then copy the data into that rather than into rm.

SET_VECTOR_ELT is giving a warning about incompatible pointer types, but the more important issue is that you need to return pointers to n different blocks of memory that are not going to go away after your function finishes.

-thomas

*>
**>
*

> for (i=0; i<n; i=i+1) {

*>
**> for (a=0; a<k; ++a){ // convert to gsl matrix
**> for (b=0; b<k; ++b){
**> ind = b*k + a;
**> m = REAL(VECTOR_ELT(list,i))[ind];
**> gsl_matrix_set(mat, b, a, m);
**> }
**> } // end conversion loops
**>
**> gsl_linalg_LU_decomp(mat, p, &s);
**> gsl_linalg_LU_invert(mat, p, inv);
**>
**>
**> for (a=0; a<k; ++a){ // convert from gsl matrix
**> for (b=0; b<k; ++b){
**> rm[b][a] = gsl_matrix_get(inv,b,a);
**> }
**> }
**> SET_VECTOR_ELT(r, i, rm);
**> }
**> UNPROTECT(1);
**> return (r);
**> }
**>
**> ______________________________________________
**> R-devel@r-project.org mailing list
**> https://stat.ethz.ch/mailman/listinfo/r-devel
**>
*

Thomas Lumley Assoc. Professor, Biostatistics tlumley@u.washington.edu University of Washington, Seattle ______________________________________________R-devel@r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-devel Received on Mon Dec 25 14:56:44 2006

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.1.8, at Tue 26 Dec 2006 - 09:31:03 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.
*