From wdunlap at tibco.com  Thu Apr 22 00:07:04 2010
From: wdunlap at tibco.com (William Dunlap)
Date: Wed, 21 Apr 2010 15:07:04 -0700
Subject: [Rd] suggestion how to use memcpy in duplicate.c
In-Reply-To: <4BCF52C6.8020303@r-enthusiasts.com>
References: <hqn72q$uev$1@dough.gmane.org> <4BCF52C6.8020303@r-enthusiasts.com>
Message-ID: <77EB52C6DD32BA4D87471DCD70C8D70002CE6BFA@NA-PA-VBE03.na.tibco.com>

If I were worried about the time this loop takes,
I would avoid using i%nt.  For the attached C code
compile with gcc 4.3.3 with -O2 I get 
  > # INTEGER() in loop
  > system.time( r1 <- .Call("my_rep1", 1:3, 1e7) )
     user  system elapsed
    0.060   0.012   0.071

  > # INTEGER() before loop
  > system.time( r2 <- .Call("my_rep2", 1:3, 1e7) )
     user  system elapsed
    0.076   0.008   0.086

  > # replace i%src_length in loop with j=0 before loop and
  > #    if(++j==src_length) j=0 ;
  > # in the loop.
  > system.time( r3 <- .Call("my_rep3", 1:3, 1e7) )
     user  system elapsed
    0.024   0.028   0.050
  > identical(r1,r2) && identical(r2,r3)
  [1] TRUE

The C code is:
#define USE_RINTERNALS /* pretend we are in the R kernel */
#include <R.h>
#include <Rinternals.h>


SEXP my_rep1(SEXP s_src, SEXP s_dest_length)
{
    int src_length = length(s_src) ;
    int dest_length = asInteger(s_dest_length) ;
    int i,j ;
    SEXP s_dest ;
    PROTECT(s_dest = allocVector(INTSXP, dest_length)) ;
    if(TYPEOF(s_src) != INTSXP) error("src must be integer data") ;
    for(i=0;i<dest_length;i++) {
        INTEGER(s_dest)[i] = INTEGER(s_src)[i % src_length] ;
    }
    UNPROTECT(1) ;
    return s_dest ;
}
SEXP my_rep2(SEXP s_src, SEXP s_dest_length)
{
    int src_length = length(s_src) ;
    int dest_length = asInteger(s_dest_length) ;
    int *psrc = INTEGER(s_src) ;
    int *pdest ;
    int i ;
    SEXP s_dest ;
    PROTECT(s_dest = allocVector(INTSXP, dest_length)) ;
    pdest = INTEGER(s_dest) ;
    if(TYPEOF(s_src) != INTSXP) error("src must be integer data") ;
    /* end of boilerplate */
    for(i=0;i<dest_length;i++) {
        pdest[i] = psrc[i % src_length] ;
    }
    UNPROTECT(1) ;
    return s_dest ;
}
SEXP my_rep3(SEXP s_src, SEXP s_dest_length)
{
    int src_length = length(s_src) ;
    int dest_length = asInteger(s_dest_length) ;
    int *psrc = INTEGER(s_src) ;
    int *pdest ;
    int i,j ;
    SEXP s_dest ;
    PROTECT(s_dest = allocVector(INTSXP, dest_length)) ;
    pdest = INTEGER(s_dest) ;
    if(TYPEOF(s_src) != INTSXP) error("src must be integer data") ;
    /* end of boilerplate */
    for(j=0,i=0;i<dest_length;i++) {
        *pdest++ = psrc[j++] ;
        if (j==src_length) {
            j = 0 ;
        }
    }
    UNPROTECT(1) ;
    return s_dest ;
}

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com  

> -----Original Message-----
> From: r-devel-bounces at r-project.org 
> [mailto:r-devel-bounces at r-project.org] On Behalf Of Romain Francois
> Sent: Wednesday, April 21, 2010 12:32 PM
> To: Matthew Dowle
> Cc: r-devel at stat.math.ethz.ch
> Subject: Re: [Rd] suggestion how to use memcpy in duplicate.c
> 
> Le 21/04/10 17:54, Matthew Dowle a ?crit :
> >
> >> From copyVector in duplicate.c :
> >
> > void copyVector(SEXP s, SEXP t)
> > {
> >      int i, ns, nt;
> >      nt = LENGTH(t);
> >      ns = LENGTH(s);
> >      switch (TYPEOF(s)) {
> > ...
> >      case INTSXP:
> >      for (i = 0; i<  ns; i++)
> >          INTEGER(s)[i] = INTEGER(t)[i % nt];
> >      break;
> > ...
> >
> > could that be replaced with :
> >
> >      case INTSXP:
> >      for (i=0; i<ns/nt; i++)
> >          memcpy((char *)DATAPTR(s)+i*nt*sizeof(int), (char 
> *)DATAPTR(t),
> > nt*sizeof(int));
> >      break;
> 
> or at least with something like this:
> 
> int* p_s = INTEGER(s) ;
> int* p_t = INTEGER(t) ;
> for( i=0 ; i < ns ; i++){
> 	p_s[i] = p_t[i % nt];
> }
> 
> since expanding the INTEGER macro over and over has a price.
> 
> > and similar for the other types in copyVector.  This won't 
> help regular
> > vector copies, since those seem to be done by the 
> DUPLICATE_ATOMIC_VECTOR
> > macro, see next suggestion below, but it should help 
> copyMatrix which calls
> > copyVector, scan.c which calls copyVector on three lines, 
> dcf.c (once) and
> > dounzip.c (once).
> >
> > For the DUPLICATE_ATOMIC_VECTOR macro there is already a 
> comment next to it
> > :
> >
> >      <FIXME>: surely memcpy would be faster here?
> >
> > which seems to refer to the for loop  :
> >
> >      else { \
> >      int __i__; \
> >      type *__fp__ = fun(from), *__tp__ = fun(to); \
> >      for (__i__ = 0; __i__<  __n__; __i__++) \
> >        __tp__[__i__] = __fp__[__i__]; \
> >    } \
> >
> > Could that loop be replaced by the following ?
> >
> >     else { \
> >     memcpy((char *)DATAPTR(to), (char *)DATAPTR(from), 
> __n__*sizeof(type)); \
> >     }\
> >
> > In the data.table package, dogroups.c uses this technique, 
> so the principle
> > is tested and works well so far.
> >
> > Are there any road blocks preventing this change, or is 
> anyone already
> > working on it ?  If not then I'll try and test it (on 
> Ubuntu 32bit) and
> > submit patch with timings, as before.  Comments/pointers 
> much appreciated.
> >
> > Matthew
> >
> > ______________________________________________
> > R-devel at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-devel
> >
> >
> 
> 
> -- 
> Romain Francois
> Professional R Enthusiast
> +33(0) 6 28 91 30 30
> http://romainfrancois.blog.free.fr
> |- http://bit.ly/9aKDM9 : embed images in Rd documents
> |- http://tr.im/OIXN : raster images and RImageJ
> |- http://tr.im/OcQe : Rcpp 0.7.7
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
> 



[Rd] suggestion how to use memcpy in duplicate.c
William Dunlap