master
madrocketsci 4 weeks ago
parent b722664be4
commit 4e0a826903

Binary file not shown.

@ -14,56 +14,58 @@ namespace rand
//Simple Deterministic Psuedorandom Number Generator (32 bit version) // //Simple Deterministic Psuedorandom Number Generator (32 bit version) //
//////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////
static const int32_t dpr32_mod = ( ((int32_t)1) << ((int32_t)30) ) - (int32_t)1; typedef int32_t narray_rands_t;
static const int32_t dpr32_mult1 = ( (int32_t) 1201633 );
static int32_t dpr32_rseed = 0; // global random number seed for numeric array class static const narray_rands_t dpr32_mod = ( ((narray_rands_t)1) << ((narray_rands_t)30) ) - (narray_rands_t)1;
static const narray_rands_t dpr32_mult1 = ( (narray_rands_t) 1201633 );
extern narray_rands_t dpr32_rseed; // global random number seed for numeric array class
int32_t dpr32_nextseed(int32_t seed); narray_rands_t dpr32_nextseed(narray_rands_t seed);
double dpr32_randd(int32_t *seed); double dpr32_randd(narray_rands_t *seed);
float dpr32_randf(int32_t *seed); float dpr32_randf(narray_rands_t *seed);
double dpr32_gaussian(int32_t *seed); double dpr32_gaussian(narray_rands_t *seed);
float dpr32_gaussianf(int32_t *seed); float dpr32_gaussianf(narray_rands_t *seed);
int dpr32_randint(int32_t *seed, int low, int high); int dpr32_randint(narray_rands_t *seed, int low, int high);
ams::vect2 dpr32_randuvect2(int32_t *seed); ams::vect2 dpr32_randuvect2(narray_rands_t *seed);
ams::vect3 dpr32_randuvect3(int32_t *seed); ams::vect3 dpr32_randuvect3(narray_rands_t *seed);
ams::vect4 dpr32_randuvect4(int32_t *seed); ams::vect4 dpr32_randuvect4(narray_rands_t *seed);
ams::vect2f dpr32_randuvect2f(int32_t *seed); ams::vect2f dpr32_randuvect2f(narray_rands_t *seed);
ams::vect3f dpr32_randuvect3f(int32_t *seed); ams::vect3f dpr32_randuvect3f(narray_rands_t *seed);
ams::vect4f dpr32_randuvect4f(int32_t *seed); ams::vect4f dpr32_randuvect4f(narray_rands_t *seed);
//////////////////////////////////////////// ////////////////////////////////////////////
//Numeric array threaded random generators// //Numeric array threaded random generators//
//////////////////////////////////////////// ////////////////////////////////////////////
int32_t get_rseed(); narray_rands_t get_rseed();
void set_rseed(int32_t seed); void set_rseed(narray_rands_t seed);
void set_rseed_withtimer(); void set_rseed_withtimer();
narray<double> narray_rand(narray_size_t N, narray<double> narray_rand(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<float> narray_randf(narray_size_t N, narray<float> narray_randf(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<double> narray_randgauss(narray_size_t N, narray<double> narray_randgauss(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<float> narray_randgaussf(narray_size_t N, narray<float> narray_randgaussf(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<int> narray_randint(narray_size_t N, int low, int highexcl, narray<int> narray_randint(narray_size_t N, int low, int highexcl,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<vect2> narray_randuvect2(narray_size_t N, narray<vect2> narray_randuvect2(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<vect3> narray_randuvect3(narray_size_t N, narray<vect3> narray_randuvect3(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<vect4> narray_randuvect4(narray_size_t N, narray<vect4> narray_randuvect4(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<vect2f> narray_randuvect2f(narray_size_t N, narray<vect2f> narray_randuvect2f(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<vect3f> narray_randuvect3f(narray_size_t N, narray<vect3f> narray_randuvect3f(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
narray<vect4f> narray_randuvect4f(narray_size_t N, narray<vect4f> narray_randuvect4f(narray_size_t N,
int32_t *rseed = &(ams::narray::rand::dpr32_rseed)); narray_rands_t *rseed = &(ams::narray::rand::dpr32_rseed));
void test_random1(); void test_random1();

@ -57,9 +57,19 @@ namespace narray
void narray_testmath2() void narray_testmath2()
{ {
rand::narray_rands_t rs = 1943981;
rand::set_rseed_withtimer(); rand::set_rseed_withtimer();
printf("rseed: %d\n",(int)rand::get_rseed());
printf("rseed: %d\n",(int)rand::dpr32_rseed); //< -- this is different than the line above
narray<double> unif = ams::narray::rand::narray_rand(10000); //rand::set_rseed(25);
printf("rseed: %d\n",(int)rand::get_rseed());
narray<double> unif = ams::narray::rand::narray_rand(5,&rs);
printf("rseed: %d\n",(int)rand::get_rseed());
double mn1 = narray_mean(unif); double mn1 = narray_mean(unif);
double std1 = narray_stdev(unif); double std1 = narray_stdev(unif);
double mn2,std2; double mn2,std2;
@ -69,11 +79,14 @@ namespace narray
//_intl_printarr(unif); printf("\n"); //_intl_printarr(unif); printf("\n");
unif = unif*2.0-1.0; unif = unif*2.0-1.0;
//unif.setall(3.5);
//_intl_printarr(unif); printf("\n"); //_intl_printarr(unif); printf("\n");
mn2 = narray_mean(unif); mn2 = narray_mean(unif);
std2 = narray_stdev(unif); std2 = narray_stdev(unif);
printf("mn1=%1.3f std1=%1.3f\n",mn1,std1); printf("mn1=%1.3f std1=%1.3f\n",mn1,std1);
printf("mn2=%1.3f std2=%1.3f\n",mn2,std2); printf("mn2=%1.3f std2=%1.3f\n",mn2,std2);

@ -7,22 +7,27 @@ namespace narray
namespace rand namespace rand
{ {
int32_t dpr32_nextseed(int32_t seed) //defining global random number seed for numeric array random routines
narray_rands_t dpr32_rseed = 0;
narray_rands_t dpr32_nextseed(narray_rands_t seed)
{ {
int32_t sret = seed; narray_rands_t sret = seed;
sret = ams::mod(sret*dpr32_mult1+1,dpr32_mod); sret = ams::mod(sret*dpr32_mult1+1,dpr32_mod);
return sret; return sret;
} }
double dpr32_randd(int32_t *seed) double dpr32_randd(narray_rands_t *seed)
{ {
double ret; double ret;
printf("debug: seed: %d\n",(int)*seed);
*seed = ams::narray::rand::dpr32_nextseed(*seed); *seed = ams::narray::rand::dpr32_nextseed(*seed);
printf("debug: seed: %d\n",(int)*seed);
ret = (double)*seed/(double)(dpr32_mod-1); ret = (double)*seed/(double)(dpr32_mod-1);
return ret; return ret;
} }
float dpr32_randf(int32_t *seed) float dpr32_randf(narray_rands_t *seed)
{ {
float ret; float ret;
*seed = ams::narray::rand::dpr32_nextseed(*seed); *seed = ams::narray::rand::dpr32_nextseed(*seed);
@ -30,7 +35,7 @@ namespace rand
return ret; return ret;
} }
double dpr32_gaussian(int32_t *seed) double dpr32_gaussian(narray_rands_t *seed)
{ {
double ret = 0.0; double ret = 0.0;
double u1,u2; double u1,u2;
@ -46,7 +51,7 @@ namespace rand
} }
float dpr32_gaussianf(int32_t *seed) float dpr32_gaussianf(narray_rands_t *seed)
{ {
float ret = 0.0f; float ret = 0.0f;
float u1,u2; float u1,u2;
@ -61,7 +66,7 @@ namespace rand
return ret; return ret;
} }
vect2 dpr32_randuvect2(int32_t *seed) vect2 dpr32_randuvect2(narray_rands_t *seed)
{ {
vect2 ret = vect2(0,0); vect2 ret = vect2(0,0);
double th = ams::narray::rand::dpr32_randd(seed)*2.0*pi; double th = ams::narray::rand::dpr32_randd(seed)*2.0*pi;
@ -69,7 +74,7 @@ namespace rand
return ret; return ret;
} }
vect3 dpr32_randuvect3(int32_t *seed) vect3 dpr32_randuvect3(narray_rands_t *seed)
{ {
vect3 ret; vect3 ret;
double az,el; double az,el;
@ -81,7 +86,7 @@ namespace rand
return ret; return ret;
} }
vect4 dpr32_randuvect4(int32_t *seed) vect4 dpr32_randuvect4(narray_rands_t *seed)
{ {
vect4 ret; vect4 ret;
double x,y,z,w; double x,y,z,w;
@ -96,7 +101,7 @@ namespace rand
return ret; return ret;
} }
vect2f dpr32_randuvect2f(int32_t *seed) vect2f dpr32_randuvect2f(narray_rands_t *seed)
{ {
vect2f ret = vect2f(0,0); vect2f ret = vect2f(0,0);
float th = ams::narray::rand::dpr32_randf(seed)*2.0*pi; float th = ams::narray::rand::dpr32_randf(seed)*2.0*pi;
@ -104,7 +109,7 @@ namespace rand
return ret; return ret;
} }
vect3f dpr32_randuvect3f(int32_t *seed) vect3f dpr32_randuvect3f(narray_rands_t *seed)
{ {
vect3f ret; vect3f ret;
float az,el; float az,el;
@ -116,7 +121,7 @@ namespace rand
return ret; return ret;
} }
vect4f dpr32_randuvect4f(int32_t *seed) vect4f dpr32_randuvect4f(narray_rands_t *seed)
{ {
vect4f ret; vect4f ret;
float x,y,z,w; float x,y,z,w;
@ -131,27 +136,27 @@ namespace rand
return ret; return ret;
} }
int32_t get_rseed() narray_rands_t get_rseed()
{ {
return ams::narray::rand::dpr32_rseed; return ams::narray::rand::dpr32_rseed;
} }
void set_rseed(int32_t seed) void set_rseed(narray_rands_t seed)
{ {
ams::narray::rand::dpr32_rseed = seed; ams::narray::rand::dpr32_rseed = seed;
} }
void set_rseed_withtimer() void set_rseed_withtimer()
{ {
int32_t t1 = (int32_t)time(NULL); narray_rands_t t1 = (narray_rands_t)time(NULL);
int32_t t2 = (int32_t)(::fmod((double)clock()/((double)CLOCKS_PER_SEC)*1000.0f,36000.0f)); narray_rands_t t2 = (narray_rands_t)(::fmod((double)clock()/((double)CLOCKS_PER_SEC)*1000.0f,36000.0f));
ams::narray::rand::dpr32_rseed = (int32_t)t1 + (int32_t)t2; ams::narray::rand::dpr32_rseed = (narray_rands_t)t1 + (narray_rands_t)t2;
} }
template<typename T> void narray_rand_threadf1( template<typename T> void narray_rand_threadf1(
narray<T> *out, narray<T> *out,
narray<int32_t> *rseeds, narray<narray_rands_t> *rseeds,
T (*randfunc)(int32_t *), T (*randfunc)(narray_rands_t *),
int threadnum, int threadnum,
int nthreads int nthreads
) )
@ -170,15 +175,15 @@ namespace rand
template<typename T> void narray_random_threadexec1( template<typename T> void narray_random_threadexec1(
narray<T> *out, narray<T> *out,
narray_size_t N, narray_size_t N,
T (*randfunc)(int32_t *), T (*randfunc)(narray_rands_t *),
int32_t *rseed narray_rands_t *rseed
) )
{ {
narray_size_t I; narray_size_t I;
int J; int J;
int nthreads; int nthreads;
std::vector<std::thread*> threads; std::vector<std::thread*> threads;
narray<int32_t> rseeds; narray<narray_rands_t> rseeds;
int res; int res;
res = out->resize(N); res = out->resize(N);
@ -204,6 +209,7 @@ namespace rand
nthreads = (nthreads>narray_max_threads) ? narray_max_threads : nthreads; nthreads = (nthreads>narray_max_threads) ? narray_max_threads : nthreads;
threads.resize(nthreads); threads.resize(nthreads);
rseeds.resize(nthreads); rseeds.resize(nthreads);
for(J=0;J<nthreads;J++) for(J=0;J<nthreads;J++)
{ {
*rseed = dpr32_nextseed(*rseed); *rseed = dpr32_nextseed(*rseed);
@ -239,7 +245,7 @@ namespace rand
} }
narray<double> narray_rand(narray_size_t N, narray<double> narray_rand(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<double> ret; narray<double> ret;
@ -260,7 +266,7 @@ namespace rand
} }
narray<float> narray_randf(narray_size_t N, narray<float> narray_randf(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<float> ret; narray<float> ret;
@ -281,7 +287,7 @@ namespace rand
} }
narray<double> narray_randgauss(narray_size_t N, narray<double> narray_randgauss(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<double> ret; narray<double> ret;
@ -302,7 +308,7 @@ namespace rand
} }
narray<float> narray_randgaussf(narray_size_t N, narray<float> narray_randgaussf(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<float> ret; narray<float> ret;
@ -323,7 +329,7 @@ namespace rand
} }
narray<vect2> narray_randuvect2(narray_size_t N, narray<vect2> narray_randuvect2(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<vect2> ret; narray<vect2> ret;
@ -344,7 +350,7 @@ namespace rand
} }
narray<vect3> narray_randuvect3(narray_size_t N, narray<vect3> narray_randuvect3(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<vect3> ret; narray<vect3> ret;
@ -365,7 +371,7 @@ namespace rand
} }
narray<vect4> narray_randuvect4(narray_size_t N, narray<vect4> narray_randuvect4(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<vect4> ret; narray<vect4> ret;
@ -386,7 +392,7 @@ namespace rand
} }
narray<vect2f> narray_randuvect2f(narray_size_t N, narray<vect2f> narray_randuvect2f(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<vect2f> ret; narray<vect2f> ret;
@ -407,7 +413,7 @@ namespace rand
} }
narray<vect3f> narray_randuvect3f(narray_size_t N, narray<vect3f> narray_randuvect3f(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<vect3f> ret; narray<vect3f> ret;
@ -428,7 +434,7 @@ namespace rand
} }
narray<vect4f> narray_randuvect4f(narray_size_t N, narray<vect4f> narray_randuvect4f(narray_size_t N,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<vect4f> ret; narray<vect4f> ret;
@ -449,7 +455,7 @@ namespace rand
} }
int dpr32_randint(int32_t *seed, int low, int high) int dpr32_randint(narray_rands_t *seed, int low, int high)
{ {
int ret = 0; int ret = 0;
int val; int val;
@ -465,8 +471,8 @@ namespace rand
template<typename T> void narray_rand_threadf2( template<typename T> void narray_rand_threadf2(
narray<T> *out, narray<T> *out,
narray<int32_t> *rseeds, narray<narray_rands_t> *rseeds,
T (*randfunc)(int32_t *, T, T), T (*randfunc)(narray_rands_t *, T, T),
T low, T low,
T highexcl, T highexcl,
int threadnum, int threadnum,
@ -487,17 +493,17 @@ namespace rand
template<typename T> void narray_random_threadexec2( template<typename T> void narray_random_threadexec2(
narray<T> *out, narray<T> *out,
narray_size_t N, narray_size_t N,
T (*randfunc)(int32_t *, T, T), T (*randfunc)(narray_rands_t *, T, T),
T low, T low,
T highexcl, T highexcl,
int32_t *rseed narray_rands_t *rseed
) )
{ {
narray_size_t I; narray_size_t I;
int J; int J;
int nthreads; int nthreads;
std::vector<std::thread*> threads; std::vector<std::thread*> threads;
narray<int32_t> rseeds; narray<narray_rands_t> rseeds;
int res; int res;
res = out->resize(N); res = out->resize(N);
@ -559,7 +565,7 @@ namespace rand
} }
narray<int> narray_randint(narray_size_t N, int low, int highexcl, narray<int> narray_randint(narray_size_t N, int low, int highexcl,
int32_t *rseed) narray_rands_t *rseed)
{ {
narray<int> ret; narray<int> ret;

Loading…
Cancel
Save