This is G o o g l e's cache of http://home.earthlink.net/~neilbawd/gbflip.html.
G o o g l e's cache is the snapshot that we took of the page as we crawled the web.
The page may have changed since that time. Click here for the current page without highlighting.
To link to or bookmark this page, use the following url: http://www.google.com/search?q=cache:MIEYp5EyesUC:home.earthlink.net/~neilbawd/gbflip.html+&hl=en&ie=UTF-8


Google is not affiliated with the authors of this page nor responsible for its content.

GB_FLIP

GB_FLIP

Wil Baden 1993-09-11

Knuth's recommended Subtractive RNG from The Stanford GraphBase.

This is a port of GB_FLIP, Donald E. Knuth's random number generator package from his The Stanford GraphBase, 1993, ISBN 0-201-54275-7.

Words have been renamed to agree with Forth naming conventions.

These routines give identical results whether working in C or Forth.

From the book, replacing C jargon with Forth jargon.

INIT-RAND           ( n -- )
To use the routines in this file, first call the function seed INIT-RAND.
NEXT-RAND           ( -- +n )
Subsequent use of NEXT-RAND will return pseudo-random integers between 0 and 2^31 - 1, inclusive.

GraphBase programs are designed to produce identical results on almost all existing computers and operating systems. An improved version of the portable subtractive method recommended in Seminumerical Algorithms, Section 3.6, is used to generate random numbers in the routines below. The period length is at least 2^55 - 1, and is in fact plausibly conjectured to be 2^85 - 2^30 for all but at most one choice of the seed value. The low-order bits of the generated numbers are just as random as the high-order bits.

A validation program is provided so installers can tell if it's working properly.

UNIF-RAND           ( m -- u )
Uniform integers. Here is a simple routine that produces a uniform integer between 0 and m-1 inclusive, when m is any positive integer less than 2^31. It avoids the bias toward small values that would occur if we simply calculated NEXT-RAND m MOD. (The bias is insignificant when m is small, but can be serious when m is large. For example, if m is approximately 2^32/3, the simple remainder algorithm would give an answer less than m / 2 about 2/3 of the time.)
This routine consumes fewer than two random numbers, on the average, for any fixed m.
TEXT Program Text 1
 
( Standard Forth CORE and CORE EXT with NOT )

: MOD-DIFF  S" - 2147483647 AND " EVALUATE ; IMMEDIATE

CREATE    A     56 CELLS ALLOT

VARIABLE  FPTR

: FLIP-CYCLE                   ( -- +random )
    A CELL+     32 CELLS A +              ( ii jj)
    BEGIN
        OVER >R                                   ( R: ii)
            OVER @ OVER @ MOD-DIFF        ( . . diff)
        R> !                              ( ii jj)( R: )
        SWAP CELL+ SWAP    CELL+
        DUP  55 CELLS A + U>
    UNTIL
    DROP                                  ( ii)
    A CELL+                               ( ii jj)
    BEGIN
        OVER >R                             ( R: ii)
            OVER @ OVER @ MOD-DIFF        ( . . diff)
        R> !                              ( ii jj)( R: )
        SWAP CELL+ SWAP    CELL+
        OVER  55 CELLS A + U>
    UNTIL                                 2DROP
    54 CELLS A + FPTR !
    55 CELLS A + @                        ( +random)
    ;

: INIT-RAND                      ( seed -- )
    -1 A !
    0 MOD-DIFF
    DUP 55 CELLS A + !
    DUP  1                                ( seed prev next)
    1 21 DO
        DUP I CELLS A + !
        MOD-DIFF                          ( seed next)
        >R                                ( seed)( R: next)
            DUP 1 AND IF
                2/
                [ HEX ] 40000000 [ DECIMAL ] +
            ELSE
                2/
            THEN
        R>                                ( seed next)( R: )
        OVER MOD-DIFF
        I CELLS A + @ SWAP                ( seed prev next)
    I  21 +  55 MOD  I -  +LOOP           2DROP DROP
    FLIP-CYCLE DROP
    FLIP-CYCLE DROP
    FLIP-CYCLE DROP
    FLIP-CYCLE DROP
    FLIP-CYCLE DROP
    ;

: NEXT-RAND                      ( -- +random )
    FPTR @ @                        ( +random)
    DUP 0< IF
        DROP FLIP-CYCLE
    ELSE
        -1 CELLS FPTR +!
    THEN ;

: UNIF-RAND                      ( +range -- +random )
    >R                                    ( )( R: +range)
        [ HEX ] 80000000 [ DECIMAL ]      ( t)
        DUP 0 R@ UM/MOD DROP -
        BEGIN
            NEXT-RAND                     ( t +random)
            2DUP U> NOT
        WHILE
            DROP                          ( t)
        REPEAT                            ( t +random)
        NIP                               ( +random)
    R> MOD                                ( +random)( R: )
    ;

MARKER VALIDATE 

: TEST-FLIP
    -314159 INIT-RAND
    NEXT-RAND
        119318998 <> ABORT" Failure on the first try. "

    133 0 DO  NEXT-RAND DROP  LOOP

    [ HEX ] 55555555 [ DECIMAL ] UNIF-RAND
        748103812 <> ABORT" Failure on the second try. "
    ; 
    
TEST-FLIP VALIDATE

( With one's complement arithmetic, try --
: MOD-DIFF
      - DUP 0< IF   [ HEX ] 40000000 + 40000000 + [ DECIMAL ]   THEN
;
)

Go back to Neil Bawd's home page.