[Rcpp-devel] RcppArmadillo: passing matrix columns by reference

Dear All,

I am trying to pass columns from an Armadillo matrix to a function, but
I'd like to pass just a reference to the column, not a copy of the column
and I do not seem to be able to do it "elegantly".

The code below (function f1) I think shows that passing X.col to a
function creates a copy (X.unsafe_col does too). We can pass &X as
argument, and the index of the column, as in f2. And that will not create
a copy. But I think this is not the right way of doing what I want to do
(to begin with, I'd rather not pass the column index to the function).


What am I getting wrong?

Best,





// [[Rcpp::depends(RcppArmadillo)]]
#include <Rcpp.h>
#include <RcppArmadillo.h>

using namespace Rcpp;



double f1(arma::umat Z) {
   Z(0, 0) = 111;
   Z(9, 0) = 111;
   std::cout << "f1, this is Z " << std::endl << Z << std::endl;
   return 33.3;
}

double f2(arma::umat &Z, const int c1) {
   Z(1, 0) = 222;
   return 66.6;
}


double f3(arma::umat &Z) {
   Z(1, 0) = 223;
   return 99.9;
}


// [[Rcpp::export]]
List f0(IntegerVector s1_, IntegerVector c1_){
   const int  s1 = as<int>(s1_);
   const int  c1 = as<int>(c1_);

   arma::umat X(10, s1);
   for(int j = 0; j < s1; ++j) {
     for(int i = 0; i < 10; ++i) {
       X(i, j) = i * 10 + j;
     }
   }

   // in both cases, a copy seems to be made
   //double fitness = f1(X.col(c1));
   double outf1 = f1(X.unsafe_col(c1));
   std::cout << "f0, this is X after f1" << std::endl << X << std::endl;

   double outf2 = f2(X, c1);
   std::cout << "f0, this is X after f2" << std::endl << X << std::endl;


   // double outf3 = f3(X.unsafe_col(c1)); //will not work
   //double outf3 = f3(X.col(c1)); //will not work


   return List::create(Named("X") = wrap(X),
		      Named("of1") = outf1,
		      Named("of2") = outf2);
}

the .col method gives you a subview_col :

   arma_inline       subview_col<eT> col(const uword col_num);
   arma_inline const subview_col<eT> col(const uword col_num) const;

not an umat.

The name of the class implies that it is a "view" class, so a way to
look a data from another class. hence, no data of its own, so cheap copy.

double f3(arma::subview_col<unsigned int> Z) {
   Z(1, 0) = 223;
   return 99.9;
}

Also, including RcppArmadillo.h after Rcpp.h is wrong. You should only
include RcppArmadillo.h.

I should do something so that the compiler tells you this.

Romain


Le 10/12/12 16:49, Ramon Diaz-Uriarte a ?crit :

Dear All,

I am trying to pass columns from an Armadillo matrix to a function, but
I'd like to pass just a reference to the column, not a copy of the column
and I do not seem to be able to do it "elegantly".

The code below (function f1) I think shows that passing X.col to a
function creates a copy (X.unsafe_col does too). We can pass &X as
argument, and the index of the column, as in f2. And that will not create
a copy. But I think this is not the right way of doing what I want to do
(to begin with, I'd rather not pass the column index to the function).


What am I getting wrong?

Best,





// [[Rcpp::depends(RcppArmadillo)]]
#include <Rcpp.h>
#include <RcppArmadillo.h>

using namespace Rcpp;



double f1(arma::umat Z) {
   Z(0, 0) = 111;
   Z(9, 0) = 111;
   std::cout << "f1, this is Z " << std::endl << Z << std::endl;
   return 33.3;
}

double f2(arma::umat &Z, const int c1) {
   Z(1, 0) = 222;
   return 66.6;
}


double f3(arma::umat &Z) {
   Z(1, 0) = 223;
   return 99.9;
}


// [[Rcpp::export]]
List f0(IntegerVector s1_, IntegerVector c1_){
   const int  s1 = as<int>(s1_);
   const int  c1 = as<int>(c1_);

   arma::umat X(10, s1);
   for(int j = 0; j < s1; ++j) {
     for(int i = 0; i < 10; ++i) {
       X(i, j) = i * 10 + j;
     }
   }

   // in both cases, a copy seems to be made
   //double fitness = f1(X.col(c1));
   double outf1 = f1(X.unsafe_col(c1));
   std::cout << "f0, this is X after f1" << std::endl << X << std::endl;

   double outf2 = f2(X, c1);
   std::cout << "f0, this is X after f2" << std::endl << X << std::endl;


   // double outf3 = f3(X.unsafe_col(c1)); //will not work
   //double outf3 = f3(X.col(c1)); //will not work


   return List::create(Named("X") = wrap(X),
              Named("of1") = outf1,
              Named("of2") = outf2);
}

the .col method gives you a subview_col :

   arma_inline       subview_col<eT> col(const uword col_num);
   arma_inline const subview_col<eT> col(const uword col_num) const;

not an umat.

The name of the class implies that it is a "view" class, so a way to 
look a data from another class. hence, no data of its own, so cheap copy.

double f3(arma::subview_col<unsigned int> Z) {
   Z(1, 0) = 223;
   return 99.9;
}

Also, including RcppArmadillo.h after Rcpp.h is wrong. You should only 
include RcppArmadillo.h.

I should do something so that the compiler tells you this.

Romain

Le 10/12/12 16:49, Ramon Diaz-Uriarte a ?crit :

Dear All,

I am trying to pass columns from an Armadillo matrix to a function, but
I'd like to pass just a reference to the column, not a copy of the column
and I do not seem to be able to do it "elegantly".

The code below (function f1) I think shows that passing X.col to a
function creates a copy (X.unsafe_col does too). We can pass &X as
argument, and the index of the column, as in f2. And that will not create
a copy. But I think this is not the right way of doing what I want to do
(to begin with, I'd rather not pass the column index to the function).


What am I getting wrong?

Best,





// [[Rcpp::depends(RcppArmadillo)]]
#include <Rcpp.h>
#include <RcppArmadillo.h>

using namespace Rcpp;



double f1(arma::umat Z) {
   Z(0, 0) = 111;
   Z(9, 0) = 111;
   std::cout << "f1, this is Z " << std::endl << Z << std::endl;
   return 33.3;
}

double f2(arma::umat &Z, const int c1) {
   Z(1, 0) = 222;
   return 66.6;
}


double f3(arma::umat &Z) {
   Z(1, 0) = 223;
   return 99.9;
}


// [[Rcpp::export]]
List f0(IntegerVector s1_, IntegerVector c1_){
   const int  s1 = as<int>(s1_);
   const int  c1 = as<int>(c1_);

   arma::umat X(10, s1);
   for(int j = 0; j < s1; ++j) {
     for(int i = 0; i < 10; ++i) {
       X(i, j) = i * 10 + j;
     }
   }

   // in both cases, a copy seems to be made
   //double fitness = f1(X.col(c1));
   double outf1 = f1(X.unsafe_col(c1));
   std::cout << "f0, this is X after f1" << std::endl << X << std::endl;

   double outf2 = f2(X, c1);
   std::cout << "f0, this is X after f2" << std::endl << X << std::endl;


   // double outf3 = f3(X.unsafe_col(c1)); //will not work
   //double outf3 = f3(X.col(c1)); //will not work


   return List::create(Named("X") = wrap(X),
		      Named("of1") = outf1,
		      Named("of2") = outf2);
}

-- 
Romain Francois
Professional R Enthusiast
+33(0) 6 28 91 30 30

R Graph Gallery: http://gallery.r-enthusiasts.com

blog:            http://romainfrancois.blog.free.fr
|- http://bit.ly/RE6sYH : OOP with Rcpp modules
`- http://bit.ly/Thw7IK : Rcpp modules more flexible