CLHEP 2.0.4.7 Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // $Id: Ranlux64Engine.cc,v 1.4.4.2.2.2 2009/07/05 00:02:30 garren Exp $ 00002 // -*- C++ -*- 00003 // 00004 // ----------------------------------------------------------------------- 00005 // HEP Random 00006 // --- Ranlux64Engine --- 00007 // class implementation file 00008 // ----------------------------------------------------------------------- 00009 // A double-precision implementation of the RanluxEngine generator as 00010 // decsribed by the notes of the original ranlux author (Martin Luscher) 00011 // 00012 // See the note by Martin Luscher, December 1997, entitiled 00013 // Double-precision implementation of the random number generator ranlux 00014 // 00015 // ======================================================================= 00016 // Ken Smith - Initial draft: 14th Jul 1998 00017 // - Removed pow() from flat method 14th Jul 1998 00018 // - Added conversion operators: 6th Aug 1998 00019 // 00020 // Mark Fischler The following were modified mostly to make the routine 00021 // exactly match the Luscher algorithm in generating 48-bit 00022 // randoms: 00023 // 9/9/98 - Substantial changes in what used to be flat() to match 00024 // algorithm in Luscher's ranlxd.c 00025 // - Added update() method for 12 numbers, making flat() trivial 00026 // - Added advance() method to hold the unrolled loop for update 00027 // - Distinction between three forms of seeding such that it 00028 // is impossible to get same sequence from different forms - 00029 // done by discarding some fraction of one macro cycle which 00030 // is different for the three cases 00031 // - Change the misnomer "seed_table" to the more accurate 00032 // "randoms" 00033 // - Removed the no longer needed count12, i_lag, j_lag, etc. 00034 // - Corrected seed procedure which had been filling bits past 00035 // 2^-48. This actually was very bad, invalidating the 00036 // number theory behind the proof that ranlxd is good. 00037 // - Addition of 2**(-49) to generated number to prevent zero 00038 // from being returned; this does not affect the sequence 00039 // itself. 00040 // - Corrected ecu seeding, which had been supplying only 00041 // numbers less than 1/2. This is probably moot. 00042 // 9/15/98 - Modified use of the various exponents of 2 00043 // to avoid per-instance space overhead. Note that these 00044 // are initialized in setSeed, which EVERY constructor 00045 // must invoke. 00046 // J. Marraffino - Remove dependence on hepString class 13 May 1999 00047 // M. Fischler - In restore, checkFile for file not found 03 Dec 2004 00048 // M. Fischler - put get Methods for distrib instance save/restore 12/8/04 00049 // M. Fischler - split get() into tag validation and 00050 // getState() for anonymous restores 12/27/04 00051 // M. Fischler - put/get for vectors of ulongs 3/14/05 00052 // M. Fischler - State-saving using only ints, for portability 4/12/05 00053 // 00054 // ======================================================================= 00055 00056 #include "CLHEP/Random/defs.h" 00057 #include "CLHEP/Random/Random.h" 00058 #include "CLHEP/Random/Ranlux64Engine.h" 00059 #include "CLHEP/Random/engineIDulong.h" 00060 #include "CLHEP/Random/DoubConv.hh" 00061 #include <string.h> 00062 #include <cmath> // for ldexp() and abs() 00063 #include <stdlib.h> // for abs(int) 00064 00065 using namespace std; 00066 00067 namespace CLHEP { 00068 00069 static const int MarkerLen = 64; // Enough room to hold a begin or end marker. 00070 00071 00072 // Number of instances with automatic seed selection 00073 int Ranlux64Engine::numEngines = 0; 00074 00075 // Maximum index into the seed table 00076 int Ranlux64Engine::maxIndex = 215; 00077 00078 double Ranlux64Engine::twoToMinus_32; 00079 double Ranlux64Engine::twoToMinus_48; 00080 double Ranlux64Engine::twoToMinus_49; 00081 00082 std::string Ranlux64Engine::name() const {return "Ranlux64Engine";} 00083 00084 Ranlux64Engine::Ranlux64Engine() 00085 { 00086 luxury = 1; 00087 int cycle = abs(int(numEngines/maxIndex)); 00088 int curIndex = abs(int(numEngines%maxIndex)); 00089 numEngines +=1; 00090 long mask = ((cycle & 0x007fffff) << 8); 00091 long seedlist[2]; 00092 HepRandom::getTheTableSeeds( seedlist, curIndex ); 00093 seedlist[0] ^= mask; 00094 seedlist[1] = 0; 00095 00096 setSeeds(seedlist, luxury); 00097 advance ( 8 ); // Discard some iterations and ensure that 00098 // this sequence won't match one where seeds 00099 // were provided. 00100 } 00101 00102 Ranlux64Engine::Ranlux64Engine(long seed, int lux) 00103 { 00104 luxury = lux; 00105 long seedlist[2]={seed,0}; 00106 setSeeds(seedlist, lux); 00107 advance ( 2*lux + 1 ); // Discard some iterations to use a different 00108 // point in the sequence. 00109 } 00110 00111 Ranlux64Engine::Ranlux64Engine(int rowIndex, int colIndex, int lux) 00112 { 00113 luxury = lux; 00114 int cycle = abs(int(rowIndex/maxIndex)); 00115 int row = abs(int(rowIndex%maxIndex)); 00116 // int col = abs(int(colIndex%2)); // But this is never used! 00117 long mask = (( cycle & 0x000007ff ) << 20 ); 00118 long seedlist[2]; 00119 HepRandom::getTheTableSeeds( seedlist, row ); 00120 seedlist[0] ^= mask; 00121 seedlist[1]= 0; 00122 setSeeds(seedlist, lux); 00123 } 00124 00125 Ranlux64Engine::Ranlux64Engine( std::istream& is ) 00126 { 00127 is >> *this; 00128 } 00129 00130 Ranlux64Engine::~Ranlux64Engine() {} 00131 00132 Ranlux64Engine::Ranlux64Engine(const Ranlux64Engine &p) 00133 { 00134 *this = p; 00135 } 00136 00137 Ranlux64Engine & Ranlux64Engine::operator=( const Ranlux64Engine &p ) { 00138 if (this != &p ) { 00139 theSeed = p.theSeed; 00140 theSeeds = p.theSeeds; 00141 for (int i=0; i<12; ++i) { 00142 randoms[i] = p.randoms[i]; 00143 } 00144 pDiscard = p.pDiscard; 00145 pDozens = p.pDozens; 00146 endIters = p.endIters; 00147 luxury = p.luxury; 00148 index = p.index; 00149 carry = p.carry; 00150 } 00151 return *this; 00152 } 00153 00154 00155 double Ranlux64Engine::flat() { 00156 // Luscher improves the speed by computing several numbers in a shot, 00157 // in a manner similar to that of the Tausworth in DualRand or the Hurd 00158 // engines. Thus, the real work is done in update(). Here we merely ensure 00159 // that zero, which the algorithm can produce, is never returned by flat(). 00160 00161 if (index <= 0) update(); 00162 return randoms[--index] + twoToMinus_49; 00163 } 00164 00165 void Ranlux64Engine::update() { 00166 // Update the stash of twelve random numbers. 00167 // When this routione is entered, index is always 0. The randoms 00168 // contains the last 12 numbers in the sequents: s[0] is x[a+11], 00169 // s[1] is x[a+10] ... and s[11] is x[a] for some a. Carry contains 00170 // the last carry value (c[a+11]). 00171 // 00172 // The recursion relation (3) in Luscher's note says 00173 // delta[n] = x[n-s] = x[n-r] -c[n-1] or for n=a+12, 00174 // delta[a+12] = x[a+7] - x[a] -c[a+11] where we use r=12, s=5 per eqn. (7) 00175 // This reduces to 00176 // s[11] = s[4] - s[11] - carry. 00177 // The next number similarly will be given by s[10] = s[3] - s[10] - carry, 00178 // and so forth until s[0] is filled. 00179 // 00180 // However, we need to skip 397, 202 or 109 numbers - these are not divisible 00181 // by 12 - to "fare well in the spectral test". 00182 00183 advance(pDozens); 00184 00185 // Since we wish at the end to have the 12 last numbers in the order of 00186 // s[11] first, till s[0] last, we will have to do 1, 10, or 1 iterations 00187 // and then re-arrange to place to get the oldest one in s[11]. 00188 // Generically, this will imply re-arranging the s array at the end, 00189 // but we can treat the special case of endIters = 1 separately for superior 00190 // efficiency in the cases of levels 0 and 2. 00191 00192 register double y1; 00193 00194 if ( endIters == 1 ) { // Luxury levels 0 and 2 will go here 00195 y1 = randoms[ 4] - randoms[11] - carry; 00196 if ( y1 < 0.0 ) { 00197 y1 += 1.0; 00198 carry = twoToMinus_48; 00199 } else { 00200 carry = 0.0; 00201 } 00202 randoms[11] = randoms[10]; 00203 randoms[10] = randoms[ 9]; 00204 randoms[ 9] = randoms[ 8]; 00205 randoms[ 8] = randoms[ 7]; 00206 randoms[ 7] = randoms[ 6]; 00207 randoms[ 6] = randoms[ 5]; 00208 randoms[ 5] = randoms[ 4]; 00209 randoms[ 4] = randoms[ 3]; 00210 randoms[ 3] = randoms[ 2]; 00211 randoms[ 2] = randoms[ 1]; 00212 randoms[ 1] = randoms[ 0]; 00213 randoms[ 0] = y1; 00214 00215 } else { 00216 00217 int m, nr, ns; 00218 for ( m = 0, nr = 11, ns = 4; m < endIters; ++m, --nr ) { 00219 y1 = randoms [ns] - randoms[nr] - carry; 00220 if ( y1 < 0.0 ) { 00221 y1 += 1.0; 00222 carry = twoToMinus_48; 00223 } else { 00224 carry = 0.0; 00225 } 00226 randoms[nr] = y1; 00227 --ns; 00228 if ( ns < 0 ) { 00229 ns = 11; 00230 } 00231 } // loop on m 00232 00233 double temp[12]; 00234 for (m=0; m<12; m++) { 00235 temp[m]=randoms[m]; 00236 } 00237 00238 ns = 11 - endIters; 00239 for (m=11; m>=0; --m) { 00240 randoms[m] = temp[ns]; 00241 --ns; 00242 if ( ns < 0 ) { 00243 ns = 11; 00244 } 00245 } 00246 00247 } 00248 00249 // Now when we return, there are 12 fresh usable numbers in s[11] ... s[0] 00250 00251 index = 11; 00252 00253 } // update() 00254 00255 void Ranlux64Engine::advance(int dozens) { 00256 00257 register double y1, y2, y3; 00258 register double cValue = twoToMinus_48; 00259 register double zero = 0.0; 00260 register double one = 1.0; 00261 00262 // Technical note: We use Luscher's trick to only do the 00263 // carry subtraction when we really have to. Like him, we use 00264 // three registers instead of two so that we avoid sequences 00265 // like storing y1 then immediately replacing its value: 00266 // some architectures lose time when this is done. 00267 00268 // Luscher's ranlxd.c fills the stash going 00269 // upward. We fill it downward to save a bit of time in the 00270 // flat() routine at no cost later. This means that while 00271 // Luscher's ir is jr+5, our n-r is (n-s)-5. (Note that 00272 // though ranlxd.c initializes ir and jr to 11 and 7, ir as 00273 // used is 5 more than jr because update is entered after 00274 // incrementing ir.) 00275 // 00276 00277 // I have CAREFULLY checked that the algorithms do match 00278 // in all details. 00279 00280 int k; 00281 for ( k = dozens; k > 0; --k ) { 00282 00283 y1 = randoms[ 4] - randoms[11] - carry; 00284 00285 y2 = randoms[ 3] - randoms[10]; 00286 if ( y1 < zero ) { 00287 y1 += one; 00288 y2 -= cValue; 00289 } 00290 randoms[11] = y1; 00291 00292 y3 = randoms[ 2] - randoms[ 9]; 00293 if ( y2 < zero ) { 00294 y2 += one; 00295 y3 -= cValue; 00296 } 00297 randoms[10] = y2; 00298 00299 y1 = randoms[ 1] - randoms[ 8]; 00300 if ( y3 < zero ) { 00301 y3 += one; 00302 y1 -= cValue; 00303 } 00304 randoms[ 9] = y3; 00305 00306 y2 = randoms[ 0] - randoms[ 7]; 00307 if ( y1 < zero ) { 00308 y1 += one; 00309 y2 -= cValue; 00310 } 00311 randoms[ 8] = y1; 00312 00313 y3 = randoms[11] - randoms[ 6]; 00314 if ( y2 < zero ) { 00315 y2 += one; 00316 y3 -= cValue; 00317 } 00318 randoms[ 7] = y2; 00319 00320 y1 = randoms[10] - randoms[ 5]; 00321 if ( y3 < zero ) { 00322 y3 += one; 00323 y1 -= cValue; 00324 } 00325 randoms[ 6] = y3; 00326 00327 y2 = randoms[ 9] - randoms[ 4]; 00328 if ( y1 < zero ) { 00329 y1 += one; 00330 y2 -= cValue; 00331 } 00332 randoms[ 5] = y1; 00333 00334 y3 = randoms[ 8] - randoms[ 3]; 00335 if ( y2 < zero ) { 00336 y2 += one; 00337 y3 -= cValue; 00338 } 00339 randoms[ 4] = y2; 00340 00341 y1 = randoms[ 7] - randoms[ 2]; 00342 if ( y3 < zero ) { 00343 y3 += one; 00344 y1 -= cValue; 00345 } 00346 randoms[ 3] = y3; 00347 00348 y2 = randoms[ 6] - randoms[ 1]; 00349 if ( y1 < zero ) { 00350 y1 += one; 00351 y2 -= cValue; 00352 } 00353 randoms[ 2] = y1; 00354 00355 y3 = randoms[ 5] - randoms[ 0]; 00356 if ( y2 < zero ) { 00357 y2 += one; 00358 y3 -= cValue; 00359 } 00360 randoms[ 1] = y2; 00361 00362 if ( y3 < zero ) { 00363 y3 += one; 00364 carry = cValue; 00365 } 00366 randoms[ 0] = y3; 00367 00368 } // End of major k loop doing 12 numbers at each cycle 00369 00370 } // advance(dozens) 00371 00372 void Ranlux64Engine::flatArray(const int size, double* vect) { 00373 for( int i=0; i < size; ++i ) { 00374 vect[i] = flat(); 00375 } 00376 } 00377 00378 void Ranlux64Engine::setSeed(long seed, int lux) { 00379 00380 // The initialization is carried out using a Multiplicative 00381 // Congruential generator using formula constants of L'Ecuyer 00382 // as described in "A review of pseudorandom number generators" 00383 // (Fred James) published in Computer Physics Communications 60 (1990) 00384 // pages 329-344 00385 00386 twoToMinus_32 = ldexp (1.0, -32); 00387 twoToMinus_48 = ldexp (1.0, -48); 00388 twoToMinus_49 = ldexp (1.0, -49); 00389 00390 const int ecuyer_a(53668); 00391 const int ecuyer_b(40014); 00392 const int ecuyer_c(12211); 00393 const int ecuyer_d(2147483563); 00394 00395 const int lux_levels[3] = {109, 202, 397}; 00396 theSeed = seed; 00397 00398 if( (lux > 2)||(lux < 0) ){ 00399 pDiscard = (lux >= 12) ? (lux-12) : lux_levels[1]; 00400 }else{ 00401 pDiscard = lux_levels[luxury]; 00402 } 00403 pDozens = pDiscard / 12; 00404 endIters = pDiscard % 12; 00405 00406 long init_table[24]; 00407 long next_seed = seed; 00408 long k_multiple; 00409 int i; 00410 00411 for(i = 0;i != 24;i++){ 00412 k_multiple = next_seed / ecuyer_a; 00413 next_seed = ecuyer_b * (next_seed - k_multiple * ecuyer_a) 00414 - k_multiple * ecuyer_c; 00415 if(next_seed < 0) { 00416 next_seed += ecuyer_d; 00417 } 00418 init_table[i] = next_seed & 0xffffffff; 00419 } 00420 00421 for(i = 0;i < 12; i++){ 00422 randoms[i] = (init_table[2*i ] ) * 2.0 * twoToMinus_32 + 00423 (init_table[2*i+1] >> 15) * twoToMinus_48; 00424 } 00425 00426 carry = 0.0; 00427 if ( randoms[11] == 0. ) carry = twoToMinus_48; 00428 index = 11; 00429 00430 } // setSeed() 00431 00432 void Ranlux64Engine::setSeeds(const long * seeds, int lux) { 00433 // old code only uses the first long in seeds 00434 // setSeed( *seeds ? *seeds : 32767, lux ); 00435 // theSeeds = seeds; 00436 00437 // using code from Ranlux - even those are 32bit seeds, 00438 // that is good enough to completely differentiate the sequences 00439 00440 twoToMinus_32 = ldexp (1.0, -32); 00441 twoToMinus_48 = ldexp (1.0, -48); 00442 twoToMinus_49 = ldexp (1.0, -49); 00443 00444 const int ecuyer_a = 53668; 00445 const int ecuyer_b = 40014; 00446 const int ecuyer_c = 12211; 00447 const int ecuyer_d = 2147483563; 00448 00449 const int lux_levels[3] = {109, 202, 397}; 00450 const long *seedptr; 00451 00452 theSeeds = seeds; 00453 seedptr = seeds; 00454 00455 if(seeds == 0){ 00456 setSeed(theSeed,lux); 00457 theSeeds = &theSeed; 00458 return; 00459 } 00460 00461 theSeed = *seeds; 00462 00463 // number of additional random numbers that need to be 'thrown away' 00464 // every 24 numbers is set using luxury level variable. 00465 00466 if( (lux > 2)||(lux < 0) ){ 00467 pDiscard = (lux >= 12) ? (lux-12) : lux_levels[1]; 00468 }else{ 00469 pDiscard = lux_levels[luxury]; 00470 } 00471 pDozens = pDiscard / 12; 00472 endIters = pDiscard % 12; 00473 00474 long init_table[24]; 00475 long next_seed = *seeds; 00476 long k_multiple; 00477 int i; 00478 00479 for( i = 0;(i != 24)&&(*seedptr != 0);i++){ 00480 init_table[i] = *seedptr & 0xffffffff; 00481 seedptr++; 00482 } 00483 00484 if(i != 24){ 00485 next_seed = init_table[i-1]; 00486 for(;i != 24;i++){ 00487 k_multiple = next_seed / ecuyer_a; 00488 next_seed = ecuyer_b * (next_seed - k_multiple * ecuyer_a) 00489 - k_multiple * ecuyer_c; 00490 if(next_seed < 0) { 00491 next_seed += ecuyer_d; 00492 } 00493 init_table[i] = next_seed & 0xffffffff; 00494 } 00495 } 00496 00497 for(i = 0;i < 12; i++){ 00498 randoms[i] = (init_table[2*i ] ) * 2.0 * twoToMinus_32 + 00499 (init_table[2*i+1] >> 15) * twoToMinus_48; 00500 } 00501 00502 carry = 0.0; 00503 if ( randoms[11] == 0. ) carry = twoToMinus_48; 00504 index = 11; 00505 00506 } 00507 00508 void Ranlux64Engine::saveStatus( const char filename[] ) const 00509 { 00510 std::ofstream outFile( filename, std::ios::out ) ; 00511 if (!outFile.bad()) { 00512 outFile << "Uvec\n"; 00513 std::vector<unsigned long> v = put(); 00514 #ifdef TRACE_IO 00515 std::cout << "Result of v = put() is:\n"; 00516 #endif 00517 for (unsigned int i=0; i<v.size(); ++i) { 00518 outFile << v[i] << "\n"; 00519 #ifdef TRACE_IO 00520 std::cout << v[i] << " "; 00521 if (i%6==0) std::cout << "\n"; 00522 #endif 00523 } 00524 #ifdef TRACE_IO 00525 std::cout << "\n"; 00526 #endif 00527 } 00528 #ifdef REMOVED 00529 if (!outFile.bad()) { 00530 outFile << theSeed << std::endl; 00531 for (int i=0; i<12; ++i) { 00532 outFile <<std::setprecision(20) << randoms[i] << std::endl; 00533 } 00534 outFile << std::setprecision(20) << carry << " " << index << std::endl; 00535 outFile << luxury << " " << pDiscard << std::endl; 00536 } 00537 #endif 00538 } 00539 00540 void Ranlux64Engine::restoreStatus( const char filename[] ) 00541 { 00542 std::ifstream inFile( filename, std::ios::in); 00543 if (!checkFile ( inFile, filename, engineName(), "restoreStatus" )) { 00544 std::cerr << " -- Engine state remains unchanged\n"; 00545 return; 00546 } 00547 if ( possibleKeywordInput ( inFile, "Uvec", theSeed ) ) { 00548 std::vector<unsigned long> v; 00549 unsigned long xin; 00550 for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) { 00551 inFile >> xin; 00552 #ifdef TRACE_IO 00553 std::cout << "ivec = " << ivec << " xin = " << xin << " "; 00554 if (ivec%3 == 0) std::cout << "\n"; 00555 #endif 00556 if (!inFile) { 00557 inFile.clear(std::ios::badbit | inFile.rdstate()); 00558 std::cerr << "\nJamesRandom state (vector) description improper." 00559 << "\nrestoreStatus has failed." 00560 << "\nInput stream is probably mispositioned now." << std::endl; 00561 return; 00562 } 00563 v.push_back(xin); 00564 } 00565 getState(v); 00566 return; 00567 } 00568 00569 if (!inFile.bad() && !inFile.eof()) { 00570 // inFile >> theSeed; removed -- encompased by possibleKeywordInput 00571 for (int i=0; i<12; ++i) { 00572 inFile >> randoms[i]; 00573 } 00574 inFile >> carry; inFile >> index; 00575 inFile >> luxury; inFile >> pDiscard; 00576 pDozens = pDiscard / 12; 00577 endIters = pDiscard % 12; 00578 } 00579 } 00580 00581 void Ranlux64Engine::showStatus() const 00582 { 00583 std::cout << std::endl; 00584 std::cout << "--------- Ranlux engine status ---------" << std::endl; 00585 std::cout << " Initial seed = " << theSeed << std::endl; 00586 std::cout << " randoms[] = "; 00587 for (int i=0; i<12; ++i) { 00588 std::cout << randoms[i] << std::endl; 00589 } 00590 std::cout << std::endl; 00591 std::cout << " carry = " << carry << ", index = " << index << std::endl; 00592 std::cout << " luxury = " << luxury << " pDiscard = " 00593 << pDiscard << std::endl; 00594 std::cout << "----------------------------------------" << std::endl; 00595 } 00596 00597 std::ostream & Ranlux64Engine::put( std::ostream& os ) const 00598 { 00599 char beginMarker[] = "Ranlux64Engine-begin"; 00600 os << beginMarker << "\nUvec\n"; 00601 std::vector<unsigned long> v = put(); 00602 for (unsigned int i=0; i<v.size(); ++i) { 00603 os << v[i] << "\n"; 00604 } 00605 return os; 00606 #ifdef REMOVED 00607 char endMarker[] = "Ranlux64Engine-end"; 00608 int pr = os.precision(20); 00609 os << " " << beginMarker << " "; 00610 os << theSeed << " "; 00611 for (int i=0; i<12; ++i) { 00612 os << randoms[i] << std::endl; 00613 } 00614 os << carry << " " << index << " "; 00615 os << luxury << " " << pDiscard << "\n"; 00616 os << endMarker << " "; 00617 os.precision(pr); 00618 return os; 00619 #endif 00620 } 00621 00622 std::vector<unsigned long> Ranlux64Engine::put () const { 00623 std::vector<unsigned long> v; 00624 v.push_back (engineIDulong<Ranlux64Engine>()); 00625 std::vector<unsigned long> t; 00626 for (int i=0; i<12; ++i) { 00627 t = DoubConv::dto2longs(randoms[i]); 00628 v.push_back(t[0]); v.push_back(t[1]); 00629 } 00630 t = DoubConv::dto2longs(carry); 00631 v.push_back(t[0]); v.push_back(t[1]); 00632 v.push_back(static_cast<unsigned long>(index)); 00633 v.push_back(static_cast<unsigned long>(luxury)); 00634 v.push_back(static_cast<unsigned long>(pDiscard)); 00635 return v; 00636 } 00637 00638 std::istream & Ranlux64Engine::get ( std::istream& is ) 00639 { 00640 char beginMarker [MarkerLen]; 00641 is >> std::ws; 00642 is.width(MarkerLen); // causes the next read to the char* to be <= 00643 // that many bytes, INCLUDING A TERMINATION \0 00644 // (Stroustrup, section 21.3.2) 00645 is >> beginMarker; 00646 if (strcmp(beginMarker,"Ranlux64Engine-begin")) { 00647 is.clear(std::ios::badbit | is.rdstate()); 00648 std::cerr << "\nInput stream mispositioned or" 00649 << "\nRanlux64Engine state description missing or" 00650 << "\nwrong engine type found." << std::endl; 00651 return is; 00652 } 00653 return getState(is); 00654 } 00655 00656 std::string Ranlux64Engine::beginTag ( ) { 00657 return "Ranlux64Engine-begin"; 00658 } 00659 00660 std::istream & Ranlux64Engine::getState ( std::istream& is ) 00661 { 00662 if ( possibleKeywordInput ( is, "Uvec", theSeed ) ) { 00663 std::vector<unsigned long> v; 00664 unsigned long uu; 00665 for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) { 00666 is >> uu; 00667 if (!is) { 00668 is.clear(std::ios::badbit | is.rdstate()); 00669 std::cerr << "\nRanlux64Engine state (vector) description improper." 00670 << "\ngetState() has failed." 00671 << "\nInput stream is probably mispositioned now." << std::endl; 00672 return is; 00673 } 00674 v.push_back(uu); 00675 } 00676 getState(v); 00677 return (is); 00678 } 00679 00680 // is >> theSeed; Removed, encompassed by possibleKeywordInput() 00681 00682 char endMarker [MarkerLen]; 00683 for (int i=0; i<12; ++i) { 00684 is >> randoms[i]; 00685 } 00686 is >> carry; is >> index; 00687 is >> luxury; is >> pDiscard; 00688 pDozens = pDiscard / 12; 00689 endIters = pDiscard % 12; 00690 is >> std::ws; 00691 is.width(MarkerLen); 00692 is >> endMarker; 00693 if (strcmp(endMarker,"Ranlux64Engine-end")) { 00694 is.clear(std::ios::badbit | is.rdstate()); 00695 std::cerr << "\nRanlux64Engine state description incomplete." 00696 << "\nInput stream is probably mispositioned now." << std::endl; 00697 return is; 00698 } 00699 return is; 00700 } 00701 00702 bool Ranlux64Engine::get (const std::vector<unsigned long> & v) { 00703 if ((v[0] & 0xffffffffUL) != engineIDulong<Ranlux64Engine>()) { 00704 std::cerr << 00705 "\nRanlux64Engine get:state vector has wrong ID word - state unchanged\n"; 00706 return false; 00707 } 00708 return getState(v); 00709 } 00710 00711 bool Ranlux64Engine::getState (const std::vector<unsigned long> & v) { 00712 if (v.size() != VECTOR_STATE_SIZE ) { 00713 std::cerr << 00714 "\nRanlux64Engine get:state vector has wrong length - state unchanged\n"; 00715 return false; 00716 } 00717 std::vector<unsigned long> t(2); 00718 for (int i=0; i<12; ++i) { 00719 t[0] = v[2*i+1]; t[1] = v[2*i+2]; 00720 randoms[i] = DoubConv::longs2double(t); 00721 } 00722 t[0] = v[25]; t[1] = v[26]; 00723 carry = DoubConv::longs2double(t); 00724 index = v[27]; 00725 luxury = v[28]; 00726 pDiscard = v[29]; 00727 return true; 00728 } 00729 00730 } // namespace CLHEP