Code for normal random variables Bradley Lucier 05 May 2020 00:49 UTC

The implementation's fundamental code for normal random variables is

       (lambda ()
         ;;Box-Muller
         (let ((r (sqrt (* -2 (log (rand-real-proc)))))
               (theta (* 2 PI (rand-real-proc))))
           (+ mean (* deviation r (sin theta)))))))))

Box-Muller can generate variables in pairs, I tend to write the routine
(in Gambit) as
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
(declare (flonum))

(define twopi (* 8. (atan 1.)))

(define gaussian-random-number
   (let ((state #f))
     (lambda ()
       (if state
	  (let ((result state))
	    (set! state #f)
	    result)
	  (let ((x1 (sqrt (- (* 2. (log (random-real))))))
		(x2 (random-real)))
	    (set! state (* x1 (sin (* twopi x2))))
	    (* x1 (cos (* twopi x2))))))))

(declare (generic))
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
so as to return one of the pair and save the other for later.

Brad