Rstudio аварийно завершает работу с функциями Rcpp и OpenMP

Question

Это дополнительный вопрос к dqrng с Rcpp для рисования из нормального и биномиального распределения .Я попытался реализовать ответ, но вместо того, чтобы рисовать из одного дистрибутива, я рисую из 3. Это код, который я написал:

// [[Rcpp::depends(dqrng, BH, RcppArmadillo)]]
#include <RcppArmadillo.h>
#include <boost/random/binomial_distribution.hpp>
#include <xoshiro.h>
#include <dqrng_distribution.h>
// [[Rcpp::plugins(openmp)]]
#include <omp.h>

// [[Rcpp::plugins(cpp11)]]

// [[Rcpp::export]]
arma::mat parallel_random_matrix(int n, int m, int ncores, double p=0.5) {
  dqrng::xoshiro256plus rng(42);
  arma::mat out(n*m,3);
  // ok to use rng here
#pragma omp parallel num_threads(ncores)
{
  dqrng::xoshiro256plus lrng(rng);      // make thread local copy of rng 
  lrng.jump(omp_get_thread_num() + 1);  // advance rng by 1 ... ncores jumps 
  int iter = 0;
#pragma omp for
  for (int i = 0; i < m; ++i) {
    for (int j = 0; j < n; ++j) {
      iter = i * n + j;
      // p can be a function of i and j
      boost::random::binomial_distribution<int> dist_binomial(1,p);
      auto gen_bernoulli = std::bind(dist_binomial, std::ref(lrng));
      boost::random::normal_distribution<int> dist_normal1(2.0,1.0);
      auto gen_normal1 = std::bind(dist_normal1, std::ref(lrng));
      boost::random::normal_distribution<int> dist_normal2(4.0,3.0);
      auto gen_normal2 = std::bind(dist_normal2, std::ref(lrng));
      out(iter,0) = gen_bernoulli();
      out(iter,1) = gen_normal1();
      out(iter,2) = gen_normal2();
    }
  }
}
// ok to use rng here
return out;
}

/*** R
parallel_random_matrix(5, 5, 4, 0.75)
*/

Когда я пытаюсь его запустить, Rstudio падает.Однако, когда я изменяю код следующим образом, он работает:

// [[Rcpp::depends(dqrng, BH, RcppArmadillo)]]
#include <RcppArmadillo.h>
#include <boost/random/binomial_distribution.hpp>
#include <xoshiro.h>
#include <dqrng_distribution.h>
// [[Rcpp::plugins(openmp)]]
#include <omp.h>

// [[Rcpp::plugins(cpp11)]]

// [[Rcpp::export]]
arma::mat parallel_random_matrix(int n, int m, int ncores, double p=0.5) {
  dqrng::xoshiro256plus rng(42);
  arma::mat out(n*m,3);
  // ok to use rng here
#pragma omp parallel num_threads(ncores)
{
  dqrng::xoshiro256plus lrng(rng);      // make thread local copy of rng 
  lrng.jump(omp_get_thread_num() + 1);  // advance rng by 1 ... ncores jumps 
  int iter = 0;
#pragma omp for
  for (int i = 0; i < m; ++i) {
    for (int j = 0; j < n; ++j) {
      iter = i * n + j;
      // p can be a function of i and j
      boost::random::binomial_distribution<int> dist_binomial(1,p);
      auto gen_bernoulli = std::bind(dist_binomial, std::ref(lrng));
      boost::random::normal_distribution<int> dist_normal1(2.0,1.0);
      auto gen_normal1 = std::bind(dist_normal1, std::ref(lrng));
      boost::random::normal_distribution<int> dist_normal2(4.0,3.0);
      auto gen_normal2 = std::bind(dist_normal2, std::ref(lrng));
      out(iter,0) = gen_bernoulli();
      out(iter,1) = 2.0;//gen_normal1();
      out(iter,2) = 3.0;//gen_normal2();
    }
  }
}
// ok to use rng here
return out;
}

/*** R
parallel_random_matrix(5, 5, 4, 0.75)
*/

Что я делаю не так?

Answer 1

Здесь кроется проблема:

  boost::random::normal_distribution<int> dist_normal1(2.0,1.0);
                                     ^^^ 

Это распределение предназначено для реальных типов, а не целых типов, ср https://www.boost.org/doc/libs/1_69_0/doc/html/boost/random/normal_distribution.html. Правильно будет

  boost::random::normal_distribution<double> dist_normal1(2.0,1.0);