source: project/release/4/srfi-27/trunk/well512.scm @ 34855

Last change on this file since 34855 was 34855, checked in by Kon Lovett, 4 years ago

st trk

File size: 5.2 KB
Line 
1;;;; well512.scm
2;;;; Kon Lovett, Nov '17
3
4(module well512
5
6(;export
7  make-random-source-well512)
8
9(import (except scheme <= inexact->exact exact->inexact number?))
10
11(import chicken foreign)
12
13(import
14  srfi-4
15  (only numbers <= inexact->exact exact->inexact number?)
16  random-source
17  entropy-source
18  (only srfi-27-numbers
19    check-positive-integer
20    random-large-integer random-large-real
21    native-real-precision?))
22(require-library
23  srfi-4
24  numbers
25  random-source entropy-source
26  srfi-27-numbers)
27
28(declare
29  (not usual-integrations
30    <= exact->inexact inexact->exact))
31
32#>
33#include "WELL512a.c"
34
35// uint init[ 16 ]
36static void
37init_state() {
38  InitWELLRNG512a( unsigned int *init );
39}
40
41static double
42uniformf64( uint32_t *state )
43{
44  return WELLRNG512a();
45}
46
47static uint32_t
48randomu32( uint32_t *state, uint32_t m )
49{
50  return WELLRNG512a (void) ;
51}
52
53static void
54uniformu64_ith_state( uint32_t *state, uint32_t i )
55{
56  for( ;i > 0; --i )
57    uniformu32( state );
58}
59
60static void
61uniformu64_jth_offset_state( uint32_t *state, uint32_t j )
62{
63  uint64_t x = (uint64_t) pow( A, j ) * state[ W ] + state[ C ];
64  state[ W ] = low32( x );
65  state[ C ] = high32( x );
66}
67
68<#
69
70(define init_state (foreign-lambda void "init_state" nonnull-u32vector unsigned-integer64))
71
72(define well512-random-integer (foreign-lambda unsigned-integer32 "randomu32" nonnull-u32vector unsigned-integer32))
73(define well512-random-real (foreign-lambda double "uniformf64" nonnull-u32vector))
74
75;;;
76
77;;
78
79(define-constant maximum-unsigned-integer32-flonum 4294967295.0)
80
81(cond-expand
82  (64bit
83    (define-constant maximum-unsigned-integer32 4294967295) )    ;MAX
84  (else ;32bit
85    (define-constant maximum-unsigned-integer32 1073741823) ) )  ;32bit most-positive-fixnum
86
87(define-constant fpMAX maximum-unsigned-integer32-flonum)  ;2^32 - 1
88(define-constant LOG2-PERIOD 250)
89
90(define eMAX (inexact->exact fpMAX)) ;Create a "bignum" if necessary
91
92(define-constant INITIAL-SEED maximum-unsigned-integer32-flonum)
93
94(define-constant STATE-LENGTH 10)
95
96(define INTERNAL-ID 'well512)
97(define EXTERNAL-ID 'well512)
98
99;;
100
101(define (make-state) (make-u32vector STATE-LENGTH))
102
103(define (well512-initial-state)
104  (let ((state (make-state)))
105    (init_state state INITIAL-SEED)
106    state ) )
107
108(define (well512-unpack-state state)
109  (cons EXTERNAL-ID (u32vector->list state)) )
110
111(define (well512-pack-state external-state)
112  (unless (well512-external-state? external-state)
113    (error 'well512-pack-state "malformed state" external-state) )
114  (let ((state (make-state)))
115    (do ((i 0 (fx+ i 1))
116         (ss (cdr external-state) (cdr ss)) )
117        ((null? ss) state)
118      (let ((x (car ss)))
119        (if (and (integer? x) (<= 0 x 4294967295))
120          (u32vector-set! state i x)
121          (error 'well512-pack-state "illegal value" x) ) ) ) ) )
122
123(define (well512-external-state? obj)
124  (and
125    (pair? obj)
126    (eq? EXTERNAL-ID (car obj))
127    (fx= (fx+ STATE-LENGTH 1) (length obj)) ) )
128
129(define (well512-randomize-state state entropy-source)
130  (init_state
131    state
132    (exact->inexact
133      (modulo
134        (fpabs (entropy-source-f64-integer entropy-source))
135        fpMAX)))
136  state )
137
138(define (well512-pseudo-randomize-state i j)
139  (let ((state (make-state)))
140    (init_state state (exact->inexact (+ i j)))
141    state ) )
142
143(define (well512-random-large state n)
144  (random-large-integer well512-random-integer state fpMAX eMAX n) )
145
146(define (well512-random-real-mp state prec)
147  (random-large-real well512-random-integer state fpMAX eMAX prec) )
148
149;;;
150
151(define (make-random-source-well512)
152  (let ((state (well512-initial-state)))
153    (*make-random-source
154      ;
155      make-random-source-well512
156      ;
157      EXTERNAL-ID
158      ;
159      "Well's 512-bit Generator"
160      ;
161      LOG2-PERIOD
162      ;
163      fpMAX
164      ;
165      #f
166      ;
167      (lambda ()
168        (well512-unpack-state state) )
169      ;
170      (lambda (new-state)
171        (set! state (well512-pack-state new-state)) )
172      ;
173      (lambda (entropy-source)
174        (set! state (well512-randomize-state state entropy-source)) )
175      ;
176      (lambda (i j)
177        (set! state (well512-pseudo-randomize-state i j)) )
178      ;
179      (lambda ()
180        (lambda (n)
181          (check-positive-integer INTERNAL-ID n 'range)
182          (cond-expand
183            (64bit
184              (cond
185                ((and (fixnum? n) (<= n maximum-unsigned-integer32))
186                  (well512-random-integer state n))
187                (else
188                  (well512-random-large state n) ) ) )
189            (else ;32bit
190              (cond
191                ((and (fixnum? n) (<= n maximum-unsigned-integer32))
192                  (well512-random-integer state n))
193                ;'n' maybe bignum - must be convertable to "unsigned-integer32"
194                ((<= n eMAX)
195                  (well512-random-integer state (exact->inexact n)))
196                (else
197                  (well512-random-large state n) ) ) ) ) ) )
198      ;
199      (lambda (prec)
200        (cond
201          ((native-real-precision? prec eMAX)
202            (lambda ()
203              (well512-random-real state) ) )
204          (else
205            (lambda ()
206              (well512-random-real-mp state prec) ) ) ) ) ) ) )
207
208;;;
209;;; Module Init
210;;;
211
212(register-random-source! INTERNAL-ID make-random-source-well512)
213
214) ;module well512
Note: See TracBrowser for help on using the repository browser.