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 #include <cassert> 00009 00010 #define CLEAN_OUTPUT 00011 #ifdef CLEAN_OUTPUT 00012 std::ofstream output("testInstanceRestore.cout"); 00013 #else 00014 std::ostream & output = std::cout; 00015 #endif 00016 00017 // Normally on for routine validation: 00018 00019 00020 #ifdef TURNOFF 00021 #endif 00022 00023 #define TEST_ENGINE_NAMES 00024 #define TEST_INSTANCE_METHODS 00025 00026 #define VERBOSER 00027 #define VERBOSER2 00028 00029 using namespace CLHEP; 00030 00031 // Absolutely Safe Equals Without Registers Screwing Us Up 00032 bool equals01(const std::vector<double> &ab) { 00033 return ab[1]==ab[0]; 00034 } 00035 bool equals(double a, double b) { 00036 std::vector<double> ab(2); 00037 ab[0]=a; ab[1]=b; 00038 return (equals01(ab)); 00039 } 00040 00041 std::vector<double> aSequence(int n) { 00042 std::vector<double> v; 00043 DualRand e(13542); 00044 RandFlat f(e); 00045 for (int i=0; i<n; i++) { 00046 v.push_back(f()); 00047 } 00048 return v; 00049 } 00050 00051 // ----------- Tests for instance methods ----------- 00052 00053 template <class E> 00054 int checkEngineName(const std::string & name) { 00055 int stat = 0; 00056 output << E::engineName() << "\n"; 00057 if (E::engineName() != name) { 00058 std::cout << "???? engineName mismatch for " << name << " <--> " 00059 << E::engineName() << "\n"; 00060 #ifdef CLEAN_OUTPUT 00061 output << "???? engineName mismatch for " << name << " <--> " 00062 << E::engineName() << "\n"; 00063 #endif 00064 stat |= 256; 00065 } 00066 E e(123); 00067 if (e.name() != name) { 00068 std::cout << "???? name mismatch for " << name << " <--> " 00069 << e.name() << "\n"; 00070 #ifdef CLEAN_OUTPUT 00071 output << "???? name mismatch for " << name << " <--> " 00072 << e.name() << "\n"; 00073 #endif 00074 stat |= 256; 00075 } 00076 return stat; 00077 } 00078 00079 template <class E, class D> 00080 int checkEngine() { 00081 int stat = 0; 00082 E e(1234); 00083 D d(e); 00084 if (d.engine().name() != e.name()) { 00085 std::cout << "???? Improper d.engine() \n"; 00086 #ifdef CLEAN_OUTPUT 00087 output << "???? Improper d.engine() \n"; 00088 #endif 00089 stat |= 512; 00090 } 00091 return stat; 00092 } 00093 00094 template <class E> 00095 int checkEngineInstanceSave(E & e) { 00096 int stat = 0; 00097 output << "checkEngineInstanceSave for " << e.name() << "\n"; 00098 int pr=output.precision(20); 00099 double r=0; 00100 for (int i=0; i<100; i++) r += e.flat(); 00101 {std::ofstream os ("engine.save"); os << e;} 00102 for (int i=0; i<100; i++) r += e.flat(); 00103 double keyValue1 = e.flat(); 00104 double keyValue2 = e.flat(); 00105 #ifdef VERBOSER 00106 output << keyValue1 << " " << keyValue2 << "\n"; 00107 #endif 00108 E e2; 00109 {std::ifstream is ("engine.save"); is >> e2;} 00110 for (int i=0; i<100; i++) r += e2.flat(); 00111 double k1 = e2.flat(); 00112 double k2 = e2.flat(); 00113 #ifdef VERBOSER 00114 output << k1 << " " << k2 << "\n"; 00115 #endif 00116 if ( !(equals(k1,keyValue1)) || !(equals(k2,keyValue2)) ) { 00117 std::cout << "???? checkInstanceSave failed for " << e.name() << "\n"; 00118 #ifdef CLEAN_OUTPUT 00119 output << "???? checkInstanceSave failed for " << e.name() << "\n"; 00120 #endif 00121 stat |= 1024; 00122 } 00123 output.precision(pr); 00124 return stat; 00125 } 00126 00127 template <class E, class D> 00128 int checkSaveDistribution(D & d, int nth) { 00129 // verify that engine is the expected type 00130 assert( &dynamic_cast<E &>(d.engine()) ); 00131 int stat = 0; 00132 output << "checkSaveDistribution with " << d.engine().name() 00133 << ", " << d.name() << "\n"; 00134 double r=0; 00135 r = d(); 00136 double keyValue1, keyValue2, keyValue3, keyValue4; 00137 for (int i=0; i<nth; i++) r += d(); 00138 {std::ofstream os ("distribution.save1"); os << d.engine() << d;} 00139 keyValue1 = d(); 00140 keyValue2 = d(); 00141 r += d(); 00142 // A second capture will test non-cached if first tested cached case: 00143 {std::ofstream os ("distribution.save2"); os << d.engine() << d;} 00144 keyValue3 = d(); 00145 keyValue4 = d(); 00146 int pr = output.precision(20); 00147 #ifdef VERBOSER 00148 output << "keyValue1 = " << keyValue1 << 00149 " keyValue2 = " << keyValue2 << "\n"; 00150 output << "keyValue3 = " << keyValue3 << 00151 " keyValue3 = " << keyValue4 << "\n"; 00152 #endif 00153 output.precision(pr); 00154 E e; 00155 D d2(e); 00156 { std::ifstream is ("distribution.save1"); is >> e >> d2;} 00157 double k1 = d2(); 00158 double k2 = d2(); 00159 { std::ifstream is ("distribution.save2"); is >> e >> d2;} 00160 double k3 = d2(); 00161 double k4 = d2(); 00162 #ifdef VERBOSER 00163 pr = output.precision(20); 00164 output << "k1 = " << k1 << 00165 " k2 = " << k2 << "\n"; 00166 output << "k3 = " << k3 << 00167 " k4 = " << k4 << "\n"; 00168 output.precision(pr); 00169 #endif 00170 if ( !equals(k1,keyValue1) || !equals(k2,keyValue2) || 00171 !equals(k3,keyValue3) || !equals(k4,keyValue4) ) { 00172 std::cout << "???? Incorrect restored value for distribution " 00173 << d.name() << "\n"; 00174 #ifdef CLEAN_OUTPUT 00175 output << "???? Incorrect restored value for distribution " 00176 << d.name() << "\n"; 00177 #endif 00178 stat |= 2048; 00179 } 00180 // if (stat) exit(-1); 00181 return stat; 00182 } 00183 00184 template <class E> 00185 int checkRandGeneralDistribution(RandGeneral & d, int nth) { 00186 assert( &dynamic_cast<E &>(d.engine()) ); 00187 int stat = 0; 00188 output << "checkSaveDistribution with " << d.engine().name() 00189 << ", " << d.name() << "\n"; 00190 double r=0; 00191 r = d(); 00192 double keyValue1, keyValue2, keyValue3, keyValue4; 00193 for (int i=0; i<nth; i++) r += d(); 00194 {std::ofstream os ("distribution.save1"); os << d.engine() << d;} 00195 keyValue1 = d(); 00196 keyValue2 = d(); 00197 r += d(); 00198 // A second capture will test non-cached if first tested cached case: 00199 {std::ofstream os ("distribution.save2"); os << d.engine() << d;} 00200 keyValue3 = d(); 00201 keyValue4 = d(); 00202 int pr = output.precision(20); 00203 #ifdef VERBOSER 00204 output << "keyValue1 = " << keyValue1 << 00205 " keyValue2 = " << keyValue2 << "\n"; 00206 output << "keyValue3 = " << keyValue3 << 00207 " keyValue3 = " << keyValue4 << "\n"; 00208 #endif 00209 output.precision(pr); 00210 E e; 00211 double temp = 1; 00212 RandGeneral d2(e, &temp, 1); 00213 { std::ifstream is ("distribution.save1"); is >> e >> d2;} 00214 double k1 = d2(); 00215 double k2 = d2(); 00216 { std::ifstream is ("distribution.save2"); is >> e >> d2;} 00217 double k3 = d2(); 00218 double k4 = d2(); 00219 #ifdef VERBOSER 00220 pr = output.precision(20); 00221 output << "k1 = " << k1 << 00222 " k2 = " << k2 << "\n"; 00223 output << "k3 = " << k3 << 00224 " k4 = " << k4 << "\n"; 00225 output.precision(pr); 00226 #endif 00227 if ( !equals(k1,keyValue1) || !equals(k2,keyValue2) || 00228 !equals(k3,keyValue3) || !equals(k4,keyValue4) ) { 00229 std::cout << "???? Incorrect restored value for distribution " 00230 << d.name() << "\n"; 00231 #ifdef CLEAN_OUTPUT 00232 output << "???? Incorrect restored value for distribution " 00233 << d.name() << "\n"; 00234 #endif 00235 stat |= 2048; 00236 } 00237 // if (stat) exit(-1); 00238 return stat; 00239 } 00240 00241 template <class E> 00242 int checkDistributions() { 00243 int stat = 0; 00244 {RandGauss d(new E(12561),100.0,3.0); 00245 stat |= checkSaveDistribution<E,RandGauss> (d,33); } 00246 {RandGauss d(new E(12572),100.0,3.0); 00247 stat |= checkSaveDistribution<E,RandGauss> (d,34); } 00248 {RandGaussQ d(new E(12563),10.0,4.0); 00249 stat |= checkSaveDistribution<E,RandGaussQ> (d,33); } 00250 {RandGaussT d(new E(12564),5.0,2.0); 00251 stat |= checkSaveDistribution<E,RandGaussT> (d,33); } 00252 {RandBinomial d(new E(12565),4,0.6); 00253 stat |= checkSaveDistribution<E,RandBinomial> (d,33); } 00254 {RandFlat d(new E(12576),12.5,35.0); 00255 stat |= checkSaveDistribution<E,RandFlat> (d,33); } 00256 {RandBit d(new E(12567)); 00257 stat |= checkSaveDistribution<E,RandBit> (d,31); } 00258 {RandBit d(new E(12578)); 00259 stat |= checkSaveDistribution<E,RandBit> (d,32); } 00260 {RandBit d(new E(12589)); 00261 stat |= checkSaveDistribution<E,RandBit> (d,33); } 00262 {RandBreitWigner d(new E(125611),50.0,15.0); 00263 stat |= checkSaveDistribution<E,RandBreitWigner> (d,33); } 00264 {RandChiSquare d(new E(125612),5.0); 00265 stat |= checkSaveDistribution<E,RandChiSquare> (d,33); } 00266 {RandExponential d(new E(125713),8.00); 00267 stat |= checkSaveDistribution<E,RandExponential> (d,33); } 00268 {RandGamma d(new E(125713),6.0,2.0); 00269 stat |= checkSaveDistribution<E,RandGamma> (d,33); } 00270 {RandLandau d(new E(125714)); 00271 stat |= checkSaveDistribution<E,RandLandau> (d,33); } 00272 {RandStudentT d(new E(125715),5); 00273 stat |= checkSaveDistribution<E,RandStudentT> (d,33); } 00274 00275 // Multiple tests of Poisson distributions for small desired, since 00276 // the answer in each test is a small int, and coincidental agreement 00277 // is very possible. 00278 00279 {RandPoisson d(new E(125616),2.5); 00280 stat |= checkSaveDistribution<E,RandPoisson> (d,33); } 00281 {RandPoisson d(new E(125617),105.0); 00282 stat |= checkSaveDistribution<E,RandPoisson> (d,34); } 00283 {RandPoisson d(new E(125618),2.5); 00284 stat |= checkSaveDistribution<E,RandPoisson> (d,35); } 00285 {RandPoisson d(new E(325618),2.5); 00286 stat |= checkSaveDistribution<E,RandPoisson> (d,36); } 00287 {RandPoisson d(new E(425618),2.5); 00288 stat |= checkSaveDistribution<E,RandPoisson> (d,37); } 00289 {RandPoisson d(new E(525618),2.5); 00290 stat |= checkSaveDistribution<E,RandPoisson> (d,38); } 00291 {RandPoisson d(new E(125619),110.0); 00292 stat |= checkSaveDistribution<E,RandPoisson> (d,39); } 00293 {RandPoissonQ d(new E(124616),2.5); 00294 stat |= checkSaveDistribution<E,RandPoissonQ> (d,33); } 00295 {RandPoissonQ d(new E(126616),2.5); 00296 stat |= checkSaveDistribution<E,RandPoissonQ> (d,32); } 00297 {RandPoissonQ d(new E(127616),2.5); 00298 stat |= checkSaveDistribution<E,RandPoissonQ> (d,31); } 00299 {RandPoissonQ d(new E(129616),2.5); 00300 stat |= checkSaveDistribution<E,RandPoissonQ> (d,30); } 00301 {RandPoissonQ d(new E(125616),110.0); 00302 stat |= checkSaveDistribution<E,RandPoissonQ> (d,33); } 00303 {RandPoissonQ d(new E(125616),2.5); 00304 stat |= checkSaveDistribution<E,RandPoissonQ> (d,34); } 00305 {RandPoissonQ d(new E(125616),110.0); 00306 stat |= checkSaveDistribution<E,RandPoissonQ> (d,34); } 00307 {RandPoissonT d(new E(125616),2.5); 00308 stat |= checkSaveDistribution<E,RandPoissonT> (d,33); } 00309 {RandPoissonT d(new E(125616),110.0); 00310 stat |= checkSaveDistribution<E,RandPoissonT> (d,33); } 00311 {RandPoissonT d(new E(125616),2.5); 00312 stat |= checkSaveDistribution<E,RandPoissonT> (d,34); } 00313 {RandPoissonT d(new E(125616),110.0); 00314 stat |= checkSaveDistribution<E,RandPoissonT> (d,34); } 00315 {RandPoissonT d(new E(125916),2.5); 00316 stat |= checkSaveDistribution<E,RandPoissonT> (d,10); } 00317 {RandPoissonT d(new E(125816),2.5); 00318 stat |= checkSaveDistribution<E,RandPoissonT> (d,11); } 00319 {RandPoissonT d(new E(125716),2.5); 00320 stat |= checkSaveDistribution<E,RandPoissonT> (d,12); } 00321 00322 {std::vector<double> pdf; 00323 int nbins = 20; 00324 for (int i = 0; i < nbins; ++i) 00325 pdf.push_back( 5*i + (10.5-i) * (10.5-i) ); 00326 RandGeneral d(new E(125917), &pdf[0], 20); 00327 stat |= checkRandGeneralDistribution<E> (d,33); } 00328 00329 return stat; 00330 } 00331 00332 // --------------------------------------------- 00333 // --------------------------------------------- 00334 // --------------------------------------------- 00335 00336 00337 int main() { 00338 int stat = 0; 00339 00340 #ifdef TEST_ENGINE_NAMES 00341 output << "\n=============================================\n"; 00342 output << " Part II \n"; 00343 output << "Check all engine names were entered correctly \n"; 00344 output << "=============================================\n\n"; 00345 00346 stat |= checkEngineName<DRand48Engine >("DRand48Engine"); 00347 stat |= checkEngineName<DualRand >("DualRand"); 00348 stat |= checkEngineName<Hurd160Engine >("Hurd160Engine"); 00349 stat |= checkEngineName<Hurd288Engine >("Hurd288Engine"); 00350 stat |= checkEngineName<HepJamesRandom>("HepJamesRandom"); 00351 stat |= checkEngineName<MTwistEngine >("MTwistEngine"); 00352 stat |= checkEngineName<RandEngine >("RandEngine"); 00353 stat |= checkEngineName<RanecuEngine >("RanecuEngine"); 00354 stat |= checkEngineName<Ranlux64Engine>("Ranlux64Engine"); 00355 stat |= checkEngineName<RanluxEngine >("RanluxEngine"); 00356 stat |= checkEngineName<RanshiEngine >("RanshiEngine"); 00357 stat |= checkEngineName<TripleRand >("TripleRand"); 00358 #endif 00359 00360 #ifdef TEST_INSTANCE_METHODS 00361 output << "===========================================\n\n"; 00362 output << " Part III\n"; 00363 output << "Check instance methods for specific engines \n"; 00364 output << " specific engines and distributions\n"; 00365 output << "===========================================\n\n"; 00366 00367 {DualRand e(234); stat |= checkEngineInstanceSave(e);} 00368 {Hurd160Engine e(234); stat |= checkEngineInstanceSave(e);} 00369 {Hurd288Engine e(234); stat |= checkEngineInstanceSave(e);} 00370 {HepJamesRandom e(234); stat |= checkEngineInstanceSave(e);} 00371 {MTwistEngine e(234); stat |= checkEngineInstanceSave(e);} 00372 {RandEngine e(234); stat |= checkEngineInstanceSave(e);} 00373 {RanecuEngine e(234); stat |= checkEngineInstanceSave(e);} 00374 {Ranlux64Engine e(234); stat |= checkEngineInstanceSave(e);} 00375 {RanluxEngine e(234); stat |= checkEngineInstanceSave(e);} 00376 {RanshiEngine e(234); stat |= checkEngineInstanceSave(e);} 00377 {TripleRand e(234); stat |= checkEngineInstanceSave(e);} 00378 00379 {std::vector<double> nonRand = aSequence(500); 00380 NonRandomEngine e; 00381 e.setRandomSequence(&nonRand[0], nonRand.size()); 00382 stat |= checkEngineInstanceSave(e);} 00383 00384 stat |= checkDistributions<DualRand>(); 00385 stat |= checkDistributions<Hurd160Engine>(); 00386 stat |= checkDistributions<Hurd288Engine>(); 00387 stat |= checkDistributions<HepJamesRandom>(); 00388 stat |= checkDistributions<MTwistEngine>(); 00389 stat |= checkDistributions<Ranlux64Engine>(); 00390 stat |= checkDistributions<RanluxEngine>(); 00391 stat |= checkDistributions<RanshiEngine>(); 00392 stat |= checkDistributions<TripleRand>(); 00393 00394 RandGaussQ::shoot(); // Just to verify that the static engine is OK 00395 #endif 00396 00397 output << "\n=============================================\n\n"; 00398 00399 if (stat != 0) { 00400 std::cout << "One or more problems detected: stat = " << stat << "\n"; 00401 output << "One or more problems detected: stat = " << stat << "\n"; 00402 } else { 00403 output << "testInstanceRestore passed with no problems detected.\n"; 00404 } 00405 00406 if (stat == 0) return 0; 00407 if (stat > 0) return -(stat|1); 00408 return stat|1; 00409 00410 } 00411