[Rd] Using matprod from array.c

[Prof Brian Ripley]

You need to link against -lRblas as well, since that is where the BLAS 
library entry points are.

See `Writing R Extensions' for a platform-independent way to do this via 
variables in Makevars (you need ${BLAS_LIBS) $(FLIBS)).

On Wed, 12 Oct 2005, Heather Turner wrote:

> Hi,
>
> I was wondering if I could use the matprod function from array.c in my own C routine.  I tried the following as a test
>
> /* my_matprod.c */
>
> # include <Rinternals.h> /* for REAL, SEXP etc */
> # include <R_ext/Applic.h> /* array.c says need for dgemm */
>
> /* following copied from array.c */
> static void matprod(double *x, int nrx, int ncx,
> 		    double *y, int nry, int ncy, double *z)
> {
>    char *transa = "N", *transb = "N";
>    int i,  j, k;
>    double one = 1.0, zero = 0.0, sum;
>    Rboolean have_na = FALSE;
>
>    if (nrx > 0 && ncx > 0 && nry > 0 && ncy > 0) {
> 	/* Don't trust the BLAS to handle NA/NaNs correctly: PR#4582
> 	 * The test is only O(n) here
> 	 */
> 	for (i = 0; i < nrx*ncx; i++)
> 	    if (ISNAN(x[i])) {have_na = TRUE; break;}
> 	if (!have_na)
> 	    for (i = 0; i < nry*ncy; i++)
> 		if (ISNAN(y[i])) {have_na = TRUE; break;}
> 	if (have_na) {
> 	    for (i = 0; i < nrx; i++)
> 		for (k = 0; k < ncy; k++) {
> 		    sum = 0.0;
> 		    for (j = 0; j < ncx; j++)
> 			sum += x[i + j * nrx] * y[j + k * nry];
> 		    z[i + k * nrx] = sum;
> 		}
> 	} else
> 	    F77_CALL(dgemm)(transa, transb, &nrx, &ncy, &ncx, &one,
> 			    x, &nrx, y, &nry, &zero, z, &nrx);
>    } else /* zero-extent operations should return zeroes */
> 	for(i = 0; i < nrx*ncy; i++) z[i] = 0;
> }
>
> /* test function: matrix multiplication of nr by nc matrix with nc-vector */
> SEXP my_matprod(SEXP M, SEXP v, SEXP nr, SEXP nc) {
>  R_len_t  nrm = INTEGER(nr)[0], ncm = INTEGER(nc)[0];
>  SEXP ans;
>
>  PROTECT(ans = allocMatrix(REALSXP, nrm, 1));
>  matprod(REAL(M), nrm, ncm,
> 	  REAL(v), ncm, 1, REAL(ans));
>  UNPROTECT(1);
>  return(ans);
> }
>
> When I try to make the DLL I get the following
> D:\C_routines>RCMD SHLIB my_matprod.c
> making my_matprod.d from my_matprod.c
> gcc   -IC:/R/tex/R-2.2.0/include -Wall -O2   -c my_matprod.c -o my_matprod.o
> ar cr my_matprod.a my_matprod.o
> ranlib my_matprod.a
> gcc  --shared -s  -o my_matprod.dll my_matprod.def my_matprod.a  -LC:/R/tex/R-2.
> 2.0/src/gnuwin32   -lg2c -lR
> my_matprod.a(my_matprod.o)(.text+0x19a):my_matprod.c: undefined reference to `dg
> emm_'
> collect2: ld returned 1 exit status
> make: *** [my_matprod.dll] Error 1
>
> I'm using Windows XP and I'm using the MinGW gcc. It works fine if I comment out the Fortran call so that the method for have_na is always used.
>
> Do I need to include another header file to use dgemm? Is there a better way to use matprod than just to copy the code?
>
> Any help appreciated,
>
> Heather
>
>
> Dr H Turner
> Research Assistant
> Dept. of Statistics
> The University of Warwick
> Coventry
> CV4 7AL
>
> Tel: 024 76575870
> Url: www.warwick.ac.uk/go/heatherturner
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
>

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595