[R] Call by reference or suggest workaround

[Romain Francois]

Le 19/06/10 16:32, Chidambaram Annamalai a écrit :
>
> I have written code to compute multi-indices in R [1] and due to the
> recursive nature of the computation I need to pass around the *same*
> matrix object (where each row corresponds to one multi-index). As pass
> by reference wasn't the default behavior I declared a global matrix
> (mat) and used the<<- operator to write to the global matrix. So the
> usage would be to call genMultiIndices(3,2) for side effects to
> generate all multi-indices of length 3 and sum 2. And then access the
> global matrix.
>
> However, after coding this I can't seem to export the global matrix
> object (in the NAMESPACE file) and still retain mutability since its
> binding is "locked" (R throws an error). Can I somehow unlock this?
>
> Ideally I would want to pass around the same matrix to the recursive
> function. Is that possible? If not, could someone please suggest a
> workaround to use the code in an R package?
>
> [1]: http://dpaste.com/209186/

Hi,

You can use lexical scoping and you might like ?Recall as well.

genMultiIndices <- function(N, V) {
     mat <- matrix(nrow=choose(N + V - 1, V), ncol=N)
	fillMultiIndices <- function(i, j, n, v) {
		if (n == 1) {
	        mat[i, j] <<- v
	    }
	    else if (v == 0) {
	        mat[i, j:(j + n - 1)] <<- 0L
	    }
	    else {
	        rowOffset <- 0
	        # the first element of each multi-index can be any of 0, 1, ..., v
	        for (k in v:0) {
	            times <- choose((n - 1) + (v - k) - 1, (v - k))
	            mat[(i + rowOffset):(i + rowOffset + times - 1), j] <<- k
	            Recall(i + rowOffset, j + 1, n - 1, v - k)
	            rowOffset <- rowOffset + times
	        }
	    }
	}
     fillMultiIndices(1, 1, N, V)
     mat
}



Also, you can consider writing your code in a language that supports 
references, e.g. C++. Here is a start with inline/Rcpp :

require( inline )
require( Rcpp )

genMultiIndices_internal <-  local({
inc <- '
void fillMultiIndices( Rcpp::IntegerMatrix& mat, int i, int j, int n, 
int v ){
	if( n == 1 ){
		mat( (i-1), (j-1) ) = v ;
	} else if( v == 0 ){
		for( int k=j; k < j+n; k++){
			mat( (i-1), (k-1) ) = 0 ;
		}
	} else {
		
		// using the R function
		// I leave it to you to use a C implementation
		Function choose( "choose" ) ;
		
		int rowOffset = 0 ;
		int times ;
		for( int k=v; k>=0; k--){
			times = as<int>( choose( (n-1) + (v-k) - 1, (v-k) ) );
			int start = i + rowOffset ;
			int end   = i + rowOffset + times ;
			for( int z = start; z < end; z++ ){
				mat( z-1 , j-1 ) = k ;
			}
			fillMultiIndices( mat, i + rowOffset, j+1, n-1, v-k ) ;
			rowOffset += times ;
		}
	}
}
'
code <- '
	int N  = as<int>( N_) ;
	int V  = as<int>( V_) ;
	int NR = as<int>( NR_) ;
	Rcpp::IntegerMatrix mat( NR, N ) ;
	fillMultiIndices( mat, 1, 1, N, V ) ;
	return mat ;
'
	.genMultiIndices <- cxxfunction(
		signature( N_ = "integer", V_ = "integer", NR_ = "integer" ),
		code, include = inc, plugin = "Rcpp" )
	function( N, V){
		.genMultiIndices( N, V, choose(N + V - 1, V) )
	}
} )


I've been lazy here and I am using the choose function from R, so there 
is room for some improvement.

 > ( x <- genMultiIndices( 3L , 2L )          )
      [,1] [,2] [,3]
[1,]    2    0    0
[2,]    1    1    0
[3,]    1    0    1
[4,]    0    2    0
[5,]    0    1    1
[6,]    0    0    2
 > ( y <- genMultiIndices_internal( 3L, 2L )  )
      [,1] [,2] [,3]
[1,]    2    0    0
[2,]    1    1    0
[3,]    1    0    1
[4,]    0    2    0
[5,]    0    1    1
[6,]    0    0    2
 > identical( x, y )
[1] TRUE

Romain


-- 
Romain Francois
Professional R Enthusiast
+33(0) 6 28 91 30 30
http://romainfrancois.blog.free.fr
|- http://bit.ly/98Uf7u : Rcpp 0.8.1
|- http://bit.ly/c6YnCi : graph gallery collage
`- http://bit.ly/bZ7ltC : inline 0.3.5