Not looking for advise but sharing my own experience here.
In a project using optimisation I was looking into rewriting the target function in C++ to improve performance. I have absolutely no experience in C++ but as the function was very simple I thought I'd give it a go. After some dabbeling I achieved my goal and the function was evaluating correctly and even better: It was much faster than the R counterpart. Hurray!!
However, when running the function multiple times, e.g. in microbenchmark()
or in optim()
, R would crash consistently. I couldn't for the life of me figure out what was causing the error. I let the project rest for months, because of that. At one point, with fresh motivation, I thought I'd give it a go. Again, I could not find the error but narrowed it down to a potential memory leak, something we don't have to deal with usually in R.
Reading through the C++ chapter in Advanced R by @hadley again I stumbled upon the following paragraph that gave me pause:
It could not be that easy, right? Did I ignore this simple fact, even though it was written twice and in bold? But sure enough, looking in my function I found a for loop that broke that exact rule. Fixing the error was then easy and the function worked like a charm.
What I take away from this is that
- C++ can result in different types of errors than R. (Functions can evaluate correctly and without an error message but cause problems in special cases).
- The easiest errors are sometimes the hardest to find.
- Documentation on Rcpp is still scarce for beginners (hence this post).
I've attached the devious function and the fix down below. The difference is in the for loop, where the end-point should be one less than the size (length
) of the the vector to be summed up.
Best,
Valentin
library(Rcpp)
# function that will crash R
crash_fun <-
"NumericVector cumsum_crashing(NumericVector x){
double sum = 0;
int n = x.size();
NumericVector out(n);
for (int i=0; i<= n; i++){
sum += x[i];
out[i] = sum;
}
return out;
}"
# function that wont crash R
nocrash_fun <-
"NumericVector cumsum_working(NumericVector x){
double sum = 0;
int n = x.size();
NumericVector out(n);
for (int i=0; i< n; i++){
sum += x[i];
out[i] = sum;
}
return out;
}"
Rcpp::cppFunction(crash_fun)
Rcpp::cppFunction(nocrash_fun)
Created on 2022-01-03 by the reprex package (v2.0.1)