(49G) An Inverse Congruential Pseudo Random Number Generator

[Gerald H]

Iset takes takes two positive integers S & P from the stack & sets the six parameters

Iseed

Ia

Ic

Im

Iπ

Ipow

S is the size of the prime Iπ & P the power Ipow of the prime used in the formula

Iseed := (Iseed ^ -1) * Ia + Ic MOD Iπ ^ Ipow

Iπ ^ Ipow is used to set the size of the other parameters.

Actuate IXNEXT to return the next pseudo-random integer to the stack.

ISET uses one sub-programme, see below.

Code:


ISET

::
  CK2&Dispatch
  2REAL
  ::
    DUP'
    ID Ipow
    SAFESTO
    COERCE2
    '
    ID Iseed
    '
    ID Ia
    '
    ID Ic
    5ROLL
    ID x005
    FPTR2 ^Prime+
    DUP'
    ID Iπ
    SAFESTO
    5ROLL
    FPTR2 ^PPow#
    DUP'
    ID Im
    SAFESTO
    DUP
    LENHXS
    BINT3
    ZERO_DO
    DUP
    ID x005
    3PICK
    FPTR2 ^ZMod
    4ROLL
    SAFESTO
    LOOP
    2DROP
  ;
;

IXNEXT

::
  CK0
  ::
    '
    ID Im
    '
    ID Ic
    '
    ID Ia
    '
    ID Iseed
    BINT4
    ZERO_DO
    SAFE@
    NcaseTYPEERR
    DUPTYPEZINT?
    NcaseTYPEERR
    4UNROLL
    LOOP
    4PICK
    FPTR2 ^ModDiv
    FPTR2 ^QAdd
    OVER
    FPTR2 ^ZMod
    FPTR2 ^DupZIsNeg?
    ITE
    FPTR2 ^RADDext
    SWAPDROP
    DUP'
    ID Iseed
    ?STO_HERE
  ;
;

x005

::
  DUP
  BINT12
  #/
  SWAPONE
  #MIN
  #+
  NULL$SWAP
  ZERO_DO
  %RAN
  FPTR 2 A5
  BINT9
  BINT20
  SUB$
  &$
  LOOP
  SWP1+
  1_#1-SUB$
  FPTR2 ^S>Z
;