--- capa/capa51/pProj/ranlib.c 1999/09/28 21:26:21 1.1 +++ capa/capa51/pProj/ranlib.c 2000/02/22 18:19:02 1.2 @@ -48,7 +48,7 @@ S10: S20: if(!(min(aa,bb) > 1.0)) goto S100; /* - Alborithm BB + Algorithm BB Initialize */ if(qsame) goto S30; @@ -436,6 +436,7 @@ S20: S30: return gennf; } + float gennor(float av,float sd) /* ********************************************************************** @@ -458,11 +459,51 @@ float gennor(float av,float sd) ********************************************************************** */ { -static float gennor; +float gennor; +float tmp_f; + + tmp_f = snorm(); + + gennor = sd*tmp_f+av; + return (gennor); +} + +float capa_gennor(double *num_d, float av,float sd) +/* +********************************************************************** + float gennor(float av,float sd) + GENerate random deviate from a NORmal distribution + Function + Generates a single random deviate from a normal distribution + with mean, AV, and standard deviation, SD. + Arguments + av --> Mean of the normal distribution. + sd --> Standard deviation of the normal distribution. + Method + Renames SNORM from TOMS as slightly modified by BWB to use RANF + instead of SUNIF. + For details see: + Ahrens, J.H. and Dieter, U. + Extensions of Forsythe's Method for Random + Sampling from the Normal Distribution. + Math. Comput., 27,124 (Oct. 1973), 927 - 937. +********************************************************************** +*/ +{ +float gen_num; +float tmp_f; - gennor = sd*snorm()+av; - return gennor; + tmp_f = snorm(); + + gen_num = sd*tmp_f+av; + /* printf("SNORM()=%f,GENNOR()=%f,%f*%f+%f\n",tmp_f,gen_num,sd,tmp_f,av); */ + *num_d = (double)gen_num; + + gen_num = (float)37.358341; + return (gen_num); } + + void genprm(long *iarray,int larray) /* ********************************************************************** @@ -1300,13 +1341,51 @@ float ranf(void) */ { static float ranf; +long tmp_l; +double tmp_d; /* 4.656613057E-10 is 1/M1 M1 is set in a data statement in IGNLGI and is currently 2147483563. If M1 changes, change this also. */ - ranf = ignlgi()*4.656613057E-10; + tmp_l = ignlgi(); + tmp_d = (double)tmp_l * (double)4.656613057E-10; + ranf = (float)tmp_d; + /* printf("RANF()=%f\n",ranf); */ return ranf; } +float capa_ranf(void) +/* +********************************************************************** + float ranf(void) + RANDom number generator as a Function + Returns a random floating point number from a uniform distribution + over 0 - 1 (endpoints of this interval are not returned) using the + current generator + This is a transcription from Pascal to Fortran of routine + Uniform_01 from the paper + L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package + with Splitting Facilities." ACM Transactions on Mathematical + Software, 17:98-111 (1991) +********************************************************************** +*/ +{ + float ran_f; + long my_ran; + double my_doub; +/* + 4.656613057E-10 is 1/M1 M1 is set in a data statement in IGNLGI + and is currently 2147483563. If M1 changes, change this also. +*/ + my_ran = ignlgi(); + /* printf("MY_ignlgi=%ld -- first time\n",my_ran); */ + /* ran_f = my_ran * 4.656613057E-10; */ + + my_doub = (double)my_ran * (double)4.656613057E-10; + printf("MY_ranf in double=%.15g -- first time\n",my_doub); + ran_f = (float)my_doub; + return (ran_f); +} + void setgmn(float *meanv,float *covm,long p,float *parm) /* **********************************************************************