diff options
-rw-r--r-- | sys/libn/random.c | 90 |
1 files changed, 90 insertions, 0 deletions
diff --git a/sys/libn/random.c b/sys/libn/random.c new file mode 100644 index 0000000..1e2c2d5 --- /dev/null +++ b/sys/libn/random.c @@ -0,0 +1,90 @@ +#include <u.h> + +// ---------------------------------------------------------------------------- +// Internal structure + +uint64 +rol64(uint64 x, int k) +{ + return (x << k) | (x >> (64 - k)); +} + +typedef struct Rng { + uint64 s[4]; +} Rng; + +uint64 +xoshiro256ss(Rng *state) +{ + uint64 *s = state->s; + uint64 result = rol64(s[1] * 5, 7) * 9; + uint64 t = s[1] << 17; + + s[2] ^= s[0]; + s[3] ^= s[1]; + s[1] ^= s[2]; + s[0] ^= s[3]; + + s[2] ^= t; + s[3] = rol64(s[3], 45); + + return result; +} + +typedef struct Mix +{ + uint64 s; +} Mix; + +uint64 +splitmix64(struct Mix *state) { + uint64 result = state->s; + + state->s = result + 0x9E3779B97f4A7C15; + result = (result ^ (result >> 30)) * 0xBF58476D1CE4E5B9; + result = (result ^ (result >> 27)) * 0x94D049BB133111EB; + return result ^ (result >> 31); +} + +static Rng RNG; + +// ---------------------------------------------------------------------------- +// Exported functions + +/* Initializes the global RNG */ + +error +rng·init(uint64 seed) +{ + Mix smstate = {seed}; + + for (int i=0; i < 4; i++) + RNG.s[i] = splitmix64(&smstate); + + return 0; +} + +/* Returns a random float64 between 0 and 1 */ +double +rng·random() +{ + uint64 r = xoshiro256ss(&RNG); + return (double)r / UINT64_MAX; +} + +/* Returns true or false on success of trial */ +bool +rng·bernoulli(double f) +{ + return rng·random() < f; +} + +/* Returns a random integer between 0 and max + * TODO: Modulo throws off uniformity + */ +uint64 +rng·randi(int max) +{ + uint64 r = xoshiro256ss(&RNG); + return r % max; +} |