real function urand(iy) integer iy c c urand is a uniform random number generator based on theory and c suggestions given in d.e. knuth (1969), vol 2. the integer iy c should be initialized to an arbitrary integer prior to the first call c to urand. the calling program should not alter the value of iy c between subsequent calls to urand. values of urand will be returned c in the interval (0,1). c integer ia,ic,itwo,m2,m,mic double precision halfm real s double precision datan,dsqrt data m2/0/,itwo/2/ if (m2 .ne. 0) go to 20 c c if first entry, compute machine integer word length c m = 1 10 m2 = m m = itwo*m2 if (m .gt. m2) go to 10 halfm = m2 c c compute multiplier and increment for linear congruential method c ia = 8*idint(halfm*datan(1.d0)/8.d0) + 5 ic = 2*idint(halfm*(0.5d0-dsqrt(3.d0)/6.d0)) + 1 mic = (m2 - ic) + m2 c c s is the scale factor for converting to floating point c s = 0.5/halfm c c compute next random number c 20 iy = iy*ia c c the following statement is for computers which do not allow c integer overflow on addition c if (iy .gt. mic) iy = (iy - m2) - m2 c iy = iy + ic c c the following statement is for computers where the c word length for addition is greater than for multiplication c if (iy/2 .gt. m2) iy = (iy - m2) - m2 c c the following statement is for computers where integer c overflow affects the sign bit c if (iy .lt. 0) iy = (iy + m2) + m2 urand = float(iy)*s return end