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

RanshiEngine.cc
Go to the documentation of this file.
00001 // $Id: RanshiEngine.cc,v 1.6 2010/06/16 17:24:53 garren Exp $
00002 // -*- C++ -*-
00003 //
00004 // -----------------------------------------------------------------------
00005 //                           HEP Random
00006 //                      --- RanshiEngine ---
00007 //                    class implementation file
00008 // -----------------------------------------------------------------------
00009 //
00010 // This algorithm implements the random number generator as proposed by 
00011 // "F. Gutbrod, Comp. Phys. Comm. 87 (1995) 291-306". 
00012 //
00013 // =======================================================================
00014 // Ken Smith      - Created:                                 9th June 1998
00015 //                - Removed pow() from flat method:          21st Jul 1998
00016 //                - Added conversion operators:               6th Aug 1998
00017 // J. Marraffino  - Added some explicit casts to deal with
00018 //                  machines where sizeof(int) != sizeof(long) 22 Aug 1998
00019 // M. Fischler    - Modified constructors taking seeds to not
00020 //                  depend on numEngines (same seeds should
00021 //                  produce same sequences).  Default still
00022 //                  depends on numEngines.                      16 Sep 1998
00023 //                - Modified use of the various exponents of 2
00024 //                  to avoid per-instance space overhead and
00025 //                  correct the rounding procedure              16 Sep 1998
00026 // J. Marraffino  - Remove dependence on hepString class        13 May 1999
00027 // M. Fischler    - In restore, checkFile for file not found    03 Dec 2004
00028 // M. Fischler    - Methods for instance save/restore            12/8/04    
00029 // M. Fischler    - split get() into tag validation and 
00030 //                  getState() for anonymous restores           12/27/04    
00031 // M. Fischler    - State-saving using only ints, for portability 4/12/05
00032 //
00033 // =======================================================================
00034 
00035 #include "CLHEP/Random/defs.h"
00036 #include "CLHEP/Random/RanshiEngine.h"
00037 #include "CLHEP/Random/engineIDulong.h"
00038 #include <string.h>     // for strcmp
00039 
00040 using namespace std;
00041 
00042 namespace CLHEP {
00043 
00044 static const int MarkerLen = 64; // Enough room to hold a begin or end marker. 
00045 
00046 std::string RanshiEngine::name() const {return "RanshiEngine";}
00047 
00048 // Number of instances with automatic seed selection
00049 int RanshiEngine::numEngines = 0;
00050 
00051 RanshiEngine::RanshiEngine()
00052 : HepRandomEngine(),
00053   halfBuff(0), numFlats(0) 
00054 {
00055   int i = 0;
00056   while (i < numBuff) {    
00057     buffer[i] = (unsigned int)(numEngines+19780503L*(i+1));
00058     ++i;
00059   }
00060   theSeed = numEngines+19780503L*++i;
00061   redSpin = (unsigned int)(theSeed & 0xffffffff);
00062   ++numEngines;
00063   for( i = 0; i < 10000; ++i) flat();  // Warm-up by running thorugh 10000 nums
00064 }
00065 
00066 RanshiEngine::RanshiEngine(std::istream& is)
00067 : HepRandomEngine(),
00068   halfBuff(0), numFlats(0) 
00069 {
00070   is >> *this;
00071 }
00072 
00073 RanshiEngine::RanshiEngine(long seed)
00074 : HepRandomEngine(),
00075   halfBuff(0), numFlats(0) 
00076 {
00077   for (int i = 0; i < numBuff; ++i) {
00078     buffer[i] = (unsigned int)seed&0xffffffff;
00079   }
00080   theSeed = seed;
00081   redSpin = (unsigned int)(theSeed & 0xffffffff);
00082   int j;
00083   for (j = 0; j < numBuff*20; ++j) {      // "warm-up" for engine to hit
00084     flat();                               //  every ball on average 20X.
00085   }
00086 }
00087 
00088 RanshiEngine::RanshiEngine(int rowIndex, int colIndex)
00089 : HepRandomEngine(),
00090   halfBuff(0), numFlats(0) 
00091 {
00092   int i = 0;
00093   while( i < numBuff ) {
00094     buffer[i] = (unsigned int)((rowIndex + (i+1)*(colIndex+8))&0xffffffff);
00095     ++i;
00096   }
00097   theSeed = rowIndex;
00098   redSpin = colIndex & 0xffffffff;
00099   for( i = 0; i < 100; ++i) flat();    // Warm-up by running thorugh 100 nums
00100 }
00101 
00102 RanshiEngine::~RanshiEngine() { }
00103 
00104 double RanshiEngine::flat() {
00105   unsigned int redAngle = (((numBuff/2) - 1) & redSpin) + halfBuff;
00106   unsigned int blkSpin     = buffer[redAngle] & 0xffffffff;
00107   unsigned int boostResult = blkSpin ^ redSpin;
00108 
00109   buffer[redAngle] = ((blkSpin << 17) | (blkSpin >> (32-17))) ^ redSpin;
00110   
00111   redSpin  = (blkSpin + numFlats++) & 0xffffffff;
00112   halfBuff = numBuff/2 - halfBuff;
00113   
00114   return ( blkSpin * twoToMinus_32() +            // most significant part
00115            (boostResult>>11) * twoToMinus_53() +  // fill in remaining bits
00116            nearlyTwoToMinus_54());              // non-zero
00117 }
00118 
00119 void RanshiEngine::flatArray(const int size, double* vect) {
00120   for (int i = 0; i < size; ++i) {
00121     vect[i] = flat();
00122   }
00123 }
00124 
00125 void RanshiEngine::setSeed(long seed, int) {
00126   *this = RanshiEngine(seed); 
00127 }
00128 
00129 void RanshiEngine::setSeeds(const long* seeds, int) {
00130   if (*seeds) {
00131     int i = 0;
00132     while (seeds[i] && i < numBuff) {
00133       buffer[i] = (unsigned int)seeds[i];
00134       ++i;
00135     }
00136     while (i < numBuff) {
00137       buffer[i] = buffer[i-1];
00138       ++i;
00139     }
00140     theSeed = seeds[0];
00141     redSpin = (unsigned int)theSeed;
00142   }
00143   theSeeds = seeds;
00144 }
00145      
00146 void RanshiEngine::saveStatus(const char filename[]) const {
00147   std::ofstream outFile(filename, std::ios::out);
00148   if (!outFile.bad()) {
00149     outFile << "Uvec\n";
00150     std::vector<unsigned long> v = put();
00151                      #ifdef TRACE_IO
00152                          std::cout << "Result of v = put() is:\n"; 
00153                      #endif
00154     for (unsigned int i=0; i<v.size(); ++i) {
00155       outFile << v[i] << "\n";
00156                      #ifdef TRACE_IO
00157                            std::cout << v[i] << " ";
00158                            if (i%6==0) std::cout << "\n";
00159                      #endif
00160     }
00161                      #ifdef TRACE_IO
00162                          std::cout << "\n";
00163                      #endif
00164   }
00165 #ifdef REMOVED
00166   if (!outFile.bad()) {
00167     outFile << std::setprecision(20) << theSeed << std::endl;
00168     for (int i = 0; i < numBuff; ++i) {
00169       outFile << buffer[i] << " ";
00170     }
00171     outFile << redSpin  << " " << numFlats << " " << halfBuff << std::endl;
00172   }
00173 #endif
00174 }
00175 
00176 void RanshiEngine::restoreStatus(const char filename[]) {
00177   std::ifstream inFile(filename, std::ios::in);
00178   if (!checkFile ( inFile, filename, engineName(), "restoreStatus" )) {
00179     std::cerr << "  -- Engine state remains unchanged\n";
00180     return;
00181   }
00182   if ( possibleKeywordInput ( inFile, "Uvec", theSeed ) ) {
00183     std::vector<unsigned long> v;
00184     unsigned long xin;
00185     for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
00186       inFile >> xin;
00187                #ifdef TRACE_IO
00188                std::cout << "ivec = " << ivec << "  xin = " << xin << "    ";
00189                if (ivec%3 == 0) std::cout << "\n"; 
00190                #endif
00191       if (!inFile) {
00192         inFile.clear(std::ios::badbit | inFile.rdstate());
00193         std::cerr << "\nRanshiEngine state (vector) description improper."
00194                << "\nrestoreStatus has failed."
00195                << "\nInput stream is probably mispositioned now." << std::endl;
00196         return;
00197       }
00198       v.push_back(xin);
00199     }
00200     getState(v);
00201     return;
00202   }
00203 
00204   if (!inFile.bad()) {
00205 //     inFile >> theSeed;  removed -- encompased by possibleKeywordInput
00206     for (int i = 0; i < numBuff; ++i) {
00207       inFile >> buffer[i];
00208     }
00209     inFile >> redSpin >> numFlats >> halfBuff;
00210   }
00211 }
00212 
00213 void RanshiEngine::showStatus() const {
00214   std::cout << std::setprecision(20) << std::endl;
00215   std::cout << "----------- Ranshi engine status ----------" << std::endl;
00216   std::cout << "Initial seed      = " << theSeed << std::endl;
00217   std::cout << "Current red spin  = " << redSpin << std::endl;
00218   std::cout << "Values produced   = " << numFlats << std::endl;
00219   std::cout << "Side of buffer    = " << (halfBuff ? "upper" : "lower")
00220             << std::endl;
00221   std::cout << "Current buffer    = " << std::endl;
00222   for (int i = 0; i < numBuff; i+=4) {
00223     std::cout << std::setw(10) << std::setiosflags(std::ios::right)
00224               << buffer[i]     << std::setw(11) << buffer[i+1] << std::setw(11)
00225               << buffer[i+2]   << std::setw(11) << buffer[i+3] << std::endl;
00226   }
00227   std::cout << "-------------------------------------------" << std::endl;
00228 }
00229 
00230 RanshiEngine::operator float() {
00231   unsigned int redAngle = (((numBuff/2) - 1) & redSpin) + halfBuff;
00232   unsigned int blkSpin  = buffer[redAngle] & 0xffffffff;
00233   
00234   buffer[redAngle] = ((blkSpin << 17) | (blkSpin >> (32-17))) ^ redSpin;
00235   
00236   redSpin  = (blkSpin + numFlats++) & 0xffffffff;
00237   halfBuff = numBuff/2 - halfBuff;
00238   
00239   return float(blkSpin * twoToMinus_32());
00240 }
00241 
00242 RanshiEngine::operator unsigned int() {
00243   unsigned int redAngle = (((numBuff/2) - 1) & redSpin) + halfBuff;
00244   unsigned int blkSpin  = buffer[redAngle] & 0xffffffff;
00245   
00246   buffer[redAngle] = ((blkSpin << 17) | (blkSpin >> (32-17))) ^ redSpin;
00247   
00248   redSpin  = (blkSpin + numFlats++) & 0xffffffff;
00249   halfBuff = numBuff/2 - halfBuff;
00250   
00251   return blkSpin;
00252 }
00253 
00254 std::ostream& RanshiEngine::put (std::ostream& os ) const {
00255   char beginMarker[] = "RanshiEngine-begin";
00256   os << beginMarker << "\nUvec\n";
00257   std::vector<unsigned long> v = put();
00258   for (unsigned int i=0; i<v.size(); ++i) {
00259      os <<  v[i] <<  "\n";
00260   }
00261   return os;  
00262 #ifdef REMOVED 
00263   char endMarker[]   = "RanshiEngine-end";
00264   int pr=os.precision(20);
00265   os << " " << beginMarker << " ";
00266   
00267   os << theSeed  << "\n";
00268   for (int i = 0; i < numBuff; ++i) {
00269     os << buffer[i]  << "\n";
00270   }
00271   os << redSpin  << " " << numFlats << "\n" << halfBuff; 
00272   
00273   os << " " << endMarker   << "\n";
00274   os.precision(pr);
00275   return os;
00276 #endif
00277 }
00278 
00279 std::vector<unsigned long> RanshiEngine::put () const {
00280   std::vector<unsigned long> v;
00281   v.push_back (engineIDulong<RanshiEngine>());
00282   for (int i = 0; i < numBuff; ++i) {
00283     v.push_back(static_cast<unsigned long>(buffer[i]));
00284   }
00285   v.push_back(static_cast<unsigned long>(redSpin));
00286   v.push_back(static_cast<unsigned long>(numFlats));
00287   v.push_back(static_cast<unsigned long>(halfBuff));  
00288   return v;
00289 }
00290 
00291 std::istream& RanshiEngine::get (std::istream& is) {
00292   char beginMarker [MarkerLen];
00293   is >> std::ws;
00294   is.width(MarkerLen);  // causes the next read to the char* to be <=
00295                         // that many bytes, INCLUDING A TERMINATION \0 
00296                         // (Stroustrup, section 21.3.2)
00297   is >> beginMarker;
00298   if (strcmp(beginMarker,"RanshiEngine-begin")) {
00299     is.clear(std::ios::badbit | is.rdstate());
00300     std::cerr << "\nInput mispositioned or"
00301               << "\nRanshiEngine state description missing or"
00302               << "\nwrong engine type found." << std::endl;
00303     return is;
00304   }
00305   return getState(is);
00306 }
00307 
00308 std::string RanshiEngine::beginTag ( )  { 
00309   return "RanshiEngine-begin"; 
00310 }
00311   
00312 std::istream& RanshiEngine::getState (std::istream& is) {
00313   if ( possibleKeywordInput ( is, "Uvec", theSeed ) ) {
00314     std::vector<unsigned long> v;
00315     unsigned long uu;
00316     for (unsigned int ivec=0; ivec < VECTOR_STATE_SIZE; ++ivec) {
00317       is >> uu;
00318       if (!is) {
00319         is.clear(std::ios::badbit | is.rdstate());
00320         std::cerr << "\nRanshiEngine state (vector) description improper."
00321                 << "\ngetState() has failed."
00322                << "\nInput stream is probably mispositioned now." << std::endl;
00323         return is;
00324       }
00325       v.push_back(uu);
00326     }
00327     getState(v);
00328     return (is);
00329   }
00330 
00331 //  is >> theSeed;  Removed, encompassed by possibleKeywordInput()
00332 
00333   char endMarker   [MarkerLen];
00334   for (int i = 0; i < numBuff; ++i) {
00335     is >> buffer[i];
00336   }
00337   is >> redSpin >> numFlats >> halfBuff;
00338   is >> std::ws;
00339   is.width(MarkerLen);  
00340   is >> endMarker;
00341   if (strcmp(endMarker,"RanshiEngine-end")) {
00342     is.clear(std::ios::badbit | is.rdstate());
00343     std::cerr << "\nRanshiEngine state description incomplete."
00344               << "\nInput stream is probably mispositioned now." << std::endl;
00345     return is;
00346   }
00347   return is;
00348 }
00349 
00350 bool RanshiEngine::get (const std::vector<unsigned long> & v) {
00351   if ((v[0] & 0xffffffffUL) != engineIDulong<RanshiEngine>()) {
00352     std::cerr << 
00353         "\nRanshiEngine get:state vector has wrong ID word - state unchanged\n";
00354     return false;
00355   }
00356   return getState(v);
00357 }
00358 
00359 bool RanshiEngine::getState (const std::vector<unsigned long> & v) {
00360   if (v.size() != VECTOR_STATE_SIZE ) {
00361     std::cerr << 
00362         "\nRanshiEngine get:state vector has wrong length - state unchanged\n";
00363     return false;
00364   }
00365   for (int i = 0; i < numBuff; ++i) {
00366     buffer[i] = v[i+1];
00367   }
00368   redSpin  = v[numBuff+1];
00369   numFlats = v[numBuff+2]; 
00370   halfBuff = v[numBuff+3];
00371   return true;
00372 }
00373 
00374 }  // namespace CLHEP