多次调用相同的 Rcpp 函数时返回不同的结果答案

【问题标题】：When calling same Rcpp function several times different results are returned多次调用相同的 Rcpp 函数时返回不同的结果
【发布时间】：2015-07-09 13:06:29
【问题描述】：

我已经使用 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);
}

它使用Rcpp::sourceCpp 编译。现在，当我从 R sumsingroups(1:10, rep(0:1, each = 5), 2) 多次调用它时，我得到了正确答案（15 40），然后得到了不同的答案（通常是正确答案的乘法）。跑步

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)
}

在随机迭代返回时中断

$sg
[1]  60 160

或

$sg
[1] 30 80

我是 Rcpp 和 RcppParallel 的新手，不知道什么会导致这种行为。

更新。没有帮助的事情：

已将for (std::size_t i = 0; i < sg.length(); i++) sg[i] = 0; 添加到两个构造函数。
更改了名称，使它们在 Worker 定义并在函数实现中。

【问题讨论】：

我怀疑你的线程在循环的情况下是不同步的。必须有一个命令来确保每个线程在收集结果之前已经完成，但是我从来没有使用过RcppParallel，所以我不知道这是怎么做到的。
Rcpp 有一个parallelFor() 函数来处理这种情况，它负责内部线程锁定。
我不明白在这种情况下如何使用parallelFor()。
我的意思是在你的 Rcpp 代码中。你有两个 for 循环。用parallelFor() 函数替换它们可能会更好，至少我是这样解释documentation 的。

标签： r rcpp

【解决方案1】：

试试这个。

#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::depends(RcppParallel)]]
#include <RcppParallel.h>

using namespace RcppParallel;
struct SumsInGroups5: public Worker
{
    const RVector<double> v;
    const RVector<int> g;

    std::vector<double> s;

    SumsInGroups5(const NumericVector v, const IntegerVector g): v(v), g(g),  s(*std::max_element(g.begin(), g.end()) + 1, 0.0){ }

    SumsInGroups5(const SumsInGroups5& p, Split): v(p.v), g(p.g), s(*std::max_element(g.begin(), g.end()) + 1, 0.0) {}

    void operator()(std::size_t begin, std::size_t end) {
        for (std::size_t i = begin; i < end; ++i) {
            s[g[i]]+=v[i];
        }

    }

    void join(const SumsInGroups5& rhs) {
        for(std::size_t i = 0; i < s.size(); i++) {
            s[i] += rhs.s[i];
        }
    }
};

// [[Rcpp::export]]
NumericVector sg5(NumericVector v, IntegerVector g) {
    SumsInGroups5 p(v, g);
    parallelReduce(0, v.length(), p);
    return wrap(p.s);
}

/*** R
a <- 1:10
g <- c(rep(0,5),rep(1,5))


bb <- lapply(1:10000,function(x)sg5(a,g))
cc<-do.call("rbind",bb)
unique(cc)
*/

与我的其他尝试相比，此代码在与其他代码相同的情况下没有产生奇怪的结果。不是很确定。

【讨论】：