'randx_randy_mod_congr_4.bas

    ' johnf 12 Sept 2004

    'A demo of the mod-congruent random number generator & its limitations.
    'Since this is (presumably) the mechanism Carl uses, it also shows the problems.
    '  of these pseude-random generators.
    'We English, realising this, made our first 'random generating' computer, called 'Ernie'
    '  with a non algorithmic generator based on thermal noise in gas-filled thermionic tubes.
    '  (used for what were called Premium Bonds- savings bonds with a chance of a prize)
    '  ( Electronic Random Number == ERN)

    'NB this prog version does not clear the graphic window for a second run..

    NOMAINWIN

    WindowWidth  = 530
    WindowHeight = 460
    UpperLeftX   = INT(( DisplayWidth  -WindowWidth)  /2)
    UpperLeftY   = INT(( DisplayHeight -WindowHeight) /2)

    graphicbox  #main.gb1,  10,  10, 200, 200
    graphicbox  #main.gb2, 310,  10, 200, 200
    textbox     #main.t1,   10, 220, 500,  20
    texteditor  #main.t2,   20, 250, 470, 150
    button      #main.b0, "Start",    [start], UL, 230,  10
    button      #main.b1, "Faster?",  speed, UL, 230, 160
    button      #main.b2, "Carl's",   typec, UL, 230,  60
    button      #main.b3, "johnf's",  typej, UL, 230, 110

    global r, c, n, k, d, flag     'flag selects my randoms or Carl's

    d =1000     'initially slow by calling a 1000ms delay in loop. Later toggle 0 or 1/10s.

    r =1234     'any seed between 0 and n- it becomes a succession of new
                'random values all in range 0 to n
    k = 106     'these three numbers are carefully chosen
                'chain code gives n numbers at random before repeating
                'r is replaced by (r * k +c) mod n
    c =1283
    n =6075     'Change these numbers at your peril! Look up in a good comp. math ref'ce!
                'eg 'Numerical Recipes' from Cambridge University Press

                ' Note how since successive values are calculated & therefore related
                'the results of successive pairs used as coordinates will be regular patterns.
                'It will repeat itself after n calls of the function.
                'However the distribution brings up every number equally often.
                'Such sequences are mathemagically called 'pseudorandom'.

    dim bin( 200)

    Open "rnd(x) v. modcongruent!" for window as #main
    #main "trapclose [quit]"

    #main.t2 "!font arial 8"
    #main.t2 " Left window shows result of throwing two random #s between 0 & 199"
    #main.t2 "     & using them as x,y Cartesian coordinates."
    #main.t2 " On right an (inverted) histogram shows how many of each # have been thrown."
    #main.t2 " The histogram should shows that all numbers come up on average equally often."
    #main.t2 " However the x,y plot shows that since successive values are related"
    #main.t2 "     there is a failure to 'hit' all locations in this 2-D space."
    #main.t2 " If you choose my simple one, they seem all to come out equally."
    #main.t2 "     but highly correlated in pairs"
    #main.t2 "Carl's LB function gives less correlation but also less visible equal prob'y"
    #main.t2, "In a mod-congruent pseudorandom number generator, r is replaced by (r +kc) modulo n"
    #main.t2, "The choice of k, c and n is very critical- look it up in say p276,"
    #main.t2, "http://www.library.cornell.edu/nr/bookcpdf/c7-1.pdf"
    #main.t2, "Which is a fantastic resource if you don't own 'Numerical Recipes'!"

    wait

    [start]
    for x = 1 to 10 *n/2            'choose random locations for x and y locations to draw at
                                    'By calling rand() n times we will have seen all values ( n /2 *2)
                                    'so once the initial set have been seen they are re-cycled

        a       = rand( 200)        'find a new random integer a in range 0 to 200
        bin( a) =bin( a) +1         'count number of times it has come up by incrementing its bin
        locx    = a                 'use as x coordinate
        #main.gb2, "color "; a; " "; bin( a); " "; 255-a    'colour code dot to be plotted.
        #main.gb2 "down ; set "; a ; " "; bin( a)           ' and increment histogram

        b       = rand( 200)        'find a second random integer b in same range
        bin( b) =bin( b) +1
        locy    = b                 'use as y coordinate
        #main.gb2, "color "; b; " "; bin( b); " "; 255-b
        #main.gb2 "down ; set "; b ; " "; bin( b)

        #main.gb1 "down ; color black ; set "; locx; " "; locy

        #main.t1 "Latest values were "; right$( "____" +str$( a), 3); " & "; right$( "____" +str$( b), 3); " after "; x; " throws."
  'NB ?Bug? if you use "    " rather than "0000" or "____" why does this fail????

        if bin( a) >200 or bin( b) >200 then [hold]
        call Pause d
'  NB Without the discard the repeated plot commands cause exhaustion of memory resourses.
'     This memory leak is NOT reclaimed when program ends.
'     I had hoped that doing discards here would help- it does not. Hence rem'd out.
        '#main.gb1 "discard"
        '#main.gb2 "discard"
        scan
    next x

    #main.gb1 "flush"
    #main.gb2 "flush"

    [hold]
    Wait

    [about]


    [quit]
    print #main.gb1, "getbmp drawing 1 1 530 370"
    bmpsave "drawing", "rndcong.bmp"
    close #main : END

    function rand( range)  's is a dummy & not used...
    r1  =( r *k +c)        'the value here now has to be reduced mod n
    j   =int( r1 /n)       'find r1 div n, the whole number part
    r   =r1 -j *n          'these three lines have simulated r --> (r*k +c) mod n
    if flag =1 then
        res =int( 200 *rnd( 1)) 'if you want to use Carl's random numbers.
    else
        res =int( range *r /n) 'if you want to use my random generator.
    end if
    rand =res
    end function

    sub Pause mil
    t =time$( "milliseconds")
    while time$( "milliseconds") <( t +mil)
    wend
    end sub

    sub speed h$
    if d =0 then
        d =100
        #main.b1 "Slow"
    else
        d=0
        #main.b1 "Fast"
    end if
    end sub

    sub typec h$
    flag =1
    end sub

    sub typej h$
    flag =0
    end sub