extra/random/xoshiro/xoshiro.factor

   1 ! Copyright (C) 2018 Doug Coleman.
   2 ! See http://factorcode.org/license.txt for BSD license.
   3 USING: accessors alien.c-types classes.struct kernel locals math
   4 math.bitwise random sequences slots.syntax ;
   5 IN: random.xoshiro
   6
   7 ! http://xoshiro.di.unimi.it/xoshiro256starstar.c
   8
   9 CONSTANT: JUMP-256 {
  10     0x180ec6d33cfd0aba
  11     0xd5a61266f0c9392c
  12     0xa9582618e03fc9aa
  13     0x39abdc4529b1661c
  14 }
  15
  16 CONSTANT: LONG-JUMP-256 {
  17     0x76e15d3efefdcbbf
  18     0xc5004e441c522fb3
  19     0x77710069854ee241
  20     0x39109bb02acbe635
  21 }
  22
  23 STRUCT: xoshiro-256-star-star { s0 ulonglong } { s1 ulonglong } { s2 ulonglong } { s3 ulonglong } ;
  24
  25 : <xoshiro-256-star-star> ( s0 s1 s2 s3 -- obj )
  26     xoshiro-256-star-star new
  27         swap >>s3
  28         swap >>s2
  29         swap >>s1
  30         swap >>s0 ; inline
  31
  32 : rotl-256 ( x: uint64_t k: int -- out: uint64_t )
  33     [ shift ]
  34     [ 64 swap - neg shift ] 2bi bitor 64 bits ; inline
  35
  36 :: (next-256) ( s0! s1! s2! s3! -- s0 s1 s2 s3 64-random-bits )
  37     s1 5 * 7 rotl-256 9 * 64 bits :> 64-random-bits
  38     s1 17 shift 64 bits :> t
  39     s0 s2 bitxor s2!
  40     s1 s3 bitxor s3!
  41     s2 s1 bitxor s1!
  42     s3 s0 bitxor s0!
  43     s2 t bitxor s2!
  44     s3 45 rotl-256 s3!
  45     s0 s1 s2 s3 64-random-bits ; inline
  46
  47 : next-256 ( xoshiro-256-star-star -- r64 )
  48     dup get[ s0 s1 s2 s3 ] (next-256)
  49     [ set[ s0 s1 s2 s3 ] drop ] dip ;
  50
  51 :: jump ( s0! s1! s2! s3! jump-table -- s0' s1' s2' s3' )
  52     0 0 0 0 :> ( t0! t1! t2! t3! )
  53     4 <iota> [
  54         64 <iota> [
  55         [ jump-table nth ] [ 1 swap shift ] bi* bitand 0 > [
  56             s0 t0 bitxor t0!
  57             s1 t1 bitxor t1!
  58             s2 t2 bitxor t2!
  59             s3 t3 bitxor t3!
  60         ] when
  61         s0 s1 s2 s3 (next-256) drop s3! s2! s1! s0!
  62         ] with each
  63     ] each
  64     t0 t1 t2 t3 ;
  65
  66 : jump-256 ( s0 s1 s2 s3 -- s0' s1' s2' s3' ) JUMP-256 jump ;
  67 : long-jump-256 ( s0 s1 s2 s3 -- s0' s1' s2' s3' ) LONG-JUMP-256 jump ;
  68
  69 M: xoshiro-256-star-star random-32*
  70     next-256 ;