1 /** 2 @file 3 @ingroup nr 4 */ 5 #include <bpm/bpm_nr.h> 6 7 /** 8 Random number generator as nicked from numerical recipes, c7.1, p280 9 10 @param idum random seed, note that the global seed is set by bpm_rseed 11 @return random number between 0 and 1 12 */ nr_ran1(long * idum)13double nr_ran1( long *idum ) { 14 15 int j; 16 long k; 17 static long iy = 0; 18 static long iv[RAN1_NTAB]; 19 double temp; 20 21 if (*idum <= 0 || !iy) 22 { 23 if (-(*idum) < 1) *idum=1; 24 else *idum = -(*idum); 25 for (j=RAN1_NTAB+7;j>=0;j--) 26 { 27 k=(*idum)/RAN1_IQ; 28 *idum=RAN1_IA*(*idum-k*RAN1_IQ)-RAN1_IR*k; 29 if (*idum < 0) *idum += RAN1_IM; 30 if (j < RAN1_NTAB) iv[j] = *idum; 31 } 32 iy=iv[0]; 33 } 34 k=(*idum)/RAN1_IQ; 35 *idum=RAN1_IA*(*idum-k*RAN1_IQ)-RAN1_IR*k; 36 if (*idum < 0) *idum += RAN1_IM; 37 j=iy/RAN1_NDIV; 38 iy=iv[j]; 39 iv[j]=*idum; 40 if((temp=RAN1_AM*iy) > RAN1_RNMX) return RAN1_RNMX; 41 else return temp; 42 } 43