I have written parallel implementation of sums in groups using RcppParallel.
// [[Rcpp::depends(RcppParallel)]]#include <Rcpp.h>#include <RcppParallel.h> using namespace Rcpp; using namespace RcppParallel; struct SumsG: public Worker{ const RVector<double> v; const RVector<int> gi; RVector<double> sg; SumsG(const NumericVector v, const IntegerVector gi, NumericVector sg): v(v), gi(gi), sg(sg) {} SumsG(const SumsG& p, Split): v(p.v), gi(p.gi), sg(p.sg) {} void operator()(std::size_t begin, std::size_t end) { for (std::size_t i = begin; i < end; i++) { sg[gi[i]] += v[i]; } } void join(const SumsG& p) { for(std::size_t i = 0; i < sg.length(); i++) { sg[i] += p.sg[i]; } }}; // [[Rcpp::export]] List sumsingroups(NumericVector v, IntegerVector gi, int ni) { NumericVector sg(ni); SumsG p(v, gi, sg); parallelReduce(0, v.length(), p); return List::create(_["sg"] = p.sg);} It compiles using Rcpp::sourceCpp. Now when I call it from R sumsingroups(1:10, rep(0:1, each = 5), 2) several times I get the right answer (15 40) and then something different (usually some multiplicative of the right answer). Running res <- sumsingroups(1:10, rep(0:1, each = 5), 2)for(i in 1:1000) { tmp <- sumsingroups(1:10, rep(0:1, each = 5), 2) if(res[[1]][1] != tmp[[1]][1]) break Sys.sleep(0.1)} breaks at random iteration returning $sg[1] 60 160 or $sg[1] 30 80 I am new to Rcpp and RcppParallel and do not know what could cause such behavior. Things that did not help: 1. Added for (std::size_t i = 0; i < sg.length(); i++) sg[i] = 0; to both of constructors. 2. Changed names so that they are different in Worker definition and in function implementation.
_______________________________________________ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel