Commit db8c360f authored by Lynn Garren's avatar Lynn Garren
Browse files

first pass fix for bug #73093, Ranlux64Engine and 64-bit seeds

parent bdf23f35
2010-10-21 Mark Fischler <mf@fnal.gov> and Lynn Garren <garren@fnal.gov>
* Random/src/Ranlux64Engine.cc: fix problem with random numbers when seed is
greater than 32bits
* Random/test/testBug73093.cc: test for the Ranlux64Engine problem
==============================
23.07.10 Release CLHEP-2.1.0.0
==============================
......
2010-10-21 Mark Fischler <mf@fnal.gov> and Lynn Garren <garren@fnal.gov>
* src/Ranlux64Engine.cc: fix problem with random numbers when seed is
greater than 32bits
* test/testBug73093.cc: test for the Ranlux64Engine problem
==============================
23.07.10 Release CLHEP-2.1.0.0
==============================
......
// $Id: Ranlux64Engine.cc,v 1.6 2010/06/16 17:24:53 garren Exp $
// $Id: Ranlux64Engine.cc,v 1.7 2010/10/21 21:32:02 garren Exp $
// -*- C++ -*-
//
// -----------------------------------------------------------------------
......@@ -378,6 +378,10 @@ void Ranlux64Engine::setSeed(long seed, int lux) {
long next_seed = seed;
long k_multiple;
int i;
next_seed &= 0xffffffff;
while( next_seed >= ecuyer_d ) {
next_seed -= ecuyer_d;
}
for(i = 0;i != 24;i++){
k_multiple = next_seed / ecuyer_a;
......@@ -386,12 +390,28 @@ void Ranlux64Engine::setSeed(long seed, int lux) {
if(next_seed < 0) {
next_seed += ecuyer_d;
}
init_table[i] = next_seed & 0xffffffff;
}
next_seed &= 0xffffffff;
init_table[i] = next_seed;
}
// are we on a 64bit machine?
if( sizeof(long) >= 8 ) {
long topbits1, topbits2;
topbits1 = ( seed >> 32) & 0xffff ;
topbits2 = ( seed >> 48) & 0xffff ;
init_table[0] ^= topbits1;
init_table[2] ^= topbits2;
//std::cout << " init_table[0] " << init_table[0] << " from " << topbits1 << std::endl;
//std::cout << " init_table[2] " << init_table[2] << " from " << topbits2 << std::endl;
}
for(i = 0;i < 12; i++){
randoms[i] = (init_table[2*i ] ) * 2.0 * twoToMinus_32() +
(init_table[2*i+1] >> 15) * twoToMinus_48();
//if( randoms[i] < 0. || randoms[i] > 1. ) {
//std::cout << "setSeed: init_table " << init_table[2*i ] << std::endl;
//std::cout << "setSeed: init_table " << init_table[2*i+1] << std::endl;
//std::cout << "setSeed: random " << i << " is " << randoms[i] << std::endl;
//}
}
carry = 0.0;
......@@ -457,7 +477,8 @@ void Ranlux64Engine::setSeeds(const long * seeds, int lux) {
if(next_seed < 0) {
next_seed += ecuyer_d;
}
init_table[i] = next_seed & 0xffffffff;
next_seed &= 0xffffffff;
init_table[i] = next_seed;
}
}
......
......@@ -31,7 +31,7 @@ check_PROGRAMS = \
testRandom testRandDists gaussSpeed gaussSmall \
testSaveEngineStatus testInstanceRestore testSaveSharedEngines \
testStaticStreamSave testAnonymousEngineRestore testVectorSave \
testBug58950 testEngineCopy testDistCopy testRanecuSequence
testBug58950 testBug73093 testEngineCopy testDistCopy testRanecuSequence
check_SCRIPTS = \
testRandom.sh testRandDists.sh
......@@ -43,13 +43,13 @@ TESTS = \
testRandom.sh testRandDists.sh \
testSaveEngineStatus testInstanceRestore testSaveSharedEngines \
testStaticStreamSave testAnonymousEngineRestore testVectorSave \
testEngineCopy testDistCopy
testBug73093 testEngineCopy testDistCopy
else
TESTS = \
testRandom.sh testRandDists.sh \
testSaveEngineStatus testInstanceRestore testSaveSharedEngines \
testStaticStreamSave testAnonymousEngineRestore testVectorSave \
testBug58950 testEngineCopy testDistCopy testRanecuSequence
testBug58950 testBug73093 testEngineCopy testDistCopy testRanecuSequence
endif
# Identify the test(s) for which failure is the intended outcome:
......@@ -68,6 +68,7 @@ testEngineCopy_SOURCES = testEngineCopy.cc
testDistCopy_SOURCES = testDistCopy.cc
testRanecuSequence_SOURCES = testRanecuSequence.cc
testBug58950_SOURCES = testBug58950.cc
testBug73093_SOURCES = testBug73093.cc
gaussSpeed_SOURCES = gaussSpeed.cc
gaussSmall_SOURCES = gaussSmall.cc
......
// ----------------------------------------------------------------------
//
// testBug73093 -- Test of CLHEP::Ranlux64Engine with 64 bit seeds
//
// Frank Winklmeier 2010-09-24
// L. Garren 2010-10-21 rewritten for test suite
//
// ----------------------------------------------------------------------
#include <iostream>
#include <cmath>
#include <stdlib.h>
#include "CLHEP/Random/Ranlux64Engine.h"
using namespace std;
int valid_range( )
{
std::ofstream output("testBug73093.cout");
int bad = 0;
long seed;
long mult=-235421;
// use several seeds
for( int il=0; il<100; ++il ) {
if( sizeof(long) > 4 ) {
// using atol so 32bit compilers won't complain
seed = atol("9899876543210000");
mult = mult + atol("120000000000");
} else {
seed = 987654321;
}
seed += il*mult;
CLHEP::Ranlux64Engine rng;
const long N = 20;
rng.setSeed(seed, /*lux*/ 1);
output << endl;
output << "sizeof(long) = " << sizeof(long) << endl;
output << "Generating " << N << " random numbers with seed " << seed << endl;
output << "Seed has " << logb(seed)+1 << " bits" << endl;
output << "Using seed " << seed << endl;
double sum(0);
for (long i=0; i<N; ++i) {
double r = rng.flat();
if( abs(r) > 1.0 ) ++bad;
output << r << endl;
sum += r;
}
output << "Sum: " << sum << endl;
output << "Average: " << sum / N << endl;
}
return bad;
}
int check_sequence()
{
// if the seed is less than 32bits long on a 64bit machine, the random
// number sequence should be the same as the sequence on a 32bit machine
std::ofstream output("testBug73093.seq");
int bad = 0;
long seed;
long mult=-235421;
// use several seeds
for( int il=0; il<50; ++il ) {
seed = 97654321;
seed += il*mult;
CLHEP::Ranlux64Engine rng;
const long N = 20;
rng.setSeed(seed, /*lux*/ 1);
double sum(0);
for (long i=0; i<N; ++i) {
double r = rng.flat();
if( abs(r) > 1.0 ) ++bad;
output << "[" << il << "][" << i << "] = " << r << ";" << endl;
sum += r;
}
}
return bad;
}
int main()
{
int bad = 0;
bad += valid_range( );
bad += check_sequence( );
return bad;
}
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