glennrp@BRL.MIL (Glenn Randers-Pehrson, WMB) (10/22/88)
[...]
There is an excellent article on random number generators in this
month's Communications of the ACM. Motivated by this article, and by my
unhappy experience with the Silicon Graphics random number generator,
I checked the generators for various BRL machines. One uniform
characteristic of all of them is that they aren't documented very well.
None of the user documentation gives the particular algorithm and parameters
used, so I had to do some experiments to find out the algorithms.
On the Crays, my experiments with the RANF() function provided by CFT77
on both bob.brl.mil and patton.brl.mil indicated that the function is
f(z) = mod ( 44,485,709,377,909 * z, 2**46)
If user provides an even numbered seed, it is converted to the next
higher number. i.e., the sequences beginning with seed "2" and seed "3"
are identical. The algorithm should provide one of two sequences of 2**44
random numbers before recycling.
On adm.brl.mil (Gould) and tgr.brl.mil (VAX runing BSD4.3), the function
is IRANDM().
The algorithm is more complex than the above, being a nonlinear additive
feedback scheme, whose paramaters I do not know. On thud.brl.mil (Alliant),
the same function is provided and is called IRAND().
Also on the Gould and VAXEN is IRAND(), whose algorithm is
f(z) = mod (1 103 515 245 * z + 12345, 2**31)
This is one of the algorithms mentioned in the article as being deficient,
and the on-line documentation at BRL recommends against using it.
On taurus.brl.mil (Iris 2500T),
I could not figure out the algorithm, but it has the characteristic
of collapsing in fairly short order into a cycle of 21494 numbers.
I was trying to use ran() in a dithering scheme some time ago, but
found that it introduced diagonal artifacts in the images.
There does not seem to be a way of specifying the seed. Ran(negative)
is supposed to start a new sequence, but the seed for the new sequence
seems to depend on the number of times ran(positive) was called.
It is also extremely slow, compared to the "minimal standard random
number generator" presented in the article. I recommend the latter
be used on the IRIS instead of the built-in function "ran()". Here is
the Fortran code:
function msran(setseed)
c minimal standard random number generator
c Stephen K. Park and Keith W. Miller
c Reference: Random number generators: good ones are hard to find
c CACM vol 31 no 10, Oct 88
c f(z) = mod ( 16,807 * z, 2**31 -1)z
implicit integer(a-z)
save seed
data seed / 1 /,
+ a / 16 807 /,
+ m / 2 147 483 647 /,
+ q / 127 773 /,
+ r / 2 836 /
if (setseed .ne. 0) seed = setseed
hi = seed/q
lo = mod(seed,q)
test = a*lo - r*hi
if(test.gt.0)then
seed = test
else
seed = test+m
endif
msran = seed
return
end
Glenn Randers-Pehrson <glennrp@brl.mil>gwyn@BRL.MIL (Doug Gwyn, VLD/VMB) (10/23/88)
Excuse me, but I would not call the CACM article "excellent". The authors were simply trying to push for adoption of their particular favorite as a "standard". There are considerably better random number generators available on BRL UNIX systems if you know where to look.