[Rcpp-devel] Mersenne Twister in RcppArmadillo?
Simon Zehnder
szehnder at uni-bonn.de
Fri Feb 15 08:56:38 CET 2013
Hi Dirk,
that is good to know. At this point I would like to know: Functions like R::qbinom or R::qnorm from R's C library can still be used right? They are not called, but the code is read and worked through in compilation?
Best
Simon
P.S.: How is it with code using RcppArmadillo: Is there a general interest to share reusable code (e.g. I programmed for my GMM estimation a Newey West covariance estimation procedure, reproducing results from R's package sandwich with prewhitening, bandwith bwNeweyWest and weights weightsAndrews. Also sandwich estimation is possible. It becomes much faster than the R equivalent when working with datasets bigger than 50.000 observations)?
On Feb 15, 2013, at 12:08 AM, Dirk Eddelbuettel <edd at debian.org> wrote:
>
> On 14 February 2013 at 22:24, Yan Zhou wrote:
> | It runs in parallel does not mean it works correctly.
>
> Well put :)
>
> | Simply increase seed certainly not avoid overlapping. The best algorithm that
> | use MT19937 in parallel environment without overlapping has a complexity O(log
> | n) AFAIK.
> |
> | Second, R's RNG has a single instance within the whole program, which means
> | your program does not work as you thought. There are obvious race conditions in
> | your approach.
>
> Precisely -- First rule of OpenMP with R: never ever call back into R from a
> thread. Which implies to not access Rcpp::Function etc.
>
> | I suggest you google a little bit of how RNG really works before purse how to
> | generate them in parallel. For example section 6 of the R's parallel package's
> | document.
>
> I meant to make the same suggestion about the same carefully written
> document. Prof Ripley also has an (older) bug on this, he _really_ knows this
> stuff.
>
> Regards, Dirk
>
> | Best,
> |
> | Yan Zhou
> |
> | On Feb 14, 2013, at 09:40 PM, Simon Zehnder <szehnder at uni-bonn.de> wrote:
> |
> |
> | Well, in OpenMP it seems to work, if I use inside the #pragma parallel for:
> |
> | Rcpp::Environment base("package:base");
> | Rcpp::Function SetSeed = base["set.seed"];
> | SetSeed(1 + i * nobs_intern * 3);
> |
> | I used OMP_NUMTHREADS = 4 and 8 iterations.
> |
> | Best Simon
> |
> | On Feb 14, 2013, at 10:10 PM, Simon Zehnder <szehnder at uni-bonn.de> wrote:
> |
> | > I made all my simple tests now by using Rcpp::Environment and then
> | Rcpp::Function for "set.seed" when calling R::rnorm. The next step would be
> | to parallelize the iterations via OpenMP. The suggestion of Yan makes quite
> | sense to produce a RNG via std::mt19937. This is also the way I did it when
> | I used Scythe Statistical Library.
> | > As I use now the R RNG it should be able to use it in parallel, if
> | Rcpp::Function is a kind of wrapper which contains the code of "set.seed".
> | But I would rather guess, that it is a call object, that just calls the R
> | function "set.seed" in R, which then would not be threadsafe (as also the
> | RNG would then just be called when using R::rnorm).....
> | >
> | > So I guess falling back to std::mt19937 is not a bad idea. For no
> | overlapping I just increase the seed in every iteration dependent on the
> | iteration and the number of random numbers to generate.
> | >
> | > Best Simon
> | >
> | >
> | >
> | >
> | > On Feb 11, 2013, at 2:16 PM, Chris Jefferson <chris at bubblescope.net>
> | wrote:
> | >
> | >> On 11/02/13 10:23, c s wrote:
> | >>> On Sun, Feb 10, 2013 at 12:32 AM, Yan Zhou <zhouyan at me.com> wrote:
> | >>>> To have Armadillo randn use MT19937 is not easy.
> | >>>> Since it use srand() for seed, I guess it also use C rand(),
> | >>>> whose quality is quite questionable.
> | >>> The quality of the rand() function from C depends on the
> | >>> implementation in libc, which varies from system to system.
> | >>>
> | >>> While I'd like to incorporate a Mersenne-Twister random number
> | >>> generator directly in Armadillo, it would either add a dependency on
> | >>> Boost, or on C++11. Boost might not be available on a user's system,
> | >>> and the degree of support for C++11 features varies from compiler to
> | >>> compiler. We also have to bear in mind that R folks currently
> | >>> disallow CRAN packages that use C++11.
> | >>
> | >> If you just want a mersenne-twister random number generator, I will
> | extract it from boost (which is under a fairly free licence).
> | >>
> | >> Chris
> | >> _______________________________________________
> | >> Rcpp-devel mailing list
> | >> Rcpp-devel at lists.r-forge.r-project.org
> | >> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
> | >
> | > _______________________________________________
> | > Rcpp-devel mailing list
> | > Rcpp-devel at lists.r-forge.r-project.org
> | > https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
> |
> | _______________________________________________
> | Rcpp-devel mailing list
> | Rcpp-devel at lists.r-forge.r-project.org
> | https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
> |
> |
> | ----------------------------------------------------------------------
> | _______________________________________________
> | Rcpp-devel mailing list
> | Rcpp-devel at lists.r-forge.r-project.org
> | https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
> --
> Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com
More information about the Rcpp-devel
mailing list