[Rcpp-devel] bizarre behavior of RcppMatrix
Romain Francois
romain at r-enthusiasts.com
Tue May 25 09:37:48 CEST 2010
Hello,
This has nothing to do with default arguments, wrappers, etc ... on the
R side.
The issue is the difference between :
> typeof( matrix((1:9), 3, 3) )
[1] "integer"
> typeof( matrix(seq(1,9)^2, ncol=3) )
[1] "double"
The constructor of NumericMatrix has no choice but to cast the input
matrix to a numeric matrix if it is passed an integer matrix, therefore
when "matrix" is an integer matrix, there are three R objects involved,
giving you what you call "correct output".
When you pass it a numeric matrix, no cast is required, so there is only
one R object, shared by matrix, orig and plus. which gives you what you
call "incorrect output".
Further comments:
- don't use Pairlist unless you understand what a pairlist is in R, e.g.
what the R functions pairlist, as.pairlist, etc ... do. If you really
want to return a list, then use List::create, something like:
using namespace Rcpp;
return List::create( _["original"] = orig, _["plus"] = plus ) ;
This way, you'll get a real list (i.e. what the R function "list"
produces) which is more likely to be what you want
- Rcpp objects (Rcpp::NumericMatrix) don't copy data if they don't need
to, they just alter the way the data is perceived, so that you can
manipulate the matrix using c++ syntax, operators, etc ... hence Dirk's
comments. If you really wanted two distinct objects, you can use the
clone function:
Rcpp::NumericMatrix orig( matrix ) ;
Rcpp::NumerixMatrix plus = Rcpp::clone( orig ) ;
This will give you the deep copy Dirk teased you with.
Romain
Le 25/05/10 03:58, Yi-Hsin Erica Tsai a écrit :
> Hi,
>
> I've just started using Rcpp with the hopes of speeding up some R loops
> through replacement with C++ code. I've been struggling with
> understanding the behavior of the NumericMatrix datatype and more
> generally how to pass matrices between R and C++.
>
> I seem to be getting bizarre behavior when I pass default parameters to
> my C++ functions. For example, take this C++ function and the
> corresponding R wrapper function:
>
> ______________________________________
>
> //testRcppMatrix.cpp
>
> #include <Rcpp.h>
>
> RcppExport SEXP plusonediag_RcppMatrix(SEXP matrix)
> {
> Rcpp::NumericMatrix orig(matrix);
> Rcpp::NumericMatrix plus(matrix);
>
> int n = orig.rows();
> int k = orig.cols();
>
> for (int i=0; i<n; i++)
> {
> if(i <= k)
> {
> plus(i,i) = plus(i,i) + 1;
> }
> }
>
> Rcpp::Pairlist res(Rcpp::Named( "original", orig),
> Rcpp::Named( "plus", plus)
> );
>
> return res;
> }
> ______________________________________
>
> #testRcpp.r
>
> plusonediag_RcppMatrix_wrapper <- function(mat=matrix(seq(1,9)^2, ncol=3))
> {
> ## Make the call...
> val <- .Call("plusonediag_RcppMatrix",
> mat)
> return(val)
> }
>
> ______________________________________
>
> If I pass actual values to the plusonediag_RcppMatrix_wrapper, then I
> get correct output; however, if I try to use the default parameter (mat)
> or pass a matrix equal to the default matrix, then the original matrix
> appears to be overwritten. See the session output below.
>
> Anyone have any idea what I'm doing wrong? I'd appreciate any help or
> insight into this problem.
>
> I've attached the test files if that's helpful.
>
> Thanks and best wishes,
> --Erica Tsai
>
> > # main.r
> >
> > library(Rcpp)
> Loading required package: inline
> > dyn.load("testRcppMatrix.so")
> > source("testRcppMatrix.r")
> >
> > foo <- matrix(1:9, 3, 3)
> > bar <- foo ^ 2
> >
> > plusonediag_RcppMatrix_wrapper #default mat is same as bar
> function(mat=matrix(seq(1,9)^2, ncol=3))
> {
> ## Make the call...
> val <- .Call("plusonediag_RcppMatrix",
> mat)
> return(val)
> }
> > plusonediag_RcppMatrix_wrapper() #does not produce correct output!!
> $original
> [,1] [,2] [,3]
> [1,] 2 16 49
> [2,] 4 26 64
> [3,] 9 36 82
>
> $plus
> [,1] [,2] [,3]
> [1,] 2 16 49
> [2,] 4 26 64
> [3,] 9 36 82
>
> > plusonediag_RcppMatrix_wrapper(foo) #produces correct output
> $original
> [,1] [,2] [,3]
> [1,] 1 4 7
> [2,] 2 5 8
> [3,] 3 6 9
>
> $plus
> [,1] [,2] [,3]
> [1,] 2 4 7
> [2,] 2 6 8
> [3,] 3 6 10
>
> > plusonediag_RcppMatrix_wrapper(bar) #does not produce correct output!!
> $original
> [,1] [,2] [,3]
> [1,] 2 16 49
> [2,] 4 26 64
> [3,] 9 36 82
>
> $plus
> [,1] [,2] [,3]
> [1,] 2 16 49
> [2,] 4 26 64
> [3,] 9 36 82
>
> > plusonediag_RcppMatrix_wrapper(matrix((1:9), 3, 3)) #produces correct
> output
> $original
> [,1] [,2] [,3]
> [1,] 1 4 7
> [2,] 2 5 8
> [3,] 3 6 9
>
> $plus
> [,1] [,2] [,3]
> [1,] 2 4 7
> [2,] 2 6 8
> [3,] 3 6 10
>
> > plusonediag_RcppMatrix_wrapper(matrix((1:9)^2, 3, 3)) #does not
> produce correct output!!
> $original
> [,1] [,2] [,3]
> [1,] 2 16 49
> [2,] 4 26 64
> [3,] 9 36 82
>
> $plus
> [,1] [,2] [,3]
> [1,] 2 16 49
> [2,] 4 26 64
> [3,] 9 36 82
>
> >
>
>
>
> _______________________________________________
> Rcpp-devel mailing list
> Rcpp-devel at lists.r-forge.r-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
--
Romain Francois
Professional R Enthusiast
+33(0) 6 28 91 30 30
http://romainfrancois.blog.free.fr
|- http://bit.ly/cork4b : highlight 0.1-8
|- http://bit.ly/bklUXt : RcppArmadillo 0.2.1
`- http://bit.ly/936ck2 : Rcpp 0.8.0
More information about the Rcpp-devel
mailing list