1 | /*************************************************************************** |
---|
2 | * blitz/rand-tt800.h Matsumoto and Kurita's TT800 uniform random |
---|
3 | * number generator. |
---|
4 | * |
---|
5 | * $Id: rand-tt800.h,v 1.2 2003/01/14 11:29:18 patricg Exp $ |
---|
6 | * |
---|
7 | *************************************************************************** |
---|
8 | * |
---|
9 | * The class TT800 encapsulates Makoto Matsumoto and Yoshiharu Kurita's |
---|
10 | * TT800 twisted generalized feedback shift register (TGFSR) random number |
---|
11 | * generator. The generator has period 2^800 - 1. |
---|
12 | * |
---|
13 | * Contact: M. Matsumoto <matumoto@math.keio.ac.jp> |
---|
14 | * |
---|
15 | * See: M. Matsumoto and Y. Kurita, Twisted GFSR Generators II, |
---|
16 | * ACM Transactions on Modelling and Computer Simulation, |
---|
17 | * Vol. 4, No. 3, 1994, pages 254-266. |
---|
18 | * |
---|
19 | * (c) 1994 Association for Computing Machinery. |
---|
20 | * |
---|
21 | * Distributed with consent of the authors. |
---|
22 | * |
---|
23 | ***************************************************************************/ |
---|
24 | |
---|
25 | #ifndef BZ_RAND_TT800_H |
---|
26 | #define BZ_RAND_TT800_H |
---|
27 | |
---|
28 | #ifndef BZ_RANDOM_H |
---|
29 | #include <blitz/random.h> |
---|
30 | #endif |
---|
31 | |
---|
32 | BZ_NAMESPACE(blitz) |
---|
33 | |
---|
34 | class TT800 { |
---|
35 | |
---|
36 | public: |
---|
37 | typedef double T_numtype; |
---|
38 | |
---|
39 | TT800(double low = 0.0, double high = 1.0, double = 0.0) |
---|
40 | : low_(low), length_(high-low) |
---|
41 | { |
---|
42 | // Initial 25 seeds |
---|
43 | x[0] = 0x95f24dab; x[1] = 0x0b685215; x[2] = 0xe76ccae7; |
---|
44 | x[3] = 0xaf3ec239; x[4] = 0x715fad23; x[5] = 0x24a590ad; |
---|
45 | x[6] = 0x69e4b5ef; x[7] = 0xbf456141; x[8] = 0x96bc1b7b; |
---|
46 | x[9] = 0xa7bdf825; x[10] = 0xc1de75b7; x[11] = 0x8858a9c9; |
---|
47 | x[12] = 0x2da87693; x[13] = 0xb657f9dd; x[14] = 0xffdc8a9f; |
---|
48 | x[15] = 0x8121da71; x[16] = 0x8b823ecb; x[17] = 0x885d05f5; |
---|
49 | x[18] = 0x4e20cd47; x[19] = 0x5a9ad5d9; x[20] = 0x512c0c03; |
---|
50 | x[21] = 0xea857ccd; x[22] = 0x4cc1d30f; x[23] = 0x8891a8a1; |
---|
51 | x[24] = 0xa6b7aadb; |
---|
52 | |
---|
53 | // Magic vector 'a', don't change |
---|
54 | mag01[0] = 0; |
---|
55 | mag01[1] = 0x8ebfd028; |
---|
56 | |
---|
57 | k = 0; |
---|
58 | |
---|
59 | // f = 1/(2^32); |
---|
60 | |
---|
61 | f = (1.0 / 65536) / 65536; |
---|
62 | } |
---|
63 | |
---|
64 | void randomize() |
---|
65 | { |
---|
66 | BZ_NOT_IMPLEMENTED(); // NEEDS_WORK |
---|
67 | } |
---|
68 | |
---|
69 | unsigned long randomUint32() |
---|
70 | { |
---|
71 | if (k==N) |
---|
72 | generate(); |
---|
73 | |
---|
74 | unsigned long y = x[k]; |
---|
75 | y ^= (y << 7) & 0x2b5b2500; /* s and b, magic vectors */ |
---|
76 | y ^= (y << 15) & 0xdb8b0000; /* t and c, magic vectors */ |
---|
77 | y &= 0xffffffff; /* you may delete this line if word size = 32 */ |
---|
78 | |
---|
79 | // the following line was added by Makoto Matsumoto in the 1996 version |
---|
80 | // to improve lower bit's corellation. |
---|
81 | // Delete this line to use the code published in 1994. |
---|
82 | |
---|
83 | y ^= (y >> 16); /* added to the 1994 version */ |
---|
84 | k++; |
---|
85 | } |
---|
86 | |
---|
87 | double random() |
---|
88 | { |
---|
89 | unsigned long y1 = randomUint32(); |
---|
90 | unsigned long y2 = randomUint32(); |
---|
91 | |
---|
92 | return low_ + length_ * ((y1 * f) * y2 * f); |
---|
93 | } |
---|
94 | |
---|
95 | protected: |
---|
96 | void generate() |
---|
97 | { |
---|
98 | /* generate N words at one time */ |
---|
99 | int kk; |
---|
100 | for (kk=0;kk<N-M;kk++) { |
---|
101 | x[kk] = x[kk+M] ^ (x[kk] >> 1) ^ mag01[x[kk] % 2]; |
---|
102 | } |
---|
103 | for (; kk<N;kk++) { |
---|
104 | x[kk] = x[kk+(M-N)] ^ (x[kk] >> 1) ^ mag01[x[kk] % 2]; |
---|
105 | } |
---|
106 | k=0; |
---|
107 | } |
---|
108 | |
---|
109 | private: |
---|
110 | enum { N = 25, M = 7 }; |
---|
111 | |
---|
112 | double low_, length_; |
---|
113 | double f; |
---|
114 | int k; |
---|
115 | unsigned long x[N]; |
---|
116 | unsigned long mag01[2]; |
---|
117 | }; |
---|
118 | |
---|
119 | BZ_NAMESPACE_END |
---|
120 | |
---|
121 | #endif // BZ_RAND_TT800_H |
---|
122 | |
---|