added a thread-safe and cross-plateform random number generator

git-svn-id: http://svn.code.sf.net/p/retroshare/code/trunk@3373 b45a01b8-16f6-495d-af2f-9b41ad6348cc
This commit is contained in:
csoler 2010-08-15 20:20:51 +00:00
parent f0da158dbe
commit 30e4509043
5 changed files with 302 additions and 3 deletions

View File

@ -369,6 +369,7 @@ HEADERS += util/folderiterator.h \
util/rsthreads.h \
util/rsversion.h \
util/rswin.h \
util/rsrandom.h \
SOURCES += dbase/cachestrapper.cc \
dbase/fimonitor.cc \
@ -485,4 +486,5 @@ SOURCES += util/folderiterator.cc \
util/rsprint.cc \
util/rsthreads.cc \
util/rsversion.cc \
util/rswin.cc
util/rswin.cc \
util/rsrandom.cc

View File

@ -7,14 +7,17 @@ RS_TOP_DIR = ../..
include $(RS_TOP_DIR)/tests/scripts/config.mk
###############################################################
TESTOBJ = netsetup_test.o
TESTS = netsetup_test
TESTOBJ = netsetup_test.o random_test.o
TESTS = netsetup_test random_test
all: tests
netsetup_test: netsetup_test.o
$(CC) $(CFLAGS) -o netsetup_test netsetup_test.o $(LIBS)
random_test: random_test.o
$(CC) $(CFLAGS) -o random_test random_test.o $(LIBS)
###############################################################
include $(RS_TOP_DIR)/tests/scripts/rules.mk
###############################################################

View File

@ -0,0 +1,161 @@
#ifdef LINUX
#include <fenv.h>
#endif
#include <iostream>
#include <vector>
#include <stdlib.h>
#include <math.h>
#include "util/utest.h"
#include "util/rsrandom.h"
INITTEST();
typedef double (* Ifctn)( double t);
/* Numerical integration method */
double Simpson3_8( Ifctn f, double a, double b, int N)
{
int j;
double l1;
double h = (b-a)/N;
double h1 = h/3.0;
double sum = f(a) + f(b);
for (j=3*N-1; j>0; j--) {
l1 = (j%3)? 3.0 : 2.0;
sum += l1*f(a+h1*j) ;
}
return h*sum/8.0;
}
#define A 12
double Gamma_Spouge( double z )
{
int k;
static double cspace[A];
static double *coefs = NULL;
double accum;
double a = A;
if (!coefs) {
double k1_factrl = 1.0;
coefs = cspace;
coefs[0] = sqrt(2.0*M_PI);
for(k=1; k<A; k++) {
coefs[k] = exp(a-k) * pow(a-k,k-0.5) / k1_factrl;
k1_factrl *= -k;
}
}
accum = coefs[0];
for (k=1; k<A; k++) {
accum += coefs[k]/(z+k);
}
accum *= exp(-(z+a)) * pow(z+a, z+0.5);
return accum/z;
}
double aa1;
double f0( double t)
{
return pow(t, aa1)*exp(-t);
}
double GammaIncomplete_Q( double a, double x)
{
double y, h = 1.5e-2; /* approximate integration step size */
/* this cuts off the tail of the integration to speed things up */
y = aa1 = a-1;
while((f0(y) * (x-y) > 2.0e-8) && (y < x)) y += .4;
if (y>x) y=x;
return 1.0 - Simpson3_8( &f0, 0, y, std::max(5,(int)(y/h)))/Gamma_Spouge(a);
}
double chi2Probability( int dof, double distance)
{
return GammaIncomplete_Q( 0.5*dof, 0.5*distance);
}
class myThread: public RsThread
{
public:
myThread()
{
_finished = false ;
}
virtual void run()
{
// test that random numbers are regularly disposed
//
int N = 500 ;
int B = 8 ;
std::vector<int> buckets(B,0) ;
for(int i=0;i<N;++i)
{
float f = RSRandom::random_f32() ;
int b = (int)floor(f*B*0.99999999) ;
++buckets[b] ;
}
// Chi2 test
//
float chi2 = 0.0f ;
float expected = 1.0/(float)B ;
for(int k=0;k<B;++k)
chi2 += pow( buckets[k]/float(N) - expected,2)/expected ;
double chi2prob = chi2Probability(B-1,chi2) ;
double significance = 0.05 ;
std::cerr << "Compariing chi2 uniform distance " << chi2 << " to chi2 probability " << chi2prob ;
if(chi2prob > significance)
std::cerr << ": passed" << std::endl ;
else
std::cerr << ": failed" << std::endl ;
_finished = true ;
}
bool finished() const { return _finished ; }
private:
bool _finished ;
};
int main(int argc, char **argv)
{
#ifdef LINUX
feenableexcept(FE_INVALID) ;
feenableexcept(FE_DIVBYZERO) ;
#endif
int nt = 10 ; // number of threads.
std::vector<myThread *> threads(nt,(myThread*)NULL) ;
for(int i=0;i<nt;++i)
{
threads[i] = new myThread ;
threads[i]->start() ;
}
while(true)
{
bool finished = true ;
for(int i=0;i<nt;++i)
if(!threads[i]->finished())
finished = false ;
if(finished)
break ;
}
for(int i=0;i<nt;++i)
delete threads[i] ;
FINALREPORT("random_test");
exit(TESTRESULT());
}

View File

@ -0,0 +1,73 @@
#include <stdlib.h>
#include "rsrandom.h"
uint32_t RSRandom::index = 0 ;
static bool auto_seed = RSRandom::seed(time(NULL)) ;
std::vector<uint32_t> RSRandom::MT(RSRandom::N,0u) ;
RsMutex RSRandom::rndMtx ;
bool RSRandom::seed(uint32_t s)
{
RsStackMutex mtx(rndMtx) ;
MT.resize(N,0) ; // because MT might not be already resized
uint32_t j ;
MT[0]= s & 0xffffffffUL;
for (j=1; j<N; j++)
MT[j] = (1812433253UL * (MT[j-1] ^ (MT[j-1] >> 30)) + j) & 0xffffffffUL ;
return true ;
}
void RSRandom::locked_next_state()
{
for(uint32_t i=0;i<N;++i)
{
uint32_t y = ((MT[i]) & UMASK) | ((MT[(i+1)%(int)N]) & LMASK) ;
MT[i] = MT[(i + M) % (int)N] ^ (y >> 1) ;
if((y & 1) == 1)
MT[i] = MT[i] ^ 0x9908b0df ;
}
index = 0 ;
}
uint32_t RSRandom::random_u32()
{
uint32_t y;
{
RsStackMutex mtx(rndMtx) ;
y = MT[index++] ;
if(index == N)
locked_next_state();
}
// Tempering
y ^= (y >> 11);
y ^= (y << 7 ) & 0x9d2c5680UL;
y ^= (y << 15) & 0xefc60000UL;
y ^= (y >> 18);
return y;
}
uint64_t RSRandom::random_u64()
{
return ((uint64_t)random_u32() << 32ul) + random_u32() ;
}
float RSRandom::random_f32()
{
return random_u32() / (float)(~(uint32_t)0) ;
}
double RSRandom::random_f64()
{
return random_u64() / (double)(~(uint64_t)0) ;
}

View File

@ -0,0 +1,60 @@
/****************************************************************
* RetroShare is distributed under the following license:
*
* Copyright (C) 2010 Cyril Soler <csoler@users.sourceforge.net>
*
* This program is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License
* as published by the Free Software Foundation; either version 2
* of the License, or (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin Street, Fifth Floor,
* Boston, MA 02110-1301, USA.
****************************************************************/
#pragma once
// RSRandom contains a random number generator that is
// - thread safe
// - system independant
// - fast
// - cryptographically safe
//
// The implementation is adapted from the Mersenne Twister page of Wikipedia.
//
// http://en.wikipedia.org/wiki/Mersenne_twister
#include <vector>
#include <util/rsthreads.h>
class RSRandom
{
public:
static uint32_t random_u32() ;
static uint64_t random_u64() ;
static float random_f32() ;
static double random_f64() ;
static bool seed(uint32_t s) ;
private:
static RsMutex rndMtx ;
static const uint32_t N = 624;
static const uint32_t M = 397;
static const uint32_t MATRIX_A = 0x9908b0dfUL;
static const uint32_t UMASK = 0x80000000UL;
static const uint32_t LMASK = 0x7fffffffUL;
static void locked_next_state() ;
static uint32_t index ;
static std::vector<uint32_t> MT ;
};