|
13 | 13 | /* Make statistics 'bout what we do.. whatever it is we _do_ */ |
14 | 14 |
|
15 | 15 | #include <assert.h> |
16 | | -#include <time.h> |
| 16 | +#include <stdint.h> |
| 17 | +#include <stdio.h> |
| 18 | +#include <stdlib.h> |
17 | 19 |
|
18 | 20 | #include "btech/context.h" |
19 | 21 | #include "command_handlers_api.h" |
|
26 | 28 | #include "registry_api.h" |
27 | 29 |
|
28 | 30 | void init_stat(BtechContext *context) { |
29 | | - btech_random_seed(&context->random, (unsigned long)time(nullptr)); |
| 31 | + if (!btech_random_seed_from_system(&context->random)) { |
| 32 | + perror("getrandom"); |
| 33 | + exit(EXIT_FAILURE); |
| 34 | + } |
30 | 35 | } |
31 | 36 |
|
32 | 37 | static const int chances[11] = {1, 2, 3, 4, 5, 6, 5, 4, 3, 2, 1}; |
@@ -67,62 +72,31 @@ void do_show_stat(CommandInvocation *invocation) { |
67 | 72 | /* |
68 | 73 | * Returns an integer chosen randomly from the interval [low,high]. |
69 | 74 | * |
70 | | - * To eliminate bias from rounding error, this routine repeatedly takes some |
71 | | - * number of high order bits from the Mersenne Twister, until it finds a value |
72 | | - * <= (high - low). If we take n bits, such that 2^n is the smallest power of |
73 | | - * two greater than (high - low), then this procedure should only require |
74 | | - * another iteration 50% or less of the time. (The actual value would be |
75 | | - * (2^n - (high - low)) / (high - low).) It also always terminates due to the |
76 | | - * statistical qualities of the Mersenne Twister, although possibly only after |
77 | | - * several (but generally very few) iterations. |
78 | | - * |
79 | | - * For example, computing a D6 should require a second iteration 33% (1/3rd) of |
80 | | - * the time, a third iteration 11% (1/9th) of the time, a fourth iteration 3.7% |
81 | | - * (1/27th) of the time, a fifth iteration 1.2% of the time (1/81st) of the |
82 | | - * time, a sixth iteration 0.4% (1/243rd) of the time, and so on. Or in other |
83 | | - * words, this will require fewer than six iterations 99.6% of the time, while |
84 | | - * completely eliminating rounding bias. |
| 75 | + * To eliminate modulo bias, this routine repeatedly draws from xoshiro256** |
| 76 | + * until it finds a value in the largest multiple of the interval width that |
| 77 | + * fits in a uint64_t. This requires at most one additional draw on average. |
85 | 78 | * |
86 | 79 | * This code is on the critical path, but modern processors can compute this |
87 | 80 | * stuff really fast. There's really no need to have the compiler inline it to |
88 | 81 | * perform further optimization. |
89 | 82 | */ |
90 | 83 | long btech_random_range(BtechContext *context, long low, long high) { |
91 | | - const unsigned long int range = (unsigned long int)(high - low); |
92 | | - |
93 | | - unsigned long value; |
94 | | - unsigned int nn; |
| 84 | + uint64_t width; |
| 85 | + uint64_t limit; |
| 86 | + uint64_t value; |
95 | 87 |
|
| 88 | + assert(context != nullptr); |
96 | 89 | assert(high >= low); |
97 | 90 |
|
98 | | - /* |
99 | | - * Compute n, the shift value. We're using the 32-bit version of the |
100 | | - * Mersenne Twister, so we only need shifts up to 32. (If we did need a |
101 | | - * larger value, we would also need to expand our random number size.) |
102 | | - * |
103 | | - * We can special case some of the common values (such as n = 8 for |
104 | | - * range = 5, for the D6) if this loop becomes a concern. |
105 | | - */ |
106 | | - for (nn = 0; nn < 32; nn++) { |
107 | | - if ((range >> nn) == 0) |
108 | | - break; |
| 91 | + width = (uint64_t)high - (uint64_t)low + UINT64_C(1); |
| 92 | + if (width == 0) { |
| 93 | + return (long)btech_random_u64(&context->random); |
109 | 94 | } |
110 | 95 |
|
111 | | - nn = 32 - nn; |
112 | | - |
113 | | - /* Shifts >= bit width are undefined in C. At least on x86, they |
114 | | - * apparently do nothing, which causes the following do-while loop to |
115 | | - * run until the generator returns 0. */ |
116 | | - if (nn == 32) { |
117 | | - return 0; |
118 | | - } |
119 | | - |
120 | | - assert(nn < 32); |
121 | | - |
122 | | - /* Repeatedly select random numbers until we get an acceptable one. */ |
| 96 | + limit = UINT64_MAX - UINT64_MAX % width; |
123 | 97 | do { |
124 | | - value = btech_random_u32(&context->random) >> nn; |
125 | | - } while (value > range); |
| 98 | + value = btech_random_u64(&context->random); |
| 99 | + } while (value >= limit); |
126 | 100 |
|
127 | | - return low + value; |
| 101 | + return (long)((uint64_t)low + value % width); |
128 | 102 | } |
0 commit comments