blitz  Version 0.9
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
rand-tt800.h
Go to the documentation of this file.
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 
120 
121 #endif // BZ_RAND_TT800_H
122