fork(1) download
  1. // Ranlux24. (1.01)
  2. // These are probably too clunky to be generally useful.
  3. // @todo Consolidate to a single class.
  4.  
  5. #include <stdlib.h>
  6. #include <limits.h>
  7. #include <assert.h>
  8. #include <stdio.h>
  9. #include <time.h>
  10.  
  11. // Utility.
  12.  
  13. #define unlikely(c) \
  14.   __builtin_expect(!!(c), 0)
  15.  
  16. #define BITS(x) \
  17.   (sizeof(x) * CHAR_BIT)
  18.  
  19. #define NEW_T(T, count) \
  20.   ((T *) allocate(count * sizeof(T)))
  21.  
  22. void *allocate(size_t n)
  23. {
  24. void *r = malloc(n);
  25. assert(r != 0 || n == 0);
  26. return r;
  27. }
  28.  
  29. // EngineBase.
  30.  
  31. typedef struct EngineBase {
  32. int (*next) (struct EngineBase*);
  33. void (*discard) (struct EngineBase*, int);
  34. void (*delete) (struct EngineBase*);
  35. } EngineBase;
  36.  
  37. int nextEngineBase(EngineBase *self)
  38. {
  39. return self->next(self);
  40. }
  41.  
  42. void discardEngineBase(EngineBase *self, int n)
  43. {
  44. return self->discard(self, n);
  45. }
  46.  
  47. void deleteEngineBase(EngineBase *self)
  48. {
  49. return self->delete(self);
  50. }
  51.  
  52. // SubtractWithCarry.
  53.  
  54. typedef struct {
  55. EngineBase base;
  56. int *x, c, i, m, r, s;
  57. } SubtractWithCarryEngine;
  58.  
  59. static inline int next_(int *x, int mask, int i_s, int i_r, int *carry)
  60. {
  61. int y = x[i_s] - x[i_r] - *carry;
  62. *carry = -(y >> (BITS(int)-1));
  63. return y & mask;
  64. }
  65.  
  66. int nextSubtractWithCarryEngine(EngineBase *base)
  67. {
  68. SubtractWithCarryEngine *self = (SubtractWithCarryEngine *) base;
  69.  
  70. self->i = self->i + 1;
  71. int *x = self->x, i = self->i, r = self->r;
  72.  
  73. if (unlikely(i >= r))
  74. {
  75. int *c = &self->c, m = self->m, s = self->s, t = r - s;
  76.  
  77. for (int i = 0; i < s; i++)
  78. x[i] = next_(x, m, i + t, i, c);
  79. for (int i = s; i < r; i++)
  80. x[i] = next_(x, m, i - s, i, c);
  81. self->i = i = 0;
  82. }
  83. return x[i];
  84. }
  85.  
  86. void discardSubtractWithCarryEngine(EngineBase *base, int n)
  87. {
  88. for (int i = 0; i < n; i++)
  89. nextSubtractWithCarryEngine(base);
  90. }
  91.  
  92. void deleteSubtractWithCarryEngine(EngineBase *base)
  93. {
  94. SubtractWithCarryEngine *self = (SubtractWithCarryEngine *) base;
  95.  
  96. if (self != 0)
  97. {
  98. free(self->x);
  99. free(self);
  100. }
  101. }
  102.  
  103. EngineBase *newSubtractWithCarryEngine(int w, int s, int r)
  104. {
  105. assert(0 < w && w < BITS(int));
  106. assert(0 < s && s < r);
  107.  
  108. SubtractWithCarryEngine *self = NEW_T(SubtractWithCarryEngine, 1);
  109.  
  110. self->base.next = nextSubtractWithCarryEngine;
  111. self->base.discard = discardSubtractWithCarryEngine;
  112. self->base.delete = deleteSubtractWithCarryEngine;
  113. self->x = NEW_T(int, r);
  114. self->c = 0;
  115. self->i = r-1;
  116. self->m = -1U >> (BITS(int) - w);
  117. self->r = r;
  118. self->s = s;
  119.  
  120. for (int i = 0; i < r; i++)
  121. self->x[i] = rand() & self->m;
  122. return (EngineBase *) self;
  123. }
  124.  
  125. // DiscardBlock.
  126.  
  127. typedef struct {
  128. EngineBase base, *owned;
  129. int p, r, i;
  130. } DiscardBlockEngine;
  131.  
  132. int nextDiscardBlockEngine(EngineBase *base)
  133. {
  134. DiscardBlockEngine *self = (DiscardBlockEngine *) base;
  135.  
  136. if (self->i == 0)
  137. {
  138. discardEngineBase(self->owned, self->p - self->r);
  139. self->i = self->r;
  140. }
  141. self->i -= 1;
  142. return nextEngineBase(self->owned);
  143. }
  144.  
  145. void deleteDiscardBlockEngine(EngineBase *base)
  146. {
  147. DiscardBlockEngine *self = (DiscardBlockEngine *) base;
  148.  
  149. if (self != 0)
  150. {
  151. deleteEngineBase(self->owned);
  152. free(self);
  153. }
  154. }
  155.  
  156. EngineBase *newDiscardBlockEngine(EngineBase **unique, int p, int r)
  157. {
  158. assert(0 < r <= p);
  159.  
  160. DiscardBlockEngine *self = NEW_T(DiscardBlockEngine, 1);
  161.  
  162. self->base.next = nextDiscardBlockEngine;
  163. self->base.discard = 0;
  164. self->base.delete = deleteDiscardBlockEngine;
  165. self->owned = *unique; *unique = 0; // Transfer ownership.
  166. self->p = p;
  167. self->r = r;
  168. self->i = r;
  169. return (EngineBase *) self;
  170. }
  171.  
  172. // Ranlux24.
  173.  
  174. EngineBase *newRanlux24(void)
  175. {
  176. EngineBase *ranlux24_base = newSubtractWithCarryEngine(24, 10, 24);
  177. return newDiscardBlockEngine(&ranlux24_base, 223, 23);
  178. }
  179.  
  180. int main(void)
  181. {
  182. srand(time(0));
  183. EngineBase *ranlux24 = newRanlux24();
  184.  
  185. for (int i = 0; i < 24; i++)
  186. {
  187. int r = nextEngineBase(ranlux24);
  188. printf("%d\n", r);
  189. }
  190. deleteEngineBase(ranlux24);
  191. return 0;
  192. }
Success #stdin #stdout 0s 5292KB
stdin
Standard input is empty
stdout
6112115
4254633
1046442
280222
9726444
16598228
7780349
12317195
16677295
10341016
7342500
861688
10310561
41450
3273101
1646951
8034073
5838112
13490451
2403694
223761
7054338
13529678
5942504