SHOW:
|
|
- or go back to the newest paste.
1 | #include <stdio.h> | |
2 | #include <stdint.h> | |
3 | #include <stdlib.h> | |
4 | ||
5 | /* http://marc-b-reynolds.github.io/math/2017/09/18/ModInverse.html */ | |
6 | static uint32_t modinv(uint32_t a) | |
7 | { | |
8 | uint32_t x; | |
9 | x = a; // 3 bits | |
10 | x *= 2 - a*x; // 6 | |
11 | x *= 2 - a*x; // 12 | |
12 | x *= 2 - a*x; // 24 | |
13 | x *= 2 - a*x; // 48 / retaining bottom 32 | |
14 | return x; | |
15 | } | |
16 | ||
17 | /* http://number-none.com/blow/blog/programming/2016/07/08/fabian-on-lcg-fast-forward.html */ | |
18 | uint32_t skip(uint32_t x, uint32_t a, uint32_t c, unsigned n) | |
19 | { | |
20 | while (n) { | |
21 | if (n & 1) | |
22 | x = a * x + c; | |
23 | c *= a + 1; | |
24 | a *= a; | |
25 | n >>= 1; | |
26 | } | |
27 | return x; | |
28 | } | |
29 | ||
30 | int main(int argc, char **argv) | |
31 | { | |
32 | srandom(time(NULL)); | |
33 | unsigned i, s = 1, a = random() | 1UL, c = random() | 1UL; | |
34 | ||
35 | /* move forward one step at a time */ | |
36 | for (i = 0; i < 100; i++) | |
37 | s = s * a + c; | |
38 | printf("%08x\n", s); | |
39 | /* skip forwaward */ | |
40 | printf("%08x\n", skip(1, a, c, 100)); | |
41 | /* you could loop through the full LCG cycle */ | |
42 | printf("%08x\n", skip(s, a, c, -100)); | |
43 | /* or inverse of LCG is still LCG */ | |
44 | a = modinv(a); | |
45 | c = -c * a; | |
46 | /* skip backward */ | |
47 | printf("%08x\n", skip(s, a, c, 100)); | |
48 | /* move backward one step at a time */ | |
49 | for (i = 0; i < 100; i++) | |
50 | s = s * a + c; | |
51 | printf("%08x\n", s); | |
52 | } |