CLHEP 2.0.4.7 Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // ---------------------------------------------------------------------- 00002 // 00003 // testBug58950 -- test problem with RanecuEngine on 64bit machines 00004 // 00005 // R. Weller 11/11/09 initial test from Vanderbilt 00006 // L. Garren 12/1/09 rewritten for test suite 00007 // 00008 // ---------------------------------------------------------------------- 00009 #include <iostream> 00010 #include <stdexcept> 00011 #include <cmath> 00012 #include <stdlib.h> 00013 #include "CLHEP/Random/RanecuEngine.h" 00014 #include "CLHEP/Random/Random.h" 00015 00016 bool printCheck( int & i, double & r, std::ofstream & os ) 00017 { 00018 os << i << " " << r << std::endl; 00019 if (r < 0 || r > 1.0 ) { 00020 std::cout << "Error: bad random number " << r << std::endl; 00021 return false; 00022 } 00023 return true; 00024 } 00025 00026 int main() { 00027 00028 std::ofstream output("testBug58950.cout"); 00029 00030 output << std::endl << "short " << sizeof(short) << std::endl; 00031 output << "int " << sizeof(int) << std::endl; 00032 output << "unsigned int " << sizeof(unsigned int) << std::endl; 00033 output << "long " << sizeof(long) << std::endl; 00034 output << "float " << sizeof(float) << std::endl; 00035 output << "double " << sizeof(double) << std::endl; 00036 output << "long double " << sizeof(long double) << std::endl << std::endl; 00037 00038 CLHEP::RanecuEngine *eng = new CLHEP::RanecuEngine; 00039 CLHEP::HepRandom::setTheEngine(eng); 00040 CLHEP::HepRandom *g; 00041 g=CLHEP::HepRandom::getTheGenerator(); 00042 00043 long rvals[2]; 00044 try { 00045 std::ifstream in("/dev/urandom", std::ios::in | std::ios::binary); 00046 if(in.is_open()) { 00047 in.read((char *)(&rvals), 2*sizeof(long)); 00048 in.close(); 00049 if(in.fail()) { 00050 throw std::runtime_error("File read error"); 00051 } 00052 } else throw std::runtime_error("File open error"); 00053 } catch(std::runtime_error e) { 00054 std::ostringstream dStr; 00055 dStr << "Error: " << e.what() 00056 << " processing seed from file \"" << "/dev/urandom" << "\"."; 00057 throw std::runtime_error(dStr.str().c_str()); 00058 } 00059 00060 int nNumbers = 20; 00061 int badcount = 0; 00062 00063 long seeds[3]; 00064 const long *pseeds; 00065 //*********************************************************************** 00066 // Seeds are expected to be positive. Therefore, if either seed is 00067 // negative then prior to 2.0.4.5 the generator set initial conditions 00068 // and generated the same sequence of numbers no matter what the seeds were. 00069 seeds[0]=rvals[0]; 00070 seeds[1]=rvals[1]; 00071 seeds[2]=0; 00072 if( rvals[0] > 0 ) seeds[0] = -rvals[0]; 00073 00074 double negseq[20] = { 0.154707, 0.587114, 0.702059, 0.566, 0.988325, 00075 0.525921, 0.191554, 0.269338, 0.234277, 0.358997, 00076 0.549936, 0.296877, 0.162243, 0.227732, 0.528862, 00077 0.631571, 0.176462, 0.247858, 0.170025, 0.284483 }; 00078 double eps = 1.0E-6; 00079 00080 output << std::endl << "********************" << std::endl; 00081 output << "This is the case that may or may not fail." << std::endl; 00082 output << "However, if it has values in (0,1), they are a " << std::endl 00083 << "deterministic sequence beginning with 0.154707." << std::endl; 00084 output << "seeds[0] = " << seeds[0] << "\n" 00085 << "seeds[1] = " << seeds[1] << std::endl << std::endl; 00086 00087 g->setTheSeeds(seeds); 00088 int rseq = 0; 00089 for (int i=0; i < nNumbers; ++i) { 00090 double r = g->flat(); 00091 if( ! printCheck(i,r,output) ) ++badcount; 00092 // before the change, the random number sequence was reliably the same 00093 if( std::fabs(r-negseq[i]) < eps ) { 00094 std::cout << " reproducing sequence " << i << " " 00095 << r << " " << negseq[i] << std::endl; 00096 ++rseq; 00097 } 00098 } 00099 if( rseq == 20 ) ++badcount; 00100 pseeds=g->getTheSeeds(); 00101 output << "Final seeds[0] = " << pseeds[0] << "\n" 00102 << "Final seeds[1] = " << pseeds[1] << std::endl << std::endl; 00103 00104 //*********************************************************************** 00105 // Prior to the 2.0.4.5 bug fix, 64bit seeds resulted in incorrect randoms 00106 seeds[0]=labs(rvals[0]); 00107 seeds[1]=labs(rvals[1]); 00108 seeds[2]=0; 00109 00110 output << std::endl << "********************" << std::endl; 00111 output << "This is the case that always fails." << std::endl; 00112 output << "seeds[0] = " << seeds[0] << "\n" 00113 << "seeds[1] = " << seeds[1] << std::endl << std::endl; 00114 00115 g->setTheSeeds(seeds); 00116 for (int i=0; i < nNumbers; ++i) { 00117 double r = g->flat(); 00118 if( ! printCheck(i,r,output) ) ++badcount; 00119 } 00120 pseeds=g->getTheSeeds(); 00121 output << "Final seeds[0] = " << pseeds[0] << "\n" 00122 << "Final seeds[1] = " << pseeds[1] << std::endl << std::endl; 00123 00124 //*********************************************************************** 00125 // recover and reuse seeds 00126 seeds[0]=labs(rvals[0]); 00127 seeds[1]=labs(rvals[1]); 00128 seeds[2]=0; 00129 00130 output << std::endl << "********************" << std::endl; 00131 output << "Check rolling back a random number seed." << std::endl; 00132 output << "seeds[0] = " << seeds[0] << "\n" 00133 << "seeds[1] = " << seeds[1] << std::endl << std::endl; 00134 std::vector<double> v; 00135 g->setTheSeeds(seeds); 00136 00137 for (int i=0; i < nNumbers; ++i) { 00138 double r = g->flat(); 00139 if( ! printCheck(i,r,output) ) ++badcount; 00140 } 00141 pseeds=g->getTheSeeds(); 00142 seeds[0] = pseeds[0]; 00143 seeds[1] = pseeds[1]; 00144 output << " pseeds[0] = " << pseeds[0] << "\n" 00145 << "pseeds[1] = " << pseeds[1] << std::endl; 00146 for (int i=0; i < nNumbers; ++i) { 00147 double r = g->flat(); 00148 v.push_back(r); 00149 } 00150 g->setTheSeeds(seeds); 00151 for (int i=0; i < nNumbers; ++i) { 00152 double r = g->flat(); 00153 if(v[i] != r ) { 00154 ++badcount; 00155 std::cerr << " rollback fails: i, v[i], r " 00156 << i << " " << v[i] << " " << r << std::endl; 00157 } 00158 } 00159 output << std::endl; 00160 00161 //*********************************************************************** 00162 // 4-byte positive integers generate valid sequences, which remain within bounds. 00163 seeds[0]= labs(static_cast<int>(rvals[0])); 00164 seeds[1]= labs(static_cast<int>(rvals[1])); 00165 seeds[2]=0; 00166 00167 output << std::endl << "********************" << std::endl; 00168 output << "This is the case that works." << std::endl; 00169 output << std::endl << "seeds[0] = " << seeds[0] << "\n" 00170 << "seeds[1] = " << seeds[1] << "\n" 00171 << "seeds[2] = " << seeds[2] << std::endl << std::endl; 00172 00173 g->setTheSeeds(seeds); 00174 for (int i=0; i < nNumbers; ++i) { 00175 double r = g->flat(); 00176 if( ! printCheck(i,r,output) ) ++badcount; 00177 } 00178 pseeds=g->getTheSeeds(); 00179 output << "Final seeds[0] = " << pseeds[0] << "\n" 00180 << "Final seeds[1] = " << pseeds[1] << std::endl << std::endl; 00181 00182 //*********************************************************************** 00183 // Before the fix, a bad 64bit sequence would eventually rectify itself. 00184 // This starts with seeds that would have failed before the 64bit corrections 00185 // were applied and loops until both seeds are positive 32-bit integers. 00186 // This looping should no longer occur. 00187 seeds[0]=labs(rvals[0]); 00188 seeds[1]=labs(rvals[1]); 00189 seeds[2]=0; 00190 00191 output << std::endl << "********************" << std::endl; 00192 output << "This case loops until valid short seeds occur." << std::endl; 00193 output << "seeds[0] = " << seeds[0] << "\n" 00194 << "seeds[1] = " << seeds[1] << std::endl << std::endl; 00195 00196 g->setTheSeeds(seeds); 00197 // Loop as long as the values are bad. 00198 double r; 00199 unsigned int low = ~0; 00200 unsigned long mask = (~0) << 31; 00201 unsigned long skipcount = 0; 00202 output << "low = " << low << " mask = " << mask << std::endl; 00203 do {r = g->flat(); pseeds = g->getTheSeeds(); ++skipcount;} 00204 while((pseeds[0]&mask) || (pseeds[1]&mask)); 00205 if ( skipcount > 1 ) ++badcount; 00206 00207 output << std::endl << "Loop terminates on two short seeds." << std::endl; 00208 output << "Skipcount = " << skipcount << std::endl; 00209 output << "pseeds[0]&mask = " << (pseeds[0]&mask) << std::endl; 00210 output << "pseeds[1]&mask = " << (pseeds[1]&mask) << std::endl; 00211 output << "Final seeds[0] = " << pseeds[0] << "\n" 00212 << "Final seeds[1] = " << pseeds[1] << std::endl << std::endl; 00213 00214 output << "This should be a valid sequence." << std::endl; 00215 for (int i=0; i < nNumbers; ++i) { 00216 double r = g->flat(); 00217 if( ! printCheck(i,r,output) ) ++badcount; 00218 } 00219 pseeds=g->getTheSeeds(); 00220 output << "seeds[0] = " << pseeds[0] << "\n" 00221 << "seeds[1] = " << pseeds[1] << std::endl << std::endl; 00222 00223 if( badcount > 0 ) std::cout << "Error count is " << badcount << std::endl; 00224 return badcount; 00225 }