Logo ROOT   6.07/09
Reference Guide
TRandom3.cxx
Go to the documentation of this file.
1 // @(#)root/mathcore:$Id$
2 // Author: Peter Malzacher 31/08/99
3 
4 /**
5 
6 \class TRandom3
7 
8 Random number generator class based on
9  M. Matsumoto and T. Nishimura,
10  Mersenne Twister: A 623-diminsionally equidistributed
11  uniform pseudorandom number generator
12  ACM Transactions on Modeling and Computer Simulation,
13  Vol. 8, No. 1, January 1998, pp 3--30.
14 
15 For more information see the Mersenne Twister homepage
16  [http://www.math.keio.ac.jp/~matumoto/emt.html]
17 
18 Advantage:
19 
20 - large period 2**19937 -1
21 - relativly fast (slightly slower than TRandom1 and TRandom2 but much faster than TRandom1)
22 
23 Drawback: a relative large internal state of 624 integers
24 
25 @ingroup Random
26 
27 */
28 
29 //////////////////////////////////////////////////////////////////////
30 // Aug.99 ROOT implementation based on CLHEP by P.Malzacher
31 //
32 // the original code contains the following copyright notice:
33 /* This library is free software; you can redistribute it and/or */
34 /* modify it under the terms of the GNU Library General Public */
35 /* License as published by the Free Software Foundation; either */
36 /* version 2 of the License, or (at your option) any later */
37 /* version. */
38 /* This library is distributed in the hope that it will be useful, */
39 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
40 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. */
41 /* See the GNU Library General Public License for more details. */
42 /* You should have received a copy of the GNU Library General */
43 /* Public License along with this library; if not, write to the */
44 /* Free Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA */
45 /* 02111-1307 USA */
46 /* Copyright (C) 1997 Makoto Matsumoto and Takuji Nishimura. */
47 /* When you use this, send an email to: matumoto@math.keio.ac.jp */
48 /* with an appropriate reference to your work. */
49 /////////////////////////////////////////////////////////////////////
50 
51 #include "TRandom3.h"
52 #include "TBuffer.h"
53 #include "TRandom2.h"
54 #include "TClass.h"
55 #include "TUUID.h"
56 
58 #ifdef R__COMPLETE_MEM_TERMINATION
59 namespace {
60  struct TRandomCleanup {
61  ~TRandomCleanup() { delete gRandom; gRandom = 0; }
62  };
63  static TRandomCleanup gCleanupRandom;
64 }
65 #endif
66 
68 
69 ////////////////////////////////////////////////////////////////////////////////
70 ///*-*-*-*-*-*-*-*-*-*-*default constructor*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
71 /// If seed is 0, the seed is automatically computed via a TUUID object.
72 /// In this case the seed is guaranteed to be unique in space and time.
73 
75 {
76  SetName("Random3");
77  SetTitle("Random number generator: Mersenne Twister");
78  SetSeed(seed);
79 }
80 
81 ////////////////////////////////////////////////////////////////////////////////
82 ///*-*-*-*-*-*-*-*-*-*-*default destructor*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
83 ///*-* ==================
84 
86 {
87 }
88 
89 ////////////////////////////////////////////////////////////////////////////////
90 /// Machine independent random number generator.
91 /// Produces uniformly-distributed floating points in (0,1)
92 /// Method: Mersenne Twister
93 
95 {
96  UInt_t y;
97 
98  const Int_t kM = 397;
99  const Int_t kN = 624;
100  const UInt_t kTemperingMaskB = 0x9d2c5680;
101  const UInt_t kTemperingMaskC = 0xefc60000;
102  const UInt_t kUpperMask = 0x80000000;
103  const UInt_t kLowerMask = 0x7fffffff;
104  const UInt_t kMatrixA = 0x9908b0df;
105 
106  if (fCount624 >= kN) {
107  Int_t i;
108 
109  for (i=0; i < kN-kM; i++) {
110  y = (fMt[i] & kUpperMask) | (fMt[i+1] & kLowerMask);
111  fMt[i] = fMt[i+kM] ^ (y >> 1) ^ ((y & 0x1) ? kMatrixA : 0x0);
112  }
113 
114  for ( ; i < kN-1 ; i++) {
115  y = (fMt[i] & kUpperMask) | (fMt[i+1] & kLowerMask);
116  fMt[i] = fMt[i+kM-kN] ^ (y >> 1) ^ ((y & 0x1) ? kMatrixA : 0x0);
117  }
118 
119  y = (fMt[kN-1] & kUpperMask) | (fMt[0] & kLowerMask);
120  fMt[kN-1] = fMt[kM-1] ^ (y >> 1) ^ ((y & 0x1) ? kMatrixA : 0x0);
121  fCount624 = 0;
122  }
123 
124  y = fMt[fCount624++];
125  y ^= (y >> 11);
126  y ^= ((y << 7 ) & kTemperingMaskB );
127  y ^= ((y << 15) & kTemperingMaskC );
128  y ^= (y >> 18);
129 
130  // 2.3283064365386963e-10 == 1./(max<UINt_t>+1) -> then returned value cannot be = 1.0
131  if (y) return ( (Double_t) y * 2.3283064365386963e-10); // * Power(2,-32)
132  return Rndm();
133 }
134 
135 ////////////////////////////////////////////////////////////////////////////////
136 /// Return an array of n random numbers uniformly distributed in ]0,1]
137 
139 {
140  for(Int_t i=0; i<n; i++) array[i]=(Float_t)Rndm();
141 }
142 
143 ////////////////////////////////////////////////////////////////////////////////
144 /// Return an array of n random numbers uniformly distributed in ]0,1]
145 
147 {
148  Int_t k = 0;
149 
150  UInt_t y;
151 
152  const Int_t kM = 397;
153  const Int_t kN = 624;
154  const UInt_t kTemperingMaskB = 0x9d2c5680;
155  const UInt_t kTemperingMaskC = 0xefc60000;
156  const UInt_t kUpperMask = 0x80000000;
157  const UInt_t kLowerMask = 0x7fffffff;
158  const UInt_t kMatrixA = 0x9908b0df;
159 
160  while (k < n) {
161  if (fCount624 >= kN) {
162  Int_t i;
163 
164  for (i=0; i < kN-kM; i++) {
165  y = (fMt[i] & kUpperMask) | (fMt[i+1] & kLowerMask);
166  fMt[i] = fMt[i+kM] ^ (y >> 1) ^ ((y & 0x1) ? kMatrixA : 0x0);
167  }
168 
169  for ( ; i < kN-1 ; i++) {
170  y = (fMt[i] & kUpperMask) | (fMt[i+1] & kLowerMask);
171  fMt[i] = fMt[i+kM-kN] ^ (y >> 1) ^ ((y & 0x1) ? kMatrixA : 0x0);
172  }
173 
174  y = (fMt[kN-1] & kUpperMask) | (fMt[0] & kLowerMask);
175  fMt[kN-1] = fMt[kM-1] ^ (y >> 1) ^ ((y & 0x1) ? kMatrixA : 0x0);
176  fCount624 = 0;
177  }
178 
179  y = fMt[fCount624++];
180  y ^= (y >> 11);
181  y ^= ((y << 7 ) & kTemperingMaskB );
182  y ^= ((y << 15) & kTemperingMaskC );
183  y ^= (y >> 18);
184 
185  if (y) {
186  array[k] = Double_t( y * 2.3283064365386963e-10); // * Power(2,-32)
187  k++;
188  }
189  }
190 }
191 
192 ////////////////////////////////////////////////////////////////////////////////
193 /// Set the random generator sequence
194 /// if seed is 0 (default value) a TUUID is generated and used to fill
195 /// the first 8 integers of the seed array.
196 /// In this case the seed is guaranteed to be unique in space and time.
197 /// Use upgraded seeding procedure to fix a known problem when seeding with values
198 /// with many zero in the bit pattern (like 2**28).
199 /// see http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html
200 
202 {
203  TRandom::SetSeed(seed);
204  fCount624 = 624;
205  if (seed > 0) {
206  fMt[0] = fSeed;
207 
208  // use multipliers from Knuth's "Art of Computer Programming" Vol. 2, 3rd Ed. p.106
209  for(Int_t i=1; i<624; i++) {
210  fMt[i] = (1812433253 * ( fMt[i-1] ^ ( fMt[i-1] >> 30)) + i );
211  }
212 
213  } else {
214 
215  // use TRandom2 (which is based on TUUId to generate the seed
216  // TRandom2 works fairly well and has been tested against example
217  // layout in https://savannah.cern.ch/bugs/?99516
218  TRandom2 r(0);
219  for (Int_t i = 0; i< 624; i++) {
220  fMt[i] = static_cast<UInt_t> (4294967296.*r.Rndm());
221  }
222  // warm up the generator calling it 10 times
223  for (Int_t i = 0; i < 10; ++i) Rndm();
224  }
225 
226 
227 }
228 
229 ////////////////////////////////////////////////////////////////////////////////
230 /// Stream an object of class TRandom3.
231 
232 void TRandom3::Streamer(TBuffer &R__b)
233 {
234  if (R__b.IsReading()) {
235  UInt_t R__s, R__c;
236  Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
237  if (R__v > 1) {
238  R__b.ReadClassBuffer(TRandom3::Class(), this, R__v, R__s, R__c);
239  return;
240  }
241  //====process old versions before automatic schema evolution
242  TRandom::Streamer(R__b);
243  R__b.ReadStaticArray(fMt);
244  R__b >> fCount624;
245  R__b.CheckByteCount(R__s, R__c, TRandom3::IsA());
246  //====end of old versions
247 
248  } else {
249  R__b.WriteClassBuffer(TRandom3::Class(),this);
250  }
251 }
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
Random number generator class based on M.
Definition: TRandom3.h:29
Bool_t IsReading() const
Definition: TBuffer.h:83
short Version_t
Definition: RtypesCore.h:61
virtual Double_t Rndm()
Machine independent random number generator.
Definition: TRandom3.cxx:94
float Float_t
Definition: RtypesCore.h:53
virtual ~TRandom3()
*-*-*-*-*-*-*-*-*-*-*default destructor*-*-*-*-*-*-*-*-*-*-*-*-*-*-* *-* ================== ...
Definition: TRandom3.cxx:85
Random number generator class based on the maximally quidistributed combined Tausworthe generator by ...
Definition: TRandom2.h:29
virtual void SetSeed(ULong_t seed=0)
Set the random generator sequence if seed is 0 (default value) a TUUID is generated and used to fill ...
Definition: TRandom3.cxx:201
Buffer base class used for serializing objects.
Definition: TBuffer.h:42
virtual Int_t CheckByteCount(UInt_t startpos, UInt_t bcnt, const TClass *clss)=0
int Int_t
Definition: RtypesCore.h:41
const char * Class
Definition: TXMLSetup.cxx:64
virtual Int_t ReadStaticArray(Bool_t *b)=0
This is the base class for the ROOT Random number generators.
Definition: TRandom.h:31
virtual void SetSeed(ULong_t seed=0)
Set the random generator seed.
Definition: TRandom.cxx:568
TRandom2 r(17)
TRandom * gRandom
Definition: TRandom3.cxx:57
unsigned int UInt_t
Definition: RtypesCore.h:42
virtual void RndmArray(Int_t n, Float_t *array)
Return an array of n random numbers uniformly distributed in ]0,1].
Definition: TRandom3.cxx:138
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
virtual Double_t Rndm()
TausWorth generator from L&#39;Ecuyer, uses as seed 3x32bits integers Use a mask of 0xffffffffUL to make ...
Definition: TRandom2.cxx:58
static const double x1[5]
#define ClassImp(name)
Definition: Rtypes.h:279
double Double_t
Definition: RtypesCore.h:55
unsigned long ULong_t
Definition: RtypesCore.h:51
Double_t y[n]
Definition: legend1.C:17
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
const Int_t n
Definition: legend1.C:16
gr SetName("gr")
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0