Better equidistributed and faster/less consing integer RANDOM.
authorLutz Euler <lutz.euler@freenet.de>
Tue, 1 May 2012 16:59:12 +0000 (18:59 +0200)
committerLutz Euler <lutz.euler@freenet.de>
Tue, 1 May 2012 16:59:12 +0000 (18:59 +0200)
Up to now the implementation of RANDOM with an integer argument just
generated a few more random bits than the length of the argument and
took this value MOD the argument. This led to a slightly uneven
distribution of the possible values unless the argument was a power of
two. Moreover, for bignums, the algorithm was quadratic both in time and
space dependent on the number of bits of the argument.

Instead generate random integers using an accept-reject loop and change
the bignum implementation to an algorithm that is linear in time and
space.

I took some inspiration from WHN's attempt at an accept-reject loop
implementation in commit 0a7604d54581d2c846838c26ce6a7993629586fa and
following.

Thanks to Christophe Rhodes for reviewing this patch!

Some details:

The implementation works correctly with both a random chunk size equal
to the word size and equal to half the word size. This is currently
necessary as a 32-bit pseudo random generator is used both under 32 and
under 64 bit word size.

In the generic RANDOM case, fixnum and bignum limits are differentiated:

With a fixnum limit an accept-reject loop on a masked random chunk is
always used. Under 64 bit word size two random chunks are used only if
the limit is so large that one doesn't suffice. This never conses.

With a bignum limit four cases are differentiated to minimize consing.
If just one random chunk is needed to supply a sufficient number of bits
the implementation only conses if the result is indeed a bignum:
* If the limit is a power of two, a chunk is generated and shifted to
  get the correct number of bits.
* If the limit is not a power of two an accept-reject loop with shifting
  is used.
If more than one random chunk is needed, a bignum is always consed even
if it happens to normalize to a fixnum:
* If the limit is a power of two a straightforward algorithm is used to
  fill a newly created bignum with random bits.
* If the limit is not a power of two an accept-reject loop is used that
  detects rejection early by starting from the most significant bits,
  thus generating on the average only one random chunk more than needed
  to fill the result once.
The test for power of two is non-consing, too.

In the case of a compile-time constant integer argument (of at most word
size) a DEFTRANSFORM triggers, that, in the general case, compiles an
accept-reject loop. For values of the limit where this sufficiently
reduces the rejection probability the largest multiple of the limit
fitting in one or two random chunks is used instead inside the loop.
To bring the result in the correct range a division is then necessary
(which the compiler converts into a multiplication). Powers of two are
optimized by leaving out the rejection test. In those cases where a word
has more bits than a random chunk, the generated expression uses two
chunks only if necessary.

NEWS
build-order.lisp-expr
package-data-list.lisp-expr
src/code/bignum-random.lisp [new file with mode: 0644]
src/code/random.lisp
src/code/target-random.lisp
src/compiler/float-tran.lisp
tests/random.pure.lisp

diff --git a/NEWS b/NEWS
index 86d8250..5b4ac62 100644 (file)
--- a/NEWS
+++ b/NEWS
@@ -1,5 +1,18 @@
 ;;;; -*- coding: utf-8; fill-column: 78 -*-
 changes relative to sbcl-1.0.56:
+  * RANDOM enhancements and bug fixes:
+    ** bug fix: the range and distribution of random integers could be
+       catastrophically wrong when the compiler derived the type of its
+       argument as a disjoint set of small integers.
+    ** bug fix: the distribution of random integers is now completely
+       uniform even when the specified limit is not a power of two.
+       (Previously some values could be about 0.1 % more probable than
+       others in the worst case.)
+    ** RANDOM on large integer arguments is generally faster and conses
+       less than before; this is visible for fixnums above a length of
+       about 24 bits, but extremely so for bignums: the old implementation
+       used time and space quadratical in the size of the argument there,
+       the new one is linear.
   * enhancement: redesigned protocol for quitting SBCL. SB-EXT:EXIT is the new
     main entry point, SB-EXT:QUIT is deprecated.
   * enhancement: additions to the SB-THREAD API: RETURN-FROM-THREAD,
index e2ec49c..c31d518 100644 (file)
  ("src/code/target-sap"        :not-host) ; uses SAP-INT type
  ("src/code/target-package"    :not-host) ; needs "code/package"
  ("src/code/target-random"     :not-host) ; needs "code/random"
+ ("src/code/bignum-random"     :not-host) ; needs "code/random" and
+                                          ;   "code/bignum"
  ("src/code/target-hash-table" :not-host) ; needs "code/hash-table"
  ("src/code/reader"            :not-host) ; needs "code/readtable"
  ("src/code/target-stream"     :not-host) ; needs WHITESPACEP from "code/reader"
index cbc0b4a..13ede30 100644 (file)
@@ -180,7 +180,7 @@ of SBCL which maintained the CMU-CL-style split into two packages.)"
                "FLOAT-BIGNUM-RATIO" "MAKE-SMALL-BIGNUM"
                "MULTIPLY-BIGNUM-AND-FIXNUM" "MULTIPLY-BIGNUMS"
                "MULTIPLY-FIXNUMS" "NEGATE-BIGNUM"
-               "SUBTRACT-BIGNUM" "SXHASH-BIGNUM"))
+               "%RANDOM-BIGNUM" "SUBTRACT-BIGNUM" "SXHASH-BIGNUM"))
 
    #s(sb-cold:package-data
       :name "SB!C"
@@ -1865,6 +1865,7 @@ is a good idea, but see SB-SYS re. blurring of boundaries."
                #!+long-float "%RANDOM-LONG-FLOAT"
                "%RANDOM-SINGLE-FLOAT" "STATIC-CLASSOID"
                "%FUNCALLABLE-INSTANCE-INFO" "RANDOM-CHUNK" "BIG-RANDOM-CHUNK"
+               "N-RANDOM-CHUNK-BITS"
                "LAYOUT-CLOS-HASH-LIMIT"
                "BUILT-IN-CLASSOID-DIRECT-SUPERCLASSES"
                "BUILT-IN-CLASSOID-TRANSLATION" "RANDOM-LAYOUT-CLOS-HASH"
@@ -1873,7 +1874,7 @@ is a good idea, but see SB-SYS re. blurring of boundaries."
                "%SET-FUNCALLABLE-INSTANCE-LAYOUT"
                "BASIC-STRUCTURE-CLASSOID"
                "REGISTER-LAYOUT"
-               "FUNCALLABLE-INSTANCE" "RANDOM-FIXNUM-MAX"
+               "FUNCALLABLE-INSTANCE"
                "MAKE-STATIC-CLASSOID"
                "%MAKE-SYMBOL"
                "%FUNCALLABLE-INSTANCE-FUNCTION" "SYMBOL-HASH"
diff --git a/src/code/bignum-random.lisp b/src/code/bignum-random.lisp
new file mode 100644 (file)
index 0000000..821d4a7
--- /dev/null
@@ -0,0 +1,183 @@
+;;;; generation of random bignums
+;;;;
+;;;; The implementation assumes that the random chunk size is either
+;;;; equal to the word size or equal to half the word size.
+
+;;;; This software is part of the SBCL system. See the README file for
+;;;; more information.
+;;;;
+;;;; This software is derived from the CMU CL system, which was
+;;;; written at Carnegie Mellon University and released into the
+;;;; public domain. The software is in the public domain and is
+;;;; provided with absolutely no warranty. See the COPYING and CREDITS
+;;;; files for more information.
+
+(in-package "SB!BIGNUM")
+
+;;; Return T if the least significant N-BITS bits of BIGNUM are all
+;;; zero, else NIL. If the integer-length of BIGNUM is less than N-BITS,
+;;; the result is NIL, too.
+(declaim (inline bignum-lower-bits-zero-p))
+(defun bignum-lower-bits-zero-p (bignum n-bits)
+  (declare (type bignum-type bignum)
+           (type bit-index n-bits))
+  (multiple-value-bind (n-full-digits n-bits-partial-digit)
+      (floor n-bits digit-size)
+    (declare (type bignum-index n-full-digits))
+    (when (> (%bignum-length bignum) n-full-digits)
+      (dotimes (index n-full-digits)
+        (declare (type bignum-index index))
+        (unless (zerop (%bignum-ref bignum index))
+          (return-from bignum-lower-bits-zero-p nil)))
+      (zerop (logand (1- (ash 1 n-bits-partial-digit))
+                     (%bignum-ref bignum n-full-digits))))))
+
+;;; Return a nonnegative integer of DIGIT-SIZE many pseudo random bits.
+(declaim (inline random-bignum-digit))
+(defun random-bignum-digit (state)
+  (if (= n-random-chunk-bits digit-size)
+      (random-chunk state)
+      (big-random-chunk state)))
+
+;;; Return a nonnegative integer of N-BITS many pseudo random bits.
+;;; N-BITS must be nonnegative and less than DIGIT-SIZE.
+(declaim (inline random-bignum-partial-digit))
+(defun random-bignum-partial-digit (n-bits state)
+  (declare (type (integer 0 #.(1- digit-size)) n-bits)
+           (type random-state state))
+  (logand (1- (ash 1 n-bits))
+          (if (<= n-bits n-random-chunk-bits)
+              (random-chunk state)
+              (big-random-chunk state))))
+
+;;; Create a (nonnegative) bignum by concatenating RANDOM-CHUNK and
+;;; BIT-COUNT many pseudo random bits, normalise and return it.
+;;; RANDOM-CHUNK must fit into a bignum digit.
+(declaim (inline concatenate-random-bignum))
+(defun concatenate-random-bignum (random-chunk bit-count state)
+  (declare (type bignum-element-type random-chunk)
+           (type (integer 0 #.sb!xc:most-positive-fixnum) bit-count)
+           (type random-state state))
+  (let* ((n-total-bits (+ 1 n-random-chunk-bits bit-count)) ; sign bit
+         (length (ceiling n-total-bits digit-size))
+         (bignum (%allocate-bignum length)))
+    (multiple-value-bind (n-random-digits n-random-bits)
+        (floor bit-count digit-size)
+      (declare (type bignum-index n-random-digits))
+      (dotimes (index n-random-digits)
+        (setf (%bignum-ref bignum index)
+              (random-bignum-digit state)))
+      (if (zerop n-random-bits)
+          (setf (%bignum-ref bignum n-random-digits) random-chunk)
+          (progn
+            (setf (%bignum-ref bignum n-random-digits)
+                  (%logior (random-bignum-partial-digit n-random-bits
+                                                        state)
+                           (%ashl random-chunk n-random-bits)))
+            (let ((shift (- digit-size n-random-bits)))
+              (when (< shift n-random-chunk-bits)
+                (setf (%bignum-ref bignum (1+ n-random-digits))
+                      (%digit-logical-shift-right random-chunk shift)))))))
+    (%normalize-bignum bignum length)))
+
+;;; Create and return a (nonnegative) bignum of N-BITS many pseudo
+;;; random bits. The result is normalised, so may be a fixnum, too.
+(declaim (inline make-random-bignum))
+(defun make-random-bignum (n-bits state)
+  (declare (type (and fixnum (integer 0)) n-bits)
+           (type random-state state))
+  (let* ((n-total-bits (1+ n-bits)) ; sign bit
+         (length (ceiling n-total-bits digit-size))
+         (bignum (%allocate-bignum length)))
+    (declare (type bignum-index length))
+    (multiple-value-bind (n-digits n-bits-partial-digit)
+        (floor n-bits digit-size)
+      (declare (type bignum-index n-digits))
+      (dotimes (index n-digits)
+        (setf (%bignum-ref bignum index)
+              (random-bignum-digit state)))
+      (unless (zerop n-bits-partial-digit)
+        (setf (%bignum-ref bignum n-digits)
+              (random-bignum-partial-digit n-bits-partial-digit state))))
+    (%normalize-bignum bignum length)))
+
+;;; Create and return a pseudo random bignum less than ARG. The result
+;;; is normalised, so may be a fixnum, too. We try to keep the number of
+;;; times RANDOM-CHUNK is called and the amount of storage consed to a
+;;; minimum.
+;;; Four cases are differentiated:
+;;; * If ARG is a power of two and only one random chunk is needed to
+;;;   supply a sufficient number of bits, a chunk is generated and
+;;;   shifted to get the correct number of bits. This only conses if the
+;;;   result is indeed a bignum. This case can only occur if the size of
+;;;   the random chunks is equal to the word size.
+;;; * If ARG is a power of two and multiple chunks are needed, we call
+;;;   MAKE-RANDOM-BIGNUM. Here a bignum is always consed even if it
+;;;   happens to normalize to a fixnum, which can't be avoided.
+;;; * If ARG is not a power of two but one random chunk suffices an
+;;;   accept-reject loop is used. Each time through the loop a chunk is
+;;;   generated and shifted to get the correct number of bits. This only
+;;;   conses if the final accepted result is indeed a bignum. This case
+;;;   too can only occur if the size of the random chunks is equal to the
+;;;   word size.
+;;; * If ARG is not a power of two and multiple chunks are needed an
+;;;   accept-reject loop is used that detects rejection early by
+;;;   starting the generation with a random chunk aligned to the most
+;;;   significant bits of ARG. If the random value is larger than the
+;;;   corresponding chunk of ARG we don't need to generate the full
+;;;   amount of random bits but can retry immediately. If the random
+;;;   value is smaller than the ARG chunk we know already that the
+;;;   result will be accepted independently of what the remaining random
+;;;   bits will be, so we generate them and return. Only in the rare
+;;;   case that the random value and the ARG chunk are equal we need to
+;;;   generate and compare the complete random number and risk to reject
+;;;   it.
+(defun %random-bignum (arg state)
+  (declare (type (integer #.(1+ sb!xc:most-positive-fixnum)) arg)
+           (type random-state state)
+           (inline bignum-lower-bits-zero-p))
+  (let ((n-bits (bignum-integer-length arg)))
+    (declare (type (integer #.sb!vm:n-fixnum-bits) n-bits))
+    ;; Don't use (ZEROP (LOGAND ARG (1- ARG))) to test if ARG is a power
+    ;; of two as that would cons.
+    (cond ((bignum-lower-bits-zero-p arg (1- n-bits))
+           ;; ARG is a power of two. We need one bit less than its
+           ;; INTEGER-LENGTH. Not using (DECF N-BITS) here allows the
+           ;; compiler to make optimal use of the type declaration for
+           ;; N-BITS above.
+           (let ((n-bits (1- n-bits)))
+             (if (<= n-bits n-random-chunk-bits)
+                 (%digit-logical-shift-right (random-chunk state)
+                                             (- n-random-chunk-bits n-bits))
+                 (make-random-bignum n-bits state))))
+          ((<= n-bits n-random-chunk-bits)
+           (let ((shift (- n-random-chunk-bits n-bits))
+                 (arg (%bignum-ref arg 0)))
+             (loop
+               (let ((bits (%digit-logical-shift-right (random-chunk state)
+                                                       shift)))
+                 (when (< bits arg)
+                   (return bits))))))
+          (t
+           ;; ARG is not a power of two and we need more than one random
+           ;; chunk.
+           (let* ((shift (- n-bits n-random-chunk-bits))
+                  (arg-first-chunk (ldb (byte n-random-chunk-bits shift)
+                                        arg)))
+             (loop
+               (let ((random-chunk (random-chunk state)))
+                 ;; If the random value is larger than the corresponding
+                 ;; chunk from the most significant bits of ARG we can
+                 ;; retry immediately; no need to generate the remaining
+                 ;; random bits.
+                 (unless (> random-chunk arg-first-chunk)
+                   ;; We need to generate the complete random number.
+                   (let ((bits (concatenate-random-bignum random-chunk
+                                                          shift state)))
+                     ;; While the second comparison below subsumes the
+                     ;; first, the first is faster and will nearly
+                     ;; always be true, so it's worth it to try it
+                     ;; first.
+                     (when (or (< random-chunk arg-first-chunk)
+                               (< bits arg))
+                       (return bits)))))))))))
index 373a3b2..cfcc7fc 100644 (file)
 (in-package "SB!KERNEL")
 
 ;;; the size of the chunks returned by RANDOM-CHUNK
-(def!constant random-chunk-length 32)
-
-;;; the amount that we overlap chunks by when building a large integer
-;;; to make up for the loss of randomness in the low bits
-(def!constant random-integer-overlap 3)
-
-;;; extra bits of randomness that we generate before taking the value MOD the
-;;; limit, to avoid loss of randomness near the limit
-(def!constant random-integer-extra-bits 10)
-
-;;; the largest fixnum we can compute from one chunk of bits
-(def!constant random-fixnum-max
-  (1- (ash 1 (- random-chunk-length random-integer-extra-bits))))
+(def!constant n-random-chunk-bits 32)
 
 (sb!xc:defstruct (random-state (:constructor %make-random-state)
                                (:copier nil)) ; since shallow copy is wrong
index 57c57a7..b4cb278 100644 (file)
 #!-sb-fluid (declaim (inline big-random-chunk))
 (defun big-random-chunk (state)
   (declare (type random-state state))
-  (logand (1- (expt 2 64))
-          (logior (ash (random-chunk state) 32)
-                  (random-chunk state))))
+  (logior (ash (random-chunk state) 32)
+          (random-chunk state)))
 \f
 ;;; Handle the single or double float case of RANDOM. We generate a
 ;;; float between 0.0 and 1.0 by clobbering the significand of 1.0
   (* arg
      (- (make-single-float
          (dpb (ash (random-chunk state)
-                   (- sb!vm:single-float-digits random-chunk-length))
+                   (- sb!vm:single-float-digits n-random-chunk-bits))
               sb!vm:single-float-significand-byte
               (single-float-bits 1.0)))
         1.0)))
   (* arg
      (- (sb!impl::make-double-float
          (dpb (ash (random-chunk state)
-                   (- sb!vm:double-float-digits random-chunk-length 32))
+                   (- sb!vm:double-float-digits n-random-chunk-bits 32))
               sb!vm:double-float-significand-byte
               (sb!impl::double-float-high-bits 1d0))
          (random-chunk state))
     (* arg
        (- (sb!impl::make-double-float
            (dpb (ash (sb!vm::random-mt19937 state-vector)
-                     (- sb!vm:double-float-digits random-chunk-length
+                     (- sb!vm:double-float-digits n-random-chunk-bits
                         sb!vm:n-word-bits))
                 sb!vm:double-float-significand-byte
                 (sb!impl::double-float-high-bits 1d0))
           1d0))))
 
 \f
-;;;; random integers
+;;;; random fixnums
 
-(defun %random-integer (arg state)
-  (declare (type (integer 1) arg) (type random-state state))
-  (let ((shift (- random-chunk-length random-integer-overlap)))
-    (do ((bits (random-chunk state)
-               (logxor (ash bits shift) (random-chunk state)))
-         (count (+ (integer-length arg)
-                   (- random-integer-extra-bits shift))
-                (- count shift)))
-        ((minusp count)
-         (rem bits arg))
-      (declare (fixnum count)))))
+;;; Generate and return a pseudo random fixnum less than ARG. To achieve
+;;; equidistribution an accept-reject loop is used.
+;;; No extra effort is made to detect the case of ARG being a power of
+;;; two where rejection is not possible, as the cost of checking for
+;;; this case is the same as doing the rejection test. When ARG is
+;;; larger than (expt 2 N-RANDOM-CHUNK-BITS), which can only happen if
+;;; the random chunk size is half the word size, two random chunks are
+;;; used in each loop iteration, otherwise only one. Finally, the
+;;; rejection probability could often be reduced by not masking the
+;;; chunk but rejecting only values as least as large as the largest
+;;; multiple of ARG that fits in a chunk (or two), but this is not done
+;;; as the speed gains due to needing fewer loop iterations are by far
+;;; outweighted by the cost of the two divisions required (one to find
+;;; the multiplier and one to bring the result into the correct range).
+#!-sb-fluid (declaim (inline %random-fixnum))
+(defun %random-fixnum (arg state)
+  (declare (type (integer 1 #.sb!xc:most-positive-fixnum) arg)
+           (type random-state state))
+  (if (= arg 1)
+      0
+      (let* ((n-bits (integer-length (1- arg)))
+             (mask (1- (ash 1 n-bits))))
+        (macrolet ((accept-reject-loop (generator)
+                     `(loop
+                        (let ((bits (logand mask (,generator state))))
+                          (when (< bits arg)
+                            (return bits))))))
+          (assert (<= n-bits (* 2 n-random-chunk-bits)))
+          (if (<= n-bits n-random-chunk-bits)
+              (accept-reject-loop random-chunk)
+              (accept-reject-loop big-random-chunk))))))
 
 (defun random (arg &optional (state *random-state*))
-  (declare (inline %random-single-float %random-double-float
+  (declare (inline %random-fixnum %random-single-float %random-double-float
                    #!+long-float %random-long-float))
   (cond
-    ((and (fixnump arg) (<= arg random-fixnum-max) (> arg 0))
-     (rem (random-chunk state) arg))
+    ((and (fixnump arg) (> arg 0))
+     (%random-fixnum arg state))
     ((and (typep arg 'single-float) (> arg 0.0f0))
      (%random-single-float arg state))
     ((and (typep arg 'double-float) (> arg 0.0d0))
     #!+long-float
     ((and (typep arg 'long-float) (> arg 0.0l0))
      (%random-long-float arg state))
-    ((and (integerp arg) (> arg 0))
-     (%random-integer arg state))
+    ((and (bignump arg) (> arg 0))
+     (%random-bignum arg state))
     (t
      (error 'simple-type-error
             :expected-type '(or (integer 1) (float (0))) :datum arg
index 8d2eaed..b5cd536 100644 (file)
   (frob %random-single-float single-float)
   (frob %random-double-float double-float))
 
-;;; Mersenne Twister RNG
-;;;
-;;; FIXME: It's unpleasant to have RANDOM functionality scattered
-;;; through the code this way. It would be nice to move this into the
-;;; same file as the other RANDOM definitions.
+;;; Return an expression to generate an integer of N-BITS many random
+;;; bits, using the minimal number of random chunks possible.
+(defun generate-random-expr-for-power-of-2 (n-bits state)
+  (declare (type (integer 1 #.sb!vm:n-word-bits) n-bits))
+  (multiple-value-bind (n-chunk-bits chunk-expr)
+      (cond ((<= n-bits n-random-chunk-bits)
+             (values n-random-chunk-bits `(random-chunk ,state)))
+            ((<= n-bits (* 2 n-random-chunk-bits))
+             (values (* 2 n-random-chunk-bits) `(big-random-chunk ,state)))
+            (t
+             (error "Unexpectedly small N-RANDOM-CHUNK-BITS")))
+    (if (< n-bits n-chunk-bits)
+        `(logand ,(1- (ash 1 n-bits)) ,chunk-expr)
+        chunk-expr)))
+
+;;; This transform for compile-time constant word-sized integers
+;;; generates an accept-reject loop to achieve equidistribution of the
+;;; returned values. Several optimizations are done: If NUM is a power
+;;; of two no loop is needed. If the random chunk size is half the word
+;;; size only one chunk is used where sufficient. For values of NUM
+;;; where it is possible and results in faster code, the rejection
+;;; probability is reduced by accepting all values below the largest
+;;; multiple of the limit that fits into one or two chunks and and doing
+;;; a division to get the random value into the desired range.
 (deftransform random ((num &optional state)
-                      ((integer 1 #.(expt 2 sb!vm::n-word-bits)) &optional *))
-  ;; FIXME: I almost conditionalized this as #!+sb-doc. Find some way
-  ;; of automatically finding #!+sb-doc in proximity to DEFTRANSFORM
-  ;; to let me scan for places that I made this mistake and didn't
-  ;; catch myself.
-  "use inline (UNSIGNED-BYTE 32) operations"
-  (let ((type (lvar-type num))
-        (limit (expt 2 sb!vm::n-word-bits))
-        (random-chunk (ecase sb!vm::n-word-bits
-                        (32 'random-chunk)
-                        (64 'sb!kernel::big-random-chunk))))
-    (if (numeric-type-p type)
-        (let ((num-high (numeric-type-high (lvar-type num))))
-          (aver num-high)
-          (cond ((constant-lvar-p num)
-                 ;; Check the worst case sum absolute error for the
-                 ;; random number expectations.
-                 (let ((rem (rem limit num-high)))
-                   (unless (< (/ (* 2 rem (- num-high rem))
-                                 num-high limit)
-                              (expt 2 (- sb!kernel::random-integer-extra-bits)))
-                     (give-up-ir1-transform
-                      "The random number expectations are inaccurate."))
-                   (if (= num-high limit)
-                       `(,random-chunk (or state *random-state*))
-                       #!-(or x86 x86-64)
-                       `(rem (,random-chunk (or state *random-state*)) num)
-                       #!+(or x86 x86-64)
-                       ;; Use multiplication, which is faster.
-                       `(values (sb!bignum::%multiply
-                                 (,random-chunk (or state *random-state*))
-                                 num)))))
-                ((> num-high random-fixnum-max)
-                 (give-up-ir1-transform
-                  "The range is too large to ensure an accurate result."))
-                #!+(or x86 x86-64)
-                ((< num-high limit)
-                 `(values (sb!bignum::%multiply
-                           (,random-chunk (or state *random-state*))
-                           num)))
-                (t
-                 `(rem (,random-chunk (or state *random-state*)) num))))
-        ;; KLUDGE: a relatively conservative treatment, but better
-        ;; than a bug (reported by PFD sbcl-devel towards the end of
-        ;; 2004-11.
-        (give-up-ir1-transform
-         "Argument type is too complex to optimize for."))))
+                      ((constant-arg (integer 1 #.(expt 2 sb!vm:n-word-bits)))
+                       &optional *)
+                      *
+                      :policy (and (> speed compilation-speed)
+                                   (> speed space)))
+  "optimize to inlined RANDOM-CHUNK operations"
+  (let ((num (lvar-value num)))
+    (if (= num 1)
+        0
+        (flet ((chunk-n-bits-and-expr (n-bits)
+                 (cond ((<= n-bits n-random-chunk-bits)
+                        (values n-random-chunk-bits
+                                '(random-chunk (or state *random-state*))))
+                       ((<= n-bits (* 2 n-random-chunk-bits))
+                        (values (* 2 n-random-chunk-bits)
+                                '(big-random-chunk (or state *random-state*))))
+                       (t
+                        (error "Unexpectedly small N-RANDOM-CHUNK-BITS")))))
+          (if (zerop (logand num (1- num)))
+              ;; NUM is a power of 2.
+              (let ((n-bits (integer-length (1- num))))
+                (multiple-value-bind (n-chunk-bits chunk-expr)
+                    (chunk-n-bits-and-expr n-bits)
+                  (if (< n-bits n-chunk-bits)
+                      `(logand ,(1- (ash 1 n-bits)) ,chunk-expr)
+                      chunk-expr)))
+              ;; Generate an accept-reject loop.
+              (let ((n-bits (integer-length num)))
+                (multiple-value-bind (n-chunk-bits chunk-expr)
+                    (chunk-n-bits-and-expr n-bits)
+                  (if (or (> (* num 3) (expt 2 n-chunk-bits))
+                          (logbitp (- n-bits 2) num))
+                      ;; Division can't help as the quotient is below 3,
+                      ;; or is too costly as the rejection probability
+                      ;; without it is already small (namely at most 1/4
+                      ;; with the given test, which is experimentally a
+                      ;; reasonable threshold and cheap to test for).
+                      `(loop
+                         (let ((bits ,(generate-random-expr-for-power-of-2
+                                       n-bits '(or state *random-state*))))
+                           (when (< bits num)
+                             (return bits))))
+                      (let ((d (truncate (expt 2 n-chunk-bits) num)))
+                        `(loop
+                           (let ((bits ,chunk-expr))
+                             (when (< bits ,(* num d))
+                               (return (values (truncate bits ,d)))))))))))))))
+
 \f
 ;;;; float accessors
 
index 101333e..829ff4d 100644 (file)
                                nconc (list (1- (expt 2 i))
                                            (expt 2 i)
                                            (1+ (expt 2 i))))
-                       ,@(loop for i from (1- sb-kernel::random-chunk-length)
-                               to (* sb-kernel::random-chunk-length 4)
+                       ,@(loop for i from (1- sb-kernel::n-random-chunk-bits)
+                               to (* sb-kernel::n-random-chunk-bits 4)
                                collect (* 3 (expt 2 i)))
                        ,@(loop for i from 2 to sb-vm:n-word-bits
                                for n = (expt 16 i)