Commit 67e23c84 authored by Lynn Garren's avatar Lynn Garren
Browse files

add and use a private method to avoid repeating sequences

parent e1fac8e8
// $Id: RanecuEngine.h,v 1.5 2010/06/16 17:24:53 garren Exp $
// $Id: RanecuEngine.h,v 1.6 2010/07/20 18:06:02 garren Exp $
// -*- C++ -*-
//
// -----------------------------------------------------------------------
......@@ -117,6 +117,9 @@ protected:
private:
// private method used to mitigate the effects of using a lookup table
void further_randomize (int seq, int col, int index, int modulus);
// Members defining the current state of the generator.
static const int maxSeq = 215;
......
// $Id: RanecuEngine.cc,v 1.6 2010/06/16 17:24:53 garren Exp $
// $Id: RanecuEngine.cc,v 1.7 2010/07/20 18:06:02 garren Exp $
// -*- C++ -*-
//
// -----------------------------------------------------------------------
......@@ -31,6 +31,9 @@
// getState() for anonymous restores 12/27/04
// M. Fischler - put/get for vectors of ulongs 3/14/05
// M. Fischler - State-saving using only ints, for portability 4/12/05
// M. Fischler - Modify ctor and setSeed to utilize all info provided
// and avoid coincidence of same state from different
// seeds 6/22/10
//
// =======================================================================
......@@ -50,6 +53,12 @@ static const double prec = 4.6566128E-10;
std::string RanecuEngine::name() const {return "RanecuEngine";}
void RanecuEngine::further_randomize (int seq, int col, int index, int modulus)
{
table[seq][col] -= (index&0x3FFFFFFF);
while (table[seq][col] <= 0) table[seq][col] += (modulus-1);
} // mf 6/22/10
// Number of instances with automatic seed selection
int RanecuEngine::numEngines = 0;
......@@ -83,6 +92,7 @@ RanecuEngine::RanecuEngine(int index)
table[j][1] ^= mask;
}
theSeeds = &table[seq][0];
further_randomize (seq, 0, index, shift1); // mf 6/22/10
}
RanecuEngine::RanecuEngine(std::istream& is)
......@@ -93,12 +103,14 @@ RanecuEngine::RanecuEngine(std::istream& is)
RanecuEngine::~RanecuEngine() {}
void RanecuEngine::setSeed(long index, int)
void RanecuEngine::setSeed(long index, int dum)
{
seq = abs(int(index%maxSeq));
theSeed = seq;
HepRandom::getTheTableSeeds(table[seq],seq);
theSeeds = &table[seq][0];
further_randomize (seq, 0, index, shift1); // mf 6/22/10
further_randomize (seq, 1, dum, shift2); // mf 6/22/10
}
void RanecuEngine::setSeeds(const long* seeds, int pos)
......
......@@ -31,7 +31,7 @@ check_PROGRAMS = \
testRandom testRandDists gaussSpeed gaussSmall \
testSaveEngineStatus testInstanceRestore testSaveSharedEngines \
testStaticStreamSave testAnonymousEngineRestore testVectorSave \
testBug58950 testEngineCopy testDistCopy
testBug58950 testEngineCopy testDistCopy testRanecuSequence
check_SCRIPTS = \
testRandom.sh testRandDists.sh
......@@ -49,7 +49,7 @@ TESTS = \
testRandom.sh testRandDists.sh \
testSaveEngineStatus testInstanceRestore testSaveSharedEngines \
testStaticStreamSave testAnonymousEngineRestore testVectorSave \
testBug58950 testEngineCopy testDistCopy
testBug58950 testEngineCopy testDistCopy testRanecuSequence
endif
# Identify the test(s) for which failure is the intended outcome:
......@@ -66,6 +66,7 @@ testAnonymousEngineRestore_SOURCES = testAnonymousEngineRestore.cc
testVectorSave_SOURCES = testVectorSave.cc
testEngineCopy_SOURCES = testEngineCopy.cc
testDistCopy_SOURCES = testDistCopy.cc
testRanecuSequence_SOURCES = testRanecuSequence.cc
testBug58950_SOURCES = testBug58950.cc
gaussSpeed_SOURCES = gaussSpeed.cc
gaussSmall_SOURCES = gaussSmall.cc
......@@ -82,7 +83,7 @@ CLEANFILES = testRandom.sh testRandDists.sh \
testSaveEngineStatus.cout testInstanceRestore.cout \
testSaveSharedEngines.cout testStaticStreamSave.cout \
testAnonymousEngineRestore.cout testVectorSave.cout \
anonymous.save
anonymous.save testRanecuSequence.output
# supply our own suffix rule
.cc.obj:
......
// ----------------------------------------------------------------------
//
// testRanecuSequence
// look at report that certain seeds produce identical sequences
//
// ----------------------------------------------------------------------
#include <iostream>
#include <string>
#include <vector>
#include "CLHEP/Random/RanecuEngine.h"
#include "CLHEP/Random/Random.h"
struct Sample {
int seed;
std::vector<double> case1;
std::vector<double> case2;
std::vector<double> case3;
std::vector<double> case4;
std::vector<double> case5;
};
void useSeed( int, std::vector<Sample>& );
bool compareSamples( Sample&, Sample&, std::ofstream& );
int main() {
std::ofstream output("testRanecuSequence.output");
output << " Setting CLHEP Random Engine..." << std::endl;
CLHEP::HepRandom::setTheEngine(new CLHEP::RanecuEngine);
std::vector<Sample> ranList;
// the first 4 seeds once resulted in identical sequences
useSeed(1275177715, ranList);
useSeed(1275171265, ranList);
useSeed(1275168040, ranList);
useSeed(1275168040-2048*215, ranList);
useSeed(4132429, ranList);
useSeed(-1357924, ranList);
int sd = 53208134;
for ( int i = 0; i < 1000; ++i ) {
++sd;
useSeed(sd, ranList);
}
int numbad = 0;
output << std::endl;
for ( unsigned int it=0; it < ranList.size(); ++it ) {
for ( unsigned int jt=it+1; jt < ranList.size(); ++jt ) {
if( ! compareSamples(ranList[it], ranList[jt], output ) ) {
++numbad;
output << "ERROR: Seed " << ranList[it].seed
<< " and Seed " << ranList[jt].seed
<< " are " << (ranList[jt].seed - ranList[it].seed )
<< " apart and result in identical sequences" << std::endl;
std::cerr << "Seed " << ranList[it].seed
<< " and Seed " << ranList[jt].seed
<< " are " << (ranList[jt].seed - ranList[it].seed )
<< " apart and result in identical sequences" << std::endl;
}
}
}
output << " numbad is " << numbad << std::endl;
return numbad;
}
void useSeed( int seed, std::vector<Sample>& ranList )
{
Sample v;
v.seed = seed;
// case 1 -- default constructor
CLHEP::RanecuEngine c1;
for( int i = 0; i < 15; ++i ) {
v.case1.push_back( c1.flat() );
}
// case 2
CLHEP::RanecuEngine c2(seed);
for( int i = 0; i < 15; ++i ) {
v.case2.push_back( c2.flat() );
}
// case 3 - use HepRandom
CLHEP::HepRandom::setTheSeed(seed);
for( int i = 0; i < 15; ++i ) {
v.case3.push_back( CLHEP::HepRandom::getTheEngine()->flat() );
}
// case 4
CLHEP::RanecuEngine c4(1);
c4.setSeed(seed);
for( int i = 0; i < 15; ++i ) {
v.case4.push_back( c4.flat() );
}
// case 5
CLHEP::RanecuEngine c5(1);
long seedarray[2];
seedarray[0] = seed;
seedarray[1] = 1;
c5.setSeeds(seedarray,2);
for( int i = 0; i < 15; ++i ) {
v.case5.push_back( c5.flat() );
}
//
ranList.push_back(v);
}
bool compareSamples( Sample& s1, Sample& s2, std::ofstream& output )
{
if ( s1.case1 == s2.case1 ) {
output << " case1: default constructor \n" ;
output << " comparing Seed " << s1.seed << " and Seed " << s2.seed << std::endl;
output << " case1 sequence:" ;
for( unsigned int i=0; i < s1.case1.size(); ++i ) output << " " << s1.case1[i];
output << std::endl;
return false;
}
if ( s1.case2 == s2.case2 ) {
output << " case2: construct with single seed \n" ;
output << " comparing Seed " << s1.seed << " and Seed " << s2.seed << std::endl;
output << " case2 sequence:" ;
for( unsigned int i=0; i < s1.case2.size(); ++i ) output << " " << s1.case2[i];
output << std::endl;
return false;
}
if ( s1.case3 == s2.case3 ) {
output << " case3: construct with single seed \n" ;
output << " comparing Seed " << s1.seed << " and Seed " << s2.seed << std::endl;
output << " case3 sequence:" ;
for( unsigned int i=0; i < s1.case3.size(); ++i ) output << " " << s1.case3[i];
output << std::endl;
return false;
}
if ( s1.case4 == s2.case4 ) {
output << " case4: use setSeed \n" ;
output << " comparing Seed " << s1.seed << " and Seed " << s2.seed << std::endl;
output << " case4 sequence:" ;
for( unsigned int i=0; i < s1.case4.size(); ++i ) output << " " << s1.case4[i];
output << std::endl;
return false;
}
if ( s1.case5 == s2.case5 ) {
output << " case5: use setSeeds \n" ;
output << " comparing Seed " << s1.seed << " and Seed " << s2.seed << std::endl;
output << " case5 sequence:" ;
for( unsigned int i=0; i < s1.case5.size(); ++i ) output << " " << s1.case5[i];
output << std::endl;
return false;
}
return true;
}
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment