CLHEP VERSION Reference Documentation
   
CLHEP Home Page     CLHEP Documentation     CLHEP Bug Reports

Hurd160Engine.cc
Go to the documentation of this file.
1 // $Id: Hurd160Engine.cc,v 1.7 2010/07/20 18:07:17 garren Exp $
2 // -*- C++ -*-
3 //
4 // -----------------------------------------------------------------------
5 // HEP Random
6 // Ken Smith - Initial draft started: 23rd Jul 1998
7 // - Added conversion operators: 6th Aug 1998
8 // J. Marraffino - Added some explicit casts to deal with
9 // machines where sizeof(int) != sizeof(long) 22 Aug 1998
10 // M. Fischler - Modified use of the various exponents of 2
11 // to avoid per-instance space overhead and
12 // correct the rounding procedure 15 Sep 1998
13 // - Modified various methods to avoid any
14 // possibility of predicting the next number
15 // based on the last several: We skip one
16 // 32-bit word per cycle. 15 Sep 1998
17 // - modify word[0] in two of the constructors
18 // so that no sequence can ever accidentally
19 // be produced by differnt seeds. 15 Sep 1998
20 // J. Marraffino - Remove dependence on hepStrings class 13 May 1999
21 // M. Fischler - Put endl at end of a save 10 Apr 2001
22 // M. Fischler - In restore, checkFile for file not found 03 Dec 2004
23 // M. Fischler - Methods put, get for instance save/restore 12/8/04
24 // M. Fischler - split get() into tag validation and
25 // getState() for anonymous restores 12/27/04
26 // M. Fischler - put/get for vectors of ulongs 3/14/05
27 // M. Fischler - State-saving using only ints, for portability 4/12/05
28 //
29 // =======================================================================
30 
31 #include "CLHEP/Random/defs.h"
32 #include "CLHEP/Random/Random.h"
33 #include "CLHEP/Random/Hurd160Engine.h"
34 #include "CLHEP/Random/engineIDulong.h"
35 #include <string.h> // for strcmp
36 #include <cstdlib> // for std::abs(int)
37 
38 namespace CLHEP {
39 
40 static const int MarkerLen = 64; // Enough room to hold a begin or end marker.
41 
42 std::string Hurd160Engine::name() const {return "Hurd160Engine";}
43 
44 // Number of instances with automatic seed selection
45 int Hurd160Engine::numEngines = 0;
46 
47 // Maximum index into the seed table
48 int Hurd160Engine::maxIndex = 215;
49 
50 static inline unsigned int f160(unsigned int a, unsigned int b, unsigned int c)
51 {
52  return ( ((b<<2) & 0x7c) | ((a<<2) & ~0x7c) | (a>>30) ) ^ ( (c<<1)|(c>>31) );
53 }
54 
57 {
58  int cycle = std::abs(int(numEngines/maxIndex));
59  int curIndex = std::abs(int(numEngines%maxIndex));
60  long mask = ((cycle & 0x007fffff) << 8);
61  long seedlist[2];
62  HepRandom::getTheTableSeeds( seedlist, curIndex );
63  seedlist[0] ^= mask;
64  seedlist[1] = 0;
65  setSeeds(seedlist, 0);
66  words[0] ^= 0x1324abcd; // To make unique vs long or two unsigned
67  if (words[0]==0) words[0] = 1; // ints in the constructor
68  ++numEngines;
69  for( int i=0; i < 100; ++i ) flat(); // warm-up just a bit
70 }
71 
74 {
75  is >> *this;
76 }
77 
80 {
81  long seedlist[2]={seed,0};
82  setSeeds(seedlist, 0);
83  words[0] ^= 0xa5482134; // To make unique vs default two unsigned
84  if (words[0]==0) words[0] = 1; // ints in the constructor
85  for( int i=0; i < 100; ++i ) flat(); // warm-up just a bit
86 }
87 
88 Hurd160Engine::Hurd160Engine( int rowIndex, int colIndex )
90 {
91  int cycle = std::abs(int(rowIndex/maxIndex));
92  int row = std::abs(int(rowIndex%maxIndex));
93  int col = colIndex & 0x1;
94  long mask = (( cycle & 0x000007ff ) << 20 );
95  long seedlist[2];
96  HepRandom::getTheTableSeeds( seedlist, row );
97  // NOTE: is that really the seed wanted (PGC) ??
98  seedlist[0] = (seedlist[col])^mask;
99  seedlist[1]= 0;
100  setSeeds(seedlist, 0);
101  for( int i=0; i < 100; ++i ) flat(); // warm-up just a bit
102 }
103 
105 
106 void Hurd160Engine::advance() {
107  register unsigned int W0, W1, W2, W3, W4;
108 
109  W0 = words[0];
110  W1 = words[1];
111  W2 = words[2];
112  W3 = words[3];
113  W4 = words[4];
114  W1 ^= W0; W0 = f160(W4, W3, W0);
115  W2 ^= W1; W1 = f160(W0, W4, W1);
116  W3 ^= W2; W2 = f160(W1, W0, W2);
117  W4 ^= W3; W3 = f160(W2, W1, W3);
118  W0 ^= W4; W4 = f160(W3, W2, W4);
119  words[0] = W0 & 0xffffffff;
120  words[1] = W1 & 0xffffffff;
121  words[2] = W2 & 0xffffffff;
122  words[3] = W3 & 0xffffffff;
123  words[4] = W4 & 0xffffffff;
124  wordIndex = 5;
125 
126 } // advance();
127 
129 
130  if( wordIndex <= 2 ) { // MF 9/15/98:
131  // skip word 0 and use two words per flat
132  advance();
133  }
134 
135  // LG 6/30/2010
136  // define the order of execution for --wordIndex
137  double x = words[--wordIndex] * twoToMinus_32() ; // most significant part
138  double y = (words[--wordIndex]>>11) * twoToMinus_53() + // fill in rest of bits
139  nearlyTwoToMinus_54(); // make sure non-zero
140  return x + y ;
141 }
142 
143 void Hurd160Engine::flatArray( const int size, double* vect ) {
144  for (int i = 0; i < size; ++i) {
145  vect[i] = flat();
146  }
147 }
148 
149 void Hurd160Engine::setSeed( long seed, int ) {
150  words[0] = (unsigned int)seed;
151  for (wordIndex = 1; wordIndex < 5; ++wordIndex) {
152  words[wordIndex] = 69607 * words[wordIndex-1] + 54329;
153  }
154 }
155 
156 void Hurd160Engine::setSeeds( const long* seeds, int ) {
157  setSeed( *seeds ? *seeds : 32767, 0 );
158  theSeeds = seeds;
159 }
160 
161 void Hurd160Engine::saveStatus( const char filename[] ) const {
162  std::ofstream outFile(filename, std::ios::out);
163  if( !outFile.bad() ) {
164  outFile << "Uvec\n";
165  std::vector<unsigned long> v = put();
166  #ifdef TRACE_IO
167  std::cout << "Result of v = put() is:\n";
168  #endif
169  for (unsigned int i=0; i<v.size(); ++i) {
170  outFile << v[i] << "\n";
171  #ifdef TRACE_IO
172  std::cout << v[i] << " ";
173  if (i%6==0) std::cout << "\n";
174  #endif
175  }
176  #ifdef TRACE_IO
177  std::cout << "\n";
178  #endif
179  }
180 #ifdef REMOVED
181  outFile << std::setprecision(20) << theSeed << " ";
182  outFile << wordIndex << " ";
183  for( int i = 0; i < 5; ++i ) {
184  outFile << words[i] << " ";
185  }
186  outFile << std::endl;
187 #endif
188 }
189 
190 void Hurd160Engine::restoreStatus( const char filename[] ) {
191  std::ifstream inFile(filename, std::ios::in);
192  if (!checkFile ( inFile, filename, engineName(), "restoreStatus" )) {
193  std::cerr << " -- Engine state remains unchanged\n";
194  return;
195  }
196  if ( possibleKeywordInput ( inFile, "Uvec", theSeed ) ) {
197  std::vector<unsigned long> v;
198  unsigned long xin;
199  for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
200  inFile >> xin;
201  #ifdef TRACE_IO
202  std::cout << "ivec = " << ivec << " xin = " << xin << " ";
203  if (ivec%3 == 0) std::cout << "\n";
204  #endif
205  if (!inFile) {
206  inFile.clear(std::ios::badbit | inFile.rdstate());
207  std::cerr << "\nHurd160Engine state (vector) description improper."
208  << "\nrestoreStatus has failed."
209  << "\nInput stream is probably mispositioned now." << std::endl;
210  return;
211  }
212  v.push_back(xin);
213  }
214  getState(v);
215  return;
216  }
217 
218  if( !inFile.bad() ) {
219 // inFile >> theSeed; removed -- encompased by possibleKeywordInput
220  inFile >> wordIndex;
221  for( int i = 0; i < 5; ++i ) {
222  inFile >> words[i];
223  }
224  }
225 }
226 
228  int pr = std::cout.precision(20);
229  std::cout << std::endl;
230  std::cout << "----------- Hurd engine status ----------" << std::endl;
231  std::cout << "Initial seed = " << theSeed << std::endl;
232  std::cout << "Current index = " << wordIndex << std::endl;
233  std::cout << "Current words = " << std::endl;
234  for( int i = 0; i < 5 ; ++i ) {
235  std::cout << " " << words[i] << std::endl;
236  }
237  std::cout << "------------------------------------------" << std::endl;
238  std::cout.precision(pr);
239 }
240 
241 Hurd160Engine::operator float() {
242  if( wordIndex <= 1 ) { // MF 9/15/98: skip word 0
243  advance();
244  }
245  return words[--wordIndex ] * twoToMinus_32();
246 }
247 
248 Hurd160Engine::operator unsigned int() {
249  if( wordIndex <= 1 ) { // MF 9/15/98: skip word 0
250  advance();
251  }
252  return words[--wordIndex];
253 }
254 
255 std::ostream& Hurd160Engine::put(std::ostream& os) const {
256  char beginMarker[] = "Hurd160Engine-begin";
257  os << beginMarker << "\nUvec\n";
258  std::vector<unsigned long> v = put();
259  for (unsigned int i=0; i<v.size(); ++i) {
260  os << v[i] << "\n";
261  }
262  return os;
263 #ifdef REMOVED
264  char endMarker[] = "Hurd160Engine-end";
265  int pr = os.precision(20);
266  os << " " << beginMarker << " ";
267  os << theSeed << " ";
268  os << wordIndex << " ";
269  for (int i = 0; i < 5; ++i) {
270  os << words[i] << "\n";
271  }
272  os << endMarker << "\n ";
273  os.precision(pr);
274  return os;
275 #endif
276 }
277 
278 std::vector<unsigned long> Hurd160Engine::put () const {
279  std::vector<unsigned long> v;
280  v.push_back (engineIDulong<Hurd160Engine>());
281  v.push_back(static_cast<unsigned long>(wordIndex));
282  for (int i = 0; i < 5; ++i) {
283  v.push_back(static_cast<unsigned long>(words[i]));
284  }
285  return v;
286 }
287 
288 
289 std::istream& Hurd160Engine::get(std::istream& is) {
290  char beginMarker [MarkerLen];
291  is >> std::ws;
292  is.width(MarkerLen); // causes the next read to the char* to be <=
293  // that many bytes, INCLUDING A TERMINATION \0
294  // (Stroustrup, section 21.3.2)
295  is >> beginMarker;
296  if (strcmp(beginMarker,"Hurd160Engine-begin")) {
297  is.clear(std::ios::badbit | is.rdstate());
298  std::cerr << "\nInput mispositioned or"
299  << "\nHurd160Engine state description missing or"
300  << "\nwrong engine type found." << std::endl;
301  return is;
302  }
303  return getState(is);
304 }
305 
306 std::string Hurd160Engine::beginTag ( ) {
307  return "Hurd160Engine-begin";
308 }
309 
310 std::istream& Hurd160Engine::getState(std::istream& is) {
311  if ( possibleKeywordInput ( is, "Uvec", theSeed ) ) {
312  std::vector<unsigned long> v;
313  unsigned long uu;
314  for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
315  is >> uu;
316  if (!is) {
317  is.clear(std::ios::badbit | is.rdstate());
318  std::cerr << "\nHurd160Engine state (vector) description improper."
319  << "\ngetState() has failed."
320  << "\nInput stream is probably mispositioned now." << std::endl;
321  return is;
322  }
323  v.push_back(uu);
324  }
325  getState(v);
326  return (is);
327  }
328 
329 // is >> theSeed; Removed, encompassed by possibleKeywordInput()
330 
331  char endMarker [MarkerLen];
332  is >> wordIndex;
333  for (int i = 0; i < 5; ++i) {
334  is >> words[i];
335  }
336  is >> std::ws;
337  is.width(MarkerLen);
338  is >> endMarker;
339  if (strcmp(endMarker,"Hurd160Engine-end")) {
340  is.clear(std::ios::badbit | is.rdstate());
341  std::cerr << "\nHurd160Engine state description incomplete."
342  << "\nInput stream is probably mispositioned now." << std::endl;
343  return is;
344  }
345  return is;
346 }
347 
348 
349 bool Hurd160Engine::get (const std::vector<unsigned long> & v) {
350  if ((v[0] & 0xffffffffUL) != engineIDulong<Hurd160Engine>()) {
351  std::cerr <<
352  "\nHurd160Engine get:state vector has wrong ID word - state unchanged\n";
353  return false;
354  }
355  return getState(v);
356 }
357 
358 bool Hurd160Engine::getState (const std::vector<unsigned long> & v) {
359  if (v.size() != VECTOR_STATE_SIZE ) {
360  std::cerr <<
361  "\nHurd160Engine get:state vector has wrong length - state unchanged\n";
362  return false;
363  }
364  wordIndex = v[1];
365  for (int i = 0; i < 5; ++i) {
366  words[i] = v[i+2];
367  }
368  return true;
369 }
370 
371 
372 } // namespace CLHEP
CLHEP::Hurd160Engine::get
virtual std::istream & get(std::istream &is)
Definition: Hurd160Engine.cc:289
a
@ a
Definition: testCategories.cc:125
CLHEP::HepRandomEngine
Definition: Matrix/CLHEP/Random/RandomEngine.h:55
CLHEP::Hurd160Engine::restoreStatus
void restoreStatus(const char filename[]="Hurd160Engine.conf")
Definition: Hurd160Engine.cc:190
CLHEP::Hurd160Engine::flatArray
void flatArray(const int size, double *vect)
Definition: Hurd160Engine.cc:143
CLHEP::HepRandomEngine::theSeed
long theSeed
Definition: Matrix/CLHEP/Random/RandomEngine.h:144
CLHEP::Hurd160Engine::name
std::string name() const
Definition: Hurd160Engine.cc:42
b
@ b
Definition: testCategories.cc:125
CLHEP::Hurd160Engine::put
std::vector< unsigned long > put() const
Definition: Hurd160Engine.cc:278
CLHEP::HepRandom::getTheTableSeeds
static void getTheTableSeeds(long *seeds, int index)
Definition: Random.cc:152
CLHEP::Hurd160Engine::showStatus
void showStatus() const
Definition: Hurd160Engine.cc:227
is
HepRotation and so forth isNear() norm2() rectify() static Rotation row1 row4(To avoid bloat in the code pulled in for programs which don 't use all these features, we split the implementation .cc files. Only isNear() goes into the original Rotation.cc) --------------------------------------- HepAxisAngle and HepEulerAngles classes --------------------------------------- These classes are very useful and simple structures for holding the result of a nice intuituve decomposition of a rotation there is no longer much content in the distinct ZOOM PhysicsVectors library The only content left in the library is the object files representing the various Exception objects When we build the CLHEP classes for the ZOOM we will set up so as to use ZOOM SpaceVector is(but we can disable namespace usage and most of our users do so at this point). What I do is leave Hep3Vector in the global namespace
CLHEP::Hurd160Engine::flat
double flat()
Definition: Hurd160Engine.cc:128
CLHEP::HepRandomEngine::theSeeds
const long * theSeeds
Definition: Matrix/CLHEP/Random/RandomEngine.h:145
CLHEP::Hurd160Engine::engineName
static std::string engineName()
Definition: Matrix/CLHEP/Random/Hurd160Engine.h:82
CLHEP::HepRandomEngine::nearlyTwoToMinus_54
static double nearlyTwoToMinus_54()
size
user code seldom needs to call this function directly ZMerrno whether or not they are still recorded ZMerrno size() Return the(integer) number of ZMthrow 'n exceptions currently recorded. 5) ZMerrno.clear() Set an internal counter to zero. This counter is available(see next function) to user code to track ZMthrow 'n exceptions that have occurred during any arbitrary time interval. 6) ZMerrno.countSinceCleared() Return the(integer) number of ZMthrow 'n exceptions that have been recorded via ZMerrno.write()
CLHEP::Hurd160Engine::beginTag
static std::string beginTag()
Definition: Hurd160Engine.cc:306
CLHEP::Hurd160Engine::saveStatus
void saveStatus(const char filename[]="Hurd160Engine.conf") const
Definition: Hurd160Engine.cc:161
CLHEP
Definition: ClhepVersion.h:13
CLHEP::HepRandomEngine::twoToMinus_53
static double twoToMinus_53()
v
they are gone ZOOM Features Discontinued The following features of the ZOOM package were felt to be extreme overkill These have been after checking that no existing user code was utilizing as in SpaceVector v
Definition: keyMergeIssues.doc:324
CLHEP::Hurd160Engine::VECTOR_STATE_SIZE
static const unsigned int VECTOR_STATE_SIZE
Definition: Matrix/CLHEP/Random/Hurd160Engine.h:90
CLHEP::possibleKeywordInput
bool possibleKeywordInput(IS &is, const std::string &key, T &t)
Definition: Matrix/CLHEP/Random/RandomEngine.h:168
seeds
Technical Maintenance Note for CLHEP Random Consequences of seeding JamesRandom with positive seed values greater than In the source code JamesRandom The usual way of seeding a generator is via the default which makes use of the table of seeds(with some trickery to ensure that the values won 't repeat after the table rows are exhausted). The trickery preserves the fact that sees are never negative(because the table values are never negative
CLHEP::Hurd160Engine::~Hurd160Engine
virtual ~Hurd160Engine()
Definition: Hurd160Engine.cc:104
CLHEP::HepRandomEngine::twoToMinus_32
static double twoToMinus_32()
i
long i
Definition: JamesRandomSeeding.txt:27
CLHEP::Hurd160Engine::setSeeds
void setSeeds(const long *seeds, int)
Definition: Hurd160Engine.cc:156
CLHEP::Hurd160Engine::getState
virtual std::istream & getState(std::istream &is)
Definition: Hurd160Engine.cc:310
CLHEP::Hurd160Engine::setSeed
void setSeed(long seed, int)
Definition: Hurd160Engine.cc:149
CLHEP::HepRandomEngine::checkFile
static bool checkFile(std::istream &file, const std::string &filename, const std::string &classname, const std::string &methodname)
Definition: RandomEngine.cc:46
in
it has advantages For I leave the ZMthrows in
Definition: keyMergeIssues.doc:62
CLHEP::Hurd160Engine::Hurd160Engine
Hurd160Engine()
Definition: Hurd160Engine.cc:55
x
any side effects of that construction would occur twice The semantics of throw x
Definition: whyZMthrowRethrows.txt:37