Генерация случайного потока битов в Rcpp эффективно - PullRequest
1 голос
/ 28 марта 2019

У меня есть вспомогательная функция в пакете R, которую я сейчас создаю, с именем rbinom01. Обратите внимание, что он вызывает random(3).

int rbinom01(int size) {
  if (!size) {
    return 0;
  }

  int64_t result = 0;
  while (size >= 32) {
    result += __builtin_popcount(random());
    size -= 32;
  }

  result += __builtin_popcount(random() & ~(LONG_MAX << size));

  return result;
}

Когда R CMD check my_package, я получил следующее предупреждение:

* checking compiled code ... NOTE
File ‘ my_package/libs/my_package.so’:
  Found ‘_random’, possibly from ‘random’ (C)
    Object: ‘ my_function.o’

Compiled code should not call entry points which might terminate R nor
write to stdout/stderr instead of to the console, nor use Fortran I/O
nor system RNGs.

See ‘Writing portable packages’ in the ‘Writing R Extensions’ manual.

Я направился к Документу , и там говорится, что я могу использовать одну из функций *_rand вместе с семейством функций распределения . Ну, это круто, но моему пакету просто нужен поток случайных битов, а не случайный double. Самый простой способ сделать это - использовать random(3) или, возможно, прочитать из /dev/urandom, но это делает мой пакет «непереносимым».

В этом посте предлагается использовать sample, но, к сожалению, он не подходит для моего случая использования. Для моего приложения генерация случайных битов, по-видимому, критически важна для производительности, поэтому я не хочу, чтобы это тратило время на вызов unif_rand, умножьте результат на N и округлите его. В любом случае, я использую C ++ для использования параллелизма на уровне битов.

Конечно, я могу вручную свернуть свой собственный PRNG или скопировать и вставить код самого современного PRNG, например xoshiro256 **, но перед этим я хотел бы посмотреть, есть ли Есть более простые альтернативы.

Кстати, может кто-нибудь связать со мной хороший краткий учебник по Rcpp? Написание расширений R является всеобъемлющим и потрясающим, но мне потребовались бы недели, чтобы закончить. Я ищу более краткую версию, но желательно, чтобы она была более информативной, чем вызов Rcpp.package.skeleton.


Как следует из ответа @ Ralf Stubner , я переписал исходный код следующим образом. Тем не менее, я получаю один и тот же результат каждый раз. Как правильно его заполнить и в то же время сохранить мой код «переносимым»?

int rbinom01(int size) {
  dqrng::xoshiro256plus rng;

  if (!size) {
    return 0;
  }

  int result = 0;
  while (size >= 64) {
    result += __builtin_popcountll(rng());
    Rcout << sizeof(rng()) << std::endl;
    size -= 64;
  }

  result += __builtin_popcountll(rng() & ((1LLU << size) - 1));

  return result;
}

1 Ответ

4 голосов
/ 28 марта 2019

Существуют различные R-пакеты, которые делают PRNG доступными в виде библиотек только для заголовков C ++:

  • BH : все от boost.random
  • sitmo: различные версии Threefry
  • dqrng : семейство PCG, xoshiro256 + и xoroshiro128 +
  • ...

Вы можете использоватьлюбой из них, добавив LinkingTo к вашему пакету DECRIPTION.Обычно эти PRNG моделируются после заголовка C ++ 11 random, что означает, что вы должны контролировать их жизненный цикл и заполнять себя.В однопоточной среде мне нравится использовать анонимные пространства имен для управления жизненным циклом, например:

#include <Rcpp.h>
// [[Rcpp::depends(dqrng)]]
#include <xoshiro.h>
// [[Rcpp::plugins(cpp11)]]

namespace {
dqrng::xoshiro256plus rng{};
}

// [[Rcpp::export]]
void set_seed(int seed) {
  rng.seed(seed);
}

// [[Rcpp::export]]
int rbinom01(int size) {
  if (!size) {
    return 0;
  }

  int result = 0;
  while (size >= 64) {
    result += __builtin_popcountll(rng());
    size -= 64;
  }

  result += __builtin_popcountll(rng() & ((1LLU << size) - 1));

  return result;
}

/*** R
set_seed(42)
rbinom01(10)
rbinom01(10)
rbinom01(10)
*/

Однако использование runif не так уж плохо и, конечно, быстрее, чем доступ к /dev/urandomdqrng есть удобная обёртка для этого.

Что касается учебников: Помимо WRE виньетка Rcpp является обязательной для прочтения. R Packages от Hadley Wickham также есть глава о «скомпилированном коде», если вы хотите пойти по devtools.

...