[1] | 1 | #include <stdlib.h> |
---|
| 2 | |
---|
| 3 | static randbuf rand48buf; |
---|
| 4 | #define A_0 0xE66D |
---|
| 5 | #define A_1 0xDEEC |
---|
| 6 | #define A_2 0x5 |
---|
| 7 | #define C 0xB |
---|
| 8 | static randbuf a = { A_0, A_1, A_2 }; |
---|
| 9 | static unsigned short c = C; |
---|
| 10 | |
---|
| 11 | static void calc_next(randbuf buf) { |
---|
| 12 | randbuf tmp; |
---|
| 13 | long t; |
---|
| 14 | t = buf[0] * a[0] + c; |
---|
| 15 | tmp[0] = t & 0xffff; |
---|
| 16 | tmp[1] = (t >> 16) & 0xffff; |
---|
| 17 | t = buf[1] * a[0] + buf[0] * a[1] + tmp[1]; |
---|
| 18 | tmp[1] = t & 0xffff; |
---|
| 19 | tmp[2] = (t >> 16) & 0xffff; |
---|
| 20 | t = buf[2] * a[0] + buf[1] * a[1] + buf[0] * a[2] + tmp[2]; |
---|
| 21 | tmp[2] = t & 0xffff; |
---|
| 22 | buf[0] = tmp[0]; |
---|
| 23 | buf[1] = tmp[1]; |
---|
| 24 | buf[2] = tmp[2]; |
---|
| 25 | } |
---|
| 26 | |
---|
| 27 | |
---|
| 28 | long lrand48(void) { |
---|
| 29 | return nrand48(rand48buf); |
---|
| 30 | } |
---|
| 31 | |
---|
| 32 | long mrand48(void) { |
---|
| 33 | return jrand48(rand48buf); |
---|
| 34 | } |
---|
| 35 | |
---|
| 36 | double drand48(void) { |
---|
| 37 | return frand48(rand48buf); |
---|
| 38 | } |
---|
| 39 | |
---|
| 40 | void srand48(long seed) { |
---|
| 41 | rand48buf[1] = (seed >> 16) & 0xffff; |
---|
| 42 | rand48buf[2] = seed & 0xffff; |
---|
| 43 | rand48buf[0] = 0x330e; |
---|
| 44 | a[0] = A_0; |
---|
| 45 | a[1] = A_1; |
---|
| 46 | a[2] = A_2; |
---|
| 47 | c = C; |
---|
| 48 | } |
---|
| 49 | |
---|
| 50 | |
---|
| 51 | double frand48(randbuf buf) { |
---|
| 52 | double ret; |
---|
| 53 | static double two16m = 1.0 / (1L << 16); |
---|
| 54 | ret = (two16m * (two16m * (two16m * buf[0] + buf[1]) + buf[2])); |
---|
| 55 | calc_next(buf); |
---|
| 56 | return ret; |
---|
| 57 | } |
---|
| 58 | |
---|
| 59 | |
---|
| 60 | |
---|
| 61 | unsigned short *seed48(randbuf buf) { |
---|
| 62 | static randbuf oldx; |
---|
| 63 | int i; |
---|
| 64 | for (i = 0; i < 3; i++) { |
---|
| 65 | oldx[i] = rand48buf[i]; |
---|
| 66 | rand48buf[i] = buf[i]; |
---|
| 67 | } |
---|
| 68 | a[0] = A_0; |
---|
| 69 | a[1] = A_1; |
---|
| 70 | a[2] = A_2; |
---|
| 71 | c = C; |
---|
| 72 | return (unsigned short *)&oldx; |
---|
| 73 | } |
---|
| 74 | |
---|
| 75 | void lcong48(unsigned short param[7]) { |
---|
| 76 | int i; |
---|
| 77 | for (i = 0; i < 3; i++) { |
---|
| 78 | rand48buf[i] = param[i]; |
---|
| 79 | a[i] = param[i + 3]; |
---|
| 80 | } |
---|
| 81 | c = param[6]; |
---|
| 82 | } |
---|
| 83 | |
---|
| 84 | long jrand48(randbuf buf) { |
---|
| 85 | long ret; |
---|
| 86 | ret = buf[2] << 16 | buf[1]; |
---|
| 87 | calc_next(buf); |
---|
| 88 | return ret; |
---|
| 89 | } |
---|
| 90 | |
---|
| 91 | long nrand48(randbuf buf) { |
---|
| 92 | return jrand48(buf) & 0x7FFFFFFFL; |
---|
| 93 | } |
---|