Dynamically add rows in Rcpp

Question

I am trying to write a function in RcppArmadillo that dynamically appends rows to an array/matrix. It should work like rbind in R or pandas.concat in Python. (I am relying on C++ for efficiency.)

My specific objective is to take in a vector called foo and to produce a three-column matrix my_matrix, each row of which is determined by some condition. Because the condition needs to be checked for each triplet {i,j,k}, it involves a triple loop. This is what I have so far (words in BLOCK LETTERS are comments I include here):

/* (From my RcppArmadillo script) */
arma::mat myFunction(arma::vec foo) {
  int n = foo.size();
  // initialize first row of column names
  arma::vec my_matrix[] = {"i", "j", "k"}; 
  // loop and append rows
  for(int i = 0; i < n; i++) {
    for(int j = 0; j < n; j++) {
      for(int k = 0; k < n; k++) { 
        if (SOME CONDITION ABOUT i,j,k and foo) {
          APPEND ROW {i,j,k} TO my_matrix 
          arma::vec new_row = {i,j,k};
          my_matrix = join_vert(my_matrix, new_row);
        }
      }
    }
  }
  return my_matrix;
}

I'm facing three issues:

1. On the line arma::vec new_row = {i,j,k};, I am told "non-constant-expression cannot be narrowed from type 'int' to 'double' in initializer list"
2. On the line my_matrix = join_vert(my_matrix, new_row);, I am told "no matching function for call to 'join_vert'"
3. On the line return my_matrix;, I am told "no viable conversion from 'arma::vec [3]' to 'arma::Mat' (aka 'Mat<<>>')"

Because I'm not familiar with C++ (especially with issues 2 and 3 that involve iterative modification) I am stuck. Could someone here help to troubleshoot? Thanks in advance!

Answer 1

Data is (generally) natively stored as colums so adding rows is not that obvious. Because a matrix is generally represented a contiguous vector, you need full copies (to create the 'holes' you want to fill).

You would be better off to represent your growing data structure as a collection of columns and let those grow individually. Which is pretty much what a data.frame does.

Answer 2

Listen to Dirk for best practices -- the data being stored by column is a huge issue that you should be mindful of always for writing efficient code. I only write separately to show mechanically how to accomplish what you describe if for some reason you can't follow his advice.

Here's how I would (and did) make your code work:

1. Have new_row be a arma::rowvec
2. Have my_matrix be a arma::mat
3. Set the column names on the R side

Now let's see what that looks like. On the C++ side we have:

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

// [[Rcpp::export(.myFunction)]]
arma::mat myFunction(arma::vec foo) {
    int n = foo.size();
    arma::mat my_matrix;
    // loop and append rows
    for(int i = 0; i < n; i++) {
        for(int j = 0; j < n; j++) {
            for(int k = 0; k < n; k++) {
                if ( (foo[i] + foo[j] + foo[k]) > 10.0 ) {
                    arma::rowvec new_row = {i,j,k};
                    my_matrix = arma::join_vert(my_matrix, new_row);
                }
            }
        }
    }
    return my_matrix;
}

Then on the R side:

myFunction <- function(foo) {
    res <- .myFunction(foo)
    colnames(res) <- c("i", "j", "k")
    return(res)
}

Here's an example of it in action:

foo <- 1:4
myFunction(foo)

     i j k
[1,] 2 3 3
[2,] 3 2 3
[3,] 3 3 2
[4,] 3 3 3

As an aside, in the future I would work on making better reproducible examples. For example, here, you didn't include the lines

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

In your question code. In this case, that part's a minor issue, because likely anyone with the knowledge to answer your question knows you need those lines, but it's still bad practice.

More importantly, you don't include

• The condition on i, j, k, and foo;
• Example function input; or
• Desired function output

So, as you see, I just had to make those things up myself. It helps others help you to provide better reproducible examples, just a tip for the future!

---

PS: There's no way this is an efficient way to solve this problem!