/[cmucl]/src/code/rand.lisp
ViewVC logotype

Contents of /src/code/rand.lisp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.6 - (show annotations)
Thu Jun 24 12:57:50 1993 UTC (20 years, 10 months ago) by ram
Branch: MAIN
Changes since 1.5: +2 -2 lines
format-string -> format-control
1 ;;; -*- Mode: Lisp; Package: Kernel -*-
2 ;;;
3 ;;; **********************************************************************
4 ;;; This code was written as part of the CMU Common Lisp project at
5 ;;; Carnegie Mellon University, and has been placed in the public domain.
6 ;;; If you want to use this code or any part of CMU Common Lisp, please contact
7 ;;; Scott Fahlman or slisp-group@cs.cmu.edu.
8 ;;;
9 (ext:file-comment
10 "$Header: /tiger/var/lib/cvsroots/cmucl/src/code/rand.lisp,v 1.6 1993/06/24 12:57:50 ram Exp $")
11 ;;;
12 ;;; **********************************************************************
13 ;;;
14 ;;; Functions to random number functions for Spice Lisp
15 ;;;
16 ;;; Originally written by David Adam. Python tuning, better large integer
17 ;;; randomness and efficient IEEE float support by Rob MacLachlan.
18 ;;;
19 (in-package "LISP")
20 (export '(random-state random-state-p random *random-state*
21 make-random-state))
22
23 (in-package "KERNEL")
24 (export '(%random-single-float %random-double-float random-chunk
25 random-fixnum-max))
26
27
28 ;;;; Random state hackery:
29
30 (defconstant random-const-a 8373)
31 (defconstant random-const-c 101010101)
32 (defconstant random-max 54)
33
34 ;;; Inclusive upper bound on the size of fixnum kept in the state (and returned
35 ;;; by random-chunk.) Must be even.
36 ;;;
37 (defconstant random-upper-bound (- most-positive-fixnum 3))
38 (defconstant random-chunk-length (integer-length random-upper-bound))
39 (deftype random-chunk () `(integer 0 ,random-upper-bound))
40
41 (defvar rand-seed 0)
42 (defstruct (random-state
43 (:constructor make-random-object)
44 (:make-load-form-fun :just-dump-it-normally))
45 (j 24 :type index)
46 (k 0 :type index)
47 (seed (make-array (1+ random-max) :initial-contents
48 (do ((list-rands () (cons (rand1) list-rands))
49 (i 0 (1+ i)))
50 ((> i random-max) list-rands)
51 (declare (fixnum i))))
52 :type simple-vector))
53
54
55 ;;; Generates a random number from rand-seed.
56 (defun rand1 ()
57 (declare (optimize (inhibit-warnings 3)))
58 (setq rand-seed
59 (mod (+ (* rand-seed random-const-a) random-const-c)
60 (1+ random-upper-bound))))
61
62
63 (defvar *random-state* (make-random-object))
64
65
66 (defun copy-state (cur-state)
67 (let ((state (make-random-object
68 :seed (make-array 55)
69 :j (random-state-j cur-state)
70 :k (random-state-k cur-state))))
71 (do ((i 0 (1+ i)))
72 ((= i 55) state)
73 (declare (fixnum i))
74 (setf (aref (random-state-seed state) i)
75 (aref (random-state-seed cur-state) i)))))
76
77 (defun make-random-state (&optional state)
78 "Make a random state object. If State is not supplied, return a copy
79 of the default random state. If State is a random state, then return a
80 copy of it. If state is T then return a random state generated from
81 the universal time."
82 (cond ((not state) (copy-state *random-state*))
83 ((random-state-p state) (copy-state state))
84 ((eq state t) (setq rand-seed (get-universal-time))
85 (make-random-object))
86 (t (error "Argument is not a RANDOM-STATE, T or NIL: ~S" state))))
87
88
89 ;;;; Random entries:
90
91 (declaim (start-block random %random-single-float %random-double-float
92 random-chunk))
93
94 ;;; random-chunk -- Internal
95 ;;;
96 ;;; This function generates fixnums between 0 and random-upper-bound,
97 ;;; inclusive. For the algorithm to work random-upper-bound must be an
98 ;;; even positive fixnum. State is the random state to use.
99 ;;;
100 (declaim (ftype (function (random-state) random-chunk) random-chunk))
101 (defun random-chunk (state)
102 (let* ((seed (random-state-seed state))
103 (j (random-state-j state))
104 (k (random-state-k state))
105 (a (- (- random-upper-bound
106 (the random-chunk
107 (svref seed
108 (setf (random-state-j state)
109 (if (= j 0) random-max (1- j))))))
110 (the random-chunk
111 (svref seed
112 (setf (random-state-k state)
113 (if (= k 0) random-max (1- k))))))))
114 (declare (fixnum a))
115 (setf (svref seed k)
116 (the random-chunk (if (minusp a) (- a) (- random-upper-bound a))))))
117
118
119 ;;; %RANDOM-SINGLE-FLOAT, %RANDOM-DOUBLE-FLOAT -- Interface
120 ;;;
121 ;;; Handle the single or double float case of RANDOM. We generate a float
122 ;;; between 0.0 and 1.0 by clobbering the significand of 1.0 with random bits,
123 ;;; then subtracting 1.0. This hides the fact that we have a hidden bit.
124 ;;;
125 (declaim (inline %random-single-float %random-double-float))
126 (defun %random-single-float (arg state)
127 (declare (type (single-float (0f0)) arg)
128 (type (or random-state null) state))
129 (* arg
130 (- (make-single-float
131 (dpb (ash (random-chunk (or state *random-state*))
132 (- vm:single-float-digits random-chunk-length))
133 vm:single-float-significand-byte
134 (single-float-bits 1.0)))
135 1.0)))
136 ;;;
137 (defun %random-double-float (arg state)
138 (declare (type (double-float (0d0)) arg)
139 (type (or random-state null) state))
140 (let ((state (or state *random-state*)))
141 (* arg
142 (- (make-double-float
143 (dpb (ash (random-chunk state)
144 (- vm:double-float-digits random-chunk-length
145 vm:word-bits))
146 vm:double-float-significand-byte
147 (double-float-high-bits 1d0))
148 (logxor (ash (random-chunk state)
149 (- vm:word-bits random-chunk-length))
150 (ash (random-chunk state)
151 (- random-chunk-length vm:word-bits))))
152 1d0))))
153
154
155 ;;;; Random integers:
156
157 ;;; Amount we overlap chunks by when building a large integer to make up for
158 ;;; the loss of randomness in the low bits.
159 ;;;
160 (defconstant random-integer-overlap 3)
161
162 ;;; Extra bits of randomness that we generate before taking the value MOD the
163 ;;; limit, to avoid loss of randomness near the limit.
164 ;;;
165 (defconstant random-integer-extra-bits 10)
166
167 ;;; Largest fixnum we can compute from one chunk of bits.
168 ;;;
169 (defconstant random-fixnum-max
170 (1- (ash 1 (- random-chunk-length random-integer-extra-bits))))
171
172
173 ;;; %RANDOM-INTEGER -- Internal
174 ;;;
175 (defun %random-integer (arg state)
176 (declare (type (integer 1) arg) (type random-state state))
177 (let ((shift (- random-chunk-length random-integer-overlap)))
178 (do ((bits (random-chunk state)
179 (logxor (ash bits shift) (random-chunk state)))
180 (count (+ (integer-length arg)
181 (- random-integer-extra-bits shift))
182 (- count shift)))
183 ((minusp count)
184 (rem bits arg))
185 (declare (fixnum count)))))
186
187 (defun random (arg &optional (state *random-state*))
188 "Generate a uniformly distributed pseudo-random number between zero
189 and Arg. State, if supplied, is the random state to use."
190 (declare (inline %random-single-float %random-double-float))
191 (cond
192 ((and (fixnump arg) (<= arg random-fixnum-max))
193 (rem (random-chunk state) arg))
194 ((typep arg 'single-float)
195 (%random-single-float arg state))
196 ((typep arg 'double-float)
197 (%random-double-float arg state))
198 ((integerp arg)
199 (%random-integer arg state))
200 (t
201 (error 'simple-type-error :expected-type '(real (0)) :datum arg
202 :format-control "Argument is not a positive real number: ~S"
203 :format-arguments (list arg)))))

  ViewVC Help
Powered by ViewVC 1.1.5