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