libbpm: bpmnr/nr_ran1.c Source File

00001 
00005 #include <bpm/bpm_nr.h>
00006 
00013 double nr_ran1( long *idum ) {
00014   
00015   int j;
00016   long k;
00017   static long iy = 0;
00018   static long iv[RAN1_NTAB];
00019   double temp;
00020 
00021   if (*idum <= 0 || !iy)
00022     {
00023       if (-(*idum) < 1) *idum=1;
00024       else *idum = -(*idum);
00025       for (j=RAN1_NTAB+7;j>=0;j--)
00026         {
00027           k=(*idum)/RAN1_IQ;
00028           *idum=RAN1_IA*(*idum-k*RAN1_IQ)-RAN1_IR*k;
00029           if (*idum < 0) *idum += RAN1_IM;
00030           if (j < RAN1_NTAB) iv[j] = *idum;
00031         }
00032       iy=iv[0];
00033     }
00034   k=(*idum)/RAN1_IQ;
00035   *idum=RAN1_IA*(*idum-k*RAN1_IQ)-RAN1_IR*k;
00036   if (*idum < 0) *idum += RAN1_IM;
00037   j=iy/RAN1_NDIV;
00038   iy=iv[j];
00039   iv[j]=*idum;
00040   if((temp=RAN1_AM*iy) > RAN1_RNMX) return RAN1_RNMX;
00041   else return temp;
00042 }