redo random/

db4
Doug Coleman 2008-03-19 16:18:03 -05:00
parent 40aab45282
commit b3527a17df
11 changed files with 209 additions and 110 deletions

View File

@ -0,0 +1,36 @@
USING: kernel math sequences namespaces
math.miller-rabin combinators.cleave combinators.lib
math.functions new-slots accessors random ;
IN: random.blum-blum-shub
! TODO: take (log log M) bits instead of 1 bit
! Blum Blum Shub, M = pq
TUPLE: blum-blum-shub x n ;
C: <blum-blum-shub> blum-blum-shub
: generate-bbs-primes ( numbits -- p q )
#! two primes congruent to 3 (mod 4)
[ [ random-prime ] curry [ 4 mod 3 = ] generate ] dup bi ;
IN: crypto
: <blum-blum-shub> ( numbits -- blum-blum-shub )
#! returns a Blum-Blum-Shub tuple
generate-bbs-primes *
[ find-relative-prime ] keep
blum-blum-shub construct-boa ;
! 256 make-bbs blum-blum-shub set-global
: next-bbs-bit ( bbs -- bit )
#! x = x^2 mod n, return low bit of calculated x
[ [ x>> 2 ] [ n>> ] bi ^mod ]
[ [ >>x ] keep x>> 1 bitand ] bi ;
IN: crypto
! : random ( n -- n )
! ! #! Cryptographically secure random number using Blum-Blum-Shub 256
! [ log2 1+ random-bits ] keep dupd >= [ -1 shift ] when ;
M: blum-blum-shub random-32 ( bbs -- r )
;

View File

@ -0,0 +1,11 @@
USING: kernel random math new-slots accessors ;
IN: random.dummy
TUPLE: random-dummy i ;
C: <random-dummy> random-dummy
M: random-dummy seed-random ( seed obj -- )
(>>i) ;
M: random-dummy random-32 ( obj -- r )
[ dup 1+ ] change-i drop ;

View File

@ -1,17 +1,17 @@
USING: help.markup help.syntax math ;
IN: random
IN: random.mersenne-twister
ARTICLE: "random-numbers" "Generating random integers"
"The " { $vocab-link "random" } " vocabulary implements the ``Mersenne Twister'' pseudo-random number generator algorithm."
{ $subsection init-random }
! { $subsection init-random }
{ $subsection (random) }
{ $subsection random } ;
ABOUT: "random-numbers"
HELP: init-random
{ $values { "seed" integer } }
{ $description "Initializes the random number generator with the given seed. This word is called on startup to initialize the random number generator with the current time." } ;
! HELP: init-random
! { $values { "seed" integer } }
! { $description "Initializes the random number generator with the given seed. This word is called on startup to initialize the random number generator with the current time." } ;
HELP: (random)
{ $values { "rand" "an integer between 0 and 2^32-1" } }

View File

@ -0,0 +1,30 @@
USING: kernel math random namespaces random.mersenne-twister
sequences tools.test ;
IN: random.mersenne-twister.tests
USE: tools.walker
: check-random ( max -- ? )
dup >r random 0 r> between? ;
[ t ] [ 100 [ drop 674 check-random ] all? ] unit-test
: make-100-randoms
[ 100 [ 100 random , ] times ] { } make ;
: test-rng ( seed quot -- )
>r <mersenne-twister> r> with-random ;
[ f ] [ 1234 [ make-100-randoms make-100-randoms = ] test-rng ] unit-test
[ 1333075495 ] [
0 [ 1000 [ drop \ random get random-32 drop ] each \ random get random-32 ] test-rng
] unit-test
[ 1575309035 ] [
0 [ 10000 [ drop \ random get random-32 drop ] each \ random get random-32 ] test-rng
] unit-test
[ 3 ] [ 101 [ 3 random-bytes length ] test-rng ] unit-test
[ 33 ] [ 101 [ 33 random-bytes length ] test-rng ] unit-test
[ t ] [ 101 [ 100 random-bits log2 90 > ] test-rng ] unit-test

View File

@ -0,0 +1,80 @@
! Copyright (C) 2005, 2008 Doug Coleman.
! See http://factorcode.org/license.txt for BSD license.
! mersenne twister based on
! http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
USING: arrays kernel math namespaces sequences
system init new-slots accessors
math.ranges combinators.cleave circular random ;
IN: random.mersenne-twister
<PRIVATE
: new-nth ( seq i -- elt ) swap nth ; inline
: new-set-nth ( seq obj n -- seq ) pick set-nth ; inline
TUPLE: mersenne-twister seq i ;
: mt-n 624 ; inline
: mt-m 397 ; inline
: mt-a HEX: 9908b0df ; inline
: mt-hi HEX: 80000000 bitand ; inline
: mt-lo HEX: 7fffffff bitand ; inline
: wrap ( x n -- y ) 2dup >= [ - ] [ drop ] if ; inline
: mt-wrap ( x -- y ) mt-n wrap ; inline
: set-generated ( mt y from-elt to -- )
>r >r [ 2/ ] [ odd? mt-a 0 ? ] bi
r> bitxor bitxor r> new-set-nth drop ; inline
: calculate-y ( mt y1 y2 -- y )
>r over r>
[ new-nth mt-hi ] [ new-nth mt-lo ] 2bi* bitor ; inline
: (mt-generate) ( mt-seq n -- y to from-elt )
[ dup 1+ mt-wrap calculate-y ]
[ mt-m + mt-wrap new-nth ]
[ nip ] 2tri ;
: mt-generate ( mt -- )
[ seq>> mt-n [ dupd (mt-generate) set-generated ] with each ]
[ 0 >>i drop ] bi ;
: init-mt-first ( seed -- seq )
>r mt-n 0 <array> r>
HEX: ffffffff bitand 0 new-set-nth ;
: init-mt-formula ( seq i -- f(seq[i]) )
tuck new-nth dup -30 shift bitxor 1812433253 * +
1+ HEX: ffffffff bitand ;
: init-mt-rest ( seq -- )
mt-n 1- [0,b) [
dupd [ init-mt-formula ] keep 1+ new-set-nth drop
] with each ;
: init-mt-seq ( seed -- seq )
init-mt-first dup init-mt-rest ;
: mt-temper ( y -- yt )
dup -11 shift bitxor
dup 7 shift HEX: 9d2c5680 bitand bitxor
dup 15 shift HEX: efc60000 bitand bitxor
dup -18 shift bitxor ; inline
PRIVATE>
: <mersenne-twister> ( seed -- obj )
init-mt-seq 0 mersenne-twister construct-boa
dup mt-generate ;
M: mersenne-twister seed-random ( mt seed -- )
init-mt-seq >>seq drop ;
M: mersenne-twister random-32 ( mt -- r )
dup [ seq>> ] [ i>> ] bi
dup mt-n < [ drop 0 pick mt-generate ] unless
new-nth mt-temper
swap [ 1+ ] change-i drop ;
[ millis <mersenne-twister> \ random set-global ] "random" add-init-hook

View File

@ -1,15 +0,0 @@
USING: kernel math random namespaces sequences tools.test ;
IN: random.tests
: check-random ( max -- ? )
dup >r random 0 r> between? ;
[ t ] [ 100 [ drop 674 check-random ] all? ] unit-test
: make-100-randoms
[ 100 [ 100 random , ] times ] { } make ;
[ f ] [ make-100-randoms make-100-randoms = ] unit-test
[ 1333075495 ] [ 0 init-random 1000 [ drop (random) drop ] each (random) ] unit-test
[ 1575309035 ] [ 0 init-random 10000 [ drop (random) drop ] each (random) ] unit-test

112
extra/random/random.factor Executable file → Normal file
View File

@ -1,107 +1,39 @@
! Copyright (C) 2005, 2007 Doug Coleman.
! Copyright (C) 2008 Doug Coleman.
! See http://factorcode.org/license.txt for BSD license.
! mersenne twister based on
! http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
USING: arrays kernel math namespaces sequences
system init alien.c-types ;
USING: alien.c-types kernel math namespaces sequences
io.backend ;
IN: random
<PRIVATE
HOOK: os-crypto-random-bytes io-backend ( n -- byte-array )
HOOK: os-random-bytes io-backend ( n -- byte-array )
HOOK: os-crypto-random-32 io-backend ( -- r )
HOOK: os-random-32 io-backend ( -- r )
TUPLE: mersenne-twister seed seq i ;
GENERIC: seed-random ( tuple seed -- )
GENERIC: random-32 ( tuple -- r )
C: <mersenne-twister> mersenne-twister
: (random-bytes) ( tuple n -- byte-array )
[ drop random-32 ] with map >c-uint-array ;
: mt-n 624 ; inline
: mt-m 397 ; inline
: mt-a HEX: 9908b0df ; inline
: mt-hi HEX: 80000000 ; inline
: mt-lo HEX: 7fffffff ; inline
DEFER: random
SYMBOL: mt
: random-bytes ( n -- r )
[
4 /mod zero? [ 1+ ] unless
\ random get swap (random-bytes)
] keep head ;
: mt-seq ( -- seq )
mt get mersenne-twister-seq ; inline
: mt-nth ( n -- nth )
mt-seq nth ; inline
: mt-i ( -- i )
mt get mersenne-twister-i ; inline
: mti-inc ( -- )
mt get [ mersenne-twister-i 1+ ] keep set-mersenne-twister-i ; inline
: set-mt-ith ( y i-get i-set -- )
>r mt-nth >r
[ 2/ ] keep odd? mt-a 0 ? r> bitxor bitxor r>
mt-seq set-nth ; inline
: mt-y ( y1 y2 -- y )
mt-nth mt-lo bitand
>r mt-nth mt-hi bitand r> bitor ; inline
: mod* ( x n -- y )
#! no floating point
2dup >= [ - ] [ drop ] if ; inline
: (mt-generate) ( n -- y n n+(mt-m) )
dup [ 1+ 624 mod* mt-y ] keep [ mt-m + 624 mod* ] keep ;
: mt-generate ( -- )
mt-n [ (mt-generate) set-mt-ith ] each
0 mt get set-mersenne-twister-i ;
: init-mt-first ( seed -- seq )
>r mt-n 0 <array> r>
HEX: ffffffff bitand 0 pick set-nth ;
: init-mt-formula ( seq i -- f(seq[i]) )
dup rot nth dup -30 shift bitxor
1812433253 * + HEX: ffffffff bitand 1+ ; inline
: init-mt-rest ( seq -- )
mt-n 1 head* [
[ init-mt-formula ] 2keep 1+ swap set-nth
] with each ;
: mt-temper ( y -- yt )
dup -11 shift bitxor
dup 7 shift HEX: 9d2c5680 bitand bitxor
dup 15 shift HEX: efc60000 bitand bitxor
dup -18 shift bitxor ; inline
PRIVATE>
: init-random ( seed -- )
global [
dup init-mt-first
[ init-mt-rest ] keep
0 <mersenne-twister> mt set
mt-generate
] bind ;
: (random) ( -- rand )
global [
mt-i dup mt-n < [ drop mt-generate 0 ] unless
mt-nth mti-inc
mt-temper
] bind ;
: big-random ( n -- r )
[ drop (random) ] map >c-uint-array byte-array>bignum ;
: random-256 ( -- r ) 8 big-random ; inline
: random-bits ( n -- r ) 2^ random ;
: random ( seq -- elt )
dup empty? [
drop f
] [
[
length dup log2 31 + 32 /i big-random swap mod
length dup log2 7 + 8 /i
random-bytes byte-array>bignum swap mod
] keep nth
] if ;
[ millis init-random ] "random" add-init-hook
: with-random ( tuple quot -- )
\ random swap with-variable ; inline

View File

@ -0,0 +1,22 @@
USING: alien.c-types io io.files io.nonblocking kernel
namespaces random io.encodings.binary singleton ;
IN: random.unix
SINGLETON: unix-random
: file-read-unbuffered ( n path -- bytes )
over default-buffer-size [
binary <file-reader> [ read ] with-stream
] with-variable ;
M: unix-random os-crypto-random-bytes ( n -- byte-array )
"/dev/random" file-read-unbuffered ;
M: unix-random os-random-bytes ( n -- byte-array )
"/dev/urandom" file-read-unbuffered ;
M: unix-random os-crypto-random-32 ( -- r )
4 os-crypto-random-bytes *uint ;
M: unix-random os-random-32 ( -- r )
4 os-random-bytes *uint ;

View File

@ -0,0 +1,3 @@
IN: random.windows
! M: windows-io