[Rcpp-devel] Rcpp shortcut for R's as.matrix and as.vector

[Romain Francois]

Le 30/03/10 15:46, Douglas Bates a écrit :
> On Tue, Mar 30, 2010 at 1:44 AM, Romain Francois
> <romain at r-enthusiasts.com>  wrote:
>> Le 29/03/10 15:36, Douglas Bates a écrit :
>>>
>>> Most Lapack routines for solving linear systems from factorizations
>>> allow the right hand side to be a matrix, thereby providing for
>>> solving multiple systems in one  call.  I would like to preserve this
>>> capability while also allowing for the most common case where the rhs
>>> is a vector.  Also, I would like to return a vector solution rather
>>> than a single matrix solution when it is a single column.
>>>
>>> I can accomplish these "as.matrix" and "as.vector" translations by
>>> calling the corresponding R functions.
>>>
>>> static Rcpp::Function asMatrix("as.matrix")
>>> static Rcpp::Function asVector("as.vector")
>>>
>>> and using those within Rcpp.  Is there an alternative within Rcpp?
>>
>> Not really. You can do :
>>
>> x.attr( "dim" ) = R_NilValue ;
>>
>> or
>>
>> x.attr( "dim" ) = IntegerVector::create( 5, 5 ) ;
>>
>>
>> Here is a sketch for something we might want to add in the Matrix template:
>>
>> SEXP as_vector(){
>>         SEXP x = ::Rf_duplicate( m_sexp ) ;
>>         ::Rf_setAttrib( x, ::Rf_install( "dim" ), R_NilValue ) ;
>>         // perhaps also dimnames
>>         return x ;
>> }
>>
>> or maybe this should return a VECTOR (which is a typedef to the appropriate
>> Vector<RTYPE>  in Matrix.
>
> I'm asking all these questions because I haven't used templates in C++
> and I keep thinking of ways that I would do things in other languages.

Me too. templates are quite weird to start with, and then you can never 
look back.

>   With self-describing objects in R you can be rather flexible about
> return types and for C code written against R you tend to have all
> your objects expressed as SEXPs so again you have flexibility.

Rcpp only gives the illusion of strong-typing. It can not afford actual 
strong typing because it is too dependent on R, and R uses loose typing. 
If you create a NumericVector and then assign a "dim" attribute to 
something of dimension 2, it is an R matrix, but it is not an 
Rcpp::NumericMatrix, as there is no way to change the c++ class of an 
object dynamically.

this is also why we have indexing operators that takes two ints defined 
in Vector and not Matrix.

The only reason for the Matrix template is that it can expose a 
different set of constructors, but what we need to make explicit is that 
it does not guarantee that the object will allways be a matrix. If you 
create a NumericMatrix and then you wipe the "dim" attribute, you still 
have a NumericMatrix, but this is no longer a matrix.

> For
> the time being I will simply use R's as.matrix in the R interface code
> or as a Function object in Rcpp code.

Sure. But would you say it is a good idea or a bad idea to have these 
helper functions in Rcpp ?

Matrix<RTYPE> Vector::as_matrix() ;
Vector<RTYPE> Matrix::as_vector() ;

or maybe they should both be in Vector, or maybe we should just have
Vector::dim().

>> I don't really have time to play with this in the coming days, but if you
>> try to add the function and some unit tests, I'll take time to review them.
>> I could also sketch out Vector::as_matrix but I guess you also can.
>>
>>> I realize that I need to return another object so that it can have a
>>> different class within C++ but I guess I would like
>>>
>>> NumericMatrix mat = Rcpp::as<NumericMatrix>(myNumericVector);
>>
>> I'm mot sure this will work (actually I'm almost sure this does not, but I
>> did not test).
>>
>> You would essentially be calling the NumericMatrix( SEXP ) constructor which
>> would not be happy about the argument not being a matrix.
>>
>>> or something similar.  It is not a big deal because the calls to the R
>>> function are not a problem (or I could manipulate the SEXP internally)
>>> but I'm wondering if I am missing something simple here.



-- 
Romain Francois
Professional R Enthusiast
+33(0) 6 28 91 30 30
http://romainfrancois.blog.free.fr
|- http://tr.im/OIXN : raster images and RImageJ
|- http://tr.im/OcQe : Rcpp 0.7.7
`- http://tr.im/O1wO : highlight 0.1-5