You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

74 lines
1.6KB

  1. #include <random.hpp>
  2. #include <math.hpp>
  3. #include <time.h>
  4. #include <sys/time.h>
  5. namespace rack {
  6. namespace random {
  7. // xoroshiro128+
  8. // from http://xoroshiro.di.unimi.it/xoroshiro128plus.c
  9. thread_local uint64_t xoroshiro128plus_state[2];
  10. static uint64_t rotl(const uint64_t x, int k) {
  11. return (x << k) | (x >> (64 - k));
  12. }
  13. static uint64_t xoroshiro128plus_next(void) {
  14. const uint64_t s0 = xoroshiro128plus_state[0];
  15. uint64_t s1 = xoroshiro128plus_state[1];
  16. const uint64_t result = s0 + s1;
  17. s1 ^= s0;
  18. xoroshiro128plus_state[0] = rotl(s0, 55) ^ s1 ^ (s1 << 14); // a, b
  19. xoroshiro128plus_state[1] = rotl(s1, 36); // c
  20. return result;
  21. }
  22. void init() {
  23. struct timeval tv;
  24. gettimeofday(&tv, NULL);
  25. xoroshiro128plus_state[0] = tv.tv_sec;
  26. xoroshiro128plus_state[1] = tv.tv_usec;
  27. // Generate a few times to fix the fact that the time is not a uniform u64
  28. for (int i = 0; i < 10; i++) {
  29. xoroshiro128plus_next();
  30. }
  31. }
  32. uint32_t u32() {
  33. return xoroshiro128plus_next() >> 32;
  34. }
  35. uint64_t u64() {
  36. return xoroshiro128plus_next();
  37. }
  38. float uniform() {
  39. // 24 bits of granularity is the best that can be done with floats while ensuring that the return value lies in [0.0, 1.0).
  40. return xoroshiro128plus_next() / std::pow(2.f, 64);
  41. }
  42. float normal() {
  43. // Box-Muller transform
  44. float radius = std::sqrt(-2.f * std::log(1.f - uniform()));
  45. float theta = 2.f * M_PI * uniform();
  46. return radius * std::sin(theta);
  47. // // Central Limit Theorem
  48. // const int n = 8;
  49. // float sum = 0.0;
  50. // for (int i = 0; i < n; i++) {
  51. // sum += uniform();
  52. // }
  53. // return (sum - n / 2.f) / std::sqrt(n / 12.f);
  54. }
  55. } // namespace random
  56. } // namespace rack