Skip to content

[Rcpp-devel] fill a NumericMatrix row by row

5 messages · Conrad S, Romain Francois, Michael Love +1 more

#
Not a question here, but a note.  Earlier there was an email about filling
a NumericMatrix row by row.

http://lists.r-forge.r-project.org/pipermail/rcpp-devel/2013-March/005390.html

I noticed that for Rcpp 0.10.2 it was possible to fill a NumericMatrix with
NumericVectors, but since Rcpp 0.10.3 it does so while printing out some
lines.

As recommended, I have switched to using Armadillo matrices, which doesn't
print out these lines.

-Mike

example with NumericMatrix, fill a matrix mat0 by row with the vector vec0:

library(Rcpp)
library(inline)
code <- '
   Rcpp::NumericMatrix mat(mat0);
   Rcpp::NumericVector vec(vec0);
   for (int i = 0; i < mat.nrow(); i++) {
      mat(i,_) = vec;
   }
   return Rcpp::wrap(mat);
'
f <-
cxxfunction(signature(mat0="numeric",vec0="numeric"),code,plugin="Rcpp")
mat0 <- matrix(c(1:4),nrow=2)
vec0 <- c(5,6)
f(mat0,vec0)

produces:

MatrixRow::get_parent_index(int = 0), parent_nrow = 2 >> 0
MatrixRow::get_parent_index(int = 1), parent_nrow = 2 >> 2
MatrixRow::get_parent_index(int = 0), parent_nrow = 2 >> 0
MatrixRow::get_parent_index(int = 1), parent_nrow = 2 >> 2
     [,1] [,2]
[1,]    5    6
[2,]    5    6


Example with Armadillo:

arma.code <- '
   arma::mat mat = Rcpp::as<arma::mat>(mat0);
   arma::rowvec vec = Rcpp::as<arma::rowvec>(vec0);
   for (int i = 0; i < mat.n_rows; i++) {
      mat.row(i) = vec;
   }
   return Rcpp::wrap(mat);
'
arma.f <-
cxxfunction(signature(mat0="numeric",vec0="numeric"),arma.code,plugin="RcppArmadillo")
arma.f(mat0,vec0)

produces:

     [,1] [,2]
[1,]    5    6
[2,]    5    6
R Under development (unstable) (2013-03-20 r62335)
Platform: x86_64-unknown-linux-gnu (64-bit)

locale:
[1] C

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] RcppArmadillo_0.3.800.1 inline_0.3.11           Rcpp_0.10.3
[4] Defaults_1.1-1          BiocInstaller_1.9.8

loaded via a namespace (and not attached):
[1] tools_3.1.0
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20130327/158726bb/attachment.html>
#
There's an even simpler way when using Armadillo, so that no loop is required:
http://arma.sourceforge.net/docs.html#each_colrow

For example:

mat X(4,5);
rowvec r(5);

X.each_row() = r;


(btw, in general it's more efficient to access matrices as columns
rather than rows).


On Wed, Mar 27, 2013 at 6:23 PM, Michael Love
<michaelisaiahlove at gmail.com> wrote:
#
I removed it in the revision 4294. Sorry for the inconvenience.

Romain

Le 27/03/13 09:23, Michael Love a ?crit :

  
    
#
Thanks for the reply and great libraries Conrad and Romain.

Conrad, my example is a bit too simple maybe, I am really doing:

for i from 1 to m:
   do computation
   store a row of results into a tall m x n matrix

...where i am pretty sure that the time required for computation is more
than storing the results.

Out of curiosity, I tried to do some profiling for row filling of tall
Armadillo matrices vs. column filling. Row filling becomes worse when the
rows are long.  I guess the memory allocation could obscure some speed
differences:

library(Rcpp)
library(inline)

# fill a matrix by row with m0 rows and n0 columns
row.code <- '
   int m = Rcpp::as<int>(m0);
   int n = Rcpp::as<int>(n0);
   arma::mat mat = arma::zeros(m,n);
   arma::rowvec vec = Rcpp::as<arma::rowvec>(vec0);
   for (int i = 0; i < m; i++) {
      mat.row(i) = vec;
   }
   return Rcpp::wrap(mat);
'
row.f <-
cxxfunction(signature(m0="integer",n0="integer",vec0="numeric"),row.code,plugin="RcppArmadillo")

# so here we swap, fill by column, m0 is the number of columns and n0 is
the number of rows
col.code <- '
   int m = Rcpp::as<int>(m0);
   int n = Rcpp::as<int>(n0);
   arma::mat mat = arma::zeros(n,m);
   arma::colvec vec = Rcpp::as<arma::colvec>(vec0);
   for (int i = 0; i < m; i++) {
      mat.col(i) = vec;
   }
   return Rcpp::wrap(mat);
'
col.f <- cxxfunction(signature(m0="
integer",n0="integer",vec0="numeric"),col.code,plugin="RcppArmadillo")
user  system elapsed
  0.779   1.022   1.801
user  system elapsed
  0.957   0.966   1.922
[1] TRUE
user  system elapsed
  1.925   1.425   3.351
user  system elapsed
  1.244   1.423   2.667
[1] TRUE
On Wed, Mar 27, 2013 at 9:44 AM, Conrad S <conradsand.arma at gmail.com> wrote:

            
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/rcpp-devel/attachments/20130327/fb13618d/attachment.html>
#
On 27 March 2013 at 12:06, Michael Love wrote:
| Thanks for the reply and great libraries Conrad and Romain.
| 
| Conrad, my example is a bit too simple maybe, I am really doing:
| 
| for i from 1 to m:
| ?? do computation
| ?? store a row of results into a tall m x n matrix
| 
| ...where i am pretty sure that the time required for computation is more than
| storing the results.
| 
| Out of curiosity, I tried to do some profiling for row filling of tall
| Armadillo matrices vs. column filling. Row filling becomes worse when the rows
| are long.? I guess the memory allocation could obscure some speed differences:

You did read the line where we explained that storage is by _columns_ so
insertion by row _has to_ result in copies?

This too has been discussed before.  If performance matters, grow columns and
transose at the end.

Dirk

| 
| library(Rcpp)
| library(inline)
| 
| # fill a matrix by row with m0 rows and n0 columns
| row.code <- '
| ?? int m = Rcpp::as<int>(m0);
| ?? int n = Rcpp::as<int>(n0);
| ?? arma::mat mat = arma::zeros(m,n);
| ?? arma::rowvec vec = Rcpp::as<arma::rowvec>(vec0);
| ?? for (int i = 0; i < m; i++) {
| ????? mat.row(i) = vec;
| ?? }
| ?? return Rcpp::wrap(mat);
| '
| row.f <- cxxfunction(signature(m0="integer",n0="integer",vec0=
| "numeric"),row.code,plugin="RcppArmadillo")
| 
| # so here we swap, fill by column, m0 is the number of columns and n0 is the
| number of rows
| col.code <- '
| ?? int m = Rcpp::as<int>(m0);
| ?? int n = Rcpp::as<int>(n0);
| ?? arma::mat mat = arma::zeros(n,m);
| ?? arma::colvec vec = Rcpp::as<arma::colvec>(vec0);
| ?? for (int i = 0; i < m; i++) {
| ????? mat.col(i) = vec;
| ?? }
| ?? return Rcpp::wrap(mat);
| '
| col.f <- cxxfunction(signature(m0="
| integer",n0="integer",vec0="numeric"),col.code,plugin="RcppArmadillo")
| 
| > m <- 1e7
| > n <- 10
| > vec0 <- seq_len(n)
| > system.time({row.f(m,n,vec0)})
| ?? user? system elapsed
| ? 0.779?? 1.022?? 1.801
| > system.time({col.f(m,n,vec0)})
| ?? user? system elapsed
| ? 0.957?? 0.966?? 1.922
| > all.equal(row.f(m,n,vec0), t(col.f(m,n,vec0)))
| [1] TRUE
| 
| > m <- 1e6
| > n <- 100
| > vec0 <- seq_len(n)
| > system.time({row.f(m,n,vec0)})
| ?? user? system elapsed
| ? 1.925?? 1.425?? 3.351
| > system.time({col.f(m,n,vec0)})
| ?? user? system elapsed
| ? 1.244?? 1.423?? 2.667
| > all.equal(row.f(m,n,vec0), t(col.f(m,n,vec0)))
| [1] TRUE
| 
|
| On Wed, Mar 27, 2013 at 9:44 AM, Conrad S <conradsand.arma at gmail.com> wrote:
| 
|     There's an even simpler way when using Armadillo, so that no loop is
|     required:
|     http://arma.sourceforge.net/docs.html#each_colrow
| 
|     For example:
| 
|     mat X(4,5);
|     rowvec r(5);
| 
|     X.each_row() = r;
| 
| 
|     (btw, in general it's more efficient to access matrices as columns
|     rather than rows).
| 
| 
|     On Wed, Mar 27, 2013 at 6:23 PM, Michael Love
| <michaelisaiahlove at gmail.com> wrote:
|     > As recommended, I have switched to using Armadillo matrices, which
|     doesn't
|     > print out these lines.
|     > ...
|     > Example with Armadillo:
|     >
|     > arma.code <- '
|     > ? ?arma::mat mat = Rcpp::as<arma::mat>(mat0);
|     > ? ?arma::rowvec vec = Rcpp::as<arma::rowvec>(vec0);
|     > ? ?for (int i = 0; i < mat.n_rows; i++) {
|     > ? ? ? mat.row(i) = vec;
|     > ? ?}
|     > ? ?return Rcpp::wrap(mat);
|     > '
|     > arma.f <-
|     > cxxfunction(signature(mat0="numeric",vec0="numeric"),arma.code,plugin=
|     "RcppArmadillo")
|     > arma.f(mat0,vec0)
| 
| 
| 
| ----------------------------------------------------------------------
| _______________________________________________
| 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