CLHEP VERSION Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // ---------------------------------------------------------------------- 00002 #include "CLHEP/Random/Randomize.h" 00003 #include "CLHEP/Random/NonRandomEngine.h" 00004 #include "CLHEP/Random/defs.h" 00005 #include <iostream> 00006 #include <iomanip> 00007 #include <vector> 00008 00009 #define CLEAN_OUTPUT 00010 #ifdef CLEAN_OUTPUT 00011 std::ofstream output("testAnonymousEngineRestore.cout"); 00012 #else 00013 std::ostream & output = std::cout; 00014 #endif 00015 00016 // Normally on for routine validation: 00017 00018 #ifdef TURNOFF 00019 #endif 00020 00021 #define TEST_ANONYMOUS_ENGINE_RESTORE 00022 #define TEST_ANONYMOUS_RESTORE_STATICS 00023 00024 #define VERBOSER 00025 #define VERBOSER2 00026 00027 using namespace CLHEP; 00028 00029 template <class E1, class E2> int anonymousRestoreStatics(); 00030 00031 00032 // Absolutely Safe Equals Without Registers Screwing Us Up 00033 bool equals01(const std::vector<double> &ab) { 00034 return ab[1]==ab[0]; 00035 } 00036 bool equals(double a, double b) { 00037 std::vector<double> ab(2); 00038 ab[0]=a; ab[1]=b; 00039 return (equals01(ab)); 00040 } 00041 00042 std::vector<double> aSequence(int n) { 00043 std::vector<double> v; 00044 DualRand e(13542); 00045 RandFlat f(e); 00046 for (int i=0; i<n; i++) { 00047 v.push_back(f()); 00048 } 00049 return v; 00050 } 00051 00052 // ----------- Tests saving all statics together ----------- 00053 00054 void randomizeStatics(int n) { 00055 for (int i=0; i<n; i++) { 00056 RandGauss::shoot(); 00057 RandGaussQ::shoot(); 00058 RandGaussT::shoot(); 00059 RandFlat::shoot(); 00060 RandBit::shoot(); 00061 RandFlat::shootBit(); 00062 RandBit::shootBit(); 00063 RandPoisson::shoot(); 00064 RandPoissonQ::shoot(); 00065 RandPoissonT::shoot(); 00066 RandBinomial::shoot(); 00067 RandBreitWigner::shoot(); 00068 RandChiSquare::shoot(); 00069 RandExponential::shoot(); 00070 RandGamma::shoot(); 00071 RandLandau::shoot(); 00072 RandSkewNormal::shoot(); 00073 RandStudentT::shoot(); 00074 } 00075 } 00076 00077 std::vector<double> captureStatics() { 00078 std::vector<double> c; 00079 c.push_back( RandGauss::shoot() ); 00080 c.push_back( RandGaussQ::shoot() ); 00081 c.push_back( RandGaussT::shoot() ); 00082 c.push_back( RandFlat::shoot() ); 00083 c.push_back( RandBit::shoot() ); 00084 for (int i=0; i<20; i++) { 00085 c.push_back( RandFlat::shootBit() ); 00086 c.push_back( RandBit::shootBit() ); 00087 } 00088 c.push_back( RandPoisson::shoot() ); 00089 c.push_back( RandPoissonQ::shoot() ); 00090 c.push_back( RandPoissonT::shoot() ); 00091 c.push_back( RandBinomial::shoot() ); 00092 c.push_back( RandBreitWigner::shoot() ); 00093 c.push_back( RandChiSquare::shoot() ); 00094 c.push_back( RandExponential::shoot() ); 00095 c.push_back( RandGamma::shoot() ); 00096 c.push_back( RandLandau::shoot() ); 00097 c.push_back( RandSkewNormal::shoot() ); 00098 c.push_back( RandStudentT::shoot() ); 00099 return c; 00100 } 00101 00102 void saveStatics(std::string filename) { 00103 std::ofstream os(filename.c_str()); 00104 RandGeneral::saveStaticRandomStates(os); 00105 // It should be possible to call this from HepRandom, or any distribution. 00106 // RandGeneral, which is meaningless as a static distribution, should be the 00107 // toughest test, so we use that here. 00108 } 00109 00110 void restoreStatics(std::string filename) { 00111 std::ifstream is(filename.c_str()); 00112 RandLandau::restoreStaticRandomStates(is); 00113 } 00114 00115 00116 00117 // ----------- Anonymous restore of engines ----------- 00118 00119 template <class E> 00120 void anonymousRestore1(int n, std::vector<double> & v) { 00121 output << "Anonymous restore for " << E::engineName() << "\n"; 00122 E e(12349876); 00123 double r=0; 00124 for (int i=0; i<n; i++) r += e.flat(); 00125 std::ofstream os("anonymous.save"); 00126 os << e; 00127 for (int j=0; j<25; j++) v.push_back(e.flat()); 00128 #ifdef VERBOSER2 00129 output << "First four of v are: " 00130 << v[0] << ", " << v[1] << ", " << v[2] << ", " << v[3] << "\n"; 00131 #endif 00132 return; 00133 } 00134 00135 template <> 00136 void anonymousRestore1<NonRandomEngine> (int n, std::vector<double> & v) { 00137 #ifdef VERBOSER 00138 output << "Anonymous restore for " << NonRandomEngine::engineName() << "\n"; 00139 #endif 00140 std::vector<double> nonRand = aSequence(500); 00141 NonRandomEngine e; 00142 e.setRandomSequence(&nonRand[0], nonRand.size()); 00143 double r=0; 00144 for (int i=0; i<n; i++) r += e.flat(); 00145 std::ofstream os("anonymous.save"); 00146 os << e; 00147 for (int j=0; j<25; j++) v.push_back(e.flat()); 00148 #ifdef VERBOSER2 00149 output << "First four of v are: " 00150 << v[0] << ", " << v[1] << ", " << v[2] << ", " << v[3] << "\n"; 00151 #endif 00152 return; 00153 } 00154 00155 template <class E> 00156 int anonymousRestore2(const std::vector<double> & v) { 00157 int stat = 0; 00158 std::vector<double> k; 00159 std::ifstream is("anonymous.save"); 00160 HepRandomEngine * a; 00161 a = HepRandomEngine::newEngine(is); 00162 for (int j=0; j<25; j++) k.push_back(a->flat()); 00163 delete a; 00164 #ifdef VERBOSER2 00165 output << "First four of k are: " 00166 << k[0] << ", " << k[1] << ", " << k[2] << ", " << k[3] << "\n"; 00167 #endif 00168 for (int m=0; m<25; m++) { 00169 if ( v[m] != k[m] ) { 00170 std::cout << "???? Incorrect restored value for anonymous engine" 00171 << E::engineName() << "\n"; 00172 #ifdef CLEAN_OUTPUT 00173 output << "???? Incorrect restored value for anonymous engine" 00174 << E::engineName() << "\n"; 00175 #endif 00176 stat |= 262144; 00177 return stat; 00178 } 00179 } 00180 return stat; 00181 } 00182 00183 00184 template <class E> 00185 int anonymousRestore(int n) { 00186 std::vector<double> v; 00187 anonymousRestore1<E>(n,v); 00188 return anonymousRestore2<E>(v); 00189 } 00190 00191 // ----------- Anonymous restore of all static distributions ----------- 00192 00193 template <class E> 00194 int anonymousRestoreStatics1() { 00195 int stat = 0; 00196 HepRandomEngine *e = new E(12456); 00197 HepRandom::setTheEngine(e); 00198 randomizeStatics(15); 00199 output << "\nRandomized, with theEngine = " << e->name() << "\n"; 00200 saveStatics("distribution.save"); 00201 output << "Saved all static distributions\n"; 00202 std::vector<double> c = captureStatics(); 00203 output << "Captured output of all static distributions\n"; 00204 randomizeStatics(11); 00205 output << "Randomized all static distributions\n"; 00206 restoreStatics("distribution.save"); 00207 output << "Restored all static distributions to saved state\n"; 00208 std::vector<double> d = captureStatics(); 00209 output << "Captured output of all static distributions\n"; 00210 for (unsigned int iv=0; iv<c.size(); iv++) { 00211 if (c[iv] != d[iv]) { 00212 std::cout << "???? restoreStaticRandomStates failed at random " 00213 << iv <<"\n"; 00214 #ifdef CLEAN_OUTPUT 00215 output << "???? restoreStaticRandomStates failed at random " 00216 << iv <<"\n"; 00217 #endif 00218 stat |= 131072; 00219 } 00220 } 00221 if ( (stat & 131072) == 0) { 00222 output << "All captured output agrees with earlier values\n"; 00223 } 00224 return stat; 00225 } 00226 00227 00228 00229 template <class E1, class E2> 00230 int anonymousRestoreStatics() { 00231 int stat = 0; 00232 if ( E1::engineName() == E2::engineName() ) { 00233 return anonymousRestoreStatics1<E1>(); 00234 } 00235 HepRandomEngine *e1 = new E1(12456); 00236 HepRandom::setTheEngine(e1); 00237 randomizeStatics(15); 00238 output << "\nRandomized, with theEngine = " << e1->name() << "\n"; 00239 saveStatics("distribution.save"); 00240 #ifdef VERBOSER2 00241 output << "Saved all static distributions\n"; 00242 #endif 00243 std::vector<double> c = captureStatics(); 00244 #ifdef VERBOSER2 00245 output << "Captured output of all static distributions\n"; 00246 #endif 00247 delete e1; 00248 HepRandomEngine *e2 = new E2(24653); 00249 HepRandom::setTheEngine(e2); 00250 output << "Switched to theEngine = " << e2->name() << "\n"; 00251 randomizeStatics(19); 00252 { std::ofstream os("engine.save"); os << *e2; } 00253 double v1 = e2->flat(); 00254 double v2 = e2->flat(); 00255 { std::ifstream is("engine.save"); is >> *e2; } 00256 #ifdef VERBOSER2 00257 output << "Saved the " << e2->name() << " engine: \n" 00258 << "Next randoms to be " << v1 << " " << v2 << "\n" 00259 << "Restored the " << e2->name() << " engine to that state\n"; 00260 #endif 00261 restoreStatics("distribution.save"); 00262 #ifdef VERBOSER2 00263 output << "Restored all static distributions to saved state\n" 00264 << "This changes the engine type back to " << E1::engineName() << "\n"; 00265 #endif 00266 std::vector<double> d = captureStatics(); 00267 #ifdef VERBOSER2 00268 output << "Captured output of all static distributions\n"; 00269 #endif 00270 for (unsigned int iv=0; iv<c.size(); iv++) { 00271 if (c[iv] != d[iv]) { 00272 std::cout << "???? restoreStaticRandomStates failed at random " 00273 << iv <<"\n"; 00274 #ifdef CLEAN_OUTPUT 00275 output << "???? restoreStaticRandomStates failed at random " 00276 << iv <<"\n"; 00277 #endif 00278 stat |= 524288; 00279 } 00280 } 00281 if ((stat & 524288) == 0) { 00282 output << "All captured output agrees with earlier values\n"; 00283 } 00284 double k1 = e2->flat(); 00285 double k2 = e2->flat(); 00286 #ifdef VERBOSER2 00287 output << "The " << e2->name() << " engine should not have been affected: \n" 00288 << "Next randoms are " << k1 << " " << k2 << "\n"; 00289 #endif 00290 if ( !equals(v1,k1) || !equals(v2,k2) ) { 00291 std::cout << "???? Engine used as theEngine was affected by restoring \n" 00292 << " static distributions to use engine of a different type.\n"; 00293 #ifdef CLEAN_OUTPUT 00294 output << "???? Engine used as theEngine was affected by restoring \n" 00295 << " static distributions to use engine of a different type.\n"; 00296 #endif 00297 stat |= 1048576; 00298 } 00299 return stat; 00300 } 00301 00302 // --------------------------------------------- 00303 // --------------------------------------------- 00304 // --------------------------------------------- 00305 00306 00307 int main() { 00308 int stat = 0; 00309 00310 #ifdef TEST_ANONYMOUS_ENGINE_RESTORE 00311 output << "\n=================================\n"; 00312 output << " Part VII \n"; 00313 output << "Anonymous restore of engines \n"; 00314 output << "=================================\n\n"; 00315 00316 stat |= anonymousRestore<DualRand>(13); 00317 stat |= anonymousRestore<DRand48Engine>(14); 00318 stat |= anonymousRestore<Hurd160Engine>(15); 00319 stat |= anonymousRestore<Hurd288Engine>(16); 00320 stat |= anonymousRestore<HepJamesRandom>(17); 00321 stat |= anonymousRestore<MTwistEngine>(18); 00322 stat |= anonymousRestore<RandEngine>(29); 00323 stat |= anonymousRestore<RanecuEngine>(39); 00324 stat |= anonymousRestore<Ranlux64Engine>(19); 00325 stat |= anonymousRestore<RanluxEngine>(20); 00326 stat |= anonymousRestore<RanshiEngine>(21); 00327 stat |= anonymousRestore<TripleRand>(22); 00328 stat |= anonymousRestore<NonRandomEngine>(22); 00329 #endif 00330 00331 #ifdef TEST_ANONYMOUS_RESTORE_STATICS 00332 output << "\n======================================\n"; 00333 output << " Part VIII \n"; 00334 output << "Anonymous restore static Distributions \n"; 00335 output << "======================================\n\n"; 00336 00337 stat |= anonymousRestoreStatics<DualRand, Ranlux64Engine> ( ); 00338 stat |= anonymousRestoreStatics<DRand48Engine, TripleRand> ( ); 00339 stat |= anonymousRestoreStatics<RandEngine, Ranlux64Engine> ( ); 00340 stat |= anonymousRestoreStatics<MTwistEngine, Hurd288Engine> ( ); 00341 stat |= anonymousRestoreStatics<RanecuEngine, MTwistEngine> ( ); 00342 stat |= anonymousRestoreStatics<HepJamesRandom, RanshiEngine> ( ); 00343 stat |= anonymousRestoreStatics<RanecuEngine, RandEngine> ( ); 00344 stat |= anonymousRestoreStatics<RanshiEngine, Hurd160Engine> ( ); 00345 stat |= anonymousRestoreStatics<TripleRand, DualRand> ( ); 00346 stat |= anonymousRestoreStatics<Hurd160Engine, HepJamesRandom> ( ); 00347 stat |= anonymousRestoreStatics<Hurd288Engine, RanecuEngine> ( ); 00348 stat |= anonymousRestoreStatics<HepJamesRandom, Ranlux64Engine> ( ); 00349 stat |= anonymousRestoreStatics<TripleRand, TripleRand> ( ); 00350 stat |= anonymousRestoreStatics<HepJamesRandom, HepJamesRandom> ( ); 00351 #endif 00352 00353 00354 output << "\n=============================================\n\n"; 00355 00356 if (stat != 0) { 00357 std::cout << "One or more problems detected: stat = " << stat << "\n"; 00358 output << "One or more problems detected: stat = " << stat << "\n"; 00359 } else { 00360 output << "testAnonymousEngineRestore passed with no problems detected.\n"; 00361 } 00362 00363 if (stat == 0) return 0; 00364 if (stat > 0) return -(stat|1); 00365 return stat|1; 00366 } 00367