diff options
Diffstat (limited to 'randinex.scm')
-rw-r--r-- | randinex.scm | 42 |
1 files changed, 23 insertions, 19 deletions
diff --git a/randinex.scm b/randinex.scm index 19b9d81..717b306 100644 --- a/randinex.scm +++ b/randinex.scm @@ -8,7 +8,7 @@ ;1. Any copy made of this software must include this copyright notice ;in full. ; -;2. I have made no warrantee or representation that the operation of +;2. I have made no warranty or representation that the operation of ;this software will be error-free, and I am under no obligation to ;provide any services, by way of maintenance, update, or otherwise. ; @@ -22,19 +22,23 @@ ;;; Sphere and normal functions corrections from: Harald Hanche-Olsen +(require 'random) +(require 'inexact) + +;;@code{(require 'random-inexact)} +;;@ftindex random-inexact + ;;; Generate an inexact real between 0 and 1. -(define random:uniform1 - ; how many chunks fill an inexact? +(define random:uniform1 ; how many chunks fill an inexact? (do ((random:chunks/float 0 (+ 1 random:chunks/float)) (smidgen 1.0 (/ smidgen 256.0))) - ((or (= 1.0 (+ 1 smidgen)) (= 4 random:chunks/float)) + ((or (= 1 (+ 1 smidgen)) (= 4 random:chunks/float)) (lambda (state) (do ((cnt random:chunks/float (+ -1 cnt)) (uni (/ (random:chunk state) 256.0) (/ (+ uni (random:chunk state)) 256.0))) ((= 1 cnt) uni)))))) - ;;@args ;;@args state ;;Returns an uniformly distributed inexact real random number in the @@ -70,24 +74,24 @@ ;;; 1-exp(-r^2/2). This latter means that u=exp(-r^2/2) is uniformly ;;; distributed on [0,1], so r=sqrt(-2 log u) can be used to generate r. -(define *2pi (* 8 (atan 1))) - ;;@args vect ;;@args vect state ;;Fills @1 with inexact real random numbers which are independent ;;and standard normally distributed (i.e., with mean 0 and variance 1). -(define (random:normal-vector! vect . args) - (let ((state (if (null? args) *random-state* (car args))) - (sum2 0)) - (let ((do! (lambda (k x) - (vector-set! vect k x) - (set! sum2 (+ sum2 (* x x)))))) - (do ((n (- (vector-length vect) 1) (- n 2))) - ((negative? n) sum2) - (let ((t (* *2pi (random:uniform1 state))) - (r (sqrt (* -2 (log (random:uniform1 state)))))) - (do! n (* r (cos t))) - (if (positive? n) (do! (- n 1) (* r (sin t))))))))) +(define random:normal-vector! + (let ((*2pi (* 8 (atan 1)))) + (lambda (vect . args) + (let ((state (if (null? args) *random-state* (car args))) + (sum2 0)) + (let ((do! (lambda (k x) + (vector-set! vect k x) + (set! sum2 (+ sum2 (* x x)))))) + (do ((n (- (vector-length vect) 1) (- n 2))) + ((negative? n) sum2) + (let ((t (* *2pi (random:uniform1 state))) + (r (sqrt (* -2 (log (random:uniform1 state)))))) + (do! n (* r (cos t))) + (if (positive? n) (do! (- n 1) (* r (sin t))))))))))) ;;; For the uniform distibution on the hollow sphere, pick a normal |