fork(1) download
  1. // Ranlux24 PRNG. (1.02)
  2. // Chaotic with long period but slow due to discard.
  3.  
  4. #include <stdlib.h>
  5. #include <limits.h>
  6. #include <assert.h>
  7. #include <stdio.h>
  8. #include <time.h>
  9.  
  10. // Utility.
  11.  
  12. #define BITS(x) \
  13.   (sizeof(x) * CHAR_BIT)
  14.  
  15. #define UNLIKELY(cond) \
  16.   __builtin_expect(!!(cond), 0)
  17.  
  18. #define NEW_T(T, n) \
  19.   ((T *) allocate(n * sizeof(T)))
  20.  
  21. void *allocate(size_t n)
  22. {
  23. void *r = malloc(n);
  24. assert(r != 0 || n == 0);
  25. return r;
  26. }
  27.  
  28. // Engine Base.
  29.  
  30. typedef struct EngineBase {
  31. int (*next) (struct EngineBase*);
  32. void (*discard) (struct EngineBase*, int);
  33. void (*delete) (struct EngineBase*);
  34. } EngineBase;
  35.  
  36. int EngineBase_next(EngineBase *self)
  37. {
  38. return self->next(self);
  39. }
  40.  
  41. void EngineBase_discard(EngineBase *self, int n)
  42. {
  43. return self->discard(self, n);
  44. }
  45.  
  46. void EngineBase_delete(EngineBase *self)
  47. {
  48. return self->delete(self);
  49. }
  50.  
  51. // Subtract With Carry.
  52.  
  53. typedef struct {
  54. EngineBase base;
  55. int *x, c, i, m, r, s;
  56. } SubtractWithCarryEngine;
  57.  
  58. static inline int next_(int *x, int mask, int i_s, int i_r, int *carry)
  59. {
  60. int y = x[i_s] - x[i_r] - *carry;
  61. *carry = -(y >> (BITS(int) - 1));
  62. return y & mask;
  63. }
  64.  
  65. int SubtractWithCarryEngine_next(EngineBase *base)
  66. {
  67. SubtractWithCarryEngine *self = (SubtractWithCarryEngine *) base;
  68.  
  69. self->i += 1;
  70. int *x = self->x, i = self->i, r = self->r;
  71.  
  72. if (UNLIKELY(i >= r))
  73. {
  74. int *c = &self->c, m = self->m, s = self->s, t = r - s;
  75.  
  76. for (int i = 0; i < s; i++)
  77. x[i] = next_(x, m, i + t, i, c);
  78. for (int i = s; i < r; i++)
  79. x[i] = next_(x, m, i - s, i, c);
  80. self->i = i = 0;
  81. }
  82. return x[i];
  83. }
  84.  
  85. void SubtractWithCarryEngine_discard(EngineBase *base, int n)
  86. {
  87. for (int i = 0; i < n; i++)
  88. SubtractWithCarryEngine_next(base);
  89. }
  90.  
  91. void SubtractWithCarryEngine_delete(EngineBase *base)
  92. {
  93. SubtractWithCarryEngine *self = (SubtractWithCarryEngine *) base;
  94.  
  95. if (self != 0)
  96. {
  97. free(self->x);
  98. free(self);
  99. }
  100. }
  101.  
  102. EngineBase *SubtractWithCarryEngine_new(int w, int s, int r)
  103. {
  104. assert(0 < w && w < BITS(int));
  105. assert(0 < s && s < r);
  106.  
  107. SubtractWithCarryEngine *self = NEW_T(SubtractWithCarryEngine, 1);
  108.  
  109. self->base.next = SubtractWithCarryEngine_next;
  110. self->base.discard = SubtractWithCarryEngine_discard;
  111. self->base.delete = SubtractWithCarryEngine_delete;
  112. self->x = NEW_T(int, r);
  113. self->c = 0;
  114. self->i = r-1;
  115. self->m = -1U >> (BITS(int) - w);
  116. self->r = r;
  117. self->s = s;
  118.  
  119. for (int i = 0; i < r; i++)
  120. self->x[i] = rand() & self->m;
  121. return (EngineBase *) self;
  122. }
  123.  
  124. // Discard Block.
  125.  
  126. typedef struct {
  127. EngineBase base, *owned;
  128. int p, r, i;
  129. } DiscardBlockEngine;
  130.  
  131. int DiscardBlockEngine_next(EngineBase *base)
  132. {
  133. DiscardBlockEngine *self = (DiscardBlockEngine *) base;
  134.  
  135. if (self->i == 0)
  136. {
  137. EngineBase_discard(self->owned, self->p - self->r);
  138. self->i = self->r;
  139. }
  140. self->i -= 1;
  141. return EngineBase_next(self->owned);
  142. }
  143.  
  144. void DiscardBlockEngine_discard(EngineBase *base, int n)
  145. {
  146. for (int i = 0; i < n; i++)
  147. DiscardBlockEngine_next(base);
  148. }
  149.  
  150. void DiscardBlockEngine_delete(EngineBase *base)
  151. {
  152. DiscardBlockEngine *self = (DiscardBlockEngine *) base;
  153.  
  154. if (self != 0)
  155. {
  156. EngineBase_delete(self->owned);
  157. free(self);
  158. }
  159. }
  160.  
  161. EngineBase *DiscardBlockEngine_new(EngineBase **unique, int p, int r)
  162. {
  163. assert(0 < r <= p);
  164.  
  165. DiscardBlockEngine *self = NEW_T(DiscardBlockEngine, 1);
  166.  
  167. self->base.next = DiscardBlockEngine_next;
  168. self->base.discard = DiscardBlockEngine_discard;
  169. self->base.delete = DiscardBlockEngine_delete;
  170. self->owned = *unique; *unique = 0; // Transfer ownership.
  171. self->p = p;
  172. self->r = r;
  173. self->i = r;
  174. return (EngineBase *) self;
  175. }
  176.  
  177. // Ranlux24.
  178.  
  179. EngineBase *Ranlux24_new(void)
  180. {
  181. EngineBase *ranlux24_base = SubtractWithCarryEngine_new(24, 10, 24);
  182. return DiscardBlockEngine_new(&ranlux24_base, 223, 23);
  183. }
  184.  
  185. // Main.
  186.  
  187. double clock_now(void)
  188. {
  189. struct timespec now;
  190. clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &now);
  191. return now.tv_sec + now.tv_nsec / 1.0E+09;
  192. }
  193.  
  194. int main(void)
  195. {
  196. srand(time(0));
  197. EngineBase *ranlux24 = Ranlux24_new();
  198.  
  199. for (int i = 0; i < 24; i++)
  200. {
  201. int r = EngineBase_next(ranlux24);
  202. printf("%d\n", r);
  203. }
  204.  
  205. int n = 10000000;
  206. double t = clock_now();
  207. EngineBase_discard(ranlux24, n);
  208. printf("Elapsed: %.9fs\n", clock_now() - t);
  209.  
  210. EngineBase_delete(ranlux24);
  211. return 0;
  212. }
Success #stdin #stdout 0.33s 5284KB
stdin
Standard input is empty
stdout
6350222
14567412
1754403
5428135
16749895
5074645
16154843
13787438
4870379
11771697
7398842
3657886
9208363
10261427
4494242
4418033
13050327
2553420
2196812
8806070
1844736
9687196
11037771
9513899
Elapsed: 0.322370179s