CLHEP 2.0.4.7 Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // -*- C++ -*- 00002 // $Id: ranRestoreTest.cc,v 1.3.4.3 2005/04/15 16:32:53 garren Exp $ 00003 // ---------------------------------------------------------------------- 00004 #include "CLHEP/Random/Randomize.h" 00005 #include "CLHEP/Random/NonRandomEngine.h" 00006 #include "CLHEP/Random/defs.h" 00007 #include <iostream> 00008 #include <iomanip> 00009 #include <vector> 00010 00011 #define CLEAN_OUTPUT 00012 #ifdef CLEAN_OUTPUT 00013 std::ofstream output("ranRestoreTest.cout"); 00014 #else 00015 std::ostream & output = std::cout; 00016 #endif 00017 00018 // Normally on for routine validation: 00019 00020 #define TEST_ORIGINAL_SAVE 00021 00022 #ifdef TURNOFF 00023 #endif 00024 00025 #define TEST_ENGINE_NAMES 00026 #define TEST_INSTANCE_METHODS 00027 #define TEST_SHARED_ENGINES 00028 #define TEST_STATIC_SAVE 00029 #define TEST_SAVE_STATIC_STATES 00030 #define TEST_ANONYMOUS_ENGINE_RESTORE 00031 #define TEST_ANONYMOUS_RESTORE_STATICS 00032 #define TEST_VECTOR_ENGINE_RESTORE 00033 00034 // Normally off for routine validation: 00035 00036 #ifdef TURNOFF 00037 #define TEST_MISSING_FILES 00038 #define CREATE_OLD_SAVES 00039 #define VERIFY_OLD_SAVES 00040 #endif 00041 00042 //#define VERBOSER 00043 //#define VERBOSER2 00044 00045 using namespace CLHEP; 00046 00047 double remembered_r2; 00048 double remembered_r1005; 00049 double remembered_r1006; 00050 double remembered_r1007; 00051 00052 // Absolutely Safe Equals Without Registers Screwing Us Up 00053 bool equals01(const std::vector<double> &ab) { 00054 return ab[1]==ab[0]; 00055 } 00056 bool equals(double a, double b) { 00057 std::vector<double> ab(2); 00058 ab[0]=a; ab[1]=b; 00059 return (equals01(ab)); 00060 } 00061 00062 // ------------------- The following should all FAIL ------------ 00063 00064 int saveStepX() { 00065 double r = RandGauss::shoot(); 00066 output << "r(1) = " << r << std::endl; 00067 HepRandom::saveEngineStatus(); 00068 r = RandGauss::shoot(); 00069 output << "r(2) = " << r << std::endl; 00070 remembered_r2 = r; 00071 r = RandGauss::shoot(); 00072 output << "r(3) = " << r << std::endl; 00073 for (int i=0; i < 1001; i++) { 00074 r = RandGauss::shoot(); 00075 } 00076 r = RandGauss::shoot(); 00077 remembered_r1005 = r; 00078 output << "r1005= " << r << std::endl; 00079 r = RandGauss::shoot(); 00080 return 0; 00081 } 00082 00083 int restoreStepX() { 00084 HepRandom::restoreEngineStatus(); 00085 double r = RandGauss::shoot(); 00086 output << "restored r(2) = " << r << std::endl; 00087 if ( ! equals(r,remembered_r2) ) { 00088 output << "THIS DOES NOT MATCH REMEMBERED VALUE BUT THAT IS EXPECTED\n"; 00089 } 00090 r = RandGauss::shoot(); 00091 output << "restored r(3) = " << r << std::endl; 00092 for (int i=0; i < 1001; i++) { 00093 r = RandGauss::shoot(); 00094 } 00095 r = RandGauss::shoot(); 00096 output << "restored r1005= " << r << std::endl; 00097 if ( !equals(r,remembered_r1005) ) { 00098 output << "THIS DOES NOT MATCH REMEMBERED VALUE BUT THAT IS EXPECTED\n"; 00099 } 00100 return 0; 00101 } 00102 00103 int BsaveStepX() { 00104 int r = RandFlat::shootBit(); 00105 output << "r(1) = " << r << std::endl; 00106 HepRandom::saveEngineStatus(); 00107 r = RandFlat::shootBit(); 00108 output << "r(2) = " << r << std::endl; 00109 remembered_r2 = r; 00110 r = RandFlat::shootBit(); 00111 output << "r(3) = " << r << std::endl; 00112 double d; 00113 for (int i=0; i < 1001; i++) { 00114 d = RandFlat::shoot(); 00115 if (d > 1) output << 00116 "This line inserted so clever compilers don't warn about not using d\n"; 00117 } 00118 r = RandFlat::shootBit(); 00119 remembered_r1005 = r; 00120 output << "r1005= " << r << std::endl; 00121 r = RandFlat::shootBit(); 00122 return 0; 00123 } 00124 00125 int BrestoreStepX() { 00126 HepRandom::restoreEngineStatus(); 00127 int r = RandFlat::shootBit(); 00128 output << "restored r(2) = " << r << std::endl; 00129 if ( r != remembered_r2 ) { 00130 output << "THIS DOES NOT MATCH REMEMBERED VALUE BUT THAT IS EXPECTED\n"; 00131 } 00132 r = RandFlat::shootBit(); 00133 output << "restored r(3) = " << r << std::endl; 00134 for (int i=0; i < 1001; i++) { 00135 r = RandFlat::shootBit(); 00136 } 00137 r = RandFlat::shootBit(); 00138 output << "restored r1005= " << r << std::endl; 00139 if ( r != remembered_r1005 ) { 00140 output << "THIS DOES NOT MATCH REMEMBERED VALUE BUT THAT IS EXPECTED\n"; 00141 } 00142 return 0; 00143 } 00144 00145 // ------------------- The following should all WORK ------------ 00146 00147 int saveStep() { 00148 int stat=0; 00149 double r = RandGauss::shoot(); 00150 output << "r(1) = " << r << std::endl; 00151 RandGauss::saveEngineStatus(); 00152 r = RandGauss::shoot(); 00153 output << "r(2) = " << r << std::endl; 00154 remembered_r2 = r; 00155 r = RandGauss::shoot(); 00156 output << "r(3) = " << r << std::endl; 00157 for (int i=0; i < 1001; i++) { 00158 r = RandGauss::shoot(); 00159 } 00160 r = RandGauss::shoot(); 00161 remembered_r1005 = r; 00162 output << "r1005= " << r << std::endl; 00163 r = RandGauss::shoot(); 00164 return stat; 00165 } 00166 00167 int restoreStep() { 00168 int stat=0; 00169 RandGauss::restoreEngineStatus(); 00170 double r = RandGauss::shoot(); 00171 output << "restored r(2) = " << r << std::endl; 00172 if ( !equals(r,remembered_r2) ) { 00173 std::cout << "restored r(2) = " << r << std::endl; 00174 std::cout << "????? THIS DOES NOT MATCH REMEMBERED VALUE!\n"; 00175 stat += 1; 00176 } 00177 r = RandGauss::shoot(); 00178 output << "restored r(3) = " << r << std::endl; 00179 for (int i=0; i < 1001; i++) { 00180 r = RandGauss::shoot(); 00181 } 00182 r = RandGauss::shoot(); 00183 output << "restored r1005= " << r << std::endl; 00184 if ( !equals(r,remembered_r1005) ) { 00185 std::cout << "restored r1005= " << r << std::endl; 00186 std::cout << "????? THIS DOES NOT MATCH REMEMBERED VALUE!\n"; 00187 stat += 2; 00188 } 00189 return stat; 00190 } 00191 00192 int BsaveStep() { 00193 int stat=0; 00194 int r = RandFlat::shootBit(); 00195 output << "r(1) = " << r << std::endl; 00196 RandFlat::saveEngineStatus(); 00197 r = RandFlat::shootBit(); 00198 output << "r(2) = " << r << std::endl; 00199 remembered_r2 = r; 00200 r = RandFlat::shootBit(); 00201 output << "r(3) = " << r << std::endl; 00202 for (int i=0; i < 1001; i++) { 00203 r = RandFlat::shootBit(); 00204 } 00205 r = RandFlat::shootBit(); 00206 remembered_r1005 = r; 00207 output << "r1005 = " << r << std::endl; 00208 r = RandFlat::shootBit(); 00209 remembered_r1006 = r; 00210 output << "r1006 = " << r << std::endl; 00211 r = RandFlat::shootBit(); 00212 remembered_r1007 = r; 00213 output << "r1007 = " << r << std::endl; 00214 r = RandFlat::shootBit(); 00215 return stat; 00216 } 00217 00218 int BrestoreStep() { 00219 int stat=0; 00220 RandFlat::restoreEngineStatus(); 00221 int r = RandFlat::shootBit(); 00222 output << "restored r(2) = " << r << std::endl; 00223 if ( r != remembered_r2 ) { 00224 stat += 4; 00225 std::cout << "restored r(2) = " << r << std::endl; 00226 std::cout << "????? THIS DOES NOT MATCH REMEMBERED VALUE!\n"; 00227 } 00228 r = RandFlat::shootBit(); 00229 output << "restored r(3) = " << r << std::endl; 00230 for (int i=0; i < 1001; i++) { 00231 r = RandFlat::shootBit(); 00232 } 00233 r = RandFlat::shootBit(); 00234 output << "restored r1005= " << r << std::endl; 00235 if ( r != remembered_r1005 ) { 00236 stat += 8; 00237 std::cout << "restored r1005= " << r << std::endl; 00238 std::cout << "????? THIS DOES NOT MATCH REMEMBERED VALUE!\n"; 00239 } 00240 r = RandFlat::shootBit(); 00241 output << "restored r1006= " << r << std::endl; 00242 if ( r != remembered_r1006 ) { 00243 stat += 16; 00244 std::cout << "restored r1006= " << r << std::endl; 00245 std::cout << "????? THIS DOES NOT MATCH REMEMBERED VALUE!\n"; 00246 } 00247 r = RandFlat::shootBit(); 00248 output << "restored r1007= " << r << std::endl; 00249 if ( r != remembered_r1007 ) { 00250 stat += 32; 00251 std::cout << "restored r1007= " << r << std::endl; 00252 std::cout << "????? THIS DOES NOT MATCH REMEMBERED VALUE!\n"; 00253 } 00254 return stat; 00255 } 00256 00257 // --- The following should work, by failing in an expected way ------- 00258 00259 template <class E, class D> 00260 int fileNotThere() { 00261 int stat = 0; 00262 HepRandomEngine * old = D::getTheEngine(); 00263 E e(123); 00264 output << "File-not-found test restoring "<<D::distributionName()<<":\n"; 00265 D::setTheEngine(&e); 00266 D::restoreEngineStatus("noSuchFile"); 00267 D::setTheEngine(old); // If we don't do this, then the static engine shared 00268 // by every shoot() method reamins e -- which is about 00269 // to go out of scope and be destructed! 00270 return stat; 00271 } 00272 00273 template <class E> 00274 int fileNotThereEngine() { 00275 int stat = 0; 00276 stat |= fileNotThere <E, RandBinomial>(); 00277 stat |= fileNotThere <E, RandBit>(); 00278 stat |= fileNotThere <E, RandBreitWigner>(); 00279 stat |= fileNotThere <E, RandChiSquare>(); 00280 stat |= fileNotThere <E, RandExponential>(); 00281 stat |= fileNotThere <E, RandFlat>(); 00282 stat |= fileNotThere <E, RandGamma>(); 00283 stat |= fileNotThere <E, RandGauss>(); 00284 stat |= fileNotThere <E, RandGaussQ>(); 00285 stat |= fileNotThere <E, RandGaussT>(); 00286 stat |= fileNotThere <E, RandLandau>(); 00287 stat |= fileNotThere <E, RandPoisson>(); 00288 stat |= fileNotThere <E, RandPoissonQ>(); 00289 stat |= fileNotThere <E, RandPoissonT>(); 00290 stat |= fileNotThere <E, RandStudentT>(); 00291 return stat; 00292 } 00293 00294 int missingFile() { 00295 int stat = 0; 00296 stat |= fileNotThereEngine<DRand48Engine>(); 00297 stat |= fileNotThereEngine<DualRand>(); 00298 stat |= fileNotThereEngine<Hurd160Engine>(); 00299 stat |= fileNotThereEngine<Hurd288Engine>(); 00300 stat |= fileNotThereEngine<HepJamesRandom>(); 00301 stat |= fileNotThereEngine<MTwistEngine>(); 00302 stat |= fileNotThereEngine<RandEngine>(); 00303 stat |= fileNotThereEngine<RanecuEngine>(); 00304 stat |= fileNotThereEngine<Ranlux64Engine>(); 00305 stat |= fileNotThereEngine<RanluxEngine>(); 00306 stat |= fileNotThereEngine<RanshiEngine>(); 00307 stat |= fileNotThereEngine<TripleRand>(); 00308 return stat; 00309 } 00310 00311 // -- The following was used to capture old-form engine states (sans name) -- 00312 00313 template <class E, class D> 00314 int saveEngine(const char* filename) { 00315 int stat = 0; 00316 HepRandomEngine * old = D::getTheEngine(); 00317 E e(123); 00318 D::setTheEngine(&e); 00319 double r=0; 00320 for (int i=0; i<3; i++) r += D::shoot(); 00321 D::saveEngineStatus(filename); 00322 if (r == -99999999.1) stat = 999; // This prevents clever compilers from 00323 // deducing that r is never needed 00324 D::setTheEngine(old); // If we don't do this, then the static engine shared 00325 // by every shoot() method reamins e -- which is about 00326 // to go out of scope and be destructed! 00327 return stat; 00328 } 00329 00330 // -- The following checks on static engine restores, from old and new forms -- 00331 00332 template <class E, class D> 00333 int checkSaveEngine(const char* filename) { 00334 int stat = 0; 00335 HepRandomEngine * old = D::getTheEngine(); 00336 00337 // Generate save with current format (default file name is fine) 00338 E e(123); 00339 D::setTheEngine(&e); 00340 double r=0; 00341 for (int i=0; i<3; i++) r += D::shoot(); 00342 D::saveEngineStatus(); 00343 00344 // Figure out what the key variate value should be 00345 double keyValue = D::shoot(); 00346 00347 // Restore state based on old file, and check for proper value 00348 D::restoreEngineStatus(filename); 00349 if (!equals(D::shoot(), keyValue)) { 00350 std::cout << "???? Value mismatch from file " << filename << "\n"; 00351 stat |= 64; 00352 } 00353 00354 // Restore state based on that save, and check for proper value 00355 D::restoreEngineStatus(); 00356 if (!equals(D::shoot(),keyValue)) { 00357 std::cout << "???? Value mismatch from new-format file \n"; 00358 stat |= 128; 00359 } 00360 00361 D::setTheEngine(old); 00362 return stat; 00363 } 00364 00365 00366 // ----------- Tests for instance methods ----------- 00367 00368 template <class E> 00369 int checkEngineName(const std::string & name) { 00370 int stat = 0; 00371 output << E::engineName() << "\n"; 00372 if (E::engineName() != name) { 00373 std::cout << "???? engineName mismatch for " << name << " <--> " 00374 << E::engineName() << "\n"; 00375 stat |= 256; 00376 } 00377 E e(123); 00378 if (e.name() != name) { 00379 std::cout << "???? name mismatch for " << name << " <--> " 00380 << e.name() << "\n"; 00381 stat |= 256; 00382 } 00383 return stat; 00384 } 00385 00386 template <class E, class D> 00387 int checkEngine() { 00388 int stat = 0; 00389 E e(1234); 00390 D d(e); 00391 if (d.engine().name() != e.name()) { 00392 std::cout << "???? Improper d.engine() \n"; 00393 stat |= 512; 00394 } 00395 return stat; 00396 } 00397 00398 template <class E> 00399 int checkEngineInstanceSave(E & e) { 00400 int stat = 0; 00401 output << "checkEngineInstanceSave for " << e.name() << "\n"; 00402 int pr=output.precision(20); 00403 double r=0; 00404 for (int i=0; i<100; i++) r += e.flat(); 00405 {std::ofstream os ("engine.save"); os << e;} 00406 for (int i=0; i<100; i++) r += e.flat(); 00407 double keyValue1 = e.flat(); 00408 double keyValue2 = e.flat(); 00409 #ifdef VERBOSER 00410 output << keyValue1 << " " << keyValue2 << "\n"; 00411 #endif 00412 E e2; 00413 {std::ifstream is ("engine.save"); is >> e2;} 00414 for (int i=0; i<100; i++) r += e2.flat(); 00415 double k1 = e2.flat(); 00416 double k2 = e2.flat(); 00417 #ifdef VERBOSER 00418 output << k1 << " " << k2 << "\n"; 00419 #endif 00420 if ( !(equals(k1,keyValue1)) || !(equals(k2,keyValue2)) ) { 00421 std::cout << "???? checkInstanceSave failed for " << e.name() << "\n"; 00422 stat |= 1024; 00423 } 00424 output.precision(pr); 00425 return stat; 00426 } 00427 00428 template <class E, class D> 00429 int checkSaveDistribution(D & d, int nth) { 00430 dynamic_cast<E &>(d.engine()); 00431 int stat = 0; 00432 output << "checkSaveDistribution with " << d.engine().name() 00433 << ", " << d.name() << "\n"; 00434 double r=0; 00435 r = d(); 00436 double keyValue1, keyValue2, keyValue3, keyValue4; 00437 for (int i=0; i<nth; i++) r += d(); 00438 {std::ofstream os ("distribution.save1"); os << d.engine() << d;} 00439 keyValue1 = d(); 00440 keyValue2 = d(); 00441 r += d(); 00442 // A second capture will test non-cached if first tested cached case: 00443 {std::ofstream os ("distribution.save2"); os << d.engine() << d;} 00444 keyValue3 = d(); 00445 keyValue4 = d(); 00446 int pr = output.precision(20); 00447 #ifdef VERBOSER 00448 output << "keyValue1 = " << keyValue1 << 00449 " keyValue2 = " << keyValue2 << "\n"; 00450 output << "keyValue3 = " << keyValue3 << 00451 " keyValue3 = " << keyValue4 << "\n"; 00452 #endif 00453 output.precision(pr); 00454 E e; 00455 D d2(e); 00456 { std::ifstream is ("distribution.save1"); is >> e >> d2;} 00457 double k1 = d2(); 00458 double k2 = d2(); 00459 { std::ifstream is ("distribution.save2"); is >> e >> d2;} 00460 double k3 = d2(); 00461 double k4 = d2(); 00462 #ifdef VERBOSER 00463 pr = output.precision(20); 00464 output << "k1 = " << k1 << 00465 " k2 = " << k2 << "\n"; 00466 output << "k3 = " << k3 << 00467 " k4 = " << k4 << "\n"; 00468 output.precision(pr); 00469 #endif 00470 if ( !equals(k1,keyValue1) || !equals(k2,keyValue2) || 00471 !equals(k3,keyValue3) || !equals(k4,keyValue4) ) { 00472 std::cout << "???? Incorrect restored value for distribution " 00473 << d.name() << "\n"; 00474 stat |= 2048; 00475 } 00476 // if (stat) exit(-1); 00477 return stat; 00478 } 00479 00480 template <class E> 00481 int checkRandGeneralDistribution(RandGeneral & d, int nth) { 00482 dynamic_cast<E &>(d.engine()); 00483 int stat = 0; 00484 output << "checkSaveDistribution with " << d.engine().name() 00485 << ", " << d.name() << "\n"; 00486 double r=0; 00487 r = d(); 00488 double keyValue1, keyValue2, keyValue3, keyValue4; 00489 for (int i=0; i<nth; i++) r += d(); 00490 {std::ofstream os ("distribution.save1"); os << d.engine() << d;} 00491 keyValue1 = d(); 00492 keyValue2 = d(); 00493 r += d(); 00494 // A second capture will test non-cached if first tested cached case: 00495 {std::ofstream os ("distribution.save2"); os << d.engine() << d;} 00496 keyValue3 = d(); 00497 keyValue4 = d(); 00498 int pr = output.precision(20); 00499 #ifdef VERBOSER 00500 output << "keyValue1 = " << keyValue1 << 00501 " keyValue2 = " << keyValue2 << "\n"; 00502 output << "keyValue3 = " << keyValue3 << 00503 " keyValue3 = " << keyValue4 << "\n"; 00504 #endif 00505 output.precision(pr); 00506 E e; 00507 double temp = 1; 00508 RandGeneral d2(e, &temp, 1); 00509 { std::ifstream is ("distribution.save1"); is >> e >> d2;} 00510 double k1 = d2(); 00511 double k2 = d2(); 00512 { std::ifstream is ("distribution.save2"); is >> e >> d2;} 00513 double k3 = d2(); 00514 double k4 = d2(); 00515 #ifdef VERBOSER 00516 pr = output.precision(20); 00517 output << "k1 = " << k1 << 00518 " k2 = " << k2 << "\n"; 00519 output << "k3 = " << k3 << 00520 " k4 = " << k4 << "\n"; 00521 output.precision(pr); 00522 #endif 00523 if ( !equals(k1,keyValue1) || !equals(k2,keyValue2) || 00524 !equals(k3,keyValue3) || !equals(k4,keyValue4) ) { 00525 std::cout << "???? Incorrect restored value for distribution " 00526 << d.name() << "\n"; 00527 stat |= 2048; 00528 } 00529 // if (stat) exit(-1); 00530 return stat; 00531 } 00532 00533 template <class E> 00534 int checkDistributions() { 00535 int stat = 0; 00536 {RandGauss d(new E(12561),100.0,3.0); 00537 stat |= checkSaveDistribution<E,RandGauss> (d,33); } 00538 {RandGauss d(new E(12572),100.0,3.0); 00539 stat |= checkSaveDistribution<E,RandGauss> (d,34); } 00540 {RandGaussQ d(new E(12563),10.0,4.0); 00541 stat |= checkSaveDistribution<E,RandGaussQ> (d,33); } 00542 {RandGaussT d(new E(12564),5.0,2.0); 00543 stat |= checkSaveDistribution<E,RandGaussT> (d,33); } 00544 {RandBinomial d(new E(12565),4,0.6); 00545 stat |= checkSaveDistribution<E,RandBinomial> (d,33); } 00546 {RandFlat d(new E(12576),12.5,35.0); 00547 stat |= checkSaveDistribution<E,RandFlat> (d,33); } 00548 {RandBit d(new E(12567)); 00549 stat |= checkSaveDistribution<E,RandBit> (d,31); } 00550 {RandBit d(new E(12578)); 00551 stat |= checkSaveDistribution<E,RandBit> (d,32); } 00552 {RandBit d(new E(12589)); 00553 stat |= checkSaveDistribution<E,RandBit> (d,33); } 00554 {RandBreitWigner d(new E(125611),50.0,15.0); 00555 stat |= checkSaveDistribution<E,RandBreitWigner> (d,33); } 00556 {RandChiSquare d(new E(125612),5.0); 00557 stat |= checkSaveDistribution<E,RandChiSquare> (d,33); } 00558 {RandExponential d(new E(125713),8.00); 00559 stat |= checkSaveDistribution<E,RandExponential> (d,33); } 00560 {RandGamma d(new E(125713),6.0,2.0); 00561 stat |= checkSaveDistribution<E,RandGamma> (d,33); } 00562 {RandLandau d(new E(125714)); 00563 stat |= checkSaveDistribution<E,RandLandau> (d,33); } 00564 {RandStudentT d(new E(125715),5); 00565 stat |= checkSaveDistribution<E,RandStudentT> (d,33); } 00566 00567 // Multiple tests of Poisson distributions for small desired, since 00568 // the answer in each test is a small int, and coincidental agreement 00569 // is very possible. 00570 00571 {RandPoisson d(new E(125616),2.5); 00572 stat |= checkSaveDistribution<E,RandPoisson> (d,33); } 00573 {RandPoisson d(new E(125617),105.0); 00574 stat |= checkSaveDistribution<E,RandPoisson> (d,34); } 00575 {RandPoisson d(new E(125618),2.5); 00576 stat |= checkSaveDistribution<E,RandPoisson> (d,35); } 00577 {RandPoisson d(new E(325618),2.5); 00578 stat |= checkSaveDistribution<E,RandPoisson> (d,36); } 00579 {RandPoisson d(new E(425618),2.5); 00580 stat |= checkSaveDistribution<E,RandPoisson> (d,37); } 00581 {RandPoisson d(new E(525618),2.5); 00582 stat |= checkSaveDistribution<E,RandPoisson> (d,38); } 00583 {RandPoisson d(new E(125619),110.0); 00584 stat |= checkSaveDistribution<E,RandPoisson> (d,39); } 00585 {RandPoissonQ d(new E(124616),2.5); 00586 stat |= checkSaveDistribution<E,RandPoissonQ> (d,33); } 00587 {RandPoissonQ d(new E(126616),2.5); 00588 stat |= checkSaveDistribution<E,RandPoissonQ> (d,32); } 00589 {RandPoissonQ d(new E(127616),2.5); 00590 stat |= checkSaveDistribution<E,RandPoissonQ> (d,31); } 00591 {RandPoissonQ d(new E(129616),2.5); 00592 stat |= checkSaveDistribution<E,RandPoissonQ> (d,30); } 00593 {RandPoissonQ d(new E(125616),110.0); 00594 stat |= checkSaveDistribution<E,RandPoissonQ> (d,33); } 00595 {RandPoissonQ d(new E(125616),2.5); 00596 stat |= checkSaveDistribution<E,RandPoissonQ> (d,34); } 00597 {RandPoissonQ d(new E(125616),110.0); 00598 stat |= checkSaveDistribution<E,RandPoissonQ> (d,34); } 00599 {RandPoissonT d(new E(125616),2.5); 00600 stat |= checkSaveDistribution<E,RandPoissonT> (d,33); } 00601 {RandPoissonT d(new E(125616),110.0); 00602 stat |= checkSaveDistribution<E,RandPoissonT> (d,33); } 00603 {RandPoissonT d(new E(125616),2.5); 00604 stat |= checkSaveDistribution<E,RandPoissonT> (d,34); } 00605 {RandPoissonT d(new E(125616),110.0); 00606 stat |= checkSaveDistribution<E,RandPoissonT> (d,34); } 00607 {RandPoissonT d(new E(125916),2.5); 00608 stat |= checkSaveDistribution<E,RandPoissonT> (d,10); } 00609 {RandPoissonT d(new E(125816),2.5); 00610 stat |= checkSaveDistribution<E,RandPoissonT> (d,11); } 00611 {RandPoissonT d(new E(125716),2.5); 00612 stat |= checkSaveDistribution<E,RandPoissonT> (d,12); } 00613 00614 {std::vector<double> pdf; 00615 int nbins = 20; 00616 for (int i = 0; i < nbins; ++i) 00617 pdf.push_back( 5*i + (10.5-i) * (10.5-i) ); 00618 RandGeneral d(new E(125917), &pdf[0], 20); 00619 stat |= checkRandGeneralDistribution<E> (d,33); } 00620 00621 return stat; 00622 } 00623 00624 template <class E, class D1, class D2> 00625 int checkSharingDistributions(D1 & d1, D2 & d2, int n1, int n2) { 00626 int stat = 0; 00627 output << "checkSharingDistribution with: \n" 00628 << d1.name() << " using " << d1.engine().name() << " and\n" 00629 << d2.name() << " using " << d2.engine().name() << "\n"; 00630 double r=0; 00631 r = d1(); 00632 r += d2(); 00633 double kv11,kv12,kv13,kv14; 00634 double kv21,kv22,kv23,kv24; 00635 for (int i=0; i<n1; i++) r += d1(); 00636 for (int j=0; j<n2; j++) r += d2(); 00637 {std::ofstream os ("shared.save1"); os << d1.engine() << d1 << d2;} 00638 kv11 = d1(); 00639 kv21 = d2(); 00640 kv12 = d1(); 00641 kv22 = d2(); 00642 r += d1() + d2(); 00643 // A second capture will test non-cached if first tested cached case: 00644 {std::ofstream os ("shared.save2"); os << d1.engine() << d1 << d2;} 00645 kv13 = d1(); 00646 kv23 = d2(); 00647 kv14 = d1(); 00648 kv24 = d2(); 00649 #ifdef VERBOSER2 00650 int pr = output.precision(20); 00651 output << "kv11 = " << kv11 << 00652 " kv21 = " << kv21 << "\n"; 00653 output << "kv12 = " << kv12 << 00654 " kv22 = " << kv22 << "\n"; 00655 output << "kv13 = " << kv13 << 00656 " kv23 = " << kv23 << "\n"; 00657 output << "kv14 = " << kv14 << 00658 " kv24 = " << kv24 << "\n"; 00659 output.precision(pr); 00660 #endif 00661 E e; 00662 D1 d1r(e); 00663 D2 d2r(e); 00664 { std::ifstream is ("shared.save1"); is >> e >> d1r >> d2r;} 00665 double k11 = d1r(); 00666 double k21 = d2r(); 00667 double k12 = d1r(); 00668 double k22 = d2r(); 00669 { std::ifstream is ("shared.save2"); is >> e >> d1r >> d2r;} 00670 double k13 = d1r(); 00671 double k23 = d2r(); 00672 double k14 = d1r(); 00673 double k24 = d2r(); 00674 #ifdef VERBOSER2 00675 pr = output.precision(20); 00676 output << "k11 = " << k11 << 00677 " k21 = " << k21 << "\n"; 00678 output << "k12 = " << k12 << 00679 " k22 = " << k22 << "\n"; 00680 output << "k13 = " << k13 << 00681 " k23 = " << k23 << "\n"; 00682 output << "k14 = " << k14 << 00683 " k24 = " << k24 << "\n"; 00684 output.precision(pr); 00685 #endif 00686 if ( !equals(k11,kv11) || !equals(k21,kv21) || 00687 !equals(k12,kv12) || !equals(k22,kv22) || 00688 !equals(k13,kv13) || !equals(k23,kv23) || 00689 !equals(k14,kv14) || !equals(k24,kv24) ) { 00690 std::cout << "???? Incorrect restored value for distributions " 00691 << d1.name() << " " << d2.name() << "\n"; 00692 stat |= 4096; 00693 } 00694 // if (stat) exit(-1); 00695 return stat; 00696 } 00697 00698 00699 00700 template <class E> 00701 int checkSharing() { 00702 int stat = 0; 00703 E e1(98746); 00704 RandGauss g1(e1,50.0,4.0); 00705 RandPoissonQ p1(e1, 112.0); 00706 RandGauss g2(e1,5.0,44.0); 00707 RandPoissonQ p2(e1, 212.0); 00708 stat |= checkSharingDistributions<E, RandGauss, RandPoissonQ>(g1,p1,5,4); 00709 stat |= checkSharingDistributions<E, RandGauss, RandPoissonQ>(g1,p2,6,6); 00710 stat |= checkSharingDistributions<E, RandGauss, RandPoissonQ>(g2,p1,8,9); 00711 stat |= checkSharingDistributions<E, RandGauss, RandPoissonQ>(g1,p1,7,5); 00712 stat |= checkSharingDistributions<E, RandPoissonQ, RandGauss>(p1,g2,5,4); 00713 stat |= checkSharingDistributions<E, RandPoissonQ, RandGauss>(p2,g1,6,6); 00714 stat |= checkSharingDistributions<E, RandPoissonQ, RandGauss>(p1,g1,8,9); 00715 stat |= checkSharingDistributions<E, RandPoissonQ, RandGauss>(p2,g1,7,5); 00716 return stat; 00717 } 00718 00719 std::vector<double> aSequence(int n) { 00720 std::vector<double> v; 00721 DualRand e(13542); 00722 RandFlat f(e); 00723 for (int i=0; i<n; i++) { 00724 v.push_back(f()); 00725 } 00726 return v; 00727 } 00728 00729 // ----------- Tests for static methods ----------- 00730 00731 template <class D> 00732 int staticSave(int n) { 00733 int stat = 0; 00734 int i; 00735 output << "staticSave for distribution " << D::distributionName() << "\n"; 00736 double r = 0; 00737 double v1, v2, k1, k2; 00738 for (i=0; i < n; i++) r += D::shoot(); 00739 { 00740 std::ofstream file ("distribution.save1"); 00741 D::saveFullState(file); 00742 v1 = D::shoot(); 00743 D::saveFullState(file); 00744 v2 = D::shoot(); 00745 #ifdef VERBOSER2 00746 int pr = output.precision(20); 00747 output << "v1 = " << v1 << " v2 = " << v2 << "\n"; 00748 output.precision(pr); 00749 #endif 00750 } 00751 for (i=0; i < n; i++) r += D::shoot(); 00752 { 00753 std::ifstream file ("distribution.save1"); 00754 D::restoreFullState(file); 00755 k1 = D::shoot(); 00756 for (i=0; i < n; i++) r += D::shoot(); 00757 D::restoreFullState(file); 00758 k2 = D::shoot(); 00759 #ifdef VERBOSER2 00760 int pr = output.precision(20); 00761 output << "k1 = " << k1 << " k2 = " << k2 << "\n"; 00762 output.precision(pr); 00763 #endif 00764 } 00765 if ( (k1 != v1) || (k2 != v2) ) { 00766 std::cout << "???? restoreFullState failed for " << D::distributionName() << "\n"; 00767 stat |= 8192; 00768 } 00769 00770 for (i=0; i < n; i++) r += D::shoot(); 00771 { 00772 std::ofstream file ("distribution.save2"); 00773 D::saveDistState(file) << *D::getTheEngine(); 00774 v1 = D::shoot(); 00775 D::saveDistState(file) << *D::getTheEngine(); 00776 v2 = D::shoot(); 00777 #ifdef VERBOSER2 00778 int pr = output.precision(20); 00779 output << "v1 = " << v1 << " v2 = " << v2 << "\n"; 00780 output.precision(pr); 00781 #endif 00782 } 00783 for (i=0; i < n; i++) r += D::shoot(); 00784 { 00785 std::ifstream file ("distribution.save2"); 00786 D::restoreDistState(file) >> *D::getTheEngine(); 00787 k1 = D::shoot(); 00788 for (i=0; i < n; i++) r += D::shoot(); 00789 D::restoreDistState(file) >> *D::getTheEngine(); 00790 k2 = D::shoot(); 00791 #ifdef VERBOSER2 00792 int pr = output.precision(20); 00793 output << "k1 = " << k1 << " k2 = " << k2 << "\n"; 00794 output.precision(pr); 00795 #endif 00796 } 00797 if ( (k1 != v1) || (k2 != v2) ) { 00798 std::cout << "???? restoreDistState failed for " << D::distributionName() << "\n"; 00799 stat |= 16384; 00800 } 00801 00802 return stat; 00803 } 00804 00805 template <class D> 00806 int staticSaveShootBit(int n) { 00807 int stat = 0; 00808 int i; 00809 output << "staticSaveShootBit for " << D::distributionName() << "\n"; 00810 double r = 0; 00811 int bit = 0; 00812 int v1, v2, k1, k2; 00813 for (i=0; i < n; i++) r += D::shoot(); 00814 for (i=0; i < n; i++) bit |= D::shootBit(); 00815 { 00816 std::ofstream file ("distribution.save1"); 00817 D::saveFullState(file); 00818 v1=0; 00819 for (i=0; i<25; i++) { 00820 v1 <<=1; 00821 v1 += D::shootBit(); 00822 } 00823 for (i=0; i < n; i++) bit |= D::shootBit(); 00824 D::saveFullState(file); 00825 v2=0; 00826 for (i=0; i<25; i++) { 00827 v2 <<=1; 00828 v2 += D::shootBit(); 00829 } 00830 #ifdef VERBOSER2 00831 int pr = output.precision(20); 00832 output << std::hex << "v1 = " << v1 << " v2 = " << v2 << std::dec << "\n"; 00833 output.precision(pr); 00834 #endif 00835 } 00836 for (i=0; i < n; i++) r += D::shoot(); 00837 { 00838 std::ifstream file ("distribution.save1"); 00839 D::restoreFullState(file); 00840 k1=0; 00841 for (i=0; i<25; i++) { 00842 k1 <<=1; 00843 k1 += D::shootBit(); 00844 } 00845 for (i=0; i < n; i++) r += D::shoot(); 00846 D::restoreFullState(file); 00847 k2=0; 00848 for (i=0; i<25; i++) { 00849 k2 <<=1; 00850 k2 += D::shootBit(); 00851 } 00852 #ifdef VERBOSER2 00853 int pr = output.precision(20); 00854 output << std::hex << "k1 = " << k1 << " k2 = " << k2 << std::dec << "\n"; 00855 output.precision(pr); 00856 #endif 00857 } 00858 if ( (k1 != v1) || (k2 != v2) ) { 00859 std::cout << "???? restoreFullState failed for D shootBit()\n"; 00860 stat |= 32768; 00861 } 00862 00863 for (i=0; i < n; i++) r += D::shoot(); 00864 for (i=0; i < n; i++) bit |= D::shootBit(); 00865 { 00866 std::ofstream file ("distribution.save2"); 00867 D::saveDistState(file) << *D::getTheEngine(); 00868 v1=0; 00869 for (i=0; i<25; i++) { 00870 v1 <<=1; 00871 v1 += D::shootBit(); 00872 } 00873 for (i=0; i < n; i++) bit |= D::shootBit(); 00874 D::saveDistState(file) << *D::getTheEngine(); 00875 v2=0; 00876 for (i=0; i<25; i++) { 00877 v2 <<=1; 00878 v2 += D::shootBit(); 00879 } 00880 #ifdef VERBOSER2 00881 int pr = output.precision(20); 00882 output << std::hex << "v1 = " << v1 << " v2 = " << v2 << std::dec << "\n"; 00883 output.precision(pr); 00884 #endif 00885 } 00886 for (i=0; i < n; i++) r += D::shoot(); 00887 { 00888 std::ifstream file ("distribution.save2"); 00889 D::restoreDistState(file) >> *D::getTheEngine(); 00890 k1=0; 00891 for (i=0; i<25; i++) { 00892 k1 <<=1; 00893 k1 += D::shootBit(); 00894 } 00895 for (i=0; i < n; i++) r += D::shoot(); 00896 for (i=0; i < n; i++) r += D::shootBit(); 00897 D::restoreDistState(file) >> *D::getTheEngine(); 00898 k2=0; 00899 for (i=0; i<25; i++) { 00900 k2 <<=1; 00901 k2 += D::shootBit(); 00902 } 00903 #ifdef VERBOSER2 00904 int pr = output.precision(20); 00905 output << std::hex << "k1 = " << k1 << " k2 = " << k2 << std::dec << "\n"; 00906 output.precision(pr); 00907 #endif 00908 } 00909 if ( (k1 != v1) || (k2 != v2) ) { 00910 std::cout << "???? restoreDistState failed for RnadFlat::shootBit()\n"; 00911 stat |= 65536; 00912 } 00913 00914 return stat; 00915 } 00916 00917 // ----------- Tests saving all statics together ----------- 00918 00919 void randomizeStatics(int n) { 00920 for (int i=0; i<n; i++) { 00921 RandGauss::shoot(); 00922 RandGaussQ::shoot(); 00923 RandGaussT::shoot(); 00924 RandFlat::shoot(); 00925 RandBit::shoot(); 00926 RandFlat::shootBit(); 00927 RandBit::shootBit(); 00928 RandPoisson::shoot(); 00929 RandPoissonQ::shoot(); 00930 RandPoissonT::shoot(); 00931 RandBinomial::shoot(); 00932 RandBreitWigner::shoot(); 00933 RandChiSquare::shoot(); 00934 RandExponential::shoot(); 00935 RandGamma::shoot(); 00936 RandLandau::shoot(); 00937 RandStudentT::shoot(); 00938 } 00939 } 00940 00941 std::vector<double> captureStatics() { 00942 std::vector<double> c; 00943 c.push_back( RandGauss::shoot() ); 00944 c.push_back( RandGaussQ::shoot() ); 00945 c.push_back( RandGaussT::shoot() ); 00946 c.push_back( RandFlat::shoot() ); 00947 c.push_back( RandBit::shoot() ); 00948 for (int i=0; i<20; i++) { 00949 c.push_back( RandFlat::shootBit() ); 00950 c.push_back( RandBit::shootBit() ); 00951 } 00952 c.push_back( RandPoisson::shoot() ); 00953 c.push_back( RandPoissonQ::shoot() ); 00954 c.push_back( RandPoissonT::shoot() ); 00955 c.push_back( RandBinomial::shoot() ); 00956 c.push_back( RandBreitWigner::shoot() ); 00957 c.push_back( RandChiSquare::shoot() ); 00958 c.push_back( RandExponential::shoot() ); 00959 c.push_back( RandGamma::shoot() ); 00960 c.push_back( RandLandau::shoot() ); 00961 c.push_back( RandStudentT::shoot() ); 00962 return c; 00963 } 00964 00965 void saveStatics(std::string filename) { 00966 std::ofstream os(filename.c_str()); 00967 RandGeneral::saveStaticRandomStates(os); 00968 // It should be possible to call this from HepRandom, or any distribution. 00969 // RandGeneral, which is meaningless as a static distribution, should be the 00970 // toughest test, so we use that here. 00971 } 00972 00973 void restoreStatics(std::string filename) { 00974 std::ifstream is(filename.c_str()); 00975 RandLandau::restoreStaticRandomStates(is); 00976 } 00977 00978 // ----------- Anonymous restore of engines ----------- 00979 00980 template <class E> 00981 void anonymousRestore1(int n, std::vector<double> & v) { 00982 output << "Anonymous restore for " << E::engineName() << "\n"; 00983 E e(12349876); 00984 double r=0; 00985 for (int i=0; i<n; i++) r += e.flat(); 00986 std::ofstream os("anonymous.save"); 00987 os << e; 00988 for (int j=0; j<25; j++) v.push_back(e.flat()); 00989 #ifdef VERBOSER2 00990 output << "First four of v are: " 00991 << v[0] << ", " << v[1] << ", " << v[2] << ", " << v[3] << "\n"; 00992 #endif 00993 return; 00994 } 00995 00996 template <> 00997 void anonymousRestore1<NonRandomEngine> (int n, std::vector<double> & v) { 00998 #ifdef VERBOSER 00999 output << "Anonymous restore for " << NonRandomEngine::engineName() << "\n"; 01000 #endif 01001 std::vector<double> nonRand = aSequence(500); 01002 NonRandomEngine e; 01003 e.setRandomSequence(&nonRand[0], nonRand.size()); 01004 double r=0; 01005 for (int i=0; i<n; i++) r += e.flat(); 01006 std::ofstream os("anonymous.save"); 01007 os << e; 01008 for (int j=0; j<25; j++) v.push_back(e.flat()); 01009 #ifdef VERBOSER2 01010 output << "First four of v are: " 01011 << v[0] << ", " << v[1] << ", " << v[2] << ", " << v[3] << "\n"; 01012 #endif 01013 return; 01014 } 01015 01016 template <class E> 01017 int anonymousRestore2(const std::vector<double> & v) { 01018 int stat = 0; 01019 std::vector<double> k; 01020 std::ifstream is("anonymous.save"); 01021 HepRandomEngine * a; 01022 a = HepRandomEngine::newEngine(is); 01023 for (int j=0; j<25; j++) k.push_back(a->flat()); 01024 delete a; 01025 #ifdef VERBOSER2 01026 output << "First four of k are: " 01027 << k[0] << ", " << k[1] << ", " << k[2] << ", " << k[3] << "\n"; 01028 #endif 01029 for (int m=0; m<25; m++) { 01030 if ( v[m] != k[m] ) { 01031 std::cout << "???? Incorrect restored value for anonymous engine" 01032 << E::engineName() << "\n"; 01033 stat |= 262144; 01034 return stat; 01035 } 01036 } 01037 return stat; 01038 } 01039 01040 01041 template <class E> 01042 int anonymousRestore(int n) { 01043 std::vector<double> v; 01044 anonymousRestore1<E>(n,v); 01045 return anonymousRestore2<E>(v); 01046 } 01047 01048 // ----------- Anonymous restore of all static distributions ----------- 01049 01050 template <class E> 01051 int anonymousRestoreStatics1() { 01052 int stat = 0; 01053 HepRandomEngine *e = new E(12456); 01054 HepRandom::setTheEngine(e); 01055 randomizeStatics(15); 01056 output << "\nRandomized, with theEngine = " << e->name() << "\n"; 01057 saveStatics("distribution.save"); 01058 output << "Saved all static distributions\n"; 01059 std::vector<double> c = captureStatics(); 01060 output << "Captured output of all static distributions\n"; 01061 randomizeStatics(11); 01062 output << "Randomized all static distributions\n"; 01063 restoreStatics("distribution.save"); 01064 output << "Restored all static distributions to saved state\n"; 01065 std::vector<double> d = captureStatics(); 01066 output << "Captured output of all static distributions\n"; 01067 for (unsigned int iv=0; iv<c.size(); iv++) { 01068 if (c[iv] != d[iv]) { 01069 std::cout << "???? restoreStaticRandomStates failed at random " 01070 << iv <<"\n"; 01071 stat |= 131072; 01072 } 01073 } 01074 if (stat & 131072 == 0) { 01075 output << "All captured output agrees with earlier values\n"; 01076 } 01077 return stat; 01078 } 01079 01080 01081 01082 template <class E1, class E2> 01083 int anonymousRestoreStatics() { 01084 int stat = 0; 01085 if ( E1::engineName() == E2::engineName() ) { 01086 return anonymousRestoreStatics1<E1>(); 01087 } 01088 HepRandomEngine *e1 = new E1(12456); 01089 HepRandom::setTheEngine(e1); 01090 randomizeStatics(15); 01091 output << "\nRandomized, with theEngine = " << e1->name() << "\n"; 01092 saveStatics("distribution.save"); 01093 #ifdef VERBOSER2 01094 output << "Saved all static distributions\n"; 01095 #endif 01096 std::vector<double> c = captureStatics(); 01097 #ifdef VERBOSER2 01098 output << "Captured output of all static distributions\n"; 01099 #endif 01100 delete e1; 01101 HepRandomEngine *e2 = new E2(24653); 01102 HepRandom::setTheEngine(e2); 01103 output << "Switched to theEngine = " << e2->name() << "\n"; 01104 randomizeStatics(19); 01105 { std::ofstream os("engine.save"); os << *e2; } 01106 double v1 = e2->flat(); 01107 double v2 = e2->flat(); 01108 { std::ifstream is("engine.save"); is >> *e2; } 01109 #ifdef VERBOSER2 01110 output << "Saved the " << e2->name() << " engine: \n" 01111 << "Next randoms to be " << v1 << " " << v2 << "\n" 01112 << "Restored the " << e2->name() << " engine to that state\n"; 01113 #endif 01114 restoreStatics("distribution.save"); 01115 #ifdef VERBOSER2 01116 output << "Restored all static distributions to saved state\n" 01117 << "This changes the engine type back to " << E1::engineName() << "\n"; 01118 #endif 01119 std::vector<double> d = captureStatics(); 01120 #ifdef VERBOSER2 01121 output << "Captured output of all static distributions\n"; 01122 #endif 01123 for (unsigned int iv=0; iv<c.size(); iv++) { 01124 if (c[iv] != d[iv]) { 01125 std::cout << "???? restoreStaticRandomStates failed at random " 01126 << iv <<"\n"; 01127 stat |= 524288; 01128 } 01129 } 01130 if (stat & 524288 == 0) { 01131 output << "All captured output agrees with earlier values\n"; 01132 } 01133 double k1 = e2->flat(); 01134 double k2 = e2->flat(); 01135 #ifdef VERBOSER2 01136 output << "The " << e2->name() << " engine should not have been affected: \n" 01137 << "Next randoms are " << k1 << " " << k2 << "\n"; 01138 #endif 01139 if ( !equals(v1,k1) || !equals(v2,k2) ) { 01140 std::cout << "???? Engine used as theEngine was affected by restoring \n" 01141 << " static distributions to use engine of a different type.\n"; 01142 stat |= 1048576; 01143 } 01144 return stat; 01145 } 01146 01147 // ----------- Vector restore of engines ----------- 01148 01149 template <class E> 01150 std::vector<unsigned long> vectorRestore1(int n, std::vector<double> & v) { 01151 output << "Vector restore for " << E::engineName() << "\n"; 01152 E e(97538466); 01153 double r=0; 01154 for (int i=0; i<n; i++) r += e.flat(); 01155 std::vector<unsigned long> state = e.put(); 01156 for (int j=0; j<25; j++) v.push_back(e.flat()); 01157 #ifdef VERBOSER2 01158 output << "First four of v are: " 01159 << v[0] << ", " << v[1] << ", " << v[2] << ", " << v[3] << "\n"; 01160 #endif 01161 return state; 01162 } 01163 01164 template <> 01165 std::vector<unsigned long> 01166 vectorRestore1<NonRandomEngine> (int n, std::vector<double> & v) { 01167 #ifdef VERBOSER2 01168 output << "Vector restore for " << NonRandomEngine::engineName() << "\n"; 01169 #endif 01170 std::vector<double> nonRand = aSequence(500); 01171 NonRandomEngine e; 01172 e.setRandomSequence(&nonRand[0], nonRand.size()); 01173 double r=0; 01174 for (int i=0; i<n; i++) r += e.flat(); 01175 std::vector<unsigned long> state = e.put(); 01176 for (int j=0; j<25; j++) v.push_back(e.flat()); 01177 #ifdef VERBOSER2 01178 output << "First four of v are: " 01179 << v[0] << ", " << v[1] << ", " << v[2] << ", " << v[3] << "\n"; 01180 #endif 01181 return state; 01182 } 01183 01184 template <class E> 01185 int vectorRestore2(const std::vector<unsigned long> state, 01186 const std::vector<double> & v) { 01187 int stat = 0; 01188 std::vector<double> k; 01189 HepRandomEngine * a; 01190 a = HepRandomEngine::newEngine(state); 01191 if (!a) { 01192 std::cout << "???? could not restore engine state from vector for " 01193 << E::engineName() << "\n"; 01194 stat |= 1048576; 01195 return stat; 01196 } 01197 if (a->name() != E::engineName()) { 01198 std::cout << "???? restored engine state from vector for " 01199 << E::engineName() << "to different type of engine: " 01200 << a->name() << "\n" 01201 << "There is probably a clash in CRC hashes for these two names!\n"; 01202 stat |= 1048576; 01203 return stat; 01204 } 01205 for (int j=0; j<25; j++) k.push_back(a->flat()); 01206 delete a; 01207 #ifdef VERBOSER2 01208 output << "First four of k are: " 01209 << k[0] << ", " << k[1] << ", " << k[2] << ", " << k[3] << "\n"; 01210 #endif 01211 for (int m=0; m<25; m++) { 01212 if ( v[m] != k[m] ) { 01213 std::cout << "???? Incorrect vector restored value for anonymous engine: " 01214 << E::engineName() << "\n"; 01215 stat |= 1048576; 01216 return stat; 01217 } 01218 } 01219 return stat; 01220 } 01221 01222 01223 template <class E> 01224 int vectorRestore(int n) { 01225 std::vector<double> v; 01226 std::vector<unsigned long> state = vectorRestore1<E>(n,v); 01227 return vectorRestore2<E>(state, v); 01228 } 01229 01230 01231 01232 // --------------------------------------------- 01233 // --------------------------------------------- 01234 // --------------------------------------------- 01235 01236 01237 int main() { 01238 int stat = 0; 01239 01240 #ifdef TEST_ORIGINAL_SAVE 01241 output << "=====================================\n"; 01242 output << " Part I \n"; 01243 output << "Original tests of static save/restore\n"; 01244 output << "=====================================\n\n"; 01245 01246 output << "Using old method or HepRandom::saveEngineStatus:\n"; 01247 output << "All these tests should have a chance of failure.\n"; 01248 01249 output << RandGauss:: getTheEngine()->name(); 01250 output << RandGaussQ::getTheEngine()->name(); 01251 01252 stat |= saveStepX(); 01253 stat |= restoreStepX(); 01254 stat |= BsaveStepX(); 01255 stat |= BrestoreStepX(); 01256 01257 output << "Using the class-specific RandGauss::saveEngineStatus:\n"; 01258 output << "All these tests should work properly.\n"; 01259 01260 stat |= saveStep(); 01261 stat |= restoreStep(); 01262 stat |= BsaveStep(); 01263 stat |= BrestoreStep(); 01264 #endif 01265 01266 #ifdef TEST_MISSING_FILES 01267 output << "\n=======================================\n"; 01268 output << " Part Ia \n"; 01269 output << "Test of behavior when a file is missing \n"; 01270 output << "=======================================\n\n"; 01271 01272 output << "Testing restoreEngineStatus with missing file:\n"; 01273 output << "Expect a number of <Failure to find or open> messages!\n"; 01274 stat |= missingFile(); 01275 #endif 01276 01277 #ifdef CREATE_OLD_SAVES 01278 stat |= saveEngine<DRand48Engine, RandPoisson>("DRand48Engine.oldsav"); 01279 stat |= saveEngine<DualRand, RandPoisson>("DualRand.oldsav"); 01280 stat |= saveEngine<Hurd160Engine, RandPoisson>("Hurd160Engine.oldsav"); 01281 stat |= saveEngine<Hurd288Engine, RandPoisson>("Hurd288Engine.oldsav"); 01282 stat |= saveEngine<HepJamesRandom,RandPoisson>("HepJamesRandom.oldsav"); 01283 stat |= saveEngine<MTwistEngine, RandPoisson>("MTwistEngine.oldsav"); 01284 stat |= saveEngine<RanecuEngine, RandPoisson>("RanecuEngine.oldsav"); 01285 stat |= saveEngine<Ranlux64Engine,RandPoisson>("Ranlux64Engine.oldsav"); 01286 stat |= saveEngine<RanluxEngine, RandPoisson>("RanluxEngine.oldsav"); 01287 stat |= saveEngine<RanshiEngine, RandPoisson>("RanshiEngine.oldsav"); 01288 stat |= saveEngine<TripleRand, RandPoisson>("TripleRand.oldsav"); 01289 #endif 01290 01291 #ifdef VERIFY_OLD_SAVES 01292 output << "\n==============================================\n"; 01293 output << " Part Ib \n"; 01294 output << " Verification that changes wont invalidate \n"; 01295 output << "invalidate engine saves from previous versions \n"; 01296 output << "==============================================\n\n"; 01297 01298 stat |= checkSaveEngine<DRand48Engine, RandPoisson>("DRand48Engine.oldsav"); 01299 stat |= checkSaveEngine<DualRand, RandPoisson>("DualRand.oldsav"); 01300 stat |= checkSaveEngine<Hurd160Engine, RandPoisson>("Hurd160Engine.oldsav"); 01301 stat |= checkSaveEngine<Hurd288Engine, RandPoisson>("Hurd288Engine.oldsav"); 01302 stat |= checkSaveEngine<HepJamesRandom,RandPoisson>("HepJamesRandom.oldsav"); 01303 stat |= checkSaveEngine<MTwistEngine, RandPoisson>("MTwistEngine.oldsav"); 01304 stat |= checkSaveEngine<Ranlux64Engine,RandPoisson>("Ranlux64Engine.oldsav"); 01305 stat |= checkSaveEngine<RanluxEngine, RandPoisson>("RanluxEngine.oldsav"); 01306 stat |= checkSaveEngine<RanshiEngine, RandPoisson>("RanshiEngine.oldsav"); 01307 stat |= checkSaveEngine<TripleRand, RandPoisson>("TripleRand.oldsav"); 01308 stat |= checkSaveEngine<RanecuEngine, RandPoisson>("RanecuEngine.oldsav"); 01309 #endif 01310 01311 #ifdef TEST_ENGINE_NAMES 01312 output << "\n=============================================\n"; 01313 output << " Part II \n"; 01314 output << "Check all engine names were entered correctly \n"; 01315 output << "=============================================\n\n"; 01316 01317 stat |= checkEngineName<DRand48Engine >("DRand48Engine"); 01318 stat |= checkEngineName<DualRand >("DualRand"); 01319 stat |= checkEngineName<Hurd160Engine >("Hurd160Engine"); 01320 stat |= checkEngineName<Hurd288Engine >("Hurd288Engine"); 01321 stat |= checkEngineName<HepJamesRandom>("HepJamesRandom"); 01322 stat |= checkEngineName<MTwistEngine >("MTwistEngine"); 01323 stat |= checkEngineName<RandEngine >("RandEngine"); 01324 stat |= checkEngineName<RanecuEngine >("RanecuEngine"); 01325 stat |= checkEngineName<Ranlux64Engine>("Ranlux64Engine"); 01326 stat |= checkEngineName<RanluxEngine >("RanluxEngine"); 01327 stat |= checkEngineName<RanshiEngine >("RanshiEngine"); 01328 stat |= checkEngineName<TripleRand >("TripleRand"); 01329 #endif 01330 01331 #ifdef TEST_INSTANCE_METHODS 01332 output << "===========================================\n\n"; 01333 output << " Part III\n"; 01334 output << "Check instance methods for specific engines \n"; 01335 output << " specific engines and distributions\n"; 01336 output << "===========================================\n\n"; 01337 01338 {DualRand e(234); stat |= checkEngineInstanceSave(e);} 01339 {Hurd160Engine e(234); stat |= checkEngineInstanceSave(e);} 01340 {Hurd288Engine e(234); stat |= checkEngineInstanceSave(e);} 01341 {HepJamesRandom e(234); stat |= checkEngineInstanceSave(e);} 01342 {MTwistEngine e(234); stat |= checkEngineInstanceSave(e);} 01343 {RandEngine e(234); stat |= checkEngineInstanceSave(e);} 01344 {RanecuEngine e(234); stat |= checkEngineInstanceSave(e);} 01345 {Ranlux64Engine e(234); stat |= checkEngineInstanceSave(e);} 01346 {RanluxEngine e(234); stat |= checkEngineInstanceSave(e);} 01347 {RanshiEngine e(234); stat |= checkEngineInstanceSave(e);} 01348 {TripleRand e(234); stat |= checkEngineInstanceSave(e);} 01349 01350 {std::vector<double> nonRand = aSequence(500); 01351 NonRandomEngine e; 01352 e.setRandomSequence(&nonRand[0], nonRand.size()); 01353 stat |= checkEngineInstanceSave(e);} 01354 01355 stat |= checkDistributions<DualRand>(); 01356 stat |= checkDistributions<Hurd160Engine>(); 01357 stat |= checkDistributions<Hurd288Engine>(); 01358 stat |= checkDistributions<HepJamesRandom>(); 01359 stat |= checkDistributions<MTwistEngine>(); 01360 stat |= checkDistributions<Ranlux64Engine>(); 01361 stat |= checkDistributions<RanluxEngine>(); 01362 stat |= checkDistributions<RanshiEngine>(); 01363 stat |= checkDistributions<TripleRand>(); 01364 01365 RandGaussQ::shoot(); // Just to verify that the static engine is OK 01366 #endif 01367 01368 #ifdef TEST_SHARED_ENGINES 01369 output << "\n=============================================\n"; 01370 output << " Part IV \n"; 01371 output << "Check behavior when engines are shared \n"; 01372 output << "=============================================\n\n"; 01373 01374 stat |= checkSharing<DualRand>(); 01375 stat |= checkSharing<Hurd160Engine>(); 01376 stat |= checkSharing<Hurd288Engine>(); 01377 stat |= checkSharing<HepJamesRandom>(); 01378 stat |= checkSharing<MTwistEngine>(); 01379 stat |= checkSharing<Ranlux64Engine>(); 01380 stat |= checkSharing<RanluxEngine>(); 01381 stat |= checkSharing<RanshiEngine>(); 01382 stat |= checkSharing<TripleRand>(); 01383 #endif 01384 01385 #ifdef TEST_STATIC_SAVE 01386 output << "\n=========================================\n"; 01387 output << " Part V \n"; 01388 output << "Static Save/restore to/from streams \n"; 01389 output << "=========================================\n\n"; 01390 01391 stat |= staticSave <RandGauss>(7); 01392 stat |= staticSave <RandFlat>(7); 01393 stat |= staticSaveShootBit<RandFlat> (19); 01394 stat |= staticSaveShootBit<RandBit> (23); 01395 for (int ibinom=0; ibinom<15; ibinom++) { 01396 stat |= staticSave <RandBinomial>(7+3*ibinom); 01397 } 01398 stat |= staticSave <RandBreitWigner>(7); 01399 stat |= staticSave <RandChiSquare>(7); 01400 stat |= staticSave <RandExponential>(7); 01401 stat |= staticSave <RandGamma>(7); 01402 stat |= staticSave <RandGaussQ>(7); 01403 stat |= staticSave <RandGaussT>(7); 01404 stat |= staticSave <RandLandau>(7); 01405 stat |= staticSave <RandPoisson>(7); 01406 stat |= staticSave <RandPoissonQ>(7); 01407 stat |= staticSave <RandPoissonT>(7); 01408 stat |= staticSave <RandStudentT>(7); 01409 #endif 01410 01411 #ifdef TEST_SAVE_STATIC_STATES 01412 output << "\n==============================================\n"; 01413 output << " Part VI \n"; 01414 output << "Save/restore all static states to/from streams \n"; 01415 output << "==============================================\n\n"; 01416 01417 randomizeStatics(15); 01418 saveStatics("distribution.save"); 01419 output << "Saved all static distributions\n"; 01420 std::vector<double> c = captureStatics(); 01421 output << "Captured output of all static distributions\n"; 01422 randomizeStatics(11); 01423 output << "Randomized all static distributions\n"; 01424 restoreStatics("distribution.save"); 01425 output << "Restored all static distributions to saved state\n"; 01426 std::vector<double> d = captureStatics(); 01427 output << "Captured output of all static distributions\n"; 01428 for (unsigned int iv=0; iv<c.size(); iv++) { 01429 if (c[iv] != d[iv]) { 01430 std::cout << "???? restoreStaticRandomStates failed at random " 01431 << iv <<"\n"; 01432 stat |= 131072; 01433 } 01434 } 01435 if (stat & 131072 == 0) { 01436 output << "All captured output agrees with earlier values\n"; 01437 } 01438 #endif 01439 01440 #ifdef TEST_ANONYMOUS_ENGINE_RESTORE 01441 output << "\n=================================\n"; 01442 output << " Part VII \n"; 01443 output << "Anonymous restore of engines \n"; 01444 output << "=================================\n\n"; 01445 01446 stat |= anonymousRestore<DualRand>(13); 01447 stat |= anonymousRestore<DRand48Engine>(14); 01448 stat |= anonymousRestore<Hurd160Engine>(15); 01449 stat |= anonymousRestore<Hurd288Engine>(16); 01450 stat |= anonymousRestore<HepJamesRandom>(17); 01451 stat |= anonymousRestore<MTwistEngine>(18); 01452 stat |= anonymousRestore<RandEngine>(29); 01453 stat |= anonymousRestore<RanecuEngine>(39); 01454 stat |= anonymousRestore<Ranlux64Engine>(19); 01455 stat |= anonymousRestore<RanluxEngine>(20); 01456 stat |= anonymousRestore<RanshiEngine>(21); 01457 stat |= anonymousRestore<TripleRand>(22); 01458 stat |= anonymousRestore<NonRandomEngine>(22); 01459 #endif 01460 01461 #ifdef TEST_ANONYMOUS_RESTORE_STATICS 01462 output << "\n======================================\n"; 01463 output << " Part VIII \n"; 01464 output << "Anonymous restore static Distributions \n"; 01465 output << "======================================\n\n"; 01466 01467 stat |= anonymousRestoreStatics<DualRand, Ranlux64Engine> ( ); 01468 stat |= anonymousRestoreStatics<DRand48Engine, TripleRand> ( ); 01469 stat |= anonymousRestoreStatics<RandEngine, Ranlux64Engine> ( ); 01470 stat |= anonymousRestoreStatics<MTwistEngine, Hurd288Engine> ( ); 01471 stat |= anonymousRestoreStatics<RanecuEngine, MTwistEngine> ( ); 01472 stat |= anonymousRestoreStatics<HepJamesRandom, RanshiEngine> ( ); 01473 stat |= anonymousRestoreStatics<RanecuEngine, RandEngine> ( ); 01474 stat |= anonymousRestoreStatics<RanshiEngine, Hurd160Engine> ( ); 01475 stat |= anonymousRestoreStatics<TripleRand, DualRand> ( ); 01476 stat |= anonymousRestoreStatics<Hurd160Engine, HepJamesRandom> ( ); 01477 stat |= anonymousRestoreStatics<Hurd288Engine, RanecuEngine> ( ); 01478 stat |= anonymousRestoreStatics<HepJamesRandom, Ranlux64Engine> ( ); 01479 stat |= anonymousRestoreStatics<TripleRand, TripleRand> ( ); 01480 stat |= anonymousRestoreStatics<HepJamesRandom, HepJamesRandom> ( ); 01481 #endif 01482 01483 #ifdef TEST_VECTOR_ENGINE_RESTORE 01484 output << "\n=================================\n"; 01485 output << " Part IX \n"; 01486 output << "Save/restore of engines to vectors\n"; 01487 output << "=================================\n\n"; 01488 01489 stat |= vectorRestore<DualRand>(113); 01490 stat |= vectorRestore<DRand48Engine>(114); 01491 stat |= vectorRestore<Hurd160Engine>(115); 01492 stat |= vectorRestore<Hurd288Engine>(116); 01493 stat |= vectorRestore<HepJamesRandom>(117); 01494 stat |= vectorRestore<MTwistEngine>(118); 01495 stat |= vectorRestore<RanecuEngine>(139); 01496 stat |= vectorRestore<Ranlux64Engine>(119); 01497 stat |= vectorRestore<RanluxEngine>(120); 01498 stat |= vectorRestore<RanshiEngine>(121); 01499 stat |= vectorRestore<TripleRand>(122); 01500 stat |= vectorRestore<NonRandomEngine>(123); 01501 stat |= vectorRestore<RandEngine>(129); 01502 #endif 01503 01504 01505 01506 01507 output << "\n=============================================\n\n"; 01508 01509 if (stat != 0) { 01510 std::cout << "One or more problems detected: stat = " << stat << "\n"; 01511 } else { 01512 output << "ranRestoreTest passed with no problems detected.\n"; 01513 } 01514 01515 return stat; 01516 } 01517