Committed on the Free edition of March Hare Software CVSNT Server. Upgrade to CVS Suite for more features and support: http://march-hare.com/cvsnt/ git-svn-id: http://moon:8086/svn/vhdl/trunk@272 cc03376c-175c-47c8-b038-4cd826a8556b
87 lines
1.7 KiB
C
87 lines
1.7 KiB
C
/* random.c
|
|
|
|
Author: Numerical recipes (ran1, gaussian), Jon Hamkins (others)
|
|
Revised by: Jon Hamkins
|
|
Date: 4-16-98
|
|
*/
|
|
|
|
#include <math.h>
|
|
|
|
#define IA 16807
|
|
#define IM 2147483647
|
|
#define AM (1.0/IM)
|
|
#define IQ 127773
|
|
#define IR 2836
|
|
#define NTAB 32
|
|
#define NDIV (1+(IM-1)/NTAB)
|
|
#define EPS 1.2e-14
|
|
#define RNMX (1.0-EPS)
|
|
/* From Numerical Recipes, p. 280, modified from float to double */
|
|
/* returns a uniform deviate in (0,1) */
|
|
double ran1(long *idum)
|
|
{
|
|
int j;
|
|
long k;
|
|
static long iy=0;
|
|
static long iv[NTAB];
|
|
double temp;
|
|
|
|
if (*idum <=0 || !iy) {
|
|
if (-(*idum) < 1) *idum=1;
|
|
else *idum = -(*idum);
|
|
for (j=NTAB+7; j>=0; j--) {
|
|
k=(*idum)/IQ;
|
|
*idum=IA*(*idum-k*IQ)-IR*k;
|
|
if (*idum < 0) *idum += IM;
|
|
if (j < NTAB) iv[j] = *idum;
|
|
}
|
|
iy = iv[0];
|
|
}
|
|
k=(*idum)/IQ;
|
|
*idum=IA*(*idum-k*IQ)-IR*k;
|
|
if (*idum < 0) *idum += IM;
|
|
j=iy/NDIV;
|
|
iy=iv[j];
|
|
iv[j] = *idum;
|
|
if ((temp=AM*iy) > RNMX) return RNMX;
|
|
else return temp;
|
|
}
|
|
#undef IA
|
|
#undef IM
|
|
#undef AM
|
|
#undef IQ
|
|
#undef IR
|
|
#undef NTAB
|
|
#undef NDIV
|
|
#undef EPS
|
|
#undef RNMX
|
|
|
|
/* Generate a N(0,1) r.v. */
|
|
double gaussian(long *idum)
|
|
{
|
|
static int iset=0;
|
|
static double gset;
|
|
double fac,r,v1,v2;
|
|
double ran1();
|
|
if (iset == 0) {
|
|
do {
|
|
v1=2.0*ran1(idum)-1.0;
|
|
v2=2.0*ran1(idum)-1.0;
|
|
r=v1*v1+v2*v2;
|
|
} while (r >= 1.0 || r == 0.0);
|
|
fac=sqrt(-2.0*log(r)/r);
|
|
gset=v1*fac;
|
|
iset=1;
|
|
return v2*fac;
|
|
} else {
|
|
iset=0;
|
|
return gset;
|
|
}
|
|
}
|
|
|
|
/* generate a random bit */
|
|
int random_bit(long *idum)
|
|
{
|
|
return((ran1(idum)<0.5) ? 0 : 1);
|
|
}
|