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 
72 Hurd160Engine::Hurd160Engine( std::istream& is )
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
static double twoToMinus_32()
static double twoToMinus_53()
static double nearlyTwoToMinus_54()
static bool checkFile(std::istream &file, const std::string &filename, const std::string &classname, const std::string &methodname)
Definition: RandomEngine.cc:46
static void getTheTableSeeds(long *seeds, int index)
Definition: Random.cc:152
virtual std::istream & get(std::istream &is)
void setSeeds(const long *seeds, int)
void restoreStatus(const char filename[]="Hurd160Engine.conf")
void showStatus() const
virtual std::istream & getState(std::istream &is)
std::vector< unsigned long > put() const
void flatArray(const int size, double *vect)
void setSeed(long seed, int)
static std::string beginTag()
static const unsigned int VECTOR_STATE_SIZE
void saveStatus(const char filename[]="Hurd160Engine.conf") const
std::string name() const
bool possibleKeywordInput(IS &is, const std::string &key, T &t)
@ b
@ a