\ Pseudo random number generator in ANS Forth \ version of November 24th, 1998 - 21:21 \ \ Leaves a pseudo random number in the range (0,1) \ on fp stack. \ \ --------------------------------------------------- \ (c) Copyright 1998 Julian V. Noble. \ \ Permission is granted by the author to \ \ use this software for any application pro- \ \ vided this copyright notice is preserved. \ \ --------------------------------------------------- \ \ Based on GGUBS algorithm: s' = 16807*s mod (2^32-1) \ P. Bratley, B.L. Fox and L.E. Schrage, A guide to simulation \ (Springer, Berlin, 1983). \ \ To simplify transport to 16-bit machines the 32-bit \ modular division is performed by synthetic division: \ note that \ \ bigdiv = divis * m1 + m2 \ \ so that ( [n] means "largest integer <= n" ) \ \ s' = s*m1 - [s*m1/b]*b = m1 * (s - [s/d]*d) - m2 * [s/d] \ \ Environmental dependences: \ \ 1. assumes at least 32-bit DOUBLEs \ 2. needs FLOATING and DOUBLE wordsets \ 3. assumes separate floating point stack \ \ MARKER -rand 2VARIABLE seed 21474.83647 D>F FCONSTANT bigdiv \ 2^31-1 1277.73 D>F FCONSTANT divis 16807. D>F FCONSTANT m1 2836. D>F FCONSTANT m2 : (rand) ( adr --) ( f: -- seed') DUP 2@ D>F divis FOVER FOVER ( f: s d s d) F/ F>D 2DUP D>F ( [s/d]) ( f: s d [s/d]) F* F- ( f: s-d*[s/d]) m1 F* ( f: m1*[s mod d]) D>F m2 F* F- FDUP F>D ( adr d) ROT 2! ; : prng ( f: -- random#) seed (rand) bigdiv ( f: -- seed 2**31-1) FSWAP FDUP F0< ( -- f) ( f: -- 2**31-1 seed) IF FOVER F+ THEN FSWAP F/ ; : test 0.1 seed 2! 1000 0 DO prng FDROP LOOP seed 2@ D. ; \ TEST 522329230 ok