[functional tests] Remove dependency on the ThunderChez library.

I've just moved the relevant code into the functional-tests dir.
This commit is contained in:
Joe Thornber
2020-04-30 12:07:42 +01:00
parent ad79b627a4
commit 3e5de399a7
61 changed files with 16836 additions and 16 deletions

View File

@@ -0,0 +1,30 @@
;; Copyright (c) 2009 Derick Eddington. All rights reserved.
;; Licensed under an MIT-style license. My license is in the file
;; named LICENSE from the original collection this file is distributed
;; with. If this file is redistributed with some other collection, my
;; license must also be included.
#!r6rs
(library (srfi s27 random-bits)
(export random-integer
random-real
default-random-source
make-random-source
random-source?
random-source-state-ref
random-source-state-set!
random-source-randomize!
random-source-pseudo-randomize!
random-source-make-integers
random-source-make-reals)
(import (rnrs)
(rnrs r5rs)
(only (srfi s19 time) time-nanosecond current-time)
(srfi s23 error tricks)
(srfi private include)
)
(SRFI-23-error->R6RS "(library (srfi s27 random-bits))"
(include/resolve ("srfi" "s27") "random.ss"))
)

View File

@@ -0,0 +1,584 @@
;; R6RS port of the Scheme48 reference implementation of SRFI-27
; MODULE DEFINITION FOR SRFI-27
; =============================
;
; Sebastian.Egner@philips.com, Mar-2002, in Scheme 48 0.57
; 1. The core generator is implemented in 'mrg32k3a-a.scm'.
; 2. The generic parts of the interface are in 'mrg32k3a.scm'.
; 3. The non-generic parts (record type, time, error) are here.
; history of this file:
; SE, 22-Mar-2002: initial version
; SE, 27-Mar-2002: checked again
; JS, 06-Dec-2007: R6RS port
(define-record-type :random-source
(fields state-ref
state-set!
randomize!
pseudo-randomize!
make-integers
make-reals))
(define :random-source-make make-:random-source)
(define state-ref :random-source-state-ref)
(define state-set! :random-source-state-set!)
(define randomize! :random-source-randomize!)
(define pseudo-randomize! :random-source-pseudo-randomize!)
(define make-integers :random-source-make-integers)
(define make-reals :random-source-make-reals)
(define (:random-source-current-time)
(time-nanosecond (current-time)))
;;; mrg32k3a-a.ss
; 54-BIT INTEGER IMPLEMENTATION OF THE "MRG32K3A"-GENERATOR
; =========================================================
;
; Sebastian.Egner@philips.com, Mar-2002.
;
; This file is an implementation of Pierre L'Ecuyer's MRG32k3a
; pseudo random number generator. Please refer to 'mrg32k3a.scm'
; for more information.
;
; compliance:
; Scheme R5RS with integers covering at least {-2^53..2^53-1}.
;
; history of this file:
; SE, 18-Mar-2002: initial version
; SE, 22-Mar-2002: comments adjusted, range added
; SE, 25-Mar-2002: pack/unpack just return their argument
; the actual generator
(define (mrg32k3a-random-m1 state)
(let ((x11 (vector-ref state 0))
(x12 (vector-ref state 1))
(x13 (vector-ref state 2))
(x21 (vector-ref state 3))
(x22 (vector-ref state 4))
(x23 (vector-ref state 5)))
(let ((x10 (modulo (- (* 1403580 x12) (* 810728 x13)) 4294967087))
(x20 (modulo (- (* 527612 x21) (* 1370589 x23)) 4294944443)))
(vector-set! state 0 x10)
(vector-set! state 1 x11)
(vector-set! state 2 x12)
(vector-set! state 3 x20)
(vector-set! state 4 x21)
(vector-set! state 5 x22)
(modulo (- x10 x20) 4294967087))))
; interface to the generic parts of the generator
(define (mrg32k3a-pack-state unpacked-state)
unpacked-state)
(define (mrg32k3a-unpack-state state)
state)
(define (mrg32k3a-random-range) ; m1
4294967087)
(define (mrg32k3a-random-integer state range) ; rejection method
(let* ((q (quotient 4294967087 range))
(qn (* q range)))
(do ((x (mrg32k3a-random-m1 state) (mrg32k3a-random-m1 state)))
((< x qn) (quotient x q)))))
(define (mrg32k3a-random-real state) ; normalization is 1/(m1+1)
(* 0.0000000002328306549295728 (+ 1.0 (mrg32k3a-random-m1 state))))
;;; mrg32k3a.ss
; GENERIC PART OF MRG32k3a-GENERATOR FOR SRFI-27
; ==============================================
;
; Sebastian.Egner@philips.com, 2002.
;
; This is the generic R5RS-part of the implementation of the MRG32k3a
; generator to be used in SRFI-27. It is based on a separate implementation
; of the core generator (presumably in native code) and on code to
; provide essential functionality not available in R5RS (see below).
;
; compliance:
; Scheme R5RS with integer covering at least {-2^53..2^53-1}.
; In addition,
; SRFI-23: error
;
; history of this file:
; SE, 22-Mar-2002: refactored from earlier versions
; SE, 25-Mar-2002: pack/unpack need not allocate
; SE, 27-Mar-2002: changed interface to core generator
; SE, 10-Apr-2002: updated spec of mrg32k3a-random-integer
; Generator
; =========
;
; Pierre L'Ecuyer's MRG32k3a generator is a Combined Multiple Recursive
; Generator. It produces the sequence {(x[1,n] - x[2,n]) mod m1 : n}
; defined by the two recursive generators
;
; x[1,n] = ( a12 x[1,n-2] + a13 x[1,n-3]) mod m1,
; x[2,n] = (a21 x[2,n-1] + a23 x[2,n-3]) mod m2,
;
; where the constants are
; m1 = 4294967087 = 2^32 - 209 modulus of 1st component
; m2 = 4294944443 = 2^32 - 22853 modulus of 2nd component
; a12 = 1403580 recursion coefficients
; a13 = -810728
; a21 = 527612
; a23 = -1370589
;
; The generator passes all tests of G. Marsaglia's Diehard testsuite.
; Its period is (m1^3 - 1)(m2^3 - 1)/2 which is nearly 2^191.
; L'Ecuyer reports: "This generator is well-behaved in all dimensions
; up to at least 45: ..." [with respect to the spectral test, SE].
;
; The period is maximal for all values of the seed as long as the
; state of both recursive generators is not entirely zero.
;
; As the successor state is a linear combination of previous
; states, it is possible to advance the generator by more than one
; iteration by applying a linear transformation. The following
; publication provides detailed information on how to do that:
;
; [1] P. L'Ecuyer, R. Simard, E. J. Chen, W. D. Kelton:
; An Object-Oriented Random-Number Package With Many Long
; Streams and Substreams. 2001.
; To appear in Operations Research.
;
; Arithmetics
; ===========
;
; The MRG32k3a generator produces values in {0..2^32-209-1}. All
; subexpressions of the actual generator fit into {-2^53..2^53-1}.
; The code below assumes that Scheme's "integer" covers this range.
; In addition, it is assumed that floating point literals can be
; read and there is some arithmetics with inexact numbers.
;
; However, for advancing the state of the generator by more than
; one step at a time, the full range {0..2^32-209-1} is needed.
; Required: Backbone Generator
; ============================
;
; At this point in the code, the following procedures are assumed
; to be defined to execute the core generator:
;
; (mrg32k3a-pack-state unpacked-state) -> packed-state
; (mrg32k3a-unpack-state packed-state) -> unpacked-state
; pack/unpack a state of the generator. The core generator works
; on packed states, passed as an explicit argument, only. This
; allows native code implementations to store their state in a
; suitable form. Unpacked states are #(x10 x11 x12 x20 x21 x22)
; with integer x_ij. Pack/unpack need not allocate new objects
; in case packed and unpacked states are identical.
;
; (mrg32k3a-random-range) -> m-max
; (mrg32k3a-random-integer packed-state range) -> x in {0..range-1}
; advance the state of the generator and return the next random
; range-limited integer.
; Note that the state is not necessarily advanced by just one
; step because we use the rejection method to avoid any problems
; with distribution anomalies.
; The range argument must be an exact integer in {1..m-max}.
; It can be assumed that range is a fixnum if the Scheme system
; has such a number representation.
;
; (mrg32k3a-random-real packed-state) -> x in (0,1)
; advance the state of the generator and return the next random
; real number between zero and one (both excluded). The type of
; the result should be a flonum if possible.
; Required: Record Data Type
; ==========================
;
; At this point in the code, the following procedures are assumed
; to be defined to create and access a new record data type:
;
; (:random-source-make a0 a1 a2 a3 a4 a5) -> s
; constructs a new random source object s consisting of the
; objects a0 .. a5 in this order.
;
; (:random-source? obj) -> bool
; tests if a Scheme object is a :random-source.
;
; (:random-source-state-ref s) -> a0
; (:random-source-state-set! s) -> a1
; (:random-source-randomize! s) -> a2
; (:random-source-pseudo-randomize! s) -> a3
; (:random-source-make-integers s) -> a4
; (:random-source-make-reals s) -> a5
; retrieve the values in the fields of the object s.
; Required: Current Time as an Integer
; ====================================
;
; At this point in the code, the following procedure is assumed
; to be defined to obtain a value that is likely to be different
; for each invokation of the Scheme system:
;
; (:random-source-current-time) -> x
; an integer that depends on the system clock. It is desired
; that the integer changes as fast as possible.
; Accessing the State
; ===================
(define (mrg32k3a-state-ref packed-state)
(cons 'lecuyer-mrg32k3a
(vector->list (mrg32k3a-unpack-state packed-state))))
(define (mrg32k3a-state-set external-state)
(define (check-value x m)
(if (and (integer? x)
(exact? x)
(<= 0 x (- m 1)))
#t
(error "illegal value" x)))
(if (and (list? external-state)
(= (length external-state) 7)
(eq? (car external-state) 'lecuyer-mrg32k3a))
(let ((s (cdr external-state)))
(check-value (list-ref s 0) mrg32k3a-m1)
(check-value (list-ref s 1) mrg32k3a-m1)
(check-value (list-ref s 2) mrg32k3a-m1)
(check-value (list-ref s 3) mrg32k3a-m2)
(check-value (list-ref s 4) mrg32k3a-m2)
(check-value (list-ref s 5) mrg32k3a-m2)
(if (or (zero? (+ (list-ref s 0) (list-ref s 1) (list-ref s 2)))
(zero? (+ (list-ref s 3) (list-ref s 4) (list-ref s 5))))
(error "illegal degenerate state" external-state))
(mrg32k3a-pack-state (list->vector s)))
(error "malformed state" external-state)))
; Pseudo-Randomization
; ====================
;
; Reference [1] above shows how to obtain many long streams and
; substream from the backbone generator.
;
; The idea is that the generator is a linear operation on the state.
; Hence, we can express this operation as a 3x3-matrix acting on the
; three most recent states. Raising the matrix to the k-th power, we
; obtain the operation to advance the state by k steps at once. The
; virtual streams and substreams are now simply parts of the entire
; periodic sequence (which has period around 2^191).
;
; For the implementation it is necessary to compute with matrices in
; the ring (Z/(m1*m1)*Z)^(3x3). By the Chinese-Remainder Theorem, this
; is isomorphic to ((Z/m1*Z) x (Z/m2*Z))^(3x3). We represent such a pair
; of matrices
; [ [[x00 x01 x02],
; [x10 x11 x12],
; [x20 x21 x22]], mod m1
; [[y00 y01 y02],
; [y10 y11 y12],
; [y20 y21 y22]] mod m2]
; as a vector of length 18 of the integers as writen above:
; #(x00 x01 x02 x10 x11 x12 x20 x21 x22
; y00 y01 y02 y10 y11 y12 y20 y21 y22)
;
; As the implementation should only use the range {-2^53..2^53-1}, the
; fundamental operation (x*y) mod m, where x, y, m are nearly 2^32,
; is computed by breaking up x and y as x = x1*w + x0 and y = y1*w + y0
; where w = 2^16. In this case, all operations fit the range because
; w^2 mod m is a small number. If proper multiprecision integers are
; available this is not necessary, but pseudo-randomize! is an expected
; to be called only occasionally so we do not provide this implementation.
(define mrg32k3a-m1 4294967087) ; modulus of component 1
(define mrg32k3a-m2 4294944443) ; modulus of component 2
(define mrg32k3a-initial-state ; 0 3 6 9 12 15 of A^16, see below
'#( 1062452522
2961816100
342112271
2854655037
3321940838
3542344109))
(define mrg32k3a-generators #f) ; computed when needed
(define (mrg32k3a-pseudo-randomize-state i j)
(define (product A B) ; A*B in ((Z/m1*Z) x (Z/m2*Z))^(3x3)
(define w 65536) ; wordsize to split {0..2^32-1}
(define w-sqr1 209) ; w^2 mod m1
(define w-sqr2 22853) ; w^2 mod m2
(define (lc i0 i1 i2 j0 j1 j2 m w-sqr) ; linear combination
(let ((a0h (quotient (vector-ref A i0) w))
(a0l (modulo (vector-ref A i0) w))
(a1h (quotient (vector-ref A i1) w))
(a1l (modulo (vector-ref A i1) w))
(a2h (quotient (vector-ref A i2) w))
(a2l (modulo (vector-ref A i2) w))
(b0h (quotient (vector-ref B j0) w))
(b0l (modulo (vector-ref B j0) w))
(b1h (quotient (vector-ref B j1) w))
(b1l (modulo (vector-ref B j1) w))
(b2h (quotient (vector-ref B j2) w))
(b2l (modulo (vector-ref B j2) w)))
(modulo
(+ (* (+ (* a0h b0h)
(* a1h b1h)
(* a2h b2h))
w-sqr)
(* (+ (* a0h b0l)
(* a0l b0h)
(* a1h b1l)
(* a1l b1h)
(* a2h b2l)
(* a2l b2h))
w)
(* a0l b0l)
(* a1l b1l)
(* a2l b2l))
m)))
(vector
(lc 0 1 2 0 3 6 mrg32k3a-m1 w-sqr1) ; (A*B)_00 mod m1
(lc 0 1 2 1 4 7 mrg32k3a-m1 w-sqr1) ; (A*B)_01
(lc 0 1 2 2 5 8 mrg32k3a-m1 w-sqr1)
(lc 3 4 5 0 3 6 mrg32k3a-m1 w-sqr1) ; (A*B)_10
(lc 3 4 5 1 4 7 mrg32k3a-m1 w-sqr1)
(lc 3 4 5 2 5 8 mrg32k3a-m1 w-sqr1)
(lc 6 7 8 0 3 6 mrg32k3a-m1 w-sqr1)
(lc 6 7 8 1 4 7 mrg32k3a-m1 w-sqr1)
(lc 6 7 8 2 5 8 mrg32k3a-m1 w-sqr1)
(lc 9 10 11 9 12 15 mrg32k3a-m2 w-sqr2) ; (A*B)_00 mod m2
(lc 9 10 11 10 13 16 mrg32k3a-m2 w-sqr2)
(lc 9 10 11 11 14 17 mrg32k3a-m2 w-sqr2)
(lc 12 13 14 9 12 15 mrg32k3a-m2 w-sqr2)
(lc 12 13 14 10 13 16 mrg32k3a-m2 w-sqr2)
(lc 12 13 14 11 14 17 mrg32k3a-m2 w-sqr2)
(lc 15 16 17 9 12 15 mrg32k3a-m2 w-sqr2)
(lc 15 16 17 10 13 16 mrg32k3a-m2 w-sqr2)
(lc 15 16 17 11 14 17 mrg32k3a-m2 w-sqr2)))
(define (power A e) ; A^e
(cond
((zero? e)
'#(1 0 0 0 1 0 0 0 1 1 0 0 0 1 0 0 0 1))
((= e 1)
A)
((even? e)
(power (product A A) (quotient e 2)))
(else
(product (power A (- e 1)) A))))
(define (power-power A b) ; A^(2^b)
(if (zero? b)
A
(power-power (product A A) (- b 1))))
(define A ; the MRG32k3a recursion
'#( 0 1403580 4294156359
1 0 0
0 1 0
527612 0 4293573854
1 0 0
0 1 0))
; check arguments
(if (not (and (integer? i)
(exact? i)
(integer? j)
(exact? j)))
(error "i j must be exact integer" i j))
; precompute A^(2^127) and A^(2^76) only once
(if (not mrg32k3a-generators)
(set! mrg32k3a-generators
(list (power-power A 127)
(power-power A 76)
(power A 16))))
; compute M = A^(16 + i*2^127 + j*2^76)
(let ((M (product
(list-ref mrg32k3a-generators 2)
(product
(power (list-ref mrg32k3a-generators 0)
(modulo i (expt 2 28)))
(power (list-ref mrg32k3a-generators 1)
(modulo j (expt 2 28)))))))
(mrg32k3a-pack-state
(vector
(vector-ref M 0)
(vector-ref M 3)
(vector-ref M 6)
(vector-ref M 9)
(vector-ref M 12)
(vector-ref M 15)))))
; True Randomization
; ==================
;
; The value obtained from the system time is feed into a very
; simple pseudo random number generator. This in turn is used
; to obtain numbers to randomize the state of the MRG32k3a
; generator, avoiding period degeneration.
(define (mrg32k3a-randomize-state state)
;; G. Marsaglia's simple 16-bit generator with carry
(let* ((m 65536)
(x (modulo (:random-source-current-time) m)))
(define (random-m)
(let ((y (modulo x m)))
(set! x (+ (* 30903 y) (quotient x m)))
y))
(define (random n) ; m < n < m^2
(modulo (+ (* (random-m) m) (random-m)) n))
; modify the state
(let ((m1 mrg32k3a-m1)
(m2 mrg32k3a-m2)
(s (mrg32k3a-unpack-state state)))
(mrg32k3a-pack-state
(vector
(+ 1 (modulo (+ (vector-ref s 0) (random (- m1 1))) (- m1 1)))
(modulo (+ (vector-ref s 1) (random m1)) m1)
(modulo (+ (vector-ref s 2) (random m1)) m1)
(+ 1 (modulo (+ (vector-ref s 3) (random (- m2 1))) (- m2 1)))
(modulo (+ (vector-ref s 4) (random m2)) m2)
(modulo (+ (vector-ref s 5) (random m2)) m2))))))
; Large Integers
; ==============
;
; To produce large integer random deviates, for n > m-max, we first
; construct large random numbers in the range {0..m-max^k-1} for some
; k such that m-max^k >= n and then use the rejection method to choose
; uniformly from the range {0..n-1}.
(define mrg32k3a-m-max
(mrg32k3a-random-range))
(define (mrg32k3a-random-power state k) ; n = m-max^k, k >= 1
(if (= k 1)
(mrg32k3a-random-integer state mrg32k3a-m-max)
(+ (* (mrg32k3a-random-power state (- k 1)) mrg32k3a-m-max)
(mrg32k3a-random-integer state mrg32k3a-m-max))))
(define (mrg32k3a-random-large state n) ; n > m-max
(do ((k 2 (+ k 1))
(mk (* mrg32k3a-m-max mrg32k3a-m-max) (* mk mrg32k3a-m-max)))
((>= mk n)
(let* ((mk-by-n (quotient mk n))
(a (* mk-by-n n)))
(do ((x (mrg32k3a-random-power state k)
(mrg32k3a-random-power state k)))
((< x a) (quotient x mk-by-n)))))))
; Multiple Precision Reals
; ========================
;
; To produce multiple precision reals we produce a large integer value
; and convert it into a real value. This value is then normalized.
; The precision goal is unit <= 1/(m^k + 1), or 1/unit - 1 <= m^k.
; If you know more about the floating point number types of the
; Scheme system, this can be improved.
(define (mrg32k3a-random-real-mp state unit)
(do ((k 1 (+ k 1))
(u (- (/ 1 unit) 1) (/ u mrg32k3a-m1)))
((<= u 1)
(/ (exact->inexact (+ (mrg32k3a-random-power state k) 1))
(exact->inexact (+ (expt mrg32k3a-m-max k) 1))))))
; Provide the Interface as Specified in the SRFI
; ==============================================
;
; An object of type random-source is a record containing the procedures
; as components. The actual state of the generator is stored in the
; binding-time environment of make-random-source.
(define (make-random-source)
(let ((state (mrg32k3a-pack-state ; make a new copy
(list->vector (vector->list mrg32k3a-initial-state)))))
(:random-source-make
(lambda ()
(mrg32k3a-state-ref state))
(lambda (new-state)
(set! state (mrg32k3a-state-set new-state)))
(lambda ()
(set! state (mrg32k3a-randomize-state state)))
(lambda (i j)
(set! state (mrg32k3a-pseudo-randomize-state i j)))
(lambda ()
(lambda (n)
(cond
((not (and (integer? n) (exact? n) (positive? n)))
(error "range must be exact positive integer" n))
((<= n mrg32k3a-m-max)
(mrg32k3a-random-integer state n))
(else
(mrg32k3a-random-large state n)))))
(lambda args
(cond
((null? args)
(lambda ()
(mrg32k3a-random-real state)))
((null? (cdr args))
(let ((unit (car args)))
(cond
((not (and (real? unit) (< 0 unit 1)))
(error "unit must be real in (0,1)" unit))
((<= (- (/ 1 unit) 1) mrg32k3a-m1)
(lambda ()
(mrg32k3a-random-real state)))
(else
(lambda ()
(mrg32k3a-random-real-mp state unit))))))
(else
(error "illegal arguments" args)))))))
(define random-source?
:random-source?)
(define (random-source-state-ref s)
((:random-source-state-ref s)))
(define (random-source-state-set! s state)
((:random-source-state-set! s) state))
(define (random-source-randomize! s)
((:random-source-randomize! s)))
(define (random-source-pseudo-randomize! s i j)
((:random-source-pseudo-randomize! s) i j))
; ---
(define (random-source-make-integers s)
((:random-source-make-integers s)))
(define (random-source-make-reals s . unit)
(apply (:random-source-make-reals s) unit))
; ---
(define default-random-source
(make-random-source))
(define random-integer
(random-source-make-integers default-random-source))
(define random-real
(random-source-make-reals default-random-source))

View File

@@ -0,0 +1,67 @@
REFERENCE IMPLEMENTATIONS FOR SRFI-27 "Sources of Random Bits"
==============================================================
Sebastian.Egner@philips.com, 10-Apr-2002.
Files
-----
readme - this file
mrg32k3a.scm - generic parts of P. L' Ecuyer's MRG32k3a PRGN
mrg32k3a-a.scm - core generator in Scheme integers
mrg32k3a-b.c - core generator in C doubles for Scheme 48
mrg32k3a-c.scm - core generator in Gambit [Scheme] flonums
srfi-27-a.scm - Scheme 48 package definition for Scheme-only impl.
srfi-27-b.scm - Scheme 48 package definition for C/Scheme impl.
srfi-27-c.scm - Gambit definition for Scheme-only impl.
conftest.scm - confidence tests for the implementation
Implementations
---------------
The implementation has been factored into three parts.
One part implements the core generator, one part provides
the more generic functionality as specified in the SRFI,
and one part combines the parts and provides the interface
as specified in the SRFI.
a) A Scheme-only implementation for Scheme 48 0.57:
srfi-27-a.scm
mrg32k3a-a.scm
mrg32k3a.scm
This implementation uses 54-bit Scheme integers for all
arithmetics of the generator. The result are Scheme integers
and inexact Scheme numbers when floating point values are
requested.
The implementation is slow but tries to stay away from
unportable features as much as possible.
b) An implementation in Scheme 48 0.57 and ANSI-C:
srfi-27-b.scm
mrg32k3a-b.scm
mrg32k3a.scm
This is a more realistic implementation using C's (double)
datatype for the core generator and 54-bit Scheme integers
for the more infrequent operations on the state like the
random-source-pseudo-randomize! operation.
This implementation is meant as an example for a realistic
native code implementation of the SRFI. Performance is good.
c) A Scheme-only implementation for Gambit 3.0:
srfi-27-c.scm
mrg32k3a-c.scm
mrg32k3a.scm
This implementation uses Gambit's 64-bit flonums. It is
entirely written in Scheme but uses a few special features
of the Gambit system to tell the compiler.
This implementation is meant as an example for a realistic
Scheme implementation using flonums in Scheme and no C-code.
Performance is good when the code is used in compiled form;
the implementation has been optimized by Brad Lucier. This
has resulted in a subtantial performance gain.