播种多线程unif_rand()

内容来源于 Stack Overflow,并遵循CC BY-SA 3.0许可协议进行翻译与使用

  • 回答 (2)
  • 关注 (0)
  • 查看 (86)

unif_rand()在多线程环境中播放R的内部。下面的代码生成2个线程内的2列均匀随机数矩阵。结果很有趣。

struct mtRunif: public RcppParallel::Worker
{
  int Nrow; // number of rows in matrix.
  double *v; // point to the 0th element of the 0th column.
  void operator() (std::size_t st, std::size_t end)
  {
    // st = 0 in the 0th thread, 1 in the 1st thread. 
    double *vst = v + st * Nrow;
    for(int i = 0; i < Nrow; ++i)
    {
      vst[i] = unif_rand();
    }
  }


  mtRunif(int Nrow, double *v): Nrow(Nrow), v(v)
  {
    RcppParallel::parallelFor(0, 2, *this);
  }
};


// [[Rcpp::export]] 
NumericMatrix testSeeding(int sampleSize)
{
  NumericMatrix rst(sampleSize, 2);
  mtRunif(sampleSize, &*rst.begin());
  return rst;
}


/***R
N = 100
set.seed(42); tmp = testSeeding(N) 
set.seed(42); tmp2 = testSeeding(N) 
# see if sequences are identical
range(tmp[, 1] - tmp2[, 1]); range(tmp[, 2] - tmp2[, 2])
# [1] 0 0
# [1] 0 0


N = 1000
set.seed(42); tmp = testSeeding(N) 
set.seed(42); tmp2 = testSeeding(N) 
range(tmp[, 1] - tmp2[, 1]); range(tmp[, 2] - tmp2[, 2])
# [1] -0.9655154  0.8989870
# [1] -0.969356  0.963239
*/

结果表明set.seed()控制小样本大小的所有线程的随机性? 最初我期望set.seed()在不超过1个线程中有效。 我不想利用这个结论,因为它可能是绝对错误的。 另一方面,unif_rand()的种子函数是否类似于std :: rand()的std :: srand()?

谢谢!

提问于
用户回答回答于

简而言之:由于R内部原因,你不能用R来做这件事,并且已经被广泛记录。

还有关于RNG和流的统计问题。因此,您最有可能想要研究适合从多个线程绘制的“流式RNG”。CRAN上有:

以及不再在CRAN上的旧版sprng

用户回答回答于

在评论中宣传 dqrng后,我意识到我没有编写任何关于如何使用该软件包中的RNG进行并行使用的文档。所以我开始了如下操作,它将成为下一个版本的一部分。这是其中一个示例,与您尝试执行的操作非常相似:

#include <Rcpp.h>
// [[Rcpp::depends(dqrng)]]
#include <pcg_random.hpp>
#include <dqrng_distribution.h>
// [[Rcpp::depends(RcppParallel)]]
#include <RcppParallel.h>
// [[Rcpp::plugins(cpp11)]]

struct RandomFill : public RcppParallel::Worker {
  RcppParallel::RMatrix<double> output;
  uint64_t seed;
  dqrng::normal_distribution dist{0.0, 1.0};

  RandomFill(Rcpp::NumericMatrix output, const uint64_t seed) : output(output), seed(seed) {};

  void operator()(std::size_t begin, std::size_t end) {
    pcg64 rng(seed, end); // ctor with seed and stream id
    auto gen = std::bind(dist, rng);
    std::generate(output.begin() + begin * output.nrow(),
                  output.begin() + end * output.nrow(),
                  std::ref(gen));
  }
};

// [[Rcpp::export]]
Rcpp::NumericMatrix parallel_random_matrix(const int n, const int m, const int ncores) {
  Rcpp::NumericMatrix res(n, m);
  RandomFill randomFill(res, 42);
  RcppParallel::parallelFor(0, m, randomFill, m/ncores + 1);
  return res;
}

/*** R
res <- parallel_random_matrix(1e6, 8, 4)
head(res)
*/

结果:

> res <- parallel_random_matrix(1e6, 8, 4)

> head(res)
           [,1]        [,2]        [,3]       [,4]       [,5]       [,6]       [,7]       [,8]
[1,]  0.7114429 -0.19759808 -0.47149983  0.6046378 -0.3709571 -0.8089533  0.8185977 0.49010575
[2,]  0.8721661 -0.47654248  1.10411136 -1.6290995 -1.3276661 -0.2585322 -1.2437521 0.90325167
[3,] -1.4959624  0.61068373 -0.54343828 -0.4623555 -1.1779352 -2.8068283 -0.4341252 1.74490995
[4,]  0.5087201 -0.05175746  0.19007581 -0.7869679  0.9672267 -0.5009787 -0.5283977 1.42487290
[5,] -0.8191448 -0.77348120 -0.03458304  0.7243224  1.0594094 -0.6951184 -0.5456669 0.00894037
[6,]  1.2289518 -2.33539762  0.40222707 -2.3346460 -0.5796549 -0.3092356  2.8961294 0.16773085

顺便说一下,请不要suestd::rand()。如果你想使用标准库,那么请使用类似std::mt19937random与C ++ 11。

所属标签

可能回答问题的人

  • 天使的炫翼

    17 粉丝531 提问9 回答
  • 优惠活动秘书

    0 粉丝2 提问8 回答
  • 最爱开车啦

    8 粉丝503 提问6 回答
  • 富有想象力的人

    3 粉丝0 提问5 回答

扫码关注云+社区

领取腾讯云代金券