X-Git-Url: http://repo.macrolet.net/gitweb/?a=blobdiff_plain;f=src%2Fcode%2Fbignum.lisp;h=054e6dfd17b56e6b77c5d4f31ea528c528358896;hb=0e3c4b4db102bd204a30402d7e5a0de44aea57ce;hp=fc6c45ee590da11fa1206108841b8707df0aa6d7;hpb=4dc4761909992ceb346d003f3fb19e5c837ee985;p=sbcl.git diff --git a/src/code/bignum.lisp b/src/code/bignum.lisp index fc6c45e..054e6df 100644 --- a/src/code/bignum.lisp +++ b/src/code/bignum.lisp @@ -20,7 +20,7 @@ ;;; bignum-ashift-right bignum-ashift-left bignum-gcd ;;; bignum-to-float bignum-integer-length ;;; bignum-logical-and bignum-logical-ior bignum-logical-xor -;;; bignum-logical-not bignum-load-byte bignum-deposit-byte +;;; bignum-logical-not bignum-load-byte ;;; bignum-truncate bignum-plus-p bignum-compare make-small-bignum ;;; bignum-logbitp bignum-logcount ;;; These symbols define the interface to the compiler: @@ -28,7 +28,7 @@ ;;; %bignum-length %bignum-set-length %bignum-ref %bignum-set ;;; %digit-0-or-plusp %add-with-carry %subtract-with-borrow ;;; %multiply-and-add %multiply %lognot %logand %logior %logxor -;;; %fixnum-to-digit %floor %fixnum-digit-with-correct-sign %ashl +;;; %fixnum-to-digit %bigfloor %fixnum-digit-with-correct-sign %ashl ;;; %ashr %digit-logical-shift-right)) ;;; The following interfaces will either be assembler routines or code @@ -41,7 +41,7 @@ ;;; %BIGNUM-SET-LENGTH ;;; %FIXNUM-DIGIT-WITH-CORRECT-SIGN ;;; %SIGN-DIGIT -;;; %ASHR +;;; %ASHR ;;; %ASHL ;;; %BIGNUM-0-OR-PLUSP ;;; %DIGIT-LOGICAL-SHIFT-RIGHT @@ -67,7 +67,7 @@ ;;; LDB ;;; %FIXNUM-TO-DIGIT ;;; TRUNCATE -;;; %FLOOR +;;; %BIGFLOOR ;;; ;;; Note: The floating routines know about the float representation. ;;; @@ -87,12 +87,12 @@ ;;; fixnums ;;; logior, logxor, logand ;;; depending on relationals, < (twice) and <= (twice) -;;; or write compare thing (twice). +;;; or write compare thing (twice). ;;; LDB on fixnum with bignum result. ;;; DPB on fixnum with bignum result. ;;; TRUNCATE returns zero or one as one value and fixnum or minus fixnum -;;; for the other value when given (truncate fixnum bignum). -;;; Returns (truncate bignum fixnum) otherwise. +;;; for the other value when given (truncate fixnum bignum). +;;; Returns (truncate bignum fixnum) otherwise. ;;; addition ;;; subtraction (twice) ;;; multiply @@ -100,11 +100,11 @@ ;;; Write MASK-FIELD and DEPOSIT-FIELD in terms of logical operations. ;;; DIVIDE ;;; IF (/ x y) with bignums: -;;; do the truncate, and if rem is 0, return quotient. -;;; if rem is non-0 -;;; gcd of x and y. -;;; "truncate" each by gcd, ignoring remainder 0. -;;; form ratio of each result, bottom is positive. +;;; do the truncate, and if rem is 0, return quotient. +;;; if rem is non-0 +;;; gcd of x and y. +;;; "truncate" each by gcd, ignoring remainder 0. +;;; form ratio of each result, bottom is positive. ;;;; What's a bignum? @@ -131,12 +131,12 @@ ;;; to be able to return the digit somewhere no one looks for real objects. (defun %bignum-ref (bignum i) (declare (type bignum-type bignum) - (type bignum-index i)) + (type bignum-index i)) (%bignum-ref bignum i)) (defun %bignum-set (bignum i value) (declare (type bignum-type bignum) - (type bignum-index i) - (type bignum-element-type value)) + (type bignum-index i) + (type bignum-element-type value)) (%bignum-set bignum i value)) ;;; Return T if digit is positive, or NIL if negative. @@ -147,7 +147,7 @@ #!-sb-fluid (declaim (inline %bignum-0-or-plusp)) (defun %bignum-0-or-plusp (bignum len) (declare (type bignum-type bignum) - (type bignum-index len)) + (type bignum-index len)) (%digit-0-or-plusp (%bignum-ref bignum (1- len)))) ;;; This should be in assembler, and should not cons intermediate @@ -155,7 +155,7 @@ ;;; together a, b, and an incoming carry. (defun %add-with-carry (a b carry) (declare (type bignum-element-type a b) - (type (mod 2) carry)) + (type (mod 2) carry)) (%add-with-carry a b carry)) ;;; This should be in assembler, and should not cons intermediate @@ -165,7 +165,7 @@ ;;; We really do: a - b - 1 + borrow, where borrow is either 0 or 1. (defun %subtract-with-borrow (a b borrow) (declare (type bignum-element-type a b) - (type (mod 2) borrow)) + (type (mod 2) borrow)) (%subtract-with-borrow a b borrow)) ;;; Multiply two digit-size numbers, returning a 2*digit-size result @@ -185,7 +185,7 @@ ;;; accumulating partial results which is where the res-digit comes ;;; from. (defun %multiply-and-add (x-digit y-digit carry-in-digit - &optional (res-digit 0)) + &optional (res-digit 0)) (declare (type bignum-element-type x-digit y-digit res-digit carry-in-digit)) (%multiply-and-add x-digit y-digit carry-in-digit res-digit)) @@ -194,7 +194,7 @@ (%lognot digit)) ;;; Each of these does the digit-size unsigned op. -#!-sb-fluid (declaim (inline %logand %logior %logxor)) +(declaim (inline %logand %logior %logxor)) (defun %logand (a b) (declare (type bignum-element-type a b)) (logand a b)) @@ -216,13 +216,13 @@ ;;; dividing the first two as a 2*digit-size integer by the third. ;;; ;;; Do weird LET and SETQ stuff to bamboozle the compiler into allowing -;;; the %FLOOR transform to expand into pseudo-assembler for which the +;;; the %BIGFLOOR transform to expand into pseudo-assembler for which the ;;; compiler can later correctly allocate registers. -(defun %floor (a b c) +(defun %bigfloor (a b c) (let ((a a) (b b) (c c)) (declare (type bignum-element-type a b c)) (setq a a b b c c) - (%floor a b c))) + (%bigfloor a b c))) ;;; Convert the digit to a regular integer assuming that the digit is signed. (defun %fixnum-digit-with-correct-sign (digit) @@ -235,27 +235,27 @@ ;;; unsigned. (defun %ashr (data count) (declare (type bignum-element-type data) - (type (mod #.sb!vm:n-word-bits) count)) + (type (mod #.sb!vm:n-word-bits) count)) (%ashr data count)) ;;; This takes a digit-size quantity and shifts it to the left, ;;; returning a digit-size quantity. (defun %ashl (data count) (declare (type bignum-element-type data) - (type (mod #.sb!vm:n-word-bits) count)) + (type (mod #.sb!vm:n-word-bits) count)) (%ashl data count)) ;;; Do an unsigned (logical) right shift of a digit by Count. (defun %digit-logical-shift-right (data count) (declare (type bignum-element-type data) - (type (mod #.sb!vm:n-word-bits) count)) + (type (mod #.sb!vm:n-word-bits) count)) (%digit-logical-shift-right data count)) ;;; Change the length of bignum to be newlen. Newlen must be the same or ;;; smaller than the old length, and any elements beyond newlen must be zeroed. (defun %bignum-set-length (bignum newlen) (declare (type bignum-type bignum) - (type bignum-index newlen)) + (type bignum-index newlen)) (%bignum-set-length bignum newlen)) ;;; This returns 0 or "-1" depending on whether the bignum is positive. This @@ -265,12 +265,12 @@ #!-sb-fluid (declaim (inline %sign-digit)) (defun %sign-digit (bignum len) (declare (type bignum-type bignum) - (type bignum-index len)) + (type bignum-index len)) (%ashr (%bignum-ref bignum (1- len)) (1- digit-size))) ;;; These take two digit-size quantities and compare or contrast them ;;; without wasting time with incorrect type checking. -#!-sb-fluid (declaim (inline %digit-compare %digit-greater)) +(declaim (inline %digit-compare %digit-greater)) (defun %digit-compare (x y) (= x y)) (defun %digit-greater (x y) @@ -283,50 +283,50 @@ (defun add-bignums (a b) (declare (type bignum-type a b)) (let ((len-a (%bignum-length a)) - (len-b (%bignum-length b))) + (len-b (%bignum-length b))) (declare (type bignum-index len-a len-b)) (multiple-value-bind (a len-a b len-b) - (if (> len-a len-b) - (values a len-a b len-b) - (values b len-b a len-a)) + (if (> len-a len-b) + (values a len-a b len-b) + (values b len-b a len-a)) (declare (type bignum-type a b) - (type bignum-index len-a len-b)) + (type bignum-index len-a len-b)) (let* ((len-res (1+ len-a)) - (res (%allocate-bignum len-res)) - (carry 0)) - (declare (type bignum-index len-res) - (type bignum-type res) - (type (mod 2) carry)) - (dotimes (i len-b) - (declare (type bignum-index i)) - (multiple-value-bind (v k) - (%add-with-carry (%bignum-ref a i) (%bignum-ref b i) carry) - (declare (type bignum-element-type v) - (type (mod 2) k)) - (setf (%bignum-ref res i) v) - (setf carry k))) - (if (/= len-a len-b) - (finish-add a res carry (%sign-digit b len-b) len-b len-a) - (setf (%bignum-ref res len-a) - (%add-with-carry (%sign-digit a len-a) - (%sign-digit b len-b) - carry))) - (%normalize-bignum res len-res))))) + (res (%allocate-bignum len-res)) + (carry 0)) + (declare (type bignum-index len-res) + (type bignum-type res) + (type (mod 2) carry)) + (dotimes (i len-b) + (declare (type bignum-index i)) + (multiple-value-bind (v k) + (%add-with-carry (%bignum-ref a i) (%bignum-ref b i) carry) + (declare (type bignum-element-type v) + (type (mod 2) k)) + (setf (%bignum-ref res i) v) + (setf carry k))) + (if (/= len-a len-b) + (finish-add a res carry (%sign-digit b len-b) len-b len-a) + (setf (%bignum-ref res len-a) + (%add-with-carry (%sign-digit a len-a) + (%sign-digit b len-b) + carry))) + (%normalize-bignum res len-res))))) ;;; This takes the longer of two bignums and propagates the carry through its ;;; remaining high order digits. (defun finish-add (a res carry sign-digit-b start end) (declare (type bignum-type a res) - (type (mod 2) carry) - (type bignum-element-type sign-digit-b) - (type bignum-index start end)) + (type (mod 2) carry) + (type bignum-element-type sign-digit-b) + (type bignum-index start end)) (do ((i start (1+ i))) ((= i end) (setf (%bignum-ref res end) - (%add-with-carry (%sign-digit a end) sign-digit-b carry))) + (%add-with-carry (%sign-digit a end) sign-digit-b carry))) (declare (type bignum-index i)) (multiple-value-bind (v k) - (%add-with-carry (%bignum-ref a i) sign-digit-b carry) + (%add-with-carry (%bignum-ref a i) sign-digit-b carry) (setf (%bignum-ref res i) v) (setf carry k))) (values)) @@ -339,27 +339,20 @@ ;;; function to call that fixes up the result returning any useful values, such ;;; as the result. This macro may evaluate its arguments more than once. (sb!xc:defmacro subtract-bignum-loop (a len-a b len-b res len-res return-fun) - (let ((borrow (gensym)) - (a-digit (gensym)) - (a-sign (gensym)) - (b-digit (gensym)) - (b-sign (gensym)) - (i (gensym)) - (v (gensym)) - (k (gensym))) + (with-unique-names (borrow a-digit a-sign b-digit b-sign i v k) `(let* ((,borrow 1) - (,a-sign (%sign-digit ,a ,len-a)) - (,b-sign (%sign-digit ,b ,len-b))) + (,a-sign (%sign-digit ,a ,len-a)) + (,b-sign (%sign-digit ,b ,len-b))) (declare (type bignum-element-type ,a-sign ,b-sign)) (dotimes (,i ,len-res) - (declare (type bignum-index ,i)) - (let ((,a-digit (if (< ,i ,len-a) (%bignum-ref ,a ,i) ,a-sign)) - (,b-digit (if (< ,i ,len-b) (%bignum-ref ,b ,i) ,b-sign))) - (declare (type bignum-element-type ,a-digit ,b-digit)) - (multiple-value-bind (,v ,k) - (%subtract-with-borrow ,a-digit ,b-digit ,borrow) - (setf (%bignum-ref ,res ,i) ,v) - (setf ,borrow ,k)))) + (declare (type bignum-index ,i)) + (let ((,a-digit (if (< ,i ,len-a) (%bignum-ref ,a ,i) ,a-sign)) + (,b-digit (if (< ,i ,len-b) (%bignum-ref ,b ,i) ,b-sign))) + (declare (type bignum-element-type ,a-digit ,b-digit)) + (multiple-value-bind (,v ,k) + (%subtract-with-borrow ,a-digit ,b-digit ,borrow) + (setf (%bignum-ref ,res ,i) ,v) + (setf ,borrow ,k)))) (,return-fun ,res ,len-res)))) ) ;EVAL-WHEN @@ -367,9 +360,9 @@ (defun subtract-bignum (a b) (declare (type bignum-type a b)) (let* ((len-a (%bignum-length a)) - (len-b (%bignum-length b)) - (len-res (1+ (max len-a len-b))) - (res (%allocate-bignum len-res))) + (len-b (%bignum-length b)) + (len-res (1+ (max len-a len-b))) + (res (%allocate-bignum len-res))) (declare (type bignum-index len-a len-b len-res)) ;Test len-res for bounds? (subtract-bignum-loop a len-a b len-b res len-res %normalize-bignum))) @@ -377,71 +370,71 @@ ;;; results, such as GCD, use this. It assumes Result is big enough for the ;;; result. (defun subtract-bignum-buffers-with-len (a len-a b len-b result len-res) - (declare (type bignum-type a b) - (type bignum-index len-a len-b)) + (declare (type bignum-type a b result) + (type bignum-index len-a len-b len-res)) (subtract-bignum-loop a len-a b len-b result len-res - %normalize-bignum-buffer)) + %normalize-bignum-buffer)) (defun subtract-bignum-buffers (a len-a b len-b result) - (declare (type bignum-type a b) - (type bignum-index len-a len-b)) + (declare (type bignum-type a b result) + (type bignum-index len-a len-b)) (subtract-bignum-loop a len-a b len-b result (max len-a len-b) - %normalize-bignum-buffer)) + %normalize-bignum-buffer)) ;;;; multiplication (defun multiply-bignums (a b) (declare (type bignum-type a b)) (let* ((a-plusp (%bignum-0-or-plusp a (%bignum-length a))) - (b-plusp (%bignum-0-or-plusp b (%bignum-length b))) - (a (if a-plusp a (negate-bignum a))) - (b (if b-plusp b (negate-bignum b))) - (len-a (%bignum-length a)) - (len-b (%bignum-length b)) - (len-res (+ len-a len-b)) - (res (%allocate-bignum len-res)) - (negate-res (not (eq a-plusp b-plusp)))) + (b-plusp (%bignum-0-or-plusp b (%bignum-length b))) + (a (if a-plusp a (negate-bignum a))) + (b (if b-plusp b (negate-bignum b))) + (len-a (%bignum-length a)) + (len-b (%bignum-length b)) + (len-res (+ len-a len-b)) + (res (%allocate-bignum len-res)) + (negate-res (not (eq a-plusp b-plusp)))) (declare (type bignum-index len-a len-b len-res)) (dotimes (i len-a) (declare (type bignum-index i)) (let ((carry-digit 0) - (x (%bignum-ref a i)) - (k i)) - (declare (type bignum-index k) - (type bignum-element-type carry-digit x)) - (dotimes (j len-b) - (multiple-value-bind (big-carry res-digit) - (%multiply-and-add x - (%bignum-ref b j) - (%bignum-ref res k) - carry-digit) - (declare (type bignum-element-type big-carry res-digit)) - (setf (%bignum-ref res k) res-digit) - (setf carry-digit big-carry) - (incf k))) - (setf (%bignum-ref res k) carry-digit))) + (x (%bignum-ref a i)) + (k i)) + (declare (type bignum-index k) + (type bignum-element-type carry-digit x)) + (dotimes (j len-b) + (multiple-value-bind (big-carry res-digit) + (%multiply-and-add x + (%bignum-ref b j) + (%bignum-ref res k) + carry-digit) + (declare (type bignum-element-type big-carry res-digit)) + (setf (%bignum-ref res k) res-digit) + (setf carry-digit big-carry) + (incf k))) + (setf (%bignum-ref res k) carry-digit))) (when negate-res (negate-bignum-in-place res)) (%normalize-bignum res len-res))) (defun multiply-bignum-and-fixnum (bignum fixnum) (declare (type bignum-type bignum) (type fixnum fixnum)) (let* ((bignum-plus-p (%bignum-0-or-plusp bignum (%bignum-length bignum))) - (fixnum-plus-p (not (minusp fixnum))) - (bignum (if bignum-plus-p bignum (negate-bignum bignum))) - (bignum-len (%bignum-length bignum)) - (fixnum (if fixnum-plus-p fixnum (- fixnum))) - (result (%allocate-bignum (1+ bignum-len))) - (carry-digit 0)) + (fixnum-plus-p (not (minusp fixnum))) + (bignum (if bignum-plus-p bignum (negate-bignum bignum))) + (bignum-len (%bignum-length bignum)) + (fixnum (if fixnum-plus-p fixnum (- fixnum))) + (result (%allocate-bignum (1+ bignum-len))) + (carry-digit 0)) (declare (type bignum-type bignum result) - (type bignum-index bignum-len) - (type bignum-element-type fixnum carry-digit)) + (type bignum-index bignum-len) + (type bignum-element-type fixnum carry-digit)) (dotimes (index bignum-len) (declare (type bignum-index index)) (multiple-value-bind (next-digit low) - (%multiply-and-add (%bignum-ref bignum index) fixnum carry-digit) - (declare (type bignum-element-type next-digit low)) - (setf carry-digit next-digit) - (setf (%bignum-ref result index) low))) + (%multiply-and-add (%bignum-ref bignum index) fixnum carry-digit) + (declare (type bignum-element-type next-digit low)) + (setf carry-digit next-digit) + (setf (%bignum-ref result index) low))) (setf (%bignum-ref result bignum-len) carry-digit) (unless (eq bignum-plus-p fixnum-plus-p) (negate-bignum-in-place result)) @@ -450,76 +443,75 @@ (defun multiply-fixnums (a b) (declare (fixnum a b)) (let* ((a-minusp (minusp a)) - (b-minusp (minusp b))) + (b-minusp (minusp b))) (multiple-value-bind (high low) - (%multiply (if a-minusp (- a) a) - (if b-minusp (- b) b)) + (%multiply (if a-minusp (- a) a) + (if b-minusp (- b) b)) (declare (type bignum-element-type high low)) (if (and (zerop high) - (%digit-0-or-plusp low)) - (let ((low (sb!ext:truly-the (unsigned-byte #.(1- sb!vm:n-word-bits)) - (%fixnum-digit-with-correct-sign low)))) - (if (eq a-minusp b-minusp) - low - (- low))) - (let ((res (%allocate-bignum 2))) - (%bignum-set res 0 low) - (%bignum-set res 1 high) - (unless (eq a-minusp b-minusp) (negate-bignum-in-place res)) - (%normalize-bignum res 2)))))) + (%digit-0-or-plusp low)) + (let ((low (sb!ext:truly-the (unsigned-byte #.(1- sb!vm:n-word-bits)) + (%fixnum-digit-with-correct-sign low)))) + (if (eq a-minusp b-minusp) + low + (- low))) + (let ((res (%allocate-bignum 2))) + (%bignum-set res 0 low) + (%bignum-set res 1 high) + (unless (eq a-minusp b-minusp) (negate-bignum-in-place res)) + (%normalize-bignum res 2)))))) ;;;; BIGNUM-REPLACE and WITH-BIGNUM-BUFFERS (eval-when (:compile-toplevel :execute) (sb!xc:defmacro bignum-replace (dest - src - &key - (start1 '0) - end1 - (start2 '0) - end2 - from-end) + src + &key + (start1 '0) + end1 + (start2 '0) + end2 + from-end) (sb!int:once-only ((n-dest dest) - (n-src src)) - (let ((n-start1 (gensym)) - (n-end1 (gensym)) - (n-start2 (gensym)) - (n-end2 (gensym)) - (i1 (gensym)) - (i2 (gensym)) - (end1 (or end1 `(%bignum-length ,n-dest))) - (end2 (or end2 `(%bignum-length ,n-src)))) - (if from-end - `(let ((,n-start1 ,start1) - (,n-start2 ,start2)) - (do ((,i1 (1- ,end1) (1- ,i1)) - (,i2 (1- ,end2) (1- ,i2))) - ((or (< ,i1 ,n-start1) (< ,i2 ,n-start2))) - (declare (fixnum ,i1 ,i2)) - (%bignum-set ,n-dest ,i1 - (%bignum-ref ,n-src ,i2)))) - `(let ((,n-end1 ,end1) - (,n-end2 ,end2)) - (do ((,i1 ,start1 (1+ ,i1)) - (,i2 ,start2 (1+ ,i2))) - ((or (>= ,i1 ,n-end1) (>= ,i2 ,n-end2))) - (declare (type bignum-index ,i1 ,i2)) - (%bignum-set ,n-dest ,i1 - (%bignum-ref ,n-src ,i2)))))))) + (n-src src)) + (with-unique-names (n-start1 n-end1 n-start2 n-end2 i1 i2) + (let ((end1 (or end1 `(%bignum-length ,n-dest))) + (end2 (or end2 `(%bignum-length ,n-src)))) + (if from-end + `(let ((,n-start1 ,start1) + (,n-start2 ,start2)) + (do ((,i1 (1- ,end1) (1- ,i1)) + (,i2 (1- ,end2) (1- ,i2))) + ((or (< ,i1 ,n-start1) (< ,i2 ,n-start2))) + (declare (fixnum ,i1 ,i2)) + (%bignum-set ,n-dest ,i1 (%bignum-ref ,n-src ,i2)))) + (if (eql start1 start2) + `(let ((,n-end1 (min ,end1 ,end2))) + (do ((,i1 ,start1 (1+ ,i1))) + ((>= ,i1 ,n-end1)) + (declare (type bignum-index ,i1)) + (%bignum-set ,n-dest ,i1 (%bignum-ref ,n-src ,i1)))) + `(let ((,n-end1 ,end1) + (,n-end2 ,end2)) + (do ((,i1 ,start1 (1+ ,i1)) + (,i2 ,start2 (1+ ,i2))) + ((or (>= ,i1 ,n-end1) (>= ,i2 ,n-end2))) + (declare (type bignum-index ,i1 ,i2)) + (%bignum-set ,n-dest ,i1 (%bignum-ref ,n-src ,i2)))))))))) (sb!xc:defmacro with-bignum-buffers (specs &body body) #!+sb-doc "WITH-BIGNUM-BUFFERS ({(var size [init])}*) Form*" (sb!int:collect ((binds) - (inits)) + (inits)) (dolist (spec specs) (let ((name (first spec)) - (size (second spec))) - (binds `(,name (%allocate-bignum ,size))) - (let ((init (third spec))) - (when init - (inits `(bignum-replace ,name ,init)))))) + (size (second spec))) + (binds `(,name (%allocate-bignum ,size))) + (let ((init (third spec))) + (when init + (inits `(bignum-replace ,name ,init)))))) `(let* ,(binds) ,@(inits) ,@body))) @@ -533,7 +525,7 @@ ;; check in normal use, and are disabled here. (sb!xc:defmacro gcd-assert (&rest args) (if nil - `(assert ,@args))) + `(assert ,@args))) ;; We'll be doing a lot of modular arithmetic. (sb!xc:defmacro modularly (form) `(logand all-ones-digit ,form))) @@ -543,8 +535,8 @@ ;;; it, we pay a heavy price in BIGNUM-GCD when compiled by the ;;; cross-compiler. -- CSR, 2004-07-19 (declaim (ftype (sfunction (bignum-type bignum-index bignum-type bignum-index) - sb!vm::positive-fixnum) - bignum-factors-of-two)) + (and unsigned-byte fixnum)) + bignum-factors-of-two)) (defun bignum-factors-of-two (a len-a b len-b) (declare (type bignum-index len-a len-b) (type bignum-type a b)) (do ((i 0 (1+ i)) @@ -553,10 +545,10 @@ (declare (type bignum-index i end)) (let ((or-digits (%logior (%bignum-ref a i) (%bignum-ref b i)))) (unless (zerop or-digits) - (return (do ((j 0 (1+ j)) - (or-digits or-digits (%ashr or-digits 1))) - ((oddp or-digits) (+ (* i digit-size) j)) - (declare (type (mod #.sb!vm:n-word-bits) j)))))))) + (return (do ((j 0 (1+ j)) + (or-digits or-digits (%ashr or-digits 1))) + ((oddp or-digits) (+ (* i digit-size) j)) + (declare (type (mod #.sb!vm:n-word-bits) j)))))))) ;;; Multiply a bignum buffer with a fixnum or a digit, storing the ;;; result in another bignum buffer, and without using any @@ -565,32 +557,32 @@ ;;; MULTIPLY-BIGNUM-AND-FIXNUM. (declaim (inline multiply-bignum-buffer-and-smallnum-to-buffer)) (defun multiply-bignum-buffer-and-smallnum-to-buffer (bignum bignum-len - smallnum res) + smallnum res) (declare (type bignum-type bignum)) (let* ((bignum-plus-p (%bignum-0-or-plusp bignum bignum-len)) - (smallnum-plus-p (not (minusp smallnum))) - (smallnum (if smallnum-plus-p smallnum (- smallnum))) - (carry-digit 0)) + (smallnum-plus-p (not (minusp smallnum))) + (smallnum (if smallnum-plus-p smallnum (- smallnum))) + (carry-digit 0)) (declare (type bignum-type bignum res) - (type bignum-index bignum-len) - (type bignum-element-type smallnum carry-digit)) + (type bignum-index bignum-len) + (type bignum-element-type smallnum carry-digit)) (unless bignum-plus-p (negate-bignum-buffer-in-place bignum bignum-len)) (dotimes (index bignum-len) (declare (type bignum-index index)) (multiple-value-bind (next-digit low) - (%multiply-and-add (%bignum-ref bignum index) - smallnum - carry-digit) - (declare (type bignum-element-type next-digit low)) - (setf carry-digit next-digit) - (setf (%bignum-ref res index) low))) + (%multiply-and-add (%bignum-ref bignum index) + smallnum + carry-digit) + (declare (type bignum-element-type next-digit low)) + (setf carry-digit next-digit) + (setf (%bignum-ref res index) low))) (setf (%bignum-ref res bignum-len) carry-digit) (unless bignum-plus-p (negate-bignum-buffer-in-place bignum bignum-len)) (let ((res-len (%normalize-bignum-buffer res (1+ bignum-len)))) (unless (eq bignum-plus-p smallnum-plus-p) - (negate-bignum-buffer-in-place res res-len)) + (negate-bignum-buffer-in-place res res-len)) res-len))) ;;; Given U and V, return U / V mod 2^32. Implements the algorithm in the @@ -598,79 +590,79 @@ (declaim (inline bmod)) (defun bmod (u v) (let ((ud (%bignum-ref u 0)) - (vd (%bignum-ref v 0)) - (umask 0) - (imask 1) - (m 0)) + (vd (%bignum-ref v 0)) + (umask 0) + (imask 1) + (m 0)) (declare (type (unsigned-byte #.sb!vm:n-word-bits) ud vd umask imask m)) (dotimes (i digit-size) (setf umask (logior umask imask)) - (unless (zerop (logand ud umask)) - (setf ud (modularly (- ud vd))) - (setf m (modularly (logior m imask)))) + (when (logtest ud umask) + (setf ud (modularly (- ud vd))) + (setf m (modularly (logior m imask)))) (setf imask (modularly (ash imask 1))) (setf vd (modularly (ash vd 1)))) m)) (defun dmod (u u-len v v-len tmp1) (loop while (> (bignum-buffer-integer-length u u-len) - (+ (bignum-buffer-integer-length v v-len) - digit-size)) + (+ (bignum-buffer-integer-length v v-len) + digit-size)) do (unless (zerop (%bignum-ref u 0)) (let* ((bmod (bmod u v)) - (tmp1-len (multiply-bignum-buffer-and-smallnum-to-buffer v v-len - bmod - tmp1))) - (setf u-len (subtract-bignum-buffers u u-len - tmp1 tmp1-len - u)) - (bignum-abs-buffer u u-len))) + (tmp1-len (multiply-bignum-buffer-and-smallnum-to-buffer v v-len + bmod + tmp1))) + (setf u-len (subtract-bignum-buffers u u-len + tmp1 tmp1-len + u)) + (bignum-abs-buffer u u-len))) (gcd-assert (zerop (%bignum-ref u 0))) (setf u-len (bignum-buffer-ashift-right u u-len digit-size))) (let* ((d (+ 1 (- (bignum-buffer-integer-length u u-len) - (bignum-buffer-integer-length v v-len)))) - (n (1- (ash 1 d)))) + (bignum-buffer-integer-length v v-len)))) + (n (1- (ash 1 d)))) (declare (type (unsigned-byte #.(integer-length #.sb!vm:n-word-bits)) d) - (type (unsigned-byte #.sb!vm:n-word-bits) n)) + (type (unsigned-byte #.sb!vm:n-word-bits) n)) (gcd-assert (>= d 0)) - (unless (zerop (logand (%bignum-ref u 0) n)) + (when (logtest (%bignum-ref u 0) n) (let ((tmp1-len - (multiply-bignum-buffer-and-smallnum-to-buffer v v-len - (logand n (bmod u - v)) - tmp1))) - (setf u-len (subtract-bignum-buffers u u-len - tmp1 tmp1-len - u)) - (bignum-abs-buffer u u-len))) + (multiply-bignum-buffer-and-smallnum-to-buffer v v-len + (logand n (bmod u + v)) + tmp1))) + (setf u-len (subtract-bignum-buffers u u-len + tmp1 tmp1-len + u)) + (bignum-abs-buffer u u-len))) u-len)) (defconstant lower-ones-digit (1- (ash 1 (truncate sb!vm:n-word-bits 2)))) - + ;;; Find D and N such that (LOGAND ALL-ONES-DIGIT (- (* D X) (* N Y))) is 0, ;;; (< 0 N LOWER-ONES-DIGIT) and (< 0 (ABS D) LOWER-ONES-DIGIT). (defun reduced-ratio-mod (x y) (let* ((c (bmod x y)) - (n1 c) - (d1 1) - (n2 (modularly (1+ (modularly (lognot n1))))) - (d2 (modularly -1))) + (n1 c) + (d1 1) + (n2 (modularly (1+ (modularly (lognot n1))))) + (d2 (modularly -1))) (declare (type (unsigned-byte #.sb!vm:n-word-bits) n1 d1 n2 d2)) (loop while (> n2 (expt 2 (truncate digit-size 2))) do - (loop for i of-type (mod #.sb!vm:n-word-bits) - downfrom (- (integer-length n1) (integer-length n2)) - while (>= n1 n2) do - (when (>= n1 (modularly (ash n2 i))) - (psetf n1 (modularly (- n1 (modularly (ash n2 i)))) - d1 (modularly (- d1 (modularly (ash d2 i))))))) - (psetf n1 n2 - d1 d2 - n2 n1 - d2 d1)) + (loop for i of-type (mod #.sb!vm:n-word-bits) + downfrom (- (integer-length n1) (integer-length n2)) + while (>= n1 n2) do + (when (>= n1 (modularly (ash n2 i))) + (psetf n1 (modularly (- n1 (modularly (ash n2 i)))) + d1 (modularly (- d1 (modularly (ash d2 i))))))) + (psetf n1 n2 + d1 d2 + n2 n1 + d2 d1)) (values n2 (if (>= d2 (expt 2 (1- digit-size))) - (lognot (logand most-positive-fixnum (lognot d2))) - (logand lower-ones-digit d2))))) + (lognot (logand most-positive-fixnum (lognot d2))) + (logand lower-ones-digit d2))))) (defun copy-bignum (a &optional (len (%bignum-length a))) @@ -678,7 +670,7 @@ (bignum-replace b a) (%bignum-set-length b len) b)) - + ;;; Allocate a single word bignum that holds fixnum. This is useful when ;;; we are trying to mix fixnum and bignum operands. #!-sb-fluid (declaim (inline make-small-bignum)) @@ -688,7 +680,7 @@ res)) ;; When the larger number is less than this many bignum digits long, revert -;; to old algorithm. +;; to old algorithm. (defparameter *accelerated-gcd-cutoff* 3) ;;; Alternate between k-ary reduction with the help of @@ -701,96 +693,97 @@ (defun bignum-gcd (u0 v0) (declare (type bignum-type u0 v0)) (let* ((u1 (if (%bignum-0-or-plusp u0 (%bignum-length u0)) - u0 - (negate-bignum u0 nil))) - (v1 (if (%bignum-0-or-plusp v0 (%bignum-length v0)) - v0 - (negate-bignum v0 nil)))) + u0 + (negate-bignum u0 nil))) + (v1 (if (%bignum-0-or-plusp v0 (%bignum-length v0)) + v0 + (negate-bignum v0 nil)))) (if (zerop v1) - (return-from bignum-gcd u1)) + (return-from bignum-gcd u1)) (when (> u1 v1) (rotatef u1 v1)) (let ((n (mod v1 u1))) (setf v1 (if (fixnump n) - (make-small-bignum n) - n))) + (make-small-bignum n) + n))) (if (and (= 1 (%bignum-length v1)) - (zerop (%bignum-ref v1 0))) - (return-from bignum-gcd (%normalize-bignum u1 - (%bignum-length u1)))) + (zerop (%bignum-ref v1 0))) + (return-from bignum-gcd (%normalize-bignum u1 + (%bignum-length u1)))) (let* ((buffer-len (+ 2 (%bignum-length u1))) - (u (%allocate-bignum buffer-len)) - (u-len (%bignum-length u1)) - (v (%allocate-bignum buffer-len)) - (v-len (%bignum-length v1)) - (tmp1 (%allocate-bignum buffer-len)) - (tmp1-len 0) - (tmp2 (%allocate-bignum buffer-len)) - (tmp2-len 0) - (factors-of-two - (bignum-factors-of-two u1 (%bignum-length u1) - v1 (%bignum-length v1)))) + (u (%allocate-bignum buffer-len)) + (u-len (%bignum-length u1)) + (v (%allocate-bignum buffer-len)) + (v-len (%bignum-length v1)) + (tmp1 (%allocate-bignum buffer-len)) + (tmp1-len 0) + (tmp2 (%allocate-bignum buffer-len)) + (tmp2-len 0) + (factors-of-two + (bignum-factors-of-two u1 (%bignum-length u1) + v1 (%bignum-length v1)))) (declare (type (or null bignum-index) - buffer-len u-len v-len tmp1-len tmp2-len)) + buffer-len u-len v-len tmp1-len tmp2-len)) (bignum-replace u u1) (bignum-replace v v1) (setf u-len - (make-gcd-bignum-odd u - (bignum-buffer-ashift-right u u-len - factors-of-two))) + (make-gcd-bignum-odd u + (bignum-buffer-ashift-right u u-len + factors-of-two))) (setf v-len - (make-gcd-bignum-odd v - (bignum-buffer-ashift-right v v-len - factors-of-two))) + (make-gcd-bignum-odd v + (bignum-buffer-ashift-right v v-len + factors-of-two))) (loop until (or (< u-len *accelerated-gcd-cutoff*) - (not v-len) - (zerop v-len) - (and (= 1 v-len) - (zerop (%bignum-ref v 0)))) - do - (gcd-assert (= buffer-len (%bignum-length u) - (%bignum-length v) - (%bignum-length tmp1) - (%bignum-length tmp2))) - (if (> (bignum-buffer-integer-length u u-len) - (+ #.(truncate sb!vm:n-word-bits 4) - (bignum-buffer-integer-length v v-len))) - (setf u-len (dmod u u-len - v v-len - tmp1)) - (multiple-value-bind (n d) (reduced-ratio-mod u v) - (setf tmp1-len - (multiply-bignum-buffer-and-smallnum-to-buffer v v-len - n tmp1)) - (setf tmp2-len - (multiply-bignum-buffer-and-smallnum-to-buffer u u-len - d tmp2)) - (gcd-assert (= (copy-bignum tmp2 tmp2-len) - (* (copy-bignum u u-len) d))) - (gcd-assert (= (copy-bignum tmp1 tmp1-len) - (* (copy-bignum v v-len) n))) - (setf u-len - (subtract-bignum-buffers-with-len tmp1 tmp1-len - tmp2 tmp2-len - u - (1+ (max tmp1-len - tmp2-len)))) - (gcd-assert (or (zerop (- (copy-bignum tmp1 tmp1-len) - (copy-bignum tmp2 tmp2-len))) - (= (copy-bignum u u-len) - (- (copy-bignum tmp1 tmp1-len) - (copy-bignum tmp2 tmp2-len))))) - (bignum-abs-buffer u u-len) - (gcd-assert (zerop (modularly u))))) - (setf u-len (make-gcd-bignum-odd u u-len)) - (rotatef u v) - (rotatef u-len v-len)) + (not v-len) + (zerop v-len) + (and (= 1 v-len) + (zerop (%bignum-ref v 0)))) + do + (gcd-assert (= buffer-len (%bignum-length u) + (%bignum-length v) + (%bignum-length tmp1) + (%bignum-length tmp2))) + (if (> (bignum-buffer-integer-length u u-len) + (+ #.(truncate sb!vm:n-word-bits 4) + (bignum-buffer-integer-length v v-len))) + (setf u-len (dmod u u-len + v v-len + tmp1)) + (multiple-value-bind (n d) (reduced-ratio-mod u v) + (setf tmp1-len + (multiply-bignum-buffer-and-smallnum-to-buffer v v-len + n tmp1)) + (setf tmp2-len + (multiply-bignum-buffer-and-smallnum-to-buffer u u-len + d tmp2)) + (gcd-assert (= (copy-bignum tmp2 tmp2-len) + (* (copy-bignum u u-len) d))) + (gcd-assert (= (copy-bignum tmp1 tmp1-len) + (* (copy-bignum v v-len) n))) + (setf u-len + (subtract-bignum-buffers-with-len tmp1 tmp1-len + tmp2 tmp2-len + u + (1+ (max tmp1-len + tmp2-len)))) + (gcd-assert (or (zerop (- (copy-bignum tmp1 tmp1-len) + (copy-bignum tmp2 tmp2-len))) + (= (copy-bignum u u-len) + (- (copy-bignum tmp1 tmp1-len) + (copy-bignum tmp2 tmp2-len))))) + (bignum-abs-buffer u u-len) + (gcd-assert (zerop (modularly u))))) + (setf u-len (make-gcd-bignum-odd u u-len)) + (rotatef u v) + (rotatef u-len v-len)) + (bignum-abs-buffer u u-len) (setf u (copy-bignum u u-len)) - (let ((n (bignum-mod-gcd v1 u))) - (ash (bignum-mod-gcd u1 (if (fixnump n) - (make-small-bignum n) - n)) - factors-of-two))))) + (let ((n (bignum-mod-gcd v1 u))) + (ash (bignum-mod-gcd u1 (if (fixnump n) + (make-small-bignum n) + n)) + factors-of-two))))) (defun bignum-mod-gcd (a b) (declare (type bignum-type a b)) @@ -801,13 +794,13 @@ ;; A and B pretty quickly). After that, use the binary GCD ;; algorithm to handle the rest. (loop until (and (= (%bignum-length b) 1) (zerop (%bignum-ref b 0))) do - (when (<= (%bignum-length a) (1+ (%bignum-length b))) - (return-from bignum-mod-gcd (bignum-binary-gcd a b))) - (let ((rem (mod a b))) - (if (fixnump rem) - (setf a (make-small-bignum rem)) - (setf a rem)) - (rotatef a b))) + (when (<= (%bignum-length a) (1+ (%bignum-length b))) + (return-from bignum-mod-gcd (bignum-binary-gcd a b))) + (let ((rem (mod a b))) + (if (fixnump rem) + (setf a (make-small-bignum rem)) + (setf a rem)) + (rotatef a b))) (if (= (%bignum-length a) 1) (%normalize-bignum a 1) a)) @@ -815,84 +808,84 @@ (defun bignum-binary-gcd (a b) (declare (type bignum-type a b)) (let* ((len-a (%bignum-length a)) - (len-b (%bignum-length b))) + (len-b (%bignum-length b))) (declare (type bignum-index len-a len-b)) (with-bignum-buffers ((a-buffer len-a a) - (b-buffer len-b b) - (res-buffer (max len-a len-b))) + (b-buffer len-b b) + (res-buffer (max len-a len-b))) (let* ((factors-of-two - (bignum-factors-of-two a-buffer len-a - b-buffer len-b)) - (len-a (make-gcd-bignum-odd - a-buffer - (bignum-buffer-ashift-right a-buffer len-a - factors-of-two))) - (len-b (make-gcd-bignum-odd - b-buffer - (bignum-buffer-ashift-right b-buffer len-b - factors-of-two)))) - (declare (type bignum-index len-a len-b)) - (let ((x a-buffer) - (len-x len-a) - (y b-buffer) - (len-y len-b) - (z res-buffer)) - (loop - (multiple-value-bind (u v len-v r len-r) - (bignum-gcd-order-and-subtract x len-x y len-y z) - (declare (type bignum-index len-v len-r)) - (when (and (= len-r 1) (zerop (%bignum-ref r 0))) - (if (zerop factors-of-two) - (let ((ret (%allocate-bignum len-v))) - (dotimes (i len-v) - (setf (%bignum-ref ret i) (%bignum-ref v i))) - (return (%normalize-bignum ret len-v))) - (return (bignum-ashift-left v factors-of-two len-v)))) - (setf x v len-x len-v) - (setf y r len-y (make-gcd-bignum-odd r len-r)) - (setf z u)))))))) + (bignum-factors-of-two a-buffer len-a + b-buffer len-b)) + (len-a (make-gcd-bignum-odd + a-buffer + (bignum-buffer-ashift-right a-buffer len-a + factors-of-two))) + (len-b (make-gcd-bignum-odd + b-buffer + (bignum-buffer-ashift-right b-buffer len-b + factors-of-two)))) + (declare (type bignum-index len-a len-b)) + (let ((x a-buffer) + (len-x len-a) + (y b-buffer) + (len-y len-b) + (z res-buffer)) + (loop + (multiple-value-bind (u v len-v r len-r) + (bignum-gcd-order-and-subtract x len-x y len-y z) + (declare (type bignum-index len-v len-r)) + (when (and (= len-r 1) (zerop (%bignum-ref r 0))) + (if (zerop factors-of-two) + (let ((ret (%allocate-bignum len-v))) + (dotimes (i len-v) + (setf (%bignum-ref ret i) (%bignum-ref v i))) + (return (%normalize-bignum ret len-v))) + (return (bignum-ashift-left v factors-of-two len-v)))) + (setf x v len-x len-v) + (setf y r len-y (make-gcd-bignum-odd r len-r)) + (setf z u)))))))) (defun bignum-gcd-order-and-subtract (a len-a b len-b res) (declare (type bignum-index len-a len-b) (type bignum-type a b)) (cond ((= len-a len-b) - (do ((i (1- len-a) (1- i))) - ((= i -1) - (setf (%bignum-ref res 0) 0) - (values a b len-b res 1)) - (let ((a-digit (%bignum-ref a i)) - (b-digit (%bignum-ref b i))) - (cond ((%digit-compare a-digit b-digit)) - ((%digit-greater a-digit b-digit) - (return - (values a b len-b res - (subtract-bignum-buffers a len-a b len-b - res)))) - (t - (return - (values b a len-a res - (subtract-bignum-buffers b len-b - a len-a - res)))))))) - ((> len-a len-b) - (values a b len-b res - (subtract-bignum-buffers a len-a b len-b res))) - (t - (values b a len-a res - (subtract-bignum-buffers b len-b a len-a res))))) + (do ((i (1- len-a) (1- i))) + ((= i -1) + (setf (%bignum-ref res 0) 0) + (values a b len-b res 1)) + (let ((a-digit (%bignum-ref a i)) + (b-digit (%bignum-ref b i))) + (cond ((%digit-compare a-digit b-digit)) + ((%digit-greater a-digit b-digit) + (return + (values a b len-b res + (subtract-bignum-buffers a len-a b len-b + res)))) + (t + (return + (values b a len-a res + (subtract-bignum-buffers b len-b + a len-a + res)))))))) + ((> len-a len-b) + (values a b len-b res + (subtract-bignum-buffers a len-a b len-b res))) + (t + (values b a len-a res + (subtract-bignum-buffers b len-b a len-a res))))) (defun make-gcd-bignum-odd (a len-a) (declare (type bignum-type a) (type bignum-index len-a)) (dotimes (index len-a) (declare (type bignum-index index)) (do ((digit (%bignum-ref a index) (%ashr digit 1)) - (increment 0 (1+ increment))) - ((zerop digit)) + (increment 0 (1+ increment))) + ((zerop digit)) (declare (type (mod #.sb!vm:n-word-bits) increment)) (when (oddp digit) - (return-from make-gcd-bignum-odd - (bignum-buffer-ashift-right a len-a - (+ (* index digit-size) - increment))))))) + (return-from make-gcd-bignum-odd + (bignum-buffer-ashift-right a len-a + (+ (* index digit-size) + increment))))))) ;;;; negation @@ -901,34 +894,30 @@ ;;; This negates bignum-len digits of bignum, storing the resulting digits into ;;; result (possibly EQ to bignum) and returning whatever end-carry there is. -(sb!xc:defmacro bignum-negate-loop (bignum - bignum-len - &optional (result nil resultp)) - (let ((carry (gensym)) - (end (gensym)) - (value (gensym)) - (last (gensym))) +(sb!xc:defmacro bignum-negate-loop + (bignum bignum-len &optional (result nil resultp)) + (with-unique-names (carry end value last) `(let* (,@(if (not resultp) `(,last)) - (,carry - (multiple-value-bind (,value ,carry) - (%add-with-carry (%lognot (%bignum-ref ,bignum 0)) 1 0) - ,(if resultp - `(setf (%bignum-ref ,result 0) ,value) - `(setf ,last ,value)) - ,carry)) - (i 1) - (,end ,bignum-len)) + (,carry + (multiple-value-bind (,value ,carry) + (%add-with-carry (%lognot (%bignum-ref ,bignum 0)) 1 0) + ,(if resultp + `(setf (%bignum-ref ,result 0) ,value) + `(setf ,last ,value)) + ,carry)) + (i 1) + (,end ,bignum-len)) (declare (type bit ,carry) - (type bignum-index i ,end)) + (type bignum-index i ,end)) (loop - (when (= i ,end) (return)) - (multiple-value-bind (,value temp) - (%add-with-carry (%lognot (%bignum-ref ,bignum i)) 0 ,carry) - ,(if resultp - `(setf (%bignum-ref ,result i) ,value) - `(setf ,last ,value)) - (setf ,carry temp)) - (incf i)) + (when (= i ,end) (return)) + (multiple-value-bind (,value temp) + (%add-with-carry (%lognot (%bignum-ref ,bignum i)) 0 ,carry) + ,(if resultp + `(setf (%bignum-ref ,result i) ,value) + `(setf ,last ,value)) + (setf ,carry temp)) + (incf i)) ,(if resultp carry `(values ,carry ,last))))) ) ; EVAL-WHEN @@ -938,15 +927,15 @@ (defun negate-bignum (x &optional (fully-normalize t)) (declare (type bignum-type x)) (let* ((len-x (%bignum-length x)) - (len-res (1+ len-x)) - (res (%allocate-bignum len-res))) + (len-res (1+ len-x)) + (res (%allocate-bignum len-res))) (declare (type bignum-index len-x len-res)) ;Test len-res for range? (let ((carry (bignum-negate-loop x len-x res))) (setf (%bignum-ref res len-x) - (%add-with-carry (%lognot (%sign-digit x len-x)) 0 carry))) + (%add-with-carry (%lognot (%sign-digit x len-x)) 0 carry))) (if fully-normalize - (%normalize-bignum res len-res) - (%mostly-normalize-bignum res len-res)))) + (%normalize-bignum res len-res) + (%mostly-normalize-bignum res len-res)))) ;;; This assumes bignum is positive; that is, the result of negating it will ;;; stay in the provided allocated bignum. @@ -960,7 +949,7 @@ (defun bignum-abs-buffer (bignum len) (unless (%bignum-0-or-plusp bignum len) - (negate-bignum-in-place bignum len))) + (negate-bignum-buffer-in-place bignum len))) ;;;; shifting @@ -979,26 +968,25 @@ ;;; digit from high bits of the i'th source digit and the start-pos number of ;;; bits from the i+1'th source digit. (sb!xc:defmacro shift-right-unaligned (source - start-digit - start-pos - res-len-form - termination - &optional result) + start-digit + start-pos + res-len-form + termination + &optional result) `(let* ((high-bits-in-first-digit (- digit-size ,start-pos)) - (res-len ,res-len-form) - (res-len-1 (1- res-len)) - ,@(if result `((,result (%allocate-bignum res-len))))) + (res-len ,res-len-form) + (res-len-1 (1- res-len)) + ,@(if result `((,result (%allocate-bignum res-len))))) (declare (type bignum-index res-len res-len-1)) - (do ((i ,start-digit i+1) - (i+1 (1+ ,start-digit) (1+ i+1)) - (j 0 (1+ j))) - ,termination - (declare (type bignum-index i i+1 j)) + (do ((i ,start-digit (1+ i)) + (j 0 (1+ j))) + ,termination + (declare (type bignum-index i j)) (setf (%bignum-ref ,(if result result source) j) - (%logior (%digit-logical-shift-right (%bignum-ref ,source i) - ,start-pos) - (%ashl (%bignum-ref ,source i+1) - high-bits-in-first-digit)))))) + (%logior (%digit-logical-shift-right (%bignum-ref ,source i) + ,start-pos) + (%ashl (%bignum-ref ,source (1+ i)) + high-bits-in-first-digit)))))) ) ; EVAL-WHEN @@ -1012,40 +1000,40 @@ ;;; locals established by the macro. (defun bignum-ashift-right (bignum count) (declare (type bignum-type bignum) - (type unsigned-byte count)) + (type unsigned-byte count)) (let ((bignum-len (%bignum-length bignum))) (declare (type bignum-index bignum-len)) (cond ((fixnump count) - (multiple-value-bind (digits n-bits) (truncate count digit-size) - (declare (type bignum-index digits)) - (cond - ((>= digits bignum-len) - (if (%bignum-0-or-plusp bignum bignum-len) 0 -1)) - ((zerop n-bits) - (bignum-ashift-right-digits bignum digits)) - (t - (shift-right-unaligned bignum digits n-bits (- bignum-len digits) - ((= j res-len-1) - (setf (%bignum-ref res j) - (%ashr (%bignum-ref bignum i) n-bits)) - (%normalize-bignum res res-len)) - res))))) - ((> count bignum-len) - (if (%bignum-0-or-plusp bignum bignum-len) 0 -1)) - ;; Since a FIXNUM should be big enough to address anything in - ;; memory, including arrays of bits, and since arrays of bits - ;; take up about the same space as corresponding fixnums, there - ;; should be no way that we fall through to this case: any shift - ;; right by a bignum should give zero. But let's check anyway: - (t (error "bignum overflow: can't shift right by ~S" count))))) + (multiple-value-bind (digits n-bits) (truncate count digit-size) + (declare (type bignum-index digits)) + (cond + ((>= digits bignum-len) + (if (%bignum-0-or-plusp bignum bignum-len) 0 -1)) + ((zerop n-bits) + (bignum-ashift-right-digits bignum digits)) + (t + (shift-right-unaligned bignum digits n-bits (- bignum-len digits) + ((= j res-len-1) + (setf (%bignum-ref res j) + (%ashr (%bignum-ref bignum i) n-bits)) + (%normalize-bignum res res-len)) + res))))) + ((> count bignum-len) + (if (%bignum-0-or-plusp bignum bignum-len) 0 -1)) + ;; Since a FIXNUM should be big enough to address anything in + ;; memory, including arrays of bits, and since arrays of bits + ;; take up about the same space as corresponding fixnums, there + ;; should be no way that we fall through to this case: any shift + ;; right by a bignum should give zero. But let's check anyway: + (t (error "bignum overflow: can't shift right by ~S" count))))) (defun bignum-ashift-right-digits (bignum digits) (declare (type bignum-type bignum) - (type bignum-index digits)) + (type bignum-index digits)) (let* ((res-len (- (%bignum-length bignum) digits)) - (res (%allocate-bignum res-len))) + (res (%allocate-bignum res-len))) (declare (type bignum-index res-len) - (type bignum-type res)) + (type bignum-type res)) (bignum-replace res bignum :start2 digits) (%normalize-bignum res res-len))) @@ -1063,15 +1051,15 @@ (cond ((zerop n-bits) (let ((new-end (- bignum-len digits))) - (bignum-replace bignum bignum :end1 new-end :start2 digits - :end2 bignum-len) - (%normalize-bignum-buffer bignum new-end))) + (bignum-replace bignum bignum :end1 new-end :start2 digits + :end2 bignum-len) + (%normalize-bignum-buffer bignum new-end))) (t (shift-right-unaligned bignum digits n-bits (- bignum-len digits) - ((= j res-len-1) - (setf (%bignum-ref bignum j) - (%ashr (%bignum-ref bignum i) n-bits)) - (%normalize-bignum-buffer bignum res-len))))))) + ((= j res-len-1) + (setf (%bignum-ref bignum j) + (%ashr (%bignum-ref bignum i) n-bits)) + (%normalize-bignum-buffer bignum res-len))))))) ;;; This handles shifting a bignum buffer to provide fresh bignum data for some ;;; internal routines. We know bignum is safe when called with bignum-len. @@ -1081,17 +1069,17 @@ ;;; branch handles the general case. (defun bignum-ashift-left (bignum x &optional bignum-len) (declare (type bignum-type bignum) - (type unsigned-byte x) - (type (or null bignum-index) bignum-len)) + (type unsigned-byte x) + (type (or null bignum-index) bignum-len)) (if (fixnump x) (multiple-value-bind (digits n-bits) (truncate x digit-size) (let* ((bignum-len (or bignum-len (%bignum-length bignum))) - (res-len (+ digits bignum-len 1))) - (when (> res-len maximum-bignum-length) - (error "can't represent result of left shift")) - (if (zerop n-bits) - (bignum-ashift-left-digits bignum bignum-len digits) - (bignum-ashift-left-unaligned bignum digits n-bits res-len)))) + (res-len (+ digits bignum-len 1))) + (when (> res-len maximum-bignum-length) + (error "can't represent result of left shift")) + (if (zerop n-bits) + (bignum-ashift-left-digits bignum bignum-len digits) + (bignum-ashift-left-unaligned bignum digits n-bits res-len)))) ;; Left shift by a number too big to be represented as a fixnum ;; would exceed our memory capacity, since a fixnum is big enough ;; to index any array, including a bit array. @@ -1100,10 +1088,10 @@ (defun bignum-ashift-left-digits (bignum bignum-len digits) (declare (type bignum-index bignum-len digits)) (let* ((res-len (+ bignum-len digits)) - (res (%allocate-bignum res-len))) + (res (%allocate-bignum res-len))) (declare (type bignum-index res-len)) (bignum-replace res bignum :start1 digits :end1 res-len :end2 bignum-len - :from-end t) + :from-end t) res)) ;;; BIGNUM-TRUNCATE uses this to store into a bignum buffer by supplying res. @@ -1116,29 +1104,28 @@ ;;; first non-zero result digit, digits. We also grab some left over high ;;; bits from the last digit of bignum. (defun bignum-ashift-left-unaligned (bignum digits n-bits res-len - &optional (res nil resp)) + &optional (res nil resp)) (declare (type bignum-index digits res-len) - (type (mod #.digit-size) n-bits)) + (type (mod #.digit-size) n-bits)) (let* ((remaining-bits (- digit-size n-bits)) - (res-len-1 (1- res-len)) - (res (or res (%allocate-bignum res-len)))) + (res-len-1 (1- res-len)) + (res (or res (%allocate-bignum res-len)))) (declare (type bignum-index res-len res-len-1)) - (do ((i 0 i+1) - (i+1 1 (1+ i+1)) - (j (1+ digits) (1+ j))) - ((= j res-len-1) - (setf (%bignum-ref res digits) - (%ashl (%bignum-ref bignum 0) n-bits)) - (setf (%bignum-ref res j) - (%ashr (%bignum-ref bignum i) remaining-bits)) - (if resp - (%normalize-bignum-buffer res res-len) - (%normalize-bignum res res-len))) - (declare (type bignum-index i i+1 j)) + (do ((i 0 (1+ i)) + (j (1+ digits) (1+ j))) + ((= j res-len-1) + (setf (%bignum-ref res digits) + (%ashl (%bignum-ref bignum 0) n-bits)) + (setf (%bignum-ref res j) + (%ashr (%bignum-ref bignum i) remaining-bits)) + (if resp + (%normalize-bignum-buffer res res-len) + (%normalize-bignum res res-len))) + (declare (type bignum-index i j)) (setf (%bignum-ref res j) - (%logior (%digit-logical-shift-right (%bignum-ref bignum i) - remaining-bits) - (%ashl (%bignum-ref bignum i+1) n-bits)))))) + (%logior (%digit-logical-shift-right (%bignum-ref bignum i) + remaining-bits) + (%ashl (%bignum-ref bignum (1+ i)) n-bits)))))) ;;;; relational operators @@ -1153,27 +1140,27 @@ (defun bignum-compare (a b) (declare (type bignum-type a b)) (let* ((len-a (%bignum-length a)) - (len-b (%bignum-length b)) - (a-plusp (%bignum-0-or-plusp a len-a)) - (b-plusp (%bignum-0-or-plusp b len-b))) + (len-b (%bignum-length b)) + (a-plusp (%bignum-0-or-plusp a len-a)) + (b-plusp (%bignum-0-or-plusp b len-b))) (declare (type bignum-index len-a len-b)) (cond ((not (eq a-plusp b-plusp)) - (if a-plusp 1 -1)) - ((= len-a len-b) - (do ((i (1- len-a) (1- i))) - (()) - (declare (type bignum-index i)) - (let ((a-digit (%bignum-ref a i)) - (b-digit (%bignum-ref b i))) - (declare (type bignum-element-type a-digit b-digit)) - (when (%digit-greater a-digit b-digit) - (return 1)) - (when (%digit-greater b-digit a-digit) - (return -1))) - (when (zerop i) (return 0)))) - ((> len-a len-b) - (if a-plusp 1 -1)) - (t (if a-plusp -1 1))))) + (if a-plusp 1 -1)) + ((= len-a len-b) + (do ((i (1- len-a) (1- i))) + (()) + (declare (type bignum-index i)) + (let ((a-digit (%bignum-ref a i)) + (b-digit (%bignum-ref b i))) + (declare (type bignum-element-type a-digit b-digit)) + (when (%digit-greater a-digit b-digit) + (return 1)) + (when (%digit-greater b-digit a-digit) + (return -1))) + (when (zerop i) (return 0)))) + ((> len-a len-b) + (if a-plusp 1 -1)) + (t (if a-plusp -1 1))))) ;;;; float conversion @@ -1183,28 +1170,28 @@ (declare (fixnum exp)) (declare (optimize #-sb-xc-host (sb!ext:inhibit-warnings 3))) (let ((res (dpb exp - sb!vm:single-float-exponent-byte - (logandc2 (logand #xffffffff - (%bignum-ref bits 1)) - sb!vm:single-float-hidden-bit)))) + sb!vm:single-float-exponent-byte + (logandc2 (logand #xffffffff + (%bignum-ref bits 1)) + sb!vm:single-float-hidden-bit)))) (make-single-float (if plusp - res - (logior res (ash -1 sb!vm:float-sign-shift)))))) + res + (logior res (ash -1 sb!vm:float-sign-shift)))))) (defun double-float-from-bits (bits exp plusp) (declare (fixnum exp)) (declare (optimize #-sb-xc-host (sb!ext:inhibit-warnings 3))) (let ((hi (dpb exp - sb!vm:double-float-exponent-byte - (logandc2 (ecase sb!vm::n-word-bits - (32 (%bignum-ref bits 2)) - (64 (ash (%bignum-ref bits 1) -32))) - sb!vm:double-float-hidden-bit))) - (lo (logand #xffffffff (%bignum-ref bits 1)))) + sb!vm:double-float-exponent-byte + (logandc2 (ecase sb!vm::n-word-bits + (32 (%bignum-ref bits 2)) + (64 (ash (%bignum-ref bits 1) -32))) + sb!vm:double-float-hidden-bit))) + (lo (logand #xffffffff (%bignum-ref bits 1)))) (make-double-float (if plusp - hi - (logior hi (ash -1 sb!vm:float-sign-shift))) - lo))) + hi + (logior hi (ash -1 sb!vm:float-sign-shift))) + lo))) #!+(and long-float x86) (defun long-float-from-bits (bits exp plusp) (declare (fixnum exp)) @@ -1220,74 +1207,74 @@ ;;; approximation. (defun bignum-to-float (bignum format) (let* ((plusp (bignum-plus-p bignum)) - (x (if plusp bignum (negate-bignum bignum))) - (len (bignum-integer-length x)) - (digits (float-format-digits format)) - (keep (+ digits digit-size)) - (shift (- keep len)) - (shifted (if (minusp shift) - (bignum-ashift-right x (- shift)) - (bignum-ashift-left x shift))) - (low (%bignum-ref shifted 0)) - (round-bit (ash 1 (1- digit-size)))) + (x (if plusp bignum (negate-bignum bignum))) + (len (bignum-integer-length x)) + (digits (float-format-digits format)) + (keep (+ digits digit-size)) + (shift (- keep len)) + (shifted (if (minusp shift) + (bignum-ashift-right x (- shift)) + (bignum-ashift-left x shift))) + (low (%bignum-ref shifted 0)) + (round-bit (ash 1 (1- digit-size)))) (declare (type bignum-index len digits keep) (fixnum shift)) (labels ((round-up () - (let ((rounded (add-bignums shifted round-bit))) - (if (> (integer-length rounded) keep) - (float-from-bits (bignum-ashift-right rounded 1) - (1+ len)) - (float-from-bits rounded len)))) - (float-from-bits (bits len) - (declare (type bignum-index len)) - (ecase format - (single-float - (single-float-from-bits - bits - (check-exponent len sb!vm:single-float-bias - sb!vm:single-float-normal-exponent-max) - plusp)) - (double-float - (double-float-from-bits - bits - (check-exponent len sb!vm:double-float-bias - sb!vm:double-float-normal-exponent-max) - plusp)) - #!+long-float - (long-float - (long-float-from-bits - bits - (check-exponent len sb!vm:long-float-bias - sb!vm:long-float-normal-exponent-max) - plusp)))) - (check-exponent (exp bias max) - (declare (type bignum-index len)) - (let ((exp (+ exp bias))) - (when (> exp max) - ;; Why a SIMPLE-TYPE-ERROR? Well, this is mainly - ;; called by COERCE, which requires an error of - ;; TYPE-ERROR if the conversion can't happen - ;; (except in certain circumstances when we are - ;; coercing to a FUNCTION) -- CSR, 2002-09-18 - (error 'simple-type-error - :format-control "Too large to be represented as a ~S:~% ~S" - :format-arguments (list format x) - :expected-type format - :datum x)) - exp))) + (let ((rounded (add-bignums shifted round-bit))) + (if (> (integer-length rounded) keep) + (float-from-bits (bignum-ashift-right rounded 1) + (1+ len)) + (float-from-bits rounded len)))) + (float-from-bits (bits len) + (declare (type bignum-index len)) + (ecase format + (single-float + (single-float-from-bits + bits + (check-exponent len sb!vm:single-float-bias + sb!vm:single-float-normal-exponent-max) + plusp)) + (double-float + (double-float-from-bits + bits + (check-exponent len sb!vm:double-float-bias + sb!vm:double-float-normal-exponent-max) + plusp)) + #!+long-float + (long-float + (long-float-from-bits + bits + (check-exponent len sb!vm:long-float-bias + sb!vm:long-float-normal-exponent-max) + plusp)))) + (check-exponent (exp bias max) + (declare (type bignum-index len)) + (let ((exp (+ exp bias))) + (when (> exp max) + ;; Why a SIMPLE-TYPE-ERROR? Well, this is mainly + ;; called by COERCE, which requires an error of + ;; TYPE-ERROR if the conversion can't happen + ;; (except in certain circumstances when we are + ;; coercing to a FUNCTION) -- CSR, 2002-09-18 + (error 'simple-type-error + :format-control "Too large to be represented as a ~S:~% ~S" + :format-arguments (list format x) + :expected-type format + :datum x)) + exp))) (cond ;; Round down if round bit is 0. - ((zerop (logand round-bit low)) + ((not (logtest round-bit low)) (float-from-bits shifted len)) ;; If only round bit is set, then round to even. ((and (= low round-bit) - (dotimes (i (- (%bignum-length x) (ceiling keep digit-size)) - t) - (unless (zerop (%bignum-ref x i)) (return nil)))) + (dotimes (i (- (%bignum-length x) (ceiling keep digit-size)) + t) + (unless (zerop (%bignum-ref x i)) (return nil)))) (let ((next (%bignum-ref shifted 1))) - (if (oddp next) - (round-up) - (float-from-bits shifted len)))) + (if (oddp next) + (round-up) + (float-from-bits shifted len)))) ;; Otherwise, round up. (t (round-up)))))) @@ -1297,9 +1284,9 @@ (defun bignum-buffer-integer-length (bignum len) (declare (type bignum-type bignum)) (let* ((len-1 (1- len)) - (digit (%bignum-ref bignum len-1))) + (digit (%bignum-ref bignum len-1))) (declare (type bignum-index len len-1) - (type bignum-element-type digit)) + (type bignum-element-type digit)) (+ (integer-length (%fixnum-digit-with-correct-sign digit)) (* len-1 digit-size)))) @@ -1312,24 +1299,25 @@ (let ((len (%bignum-length bignum))) (declare (type bignum-index len)) (multiple-value-bind (word-index bit-index) - (floor index digit-size) + (floor index digit-size) (if (>= word-index len) - (not (bignum-plus-p bignum)) - (not (zerop (logand (%bignum-ref bignum word-index) - (ash 1 bit-index)))))))) + (not (bignum-plus-p bignum)) + (logbitp bit-index (%bignum-ref bignum word-index)))))) (defun bignum-logcount (bignum) (declare (type bignum-type bignum)) - (let* ((length (%bignum-length bignum)) - (plusp (%bignum-0-or-plusp bignum length)) - (result 0)) + (let ((length (%bignum-length bignum)) + (result 0)) (declare (type bignum-index length) - (fixnum result)) + (fixnum result)) (do ((index 0 (1+ index))) - ((= index length) result) + ((= index length) + (if (%bignum-0-or-plusp bignum length) + result + (- (* length digit-size) result))) (let ((digit (%bignum-ref bignum index))) - (declare (type bignum-element-type digit)) - (incf result (logcount (if plusp digit (%lognot digit)))))))) + (declare (type bignum-element-type digit)) + (incf result (logcount digit)))))) ;;;; logical operations @@ -1338,7 +1326,7 @@ (defun bignum-logical-not (a) (declare (type bignum-type a)) (let* ((len (%bignum-length a)) - (res (%allocate-bignum len))) + (res (%allocate-bignum len))) (declare (type bignum-index len)) (dotimes (i len res) (declare (type bignum-index i)) @@ -1349,19 +1337,19 @@ (defun bignum-logical-and (a b) (declare (type bignum-type a b)) (let* ((len-a (%bignum-length a)) - (len-b (%bignum-length b)) - (a-plusp (%bignum-0-or-plusp a len-a)) - (b-plusp (%bignum-0-or-plusp b len-b))) + (len-b (%bignum-length b)) + (a-plusp (%bignum-0-or-plusp a len-a)) + (b-plusp (%bignum-0-or-plusp b len-b))) (declare (type bignum-index len-a len-b)) (cond ((< len-a len-b) (if a-plusp - (logand-shorter-positive a len-a b (%allocate-bignum len-a)) - (logand-shorter-negative a len-a b len-b (%allocate-bignum len-b)))) + (logand-shorter-positive a len-a b (%allocate-bignum len-a)) + (logand-shorter-negative a len-a b len-b (%allocate-bignum len-b)))) ((< len-b len-a) (if b-plusp - (logand-shorter-positive b len-b a (%allocate-bignum len-b)) - (logand-shorter-negative b len-b a len-a (%allocate-bignum len-a)))) + (logand-shorter-positive b len-b a (%allocate-bignum len-b)) + (logand-shorter-negative b len-b a len-a (%allocate-bignum len-a)))) (t (logand-shorter-positive a len-a b (%allocate-bignum len-a)))))) ;;; This takes a shorter bignum, a and len-a, that is positive. Because this @@ -1369,11 +1357,11 @@ ;;; sign bits will mask the other bits out of b. The result is len-a big. (defun logand-shorter-positive (a len-a b res) (declare (type bignum-type a b res) - (type bignum-index len-a)) + (type bignum-index len-a)) (dotimes (i len-a) (declare (type bignum-index i)) (setf (%bignum-ref res i) - (%logand (%bignum-ref a i) (%bignum-ref b i)))) + (%logand (%bignum-ref a i) (%bignum-ref b i)))) (%normalize-bignum res len-a)) ;;; This takes a shorter bignum, a and len-a, that is negative. Because this @@ -1381,11 +1369,11 @@ ;;; bits will include any bits from b. The result is len-b big. (defun logand-shorter-negative (a len-a b len-b res) (declare (type bignum-type a b res) - (type bignum-index len-a len-b)) + (type bignum-index len-a len-b)) (dotimes (i len-a) (declare (type bignum-index i)) (setf (%bignum-ref res i) - (%logand (%bignum-ref a i) (%bignum-ref b i)))) + (%logand (%bignum-ref a i) (%bignum-ref b i)))) (do ((i len-a (1+ i))) ((= i len-b)) (declare (type bignum-index i)) @@ -1397,19 +1385,19 @@ (defun bignum-logical-ior (a b) (declare (type bignum-type a b)) (let* ((len-a (%bignum-length a)) - (len-b (%bignum-length b)) - (a-plusp (%bignum-0-or-plusp a len-a)) - (b-plusp (%bignum-0-or-plusp b len-b))) + (len-b (%bignum-length b)) + (a-plusp (%bignum-0-or-plusp a len-a)) + (b-plusp (%bignum-0-or-plusp b len-b))) (declare (type bignum-index len-a len-b)) (cond ((< len-a len-b) (if a-plusp - (logior-shorter-positive a len-a b len-b (%allocate-bignum len-b)) - (logior-shorter-negative a len-a b len-b (%allocate-bignum len-b)))) + (logior-shorter-positive a len-a b len-b (%allocate-bignum len-b)) + (logior-shorter-negative a len-a b len-b (%allocate-bignum len-b)))) ((< len-b len-a) (if b-plusp - (logior-shorter-positive b len-b a len-a (%allocate-bignum len-a)) - (logior-shorter-negative b len-b a len-a (%allocate-bignum len-a)))) + (logior-shorter-positive b len-b a len-a (%allocate-bignum len-a)) + (logior-shorter-negative b len-b a len-a (%allocate-bignum len-a)))) (t (logior-shorter-positive a len-a b len-b (%allocate-bignum len-a)))))) ;;; This takes a shorter bignum, a and len-a, that is positive. Because this @@ -1418,11 +1406,11 @@ ;;; is len-b long. (defun logior-shorter-positive (a len-a b len-b res) (declare (type bignum-type a b res) - (type bignum-index len-a len-b)) + (type bignum-index len-a len-b)) (dotimes (i len-a) (declare (type bignum-index i)) (setf (%bignum-ref res i) - (%logior (%bignum-ref a i) (%bignum-ref b i)))) + (%logior (%bignum-ref a i) (%bignum-ref b i)))) (do ((i len-a (1+ i))) ((= i len-b)) (declare (type bignum-index i)) @@ -1434,11 +1422,11 @@ ;;; bits will include any bits from b. The result is len-b long. (defun logior-shorter-negative (a len-a b len-b res) (declare (type bignum-type a b res) - (type bignum-index len-a len-b)) + (type bignum-index len-a len-b)) (dotimes (i len-a) (declare (type bignum-index i)) (setf (%bignum-ref res i) - (%logior (%bignum-ref a i) (%bignum-ref b i)))) + (%logior (%bignum-ref a i) (%bignum-ref b i)))) (do ((i len-a (1+ i)) (sign (%sign-digit a len-a))) ((= i len-b)) @@ -1451,21 +1439,21 @@ (defun bignum-logical-xor (a b) (declare (type bignum-type a b)) (let ((len-a (%bignum-length a)) - (len-b (%bignum-length b))) + (len-b (%bignum-length b))) (declare (type bignum-index len-a len-b)) (if (< len-a len-b) - (bignum-logical-xor-aux a len-a b len-b (%allocate-bignum len-b)) - (bignum-logical-xor-aux b len-b a len-a (%allocate-bignum len-a))))) + (bignum-logical-xor-aux a len-a b len-b (%allocate-bignum len-b)) + (bignum-logical-xor-aux b len-b a len-a (%allocate-bignum len-a))))) ;;; This takes the shorter of two bignums in a and len-a. Res is len-b ;;; long. Do the XOR. (defun bignum-logical-xor-aux (a len-a b len-b res) (declare (type bignum-type a b res) - (type bignum-index len-a len-b)) + (type bignum-index len-a len-b)) (dotimes (i len-a) (declare (type bignum-index i)) (setf (%bignum-ref res i) - (%logxor (%bignum-ref a i) (%bignum-ref b i)))) + (%logxor (%bignum-ref a i) (%bignum-ref b i)))) (do ((i len-a (1+ i)) (sign (%sign-digit a len-a))) ((= i len-b)) @@ -1473,463 +1461,9 @@ (setf (%bignum-ref res i) (%logxor sign (%bignum-ref b i)))) (%normalize-bignum res len-b)) -;;;; LDB (load byte) - -#| -FOR NOW WE DON'T USE LDB OR DPB. WE USE SHIFTS AND MASKS IN NUMBERS.LISP WHICH -IS LESS EFFICIENT BUT EASIER TO MAINTAIN. BILL SAYS THIS CODE CERTAINLY WORKS! - -(defconstant maximum-fixnum-bits (- sb!vm:n-word-bits sb!vm:n-lowtag-bits)) - -(defun bignum-load-byte (byte bignum) - (declare (type bignum-type bignum)) - (let ((byte-len (byte-size byte)) - (byte-pos (byte-position byte))) - (if (< byte-len maximum-fixnum-bits) - (bignum-ldb-fixnum-res bignum byte-len byte-pos) - (bignum-ldb-bignum-res bignum byte-len byte-pos)))) - -;;; This returns a fixnum result of loading a byte from a bignum. In order, we -;;; check for the following conditions: -;;; Insufficient bignum digits to start loading a byte -- -;;; Return 0 or byte-len 1's depending on sign of bignum. -;;; One bignum digit containing the whole byte spec -- -;;; Grab 'em, shift 'em, and mask out what we don't want. -;;; Insufficient bignum digits to cover crossing a digit boundary -- -;;; Grab the available bits in the last digit, and or in whatever -;;; virtual sign bits we need to return a full byte spec. -;;; Else (we cross a digit boundary with all bits available) -- -;;; Make a couple masks, grab what we want, shift it around, and -;;; LOGIOR it all together. -;;; Because (< maximum-fixnum-bits digit-size) and -;;; (< byte-len maximum-fixnum-bits), -;;; we only cross one digit boundary if any. -(defun bignum-ldb-fixnum-res (bignum byte-len byte-pos) - (multiple-value-bind (skipped-digits pos) (truncate byte-pos digit-size) - (let ((bignum-len (%bignum-length bignum)) - (s-digits+1 (1+ skipped-digits))) - (declare (type bignum-index bignum-len s-digits+1)) - (if (>= skipped-digits bignum-len) - (if (%bignum-0-or-plusp bignum bignum-len) - 0 - (%make-ones byte-len)) - (let ((end (+ pos byte-len))) - (cond ((<= end digit-size) - (logand (ash (%bignum-ref bignum skipped-digits) (- pos)) - ;; Must LOGAND after shift here. - (%make-ones byte-len))) - ((>= s-digits+1 bignum-len) - (let* ((available-bits (- digit-size pos)) - (res (logand (ash (%bignum-ref bignum skipped-digits) - (- pos)) - ;; LOGAND should be unnecessary here - ;; with a logical right shift or a - ;; correct digit-sized one. - (%make-ones available-bits)))) - (if (%bignum-0-or-plusp bignum bignum-len) - res - (logior (%ashl (%make-ones (- end digit-size)) - available-bits) - res)))) - (t - (let* ((high-bits-in-first-digit (- digit-size pos)) - (high-mask (%make-ones high-bits-in-first-digit)) - (low-bits-in-next-digit (- end digit-size)) - (low-mask (%make-ones low-bits-in-next-digit))) - (declare (type bignum-element-type high-mask low-mask)) - (logior (%ashl (logand (%bignum-ref bignum s-digits+1) - low-mask) - high-bits-in-first-digit) - (logand (ash (%bignum-ref bignum skipped-digits) - (- pos)) - ;; LOGAND should be unnecessary here with - ;; a logical right shift or a correct - ;; digit-sized one. - high-mask)))))))))) - -;;; This returns a bignum result of loading a byte from a bignum. In order, we -;;; check for the following conditions: -;;; Insufficient bignum digits to start loading a byte -- -;;; Byte-pos starting on a digit boundary -- -;;; Byte spec contained in one bignum digit -- -;;; Grab the bits we want and stick them in a single digit result. -;;; Since we know byte-pos is non-zero here, we know our single digit -;;; will have a zero high sign bit. -;;; Else (unaligned multiple digits) -- -;;; This is like doing a shift right combined with either masking -;;; out unwanted high bits from bignum or filling in virtual sign -;;; bits if bignum had insufficient bits. We use SHIFT-RIGHT-ALIGNED -;;; and reference lots of local variables this macro establishes. -(defun bignum-ldb-bignum-res (bignum byte-len byte-pos) - (multiple-value-bind (skipped-digits pos) (truncate byte-pos digit-size) - (let ((bignum-len (%bignum-length bignum))) - (declare (type bignum-index bignum-len)) - (cond - ((>= skipped-digits bignum-len) - (make-bignum-virtual-ldb-bits bignum bignum-len byte-len)) - ((zerop pos) - (make-aligned-ldb-bignum bignum bignum-len byte-len skipped-digits)) - ((< (+ pos byte-len) digit-size) - (let ((res (%allocate-bignum 1))) - (setf (%bignum-ref res 0) - (logand (%ashr (%bignum-ref bignum skipped-digits) pos) - (%make-ones byte-len))) - res)) - (t - (make-unaligned-ldb-bignum bignum bignum-len - byte-len skipped-digits pos)))))) - -;;; This returns bits from bignum that don't physically exist. These are -;;; all zero or one depending on the sign of the bignum. -(defun make-bignum-virtual-ldb-bits (bignum bignum-len byte-len) - (if (%bignum-0-or-plusp bignum bignum-len) - 0 - (multiple-value-bind (res-len-1 extra) (truncate byte-len digit-size) - (declare (type bignum-index res-len-1)) - (let* ((res-len (1+ res-len-1)) - (res (%allocate-bignum res-len))) - (declare (type bignum-index res-len)) - (do ((j 0 (1+ j))) - ((= j res-len-1) - (setf (%bignum-ref res j) (%make-ones extra)) - (%normalize-bignum res res-len)) - (declare (type bignum-index j)) - (setf (%bignum-ref res j) all-ones-digit)))))) - -;;; Since we are picking up aligned digits, we just copy the whole digits -;;; we want and fill in extra bits. We might have a byte-len that extends -;;; off the end of the bignum, so we may have to fill in extra 1's if the -;;; bignum is negative. -(defun make-aligned-ldb-bignum (bignum bignum-len byte-len skipped-digits) - (multiple-value-bind (res-len-1 extra) (truncate byte-len digit-size) - (declare (type bignum-index res-len-1)) - (let* ((res-len (1+ res-len-1)) - (res (%allocate-bignum res-len))) - (declare (type bignum-index res-len)) - (do ((i skipped-digits (1+ i)) - (j 0 (1+ j))) - ((or (= j res-len-1) (= i bignum-len)) - (cond ((< i bignum-len) - (setf (%bignum-ref res j) - (logand (%bignum-ref bignum i) - (the bignum-element-type (%make-ones extra))))) - ((%bignum-0-or-plusp bignum bignum-len)) - (t - (do ((j j (1+ j))) - ((= j res-len-1) - (setf (%bignum-ref res j) (%make-ones extra))) - (setf (%bignum-ref res j) all-ones-digit)))) - (%normalize-bignum res res-len)) - (declare (type bignum-index i j)) - (setf (%bignum-ref res j) (%bignum-ref bignum i)))))) - -;;; This grabs unaligned bignum bits from bignum assuming byte-len causes at -;;; least one digit boundary crossing. We use SHIFT-RIGHT-UNALIGNED referencing -;;; lots of local variables established by it. -(defun make-unaligned-ldb-bignum (bignum - bignum-len - byte-len - skipped-digits - pos) - (multiple-value-bind (res-len-1 extra) (truncate byte-len digit-size) - (shift-right-unaligned - bignum skipped-digits pos (1+ res-len-1) - ((or (= j res-len-1) (= i+1 bignum-len)) - (cond ((= j res-len-1) - (cond - ((< extra high-bits-in-first-digit) - (setf (%bignum-ref res j) - (logand (ash (%bignum-ref bignum i) minus-start-pos) - ;; Must LOGAND after shift here. - (%make-ones extra)))) - (t - (setf (%bignum-ref res j) - (logand (ash (%bignum-ref bignum i) minus-start-pos) - ;; LOGAND should be unnecessary here with a logical - ;; right shift or a correct digit-sized one. - high-mask)) - (when (%bignum-0-or-plusp bignum bignum-len) - (setf (%bignum-ref res j) - (logior (%bignum-ref res j) - (%ashl (%make-ones - (- extra high-bits-in-first-digit)) - high-bits-in-first-digit))))))) - (t - (setf (%bignum-ref res j) - (logand (ash (%bignum-ref bignum i) minus-start-pos) - ;; LOGAND should be unnecessary here with a logical - ;; right shift or a correct digit-sized one. - high-mask)) - (unless (%bignum-0-or-plusp bignum bignum-len) - ;; Fill in upper half of this result digit with 1's. - (setf (%bignum-ref res j) - (logior (%bignum-ref res j) - (%ashl low-mask high-bits-in-first-digit))) - ;; Fill in any extra 1's we need to be byte-len long. - (do ((j (1+ j) (1+ j))) - ((>= j res-len-1) - (setf (%bignum-ref res j) (%make-ones extra))) - (setf (%bignum-ref res j) all-ones-digit))))) - (%normalize-bignum res res-len)) - res))) - -;;;; DPB (deposit byte) - -(defun bignum-deposit-byte (new-byte byte-spec bignum) - (declare (type bignum-type bignum)) - (let* ((byte-len (byte-size byte-spec)) - (byte-pos (byte-position byte-spec)) - (bignum-len (%bignum-length bignum)) - (bignum-plusp (%bignum-0-or-plusp bignum bignum-len)) - (byte-end (+ byte-pos byte-len)) - (res-len (1+ (max (ceiling byte-end digit-size) bignum-len))) - (res (%allocate-bignum res-len))) - (declare (type bignum-index bignum-len res-len)) - ;; Fill in an extra sign digit in case we set what would otherwise be the - ;; last digit's last bit. Normalize at the end in case this was - ;; unnecessary. - (unless bignum-plusp - (setf (%bignum-ref res (1- res-len)) all-ones-digit)) - (multiple-value-bind (end-digit end-bits) (truncate byte-end digit-size) - (declare (type bignum-index end-digit)) - ;; Fill in bits from bignum up to byte-pos. - (multiple-value-bind (pos-digit pos-bits) (truncate byte-pos digit-size) - (declare (type bignum-index pos-digit)) - (do ((i 0 (1+ i)) - (end (min pos-digit bignum-len))) - ((= i end) - (cond ((< i bignum-len) - (unless (zerop pos-bits) - (setf (%bignum-ref res i) - (logand (%bignum-ref bignum i) - (%make-ones pos-bits))))) - (bignum-plusp) - (t - (do ((i i (1+ i))) - ((= i pos-digit) - (unless (zerop pos-bits) - (setf (%bignum-ref res i) (%make-ones pos-bits)))) - (setf (%bignum-ref res i) all-ones-digit))))) - (setf (%bignum-ref res i) (%bignum-ref bignum i))) - ;; Fill in bits from new-byte. - (if (typep new-byte 'fixnum) - (deposit-fixnum-bits new-byte byte-len pos-digit pos-bits - end-digit end-bits res) - (deposit-bignum-bits new-byte byte-len pos-digit pos-bits - end-digit end-bits res))) - ;; Fill in remaining bits from bignum after byte-spec. - (when (< end-digit bignum-len) - (setf (%bignum-ref res end-digit) - (logior (logand (%bignum-ref bignum end-digit) - (%ashl (%make-ones (- digit-size end-bits)) - end-bits)) - ;; DEPOSIT-FIXNUM-BITS and DEPOSIT-BIGNUM-BITS only store - ;; bits from new-byte into res's end-digit element, so - ;; we don't need to mask out unwanted high bits. - (%bignum-ref res end-digit))) - (do ((i (1+ end-digit) (1+ i))) - ((= i bignum-len)) - (setf (%bignum-ref res i) (%bignum-ref bignum i))))) - (%normalize-bignum res res-len))) - -;;; This starts at result's pos-digit skipping pos-bits, and it stores bits -;;; from new-byte, a fixnum, into result. It effectively stores byte-len -;;; number of bits, but never stores past end-digit and end-bits in result. -;;; The first branch fires when all the bits we want from new-byte are present; -;;; if byte-len crosses from the current result digit into the next, the last -;;; argument to DEPOSIT-FIXNUM-DIGIT is a mask for those bits. The second -;;; branch handles the need to grab more bits than the fixnum new-byte has, but -;;; new-byte is positive; therefore, any virtual bits are zero. The mask for -;;; bits that don't fit in the current result digit is simply the remaining -;;; bits in the bignum digit containing new-byte; we don't care if we store -;;; some extra in the next result digit since they will be zeros. The last -;;; branch handles the need to grab more bits than the fixnum new-byte has, but -;;; new-byte is negative; therefore, any virtual bits must be explicitly filled -;;; in as ones. We call DEPOSIT-FIXNUM-DIGIT to grab what bits actually exist -;;; and to fill in the current result digit. -(defun deposit-fixnum-bits (new-byte byte-len pos-digit pos-bits - end-digit end-bits result) - (declare (type bignum-index pos-digit end-digit)) - (let ((other-bits (- digit-size pos-bits)) - (new-byte-digit (%fixnum-to-digit new-byte))) - (declare (type bignum-element-type new-byte-digit)) - (cond ((< byte-len maximum-fixnum-bits) - (deposit-fixnum-digit new-byte-digit byte-len pos-digit pos-bits - other-bits result - (- byte-len other-bits))) - ((or (plusp new-byte) (zerop new-byte)) - (deposit-fixnum-digit new-byte-digit byte-len pos-digit pos-bits - other-bits result pos-bits)) - (t - (multiple-value-bind (digit bits) - (deposit-fixnum-digit new-byte-digit byte-len pos-digit pos-bits - other-bits result - (if (< (- byte-len other-bits) digit-size) - (- byte-len other-bits) - digit-size)) - (declare (type bignum-index digit)) - (cond ((< digit end-digit) - (setf (%bignum-ref result digit) - (logior (%bignum-ref result digit) - (%ashl (%make-ones (- digit-size bits)) bits))) - (do ((i (1+ digit) (1+ i))) - ((= i end-digit) - (setf (%bignum-ref result i) (%make-ones end-bits))) - (setf (%bignum-ref result i) all-ones-digit))) - ((> digit end-digit)) - ((< bits end-bits) - (setf (%bignum-ref result digit) - (logior (%bignum-ref result digit) - (%ashl (%make-ones (- end-bits bits)) - bits)))))))))) - -;;; This fills in the current result digit from new-byte-digit. The first case -;;; handles everything we want fitting in the current digit, and other-bits is -;;; the number of bits remaining to be filled in result's current digit. This -;;; number is digit-size minus pos-bits. The second branch handles filling in -;;; result's current digit, and it shoves the unused bits of new-byte-digit -;;; into the next result digit. This is correct regardless of new-byte-digit's -;;; sign. It returns the new current result digit and how many bits already -;;; filled in the result digit. -(defun deposit-fixnum-digit (new-byte-digit byte-len pos-digit pos-bits - other-bits result next-digit-bits-needed) - (declare (type bignum-index pos-digit) - (type bignum-element-type new-byte-digit next-digit-mask)) - (cond ((<= byte-len other-bits) - ;; Bits from new-byte fit in the current result digit. - (setf (%bignum-ref result pos-digit) - (logior (%bignum-ref result pos-digit) - (%ashl (logand new-byte-digit (%make-ones byte-len)) - pos-bits))) - (if (= byte-len other-bits) - (values (1+ pos-digit) 0) - (values pos-digit (+ byte-len pos-bits)))) - (t - ;; Some of new-byte's bits go in current result digit. - (setf (%bignum-ref result pos-digit) - (logior (%bignum-ref result pos-digit) - (%ashl (logand new-byte-digit (%make-ones other-bits)) - pos-bits))) - (let ((pos-digit+1 (1+ pos-digit))) - ;; The rest of new-byte's bits go in the next result digit. - (setf (%bignum-ref result pos-digit+1) - (logand (ash new-byte-digit (- other-bits)) - ;; Must LOGAND after shift here. - (%make-ones next-digit-bits-needed))) - (if (= next-digit-bits-needed digit-size) - (values (1+ pos-digit+1) 0) - (values pos-digit+1 next-digit-bits-needed)))))) - -;;; This starts at result's pos-digit skipping pos-bits, and it stores bits -;;; from new-byte, a bignum, into result. It effectively stores byte-len -;;; number of bits, but never stores past end-digit and end-bits in result. -;;; When handling a starting bit unaligned with a digit boundary, we check -;;; in the second branch for the byte spec fitting into the pos-digit element -;;; after after pos-bits; DEPOSIT-UNALIGNED-BIGNUM-BITS expects at least one -;;; digit boundary crossing. -(defun deposit-bignum-bits (bignum-byte byte-len pos-digit pos-bits - end-digit end-bits result) - (declare (type bignum-index pos-digit end-digit)) - (cond ((zerop pos-bits) - (deposit-aligned-bignum-bits bignum-byte pos-digit end-digit end-bits - result)) - ((or (= end-digit pos-digit) - (and (= end-digit (1+ pos-digit)) - (zerop end-bits))) - (setf (%bignum-ref result pos-digit) - (logior (%bignum-ref result pos-digit) - (%ashl (logand (%bignum-ref bignum-byte 0) - (%make-ones byte-len)) - pos-bits)))) - (t (deposit-unaligned-bignum-bits bignum-byte pos-digit pos-bits - end-digit end-bits result)))) - -;;; This deposits bits from bignum-byte into result starting at pos-digit and -;;; the zero'th bit. It effectively only stores bits to end-bits in the -;;; end-digit element of result. The loop termination code takes care of -;;; picking up the last digit's bits or filling in virtual negative sign bits. -(defun deposit-aligned-bignum-bits (bignum-byte pos-digit end-digit end-bits - result) - (declare (type bignum-index pos-digit end-digit)) - (let* ((bignum-len (%bignum-length bignum-byte)) - (bignum-plusp (%bignum-0-or-plusp bignum-byte bignum-len))) - (declare (type bignum-index bignum-len)) - (do ((i 0 (1+ i )) - (j pos-digit (1+ j))) - ((or (= j end-digit) (= i bignum-len)) - (cond ((= j end-digit) - (cond ((< i bignum-len) - (setf (%bignum-ref result j) - (logand (%bignum-ref bignum-byte i) - (%make-ones end-bits)))) - (bignum-plusp) - (t - (setf (%bignum-ref result j) (%make-ones end-bits))))) - (bignum-plusp) - (t - (do ((j j (1+ j))) - ((= j end-digit) - (setf (%bignum-ref result j) (%make-ones end-bits))) - (setf (%bignum-ref result j) all-ones-digit))))) - (setf (%bignum-ref result j) (%bignum-ref bignum-byte i))))) - -;;; This assumes at least one digit crossing. -(defun deposit-unaligned-bignum-bits (bignum-byte pos-digit pos-bits - end-digit end-bits result) - (declare (type bignum-index pos-digit end-digit)) - (let* ((bignum-len (%bignum-length bignum-byte)) - (bignum-plusp (%bignum-0-or-plusp bignum-byte bignum-len)) - (low-mask (%make-ones pos-bits)) - (bits-past-pos-bits (- digit-size pos-bits)) - (high-mask (%make-ones bits-past-pos-bits)) - (minus-high-bits (- bits-past-pos-bits))) - (declare (type bignum-element-type low-mask high-mask) - (type bignum-index bignum-len)) - (do ((i 0 (1+ i)) - (j pos-digit j+1) - (j+1 (1+ pos-digit) (1+ j+1))) - ((or (= j end-digit) (= i bignum-len)) - (cond - ((= j end-digit) - (setf (%bignum-ref result j) - (cond - ((>= pos-bits end-bits) - (logand (%bignum-ref result j) (%make-ones end-bits))) - ((< i bignum-len) - (logior (%bignum-ref result j) - (%ashl (logand (%bignum-ref bignum-byte i) - (%make-ones (- end-bits pos-bits))) - pos-bits))) - (bignum-plusp - (logand (%bignum-ref result j) - ;; 0's between pos-bits and end-bits positions. - (logior (%ashl (%make-ones (- digit-size end-bits)) - end-bits) - low-mask))) - (t (logior (%bignum-ref result j) - (%ashl (%make-ones (- end-bits pos-bits)) - pos-bits)))))) - (bignum-plusp) - (t - (setf (%bignum-ref result j) - (%ashl (%make-ones bits-past-pos-bits) pos-bits)) - (do ((j j+1 (1+ j))) - ((= j end-digit) - (setf (%bignum-ref result j) (%make-ones end-bits))) - (declare (type bignum-index j)) - (setf (%bignum-ref result j) all-ones-digit))))) - (declare (type bignum-index i j j+1)) - (let ((digit (%bignum-ref bignum-byte i))) - (declare (type bignum-element-type digit)) - (setf (%bignum-ref result j) - (logior (%bignum-ref result j) - (%ashl (logand digit high-mask) pos-bits))) - (setf (%bignum-ref result j+1) - (logand (ash digit minus-high-bits) - ;; LOGAND should be unnecessary here with a logical right - ;; shift or a correct digit-sized one. - low-mask)))))) -|# +;;;; There used to be a bunch of code to implement "efficient" versions of LDB +;;;; and DPB here. But it apparently was never used, so it's been deleted. +;;;; --njf, 2007-02-04 ;;;; TRUNCATE @@ -2017,228 +1551,252 @@ IS LESS EFFICIENT BUT EASIER TO MAINTAIN. BILL SAYS THIS CODE CERTAINLY WORKS! (defun bignum-truncate (x y) (declare (type bignum-type x y)) (let (truncate-x truncate-y) - (labels + (labels ;;; Divide X by Y when Y is a single bignum digit. BIGNUM-TRUNCATE ;;; fixes up the quotient and remainder with respect to sign and ;;; normalization. - ;;; - ;;; We don't have to worry about shifting Y to make its most - ;;; significant digit sufficiently large for %FLOOR to return - ;;; digit-size quantities for the q-digit and r-digit. If Y is - ;;; a single digit bignum, it is already large enough for - ;;; %FLOOR. That is, it has some bits on pretty high in the - ;;; digit. - ((bignum-truncate-single-digit (x len-x y) - (declare (type bignum-index len-x)) - (let ((q (%allocate-bignum len-x)) - (r 0) - (y (%bignum-ref y 0))) - (declare (type bignum-element-type r y)) - (do ((i (1- len-x) (1- i))) - ((minusp i)) - (multiple-value-bind (q-digit r-digit) - (%floor r (%bignum-ref x i) y) - (declare (type bignum-element-type q-digit r-digit)) - (setf (%bignum-ref q i) q-digit) - (setf r r-digit))) - (let ((rem (%allocate-bignum 1))) - (setf (%bignum-ref rem 0) r) - (values q rem)))) - ;;; This returns a guess for the next division step. Y1 is the - ;;; highest y digit, and y2 is the second to highest y - ;;; digit. The x... variables are the three highest x digits - ;;; for the next division step. - ;;; - ;;; From Knuth, our guess is either all ones or x-i and x-i-1 - ;;; divided by y1, depending on whether x-i and y1 are the - ;;; same. We test this guess by determining whether guess*y2 - ;;; is greater than the three high digits of x minus guess*y1 - ;;; shifted left one digit: - ;;; ------------------------------ - ;;; | x-i | x-i-1 | x-i-2 | - ;;; ------------------------------ - ;;; ------------------------------ - ;;; - | g*y1 high | g*y1 low | 0 | - ;;; ------------------------------ - ;;; ... < guess*y2 ??? - ;;; If guess*y2 is greater, then we decrement our guess by one - ;;; and try again. This returns a guess that is either - ;;; correct or one too large. - (bignum-truncate-guess (y1 y2 x-i x-i-1 x-i-2) - (declare (type bignum-element-type y1 y2 x-i x-i-1 x-i-2)) - (let ((guess (if (%digit-compare x-i y1) - all-ones-digit - (%floor x-i x-i-1 y1)))) - (declare (type bignum-element-type guess)) - (loop - (multiple-value-bind (high-guess*y1 low-guess*y1) - (%multiply guess y1) - (declare (type bignum-element-type low-guess*y1 - high-guess*y1)) - (multiple-value-bind (high-guess*y2 low-guess*y2) - (%multiply guess y2) - (declare (type bignum-element-type high-guess*y2 - low-guess*y2)) - (multiple-value-bind (middle-digit borrow) - (%subtract-with-borrow x-i-1 low-guess*y1 1) - (declare (type bignum-element-type middle-digit) - (fixnum borrow)) - ;; Supplying borrow of 1 means there was no - ;; borrow, and we know x-i-2 minus 0 requires - ;; no borrow. - (let ((high-digit (%subtract-with-borrow x-i - high-guess*y1 - borrow))) - (declare (type bignum-element-type high-digit)) - (if (and (%digit-compare high-digit 0) - (or (%digit-greater high-guess*y2 - middle-digit) - (and (%digit-compare middle-digit - high-guess*y2) - (%digit-greater low-guess*y2 - x-i-2)))) - (setf guess (%subtract-with-borrow guess 1 1)) - (return guess))))))))) - ;;; Divide TRUNCATE-X by TRUNCATE-Y, returning the quotient - ;;; and destructively modifying TRUNCATE-X so that it holds - ;;; the remainder. - ;;; - ;;; LEN-X and LEN-Y tell us how much of the buffers we care about. - ;;; - ;;; TRUNCATE-X definitely has at least three digits, and it has one - ;;; more than TRUNCATE-Y. This keeps i, i-1, i-2, and low-x-digit - ;;; happy. Thanks to SHIFT-AND-STORE-TRUNCATE-BUFFERS. - (return-quotient-leaving-remainder (len-x len-y) - (declare (type bignum-index len-x len-y)) - (let* ((len-q (- len-x len-y)) - ;; Add one for extra sign digit in case high bit is on. - (q (%allocate-bignum (1+ len-q))) - (k (1- len-q)) - (y1 (%bignum-ref truncate-y (1- len-y))) - (y2 (%bignum-ref truncate-y (- len-y 2))) - (i (1- len-x)) - (i-1 (1- i)) - (i-2 (1- i-1)) - (low-x-digit (- i len-y))) - (declare (type bignum-index len-q k i i-1 i-2 low-x-digit) - (type bignum-element-type y1 y2)) - (loop - (setf (%bignum-ref q k) - (try-bignum-truncate-guess - ;; This modifies TRUNCATE-X. Must access - ;; elements each pass. - (bignum-truncate-guess y1 y2 - (%bignum-ref truncate-x i) - (%bignum-ref truncate-x i-1) - (%bignum-ref truncate-x i-2)) - len-y low-x-digit)) - (cond ((zerop k) (return)) - (t (decf k) - (decf low-x-digit) - (shiftf i i-1 i-2 (1- i-2))))) - q)) - ;;; This takes a digit guess, multiplies it by TRUNCATE-Y for a - ;;; result one greater in length than LEN-Y, and subtracts this result - ;;; from TRUNCATE-X. LOW-X-DIGIT is the first digit of X to start - ;;; the subtraction, and we know X is long enough to subtract a LEN-Y - ;;; plus one length bignum from it. Next we check the result of the - ;;; subtraction, and if the high digit in X became negative, then our - ;;; guess was one too big. In this case, return one less than GUESS - ;;; passed in, and add one value of Y back into X to account for - ;;; subtracting one too many. Knuth shows that the guess is wrong on - ;;; the order of 3/b, where b is the base (2 to the digit-size power) - ;;; -- pretty rarely. - (try-bignum-truncate-guess (guess len-y low-x-digit) - (declare (type bignum-index low-x-digit len-y) - (type bignum-element-type guess)) - (let ((carry-digit 0) - (borrow 1) - (i low-x-digit)) - (declare (type bignum-element-type carry-digit) - (type bignum-index i) - (fixnum borrow)) - ;; Multiply guess and divisor, subtracting from dividend - ;; simultaneously. - (dotimes (j len-y) - (multiple-value-bind (high-digit low-digit) - (%multiply-and-add guess - (%bignum-ref truncate-y j) - carry-digit) - (declare (type bignum-element-type high-digit low-digit)) - (setf carry-digit high-digit) - (multiple-value-bind (x temp-borrow) - (%subtract-with-borrow (%bignum-ref truncate-x i) - low-digit - borrow) - (declare (type bignum-element-type x) - (fixnum temp-borrow)) - (setf (%bignum-ref truncate-x i) x) - (setf borrow temp-borrow))) - (incf i)) - (setf (%bignum-ref truncate-x i) - (%subtract-with-borrow (%bignum-ref truncate-x i) - carry-digit borrow)) - ;; See whether guess is off by one, adding one - ;; Y back in if necessary. - (cond ((%digit-0-or-plusp (%bignum-ref truncate-x i)) - guess) - (t - ;; If subtraction has negative result, add one - ;; divisor value back in. The guess was one too - ;; large in magnitude. - (let ((i low-x-digit) - (carry 0)) - (dotimes (j len-y) - (multiple-value-bind (v k) - (%add-with-carry (%bignum-ref truncate-y j) - (%bignum-ref truncate-x i) - carry) - (declare (type bignum-element-type v)) - (setf (%bignum-ref truncate-x i) v) - (setf carry k)) - (incf i)) - (setf (%bignum-ref truncate-x i) - (%add-with-carry (%bignum-ref truncate-x i) - 0 carry))) - (%subtract-with-borrow guess 1 1))))) - ;;; This returns the amount to shift y to place a one in the - ;;; second highest bit. Y must be positive. If the last digit - ;;; of y is zero, then y has a one in the previous digit's - ;;; sign bit, so we know it will take one less than digit-size - ;;; to get a one where we want. Otherwise, we count how many - ;;; right shifts it takes to get zero; subtracting this value - ;;; from digit-size tells us how many high zeros there are - ;;; which is one more than the shift amount sought. - ;;; - ;;; Note: This is exactly the same as one less than the - ;;; integer-length of the last digit subtracted from the - ;;; digit-size. - ;;; - ;;; We shift y to make it sufficiently large that doing the - ;;; 2*digit-size by digit-size %FLOOR calls ensures the quotient and - ;;; remainder fit in digit-size. - (shift-y-for-truncate (y) - (let* ((len (%bignum-length y)) - (last (%bignum-ref y (1- len)))) - (declare (type bignum-index len) - (type bignum-element-type last)) - (- digit-size (integer-length last) 1))) - ;;; Stores two bignums into the truncation bignum buffers, - ;;; shifting them on the way in. This assumes x and y are - ;;; positive and at least two in length, and it assumes - ;;; truncate-x and truncate-y are one digit longer than x and - ;;; y. - (shift-and-store-truncate-buffers (x len-x y len-y shift) - (declare (type bignum-index len-x len-y) - (type (integer 0 (#.digit-size)) shift)) - (cond ((zerop shift) - (bignum-replace truncate-x x :end1 len-x) - (bignum-replace truncate-y y :end1 len-y)) - (t - (bignum-ashift-left-unaligned x 0 shift (1+ len-x) - truncate-x) - (bignum-ashift-left-unaligned y 0 shift (1+ len-y) - truncate-y))))) ;; LABELS + ;;; + ;;; We don't have to worry about shifting Y to make its most + ;;; significant digit sufficiently large for %BIGFLOOR to return + ;;; digit-size quantities for the q-digit and r-digit. If Y is + ;;; a single digit bignum, it is already large enough for + ;;; %BIGFLOOR. That is, it has some bits on pretty high in the + ;;; digit. + ((bignum-truncate-single-digit (x len-x y) + (declare (type bignum-index len-x)) + (let ((y (%bignum-ref y 0))) + (declare (type bignum-element-type y)) + (if (not (logtest y (1- y))) + ;; Y is a power of two. + ;; SHIFT-RIGHT-UNALIGNED won't do the right thing + ;; with a shift count of 0 or -1, so special case this. + (cond ((= y 0) + (error 'division-by-zero)) + ((= y 1) + ;; We could probably get away with (VALUES X 0) + ;; here, but it's not clear that some of the + ;; normalization logic further down would avoid + ;; mutilating X. Just go ahead and cons, consing's + ;; cheap. + (values (copy-bignum x len-x) 0)) + (t + (let ((n-bits (1- (integer-length y)))) + (values + (shift-right-unaligned x 0 n-bits len-x + ((= j res-len-1) + (setf (%bignum-ref res j) + (%ashr (%bignum-ref x i) n-bits)) + res) + res) + (logand (%bignum-ref x 0) (1- y)))))) + (do ((i (1- len-x) (1- i)) + (q (%allocate-bignum len-x)) + (r 0)) + ((minusp i) + (let ((rem (%allocate-bignum 1))) + (setf (%bignum-ref rem 0) r) + (values q rem))) + (declare (type bignum-element-type r)) + (multiple-value-bind (q-digit r-digit) + (%bigfloor r (%bignum-ref x i) y) + (declare (type bignum-element-type q-digit r-digit)) + (setf (%bignum-ref q i) q-digit) + (setf r r-digit)))))) + ;;; This returns a guess for the next division step. Y1 is the + ;;; highest y digit, and y2 is the second to highest y + ;;; digit. The x... variables are the three highest x digits + ;;; for the next division step. + ;;; + ;;; From Knuth, our guess is either all ones or x-i and x-i-1 + ;;; divided by y1, depending on whether x-i and y1 are the + ;;; same. We test this guess by determining whether guess*y2 + ;;; is greater than the three high digits of x minus guess*y1 + ;;; shifted left one digit: + ;;; ------------------------------ + ;;; | x-i | x-i-1 | x-i-2 | + ;;; ------------------------------ + ;;; ------------------------------ + ;;; - | g*y1 high | g*y1 low | 0 | + ;;; ------------------------------ + ;;; ... < guess*y2 ??? + ;;; If guess*y2 is greater, then we decrement our guess by one + ;;; and try again. This returns a guess that is either + ;;; correct or one too large. + (bignum-truncate-guess (y1 y2 x-i x-i-1 x-i-2) + (declare (type bignum-element-type y1 y2 x-i x-i-1 x-i-2)) + (let ((guess (if (%digit-compare x-i y1) + all-ones-digit + (%bigfloor x-i x-i-1 y1)))) + (declare (type bignum-element-type guess)) + (loop + (multiple-value-bind (high-guess*y1 low-guess*y1) + (%multiply guess y1) + (declare (type bignum-element-type low-guess*y1 + high-guess*y1)) + (multiple-value-bind (high-guess*y2 low-guess*y2) + (%multiply guess y2) + (declare (type bignum-element-type high-guess*y2 + low-guess*y2)) + (multiple-value-bind (middle-digit borrow) + (%subtract-with-borrow x-i-1 low-guess*y1 1) + (declare (type bignum-element-type middle-digit) + (fixnum borrow)) + ;; Supplying borrow of 1 means there was no + ;; borrow, and we know x-i-2 minus 0 requires + ;; no borrow. + (let ((high-digit (%subtract-with-borrow x-i + high-guess*y1 + borrow))) + (declare (type bignum-element-type high-digit)) + (if (and (%digit-compare high-digit 0) + (or (%digit-greater high-guess*y2 + middle-digit) + (and (%digit-compare middle-digit + high-guess*y2) + (%digit-greater low-guess*y2 + x-i-2)))) + (setf guess (%subtract-with-borrow guess 1 1)) + (return guess))))))))) + ;;; Divide TRUNCATE-X by TRUNCATE-Y, returning the quotient + ;;; and destructively modifying TRUNCATE-X so that it holds + ;;; the remainder. + ;;; + ;;; LEN-X and LEN-Y tell us how much of the buffers we care about. + ;;; + ;;; TRUNCATE-X definitely has at least three digits, and it has one + ;;; more than TRUNCATE-Y. This keeps i, i-1, i-2, and low-x-digit + ;;; happy. Thanks to SHIFT-AND-STORE-TRUNCATE-BUFFERS. + (return-quotient-leaving-remainder (len-x len-y) + (declare (type bignum-index len-x len-y)) + (let* ((len-q (- len-x len-y)) + ;; Add one for extra sign digit in case high bit is on. + (q (%allocate-bignum (1+ len-q))) + (k (1- len-q)) + (y1 (%bignum-ref truncate-y (1- len-y))) + (y2 (%bignum-ref truncate-y (- len-y 2))) + (i (1- len-x)) + (i-1 (1- i)) + (i-2 (1- i-1)) + (low-x-digit (- i len-y))) + (declare (type bignum-index len-q k i i-1 i-2 low-x-digit) + (type bignum-element-type y1 y2)) + (loop + (setf (%bignum-ref q k) + (try-bignum-truncate-guess + ;; This modifies TRUNCATE-X. Must access + ;; elements each pass. + (bignum-truncate-guess y1 y2 + (%bignum-ref truncate-x i) + (%bignum-ref truncate-x i-1) + (%bignum-ref truncate-x i-2)) + len-y low-x-digit)) + (cond ((zerop k) (return)) + (t (decf k) + (decf low-x-digit) + (shiftf i i-1 i-2 (1- i-2))))) + q)) + ;;; This takes a digit guess, multiplies it by TRUNCATE-Y for a + ;;; result one greater in length than LEN-Y, and subtracts this result + ;;; from TRUNCATE-X. LOW-X-DIGIT is the first digit of X to start + ;;; the subtraction, and we know X is long enough to subtract a LEN-Y + ;;; plus one length bignum from it. Next we check the result of the + ;;; subtraction, and if the high digit in X became negative, then our + ;;; guess was one too big. In this case, return one less than GUESS + ;;; passed in, and add one value of Y back into X to account for + ;;; subtracting one too many. Knuth shows that the guess is wrong on + ;;; the order of 3/b, where b is the base (2 to the digit-size power) + ;;; -- pretty rarely. + (try-bignum-truncate-guess (guess len-y low-x-digit) + (declare (type bignum-index low-x-digit len-y) + (type bignum-element-type guess)) + (let ((carry-digit 0) + (borrow 1) + (i low-x-digit)) + (declare (type bignum-element-type carry-digit) + (type bignum-index i) + (fixnum borrow)) + ;; Multiply guess and divisor, subtracting from dividend + ;; simultaneously. + (dotimes (j len-y) + (multiple-value-bind (high-digit low-digit) + (%multiply-and-add guess + (%bignum-ref truncate-y j) + carry-digit) + (declare (type bignum-element-type high-digit low-digit)) + (setf carry-digit high-digit) + (multiple-value-bind (x temp-borrow) + (%subtract-with-borrow (%bignum-ref truncate-x i) + low-digit + borrow) + (declare (type bignum-element-type x) + (fixnum temp-borrow)) + (setf (%bignum-ref truncate-x i) x) + (setf borrow temp-borrow))) + (incf i)) + (setf (%bignum-ref truncate-x i) + (%subtract-with-borrow (%bignum-ref truncate-x i) + carry-digit borrow)) + ;; See whether guess is off by one, adding one + ;; Y back in if necessary. + (cond ((%digit-0-or-plusp (%bignum-ref truncate-x i)) + guess) + (t + ;; If subtraction has negative result, add one + ;; divisor value back in. The guess was one too + ;; large in magnitude. + (let ((i low-x-digit) + (carry 0)) + (dotimes (j len-y) + (multiple-value-bind (v k) + (%add-with-carry (%bignum-ref truncate-y j) + (%bignum-ref truncate-x i) + carry) + (declare (type bignum-element-type v)) + (setf (%bignum-ref truncate-x i) v) + (setf carry k)) + (incf i)) + (setf (%bignum-ref truncate-x i) + (%add-with-carry (%bignum-ref truncate-x i) + 0 carry))) + (%subtract-with-borrow guess 1 1))))) + ;;; This returns the amount to shift y to place a one in the + ;;; second highest bit. Y must be positive. If the last digit + ;;; of y is zero, then y has a one in the previous digit's + ;;; sign bit, so we know it will take one less than digit-size + ;;; to get a one where we want. Otherwise, we count how many + ;;; right shifts it takes to get zero; subtracting this value + ;;; from digit-size tells us how many high zeros there are + ;;; which is one more than the shift amount sought. + ;;; + ;;; Note: This is exactly the same as one less than the + ;;; integer-length of the last digit subtracted from the + ;;; digit-size. + ;;; + ;;; We shift y to make it sufficiently large that doing the + ;;; 2*digit-size by digit-size %BIGFLOOR calls ensures the quotient and + ;;; remainder fit in digit-size. + (shift-y-for-truncate (y) + (let* ((len (%bignum-length y)) + (last (%bignum-ref y (1- len)))) + (declare (type bignum-index len) + (type bignum-element-type last)) + (- digit-size (integer-length last) 1))) + ;;; Stores two bignums into the truncation bignum buffers, + ;;; shifting them on the way in. This assumes x and y are + ;;; positive and at least two in length, and it assumes + ;;; truncate-x and truncate-y are one digit longer than x and + ;;; y. + (shift-and-store-truncate-buffers (x len-x y len-y shift) + (declare (type bignum-index len-x len-y) + (type (integer 0 (#.digit-size)) shift)) + (cond ((zerop shift) + (bignum-replace truncate-x x :end1 len-x) + (bignum-replace truncate-y y :end1 len-y)) + (t + (bignum-ashift-left-unaligned x 0 shift (1+ len-x) + truncate-x) + (bignum-ashift-left-unaligned y 0 shift (1+ len-y) + truncate-y))))) ;; LABELS ;;; Divide X by Y returning the quotient and remainder. In the ;;; general case, we shift Y to set up for the algorithm, and we ;;; use two buffers to save consing intermediate values. X gets @@ -2247,300 +1805,65 @@ IS LESS EFFICIENT BUT EASIER TO MAINTAIN. BILL SAYS THIS CODE CERTAINLY WORKS! ;;; multiple bind q and r, we first fix up the signs and then ;;; return the normalized results. (let* ((x-plusp (%bignum-0-or-plusp x (%bignum-length x))) - (y-plusp (%bignum-0-or-plusp y (%bignum-length y))) - (x (if x-plusp x (negate-bignum x nil))) - (y (if y-plusp y (negate-bignum y nil))) - (len-x (%bignum-length x)) - (len-y (%bignum-length y))) - (multiple-value-bind (q r) - (cond ((< len-y 2) - (bignum-truncate-single-digit x len-x y)) - ((plusp (bignum-compare y x)) - (let ((res (%allocate-bignum len-x))) - (dotimes (i len-x) - (setf (%bignum-ref res i) (%bignum-ref x i))) - (values 0 res))) - (t - (let ((len-x+1 (1+ len-x))) - (setf truncate-x (%allocate-bignum len-x+1)) - (setf truncate-y (%allocate-bignum (1+ len-y))) - (let ((y-shift (shift-y-for-truncate y))) - (shift-and-store-truncate-buffers x len-x y - len-y y-shift) - (values (return-quotient-leaving-remainder len-x+1 - len-y) - ;; Now that RETURN-QUOTIENT-LEAVING-REMAINDER - ;; has executed, we just tidy up the remainder - ;; (in TRUNCATE-X) and return it. - (cond - ((zerop y-shift) - (let ((res (%allocate-bignum len-y))) - (declare (type bignum-type res)) - (bignum-replace res truncate-x :end2 len-y) - (%normalize-bignum res len-y))) - (t - (shift-right-unaligned - truncate-x 0 y-shift len-y - ((= j res-len-1) - (setf (%bignum-ref res j) - (%ashr (%bignum-ref truncate-x i) - y-shift)) - (%normalize-bignum res res-len)) - res)))))))) - (let ((quotient (cond ((eq x-plusp y-plusp) q) - ((typep q 'fixnum) (the fixnum (- q))) - (t (negate-bignum-in-place q)))) - (rem (cond (x-plusp r) - ((typep r 'fixnum) (the fixnum (- r))) - (t (negate-bignum-in-place r))))) - (values (if (typep quotient 'fixnum) - quotient - (%normalize-bignum quotient (%bignum-length quotient))) - (if (typep rem 'fixnum) - rem - (%normalize-bignum rem (%bignum-length rem)))))))))) + (y-plusp (%bignum-0-or-plusp y (%bignum-length y))) + (x (if x-plusp x (negate-bignum x nil))) + (y (if y-plusp y (negate-bignum y nil))) + (len-x (%bignum-length x)) + (len-y (%bignum-length y))) + (multiple-value-bind (q r) + (cond ((< len-y 2) + (bignum-truncate-single-digit x len-x y)) + ((plusp (bignum-compare y x)) + (let ((res (%allocate-bignum len-x))) + (dotimes (i len-x) + (setf (%bignum-ref res i) (%bignum-ref x i))) + (values 0 res))) + (t + (let ((len-x+1 (1+ len-x))) + (setf truncate-x (%allocate-bignum len-x+1)) + (setf truncate-y (%allocate-bignum (1+ len-y))) + (let ((y-shift (shift-y-for-truncate y))) + (shift-and-store-truncate-buffers x len-x y + len-y y-shift) + (values (return-quotient-leaving-remainder len-x+1 + len-y) + ;; Now that RETURN-QUOTIENT-LEAVING-REMAINDER + ;; has executed, we just tidy up the remainder + ;; (in TRUNCATE-X) and return it. + (cond + ((zerop y-shift) + (let ((res (%allocate-bignum len-y))) + (declare (type bignum-type res)) + (bignum-replace res truncate-x :end2 len-y) + (%normalize-bignum res len-y))) + (t + (shift-right-unaligned + truncate-x 0 y-shift len-y + ((= j res-len-1) + (setf (%bignum-ref res j) + (%ashr (%bignum-ref truncate-x i) + y-shift)) + (%normalize-bignum res res-len)) + res)))))))) + (let ((quotient (cond ((eq x-plusp y-plusp) q) + ((typep q 'fixnum) (the fixnum (- q))) + (t (negate-bignum-in-place q)))) + (rem (cond (x-plusp r) + ((typep r 'fixnum) (the fixnum (- r))) + (t (negate-bignum-in-place r))))) + (values (if (typep quotient 'fixnum) + quotient + (%normalize-bignum quotient (%bignum-length quotient))) + (if (typep rem 'fixnum) + rem + (%normalize-bignum rem (%bignum-length rem)))))))))) -;;;; %FLOOR primitive for BIGNUM-TRUNCATE - -;;; When a machine leaves out a 2*digit-size by digit-size divide -;;; instruction (that is, two bignum-digits divided by one), we have to -;;; roll our own (the hard way). Basically, we treat the operation as -;;; four digit-size/2 digits divided by two digit-size/2 digits. This -;;; means we have duplicated most of the code above to do this nearly -;;; general digit-size/2 digit bignum divide, but we've unrolled loops -;;; and made use of other properties of this specific divide situation. - -;;;; %FLOOR for machines with a 32x32 divider. - -#!-sb-fluid -(declaim (inline 32x16-subtract-with-borrow 32x16-add-with-carry - 32x16-divide 32x16-multiply 32x16-multiply-split)) - -#!+32x16-divide -(defconstant 32x16-base-1 (1- (ash 1 (/ sb!vm:n-word-bits 2)))) - -#!+32x16-divide -(deftype bignum-half-element-type () `(unsigned-byte ,(/ sb!vm:n-word-bits 2))) -#!+32x16-divide -(defconstant half-digit-size (/ digit-size 2)) - -;;; This is similar to %SUBTRACT-WITH-BORROW. It returns a -;;; half-digit-size difference and a borrow. Returning a 1 for the -;;; borrow means there was no borrow, and 0 means there was one. -#!+32x16-divide -(defun 32x16-subtract-with-borrow (a b borrow) - (declare (type bignum-half-element-type a b) - (type (integer 0 1) borrow)) - (let ((diff (+ (- a b) borrow 32x16-base-1))) - (declare (type (unsigned-byte #.(1+ half-digit-size)) diff)) - (values (logand diff (1- (ash 1 half-digit-size))) - (ash diff (- half-digit-size))))) - -;;; This adds a and b, half-digit-size quantities, with the carry k. It -;;; returns a half-digit-size sum and a second value, 0 or 1, indicating -;;; whether there was a carry. -#!+32x16-divide -(defun 32x16-add-with-carry (a b k) - (declare (type bignum-half-element-type a b) - (type (integer 0 1) k)) - (let ((res (the fixnum (+ a b k)))) - (declare (type (unsigned-byte #.(1+ half-digit-size)) res)) - (if (zerop (the fixnum (logand (ash 1 half-digit-size) res))) - (values res 0) - (values (the bignum-half-element-type (logand (1- (ash 1 half-digit-size)) res)) - 1)))) - -;;; This is probably a digit-size by digit-size divide instruction. -#!+32x16-divide -(defun 32x16-divide (a b c) - (declare (type bignum-half-element-type a b c)) - (floor (the bignum-element-type - (logior (the bignum-element-type (ash a 16)) - b)) - c)) - -;;; This basically exists since we know the answer won't overflow -;;; bignum-element-type. It's probably just a basic multiply instruction, but -;;; it can't cons an intermediate bignum. The result goes in a non-descriptor -;;; register. -#!+32x16-divide -(defun 32x16-multiply (a b) - (declare (type bignum-half-element-type a b)) - (the bignum-element-type (* a b))) - -;;; This multiplies a and b, half-digit-size quantities, and returns the -;;; result as two half-digit-size quantities, high and low. -#!+32x16-divide -(defun 32x16-multiply-split (a b) - (let ((res (32x16-multiply a b))) - (declare (the bignum-element-type res)) - (values (the bignum-half-element-type (logand (1- (ash 1 half-digit-size)) (ash res (- half-digit-size)))) - (the bignum-half-element-type (logand (1- (ash 1 half-digit-size)) res))))) - -;;; The %FLOOR below uses this buffer the same way BIGNUM-TRUNCATE uses -;;; *truncate-x*. There's no y buffer since we pass around the two -;;; half-digit-size digits and use them slightly differently than the -;;; general truncation algorithm above. -#!+32x16-divide -(defvar *32x16-truncate-x* (make-array 4 :element-type 'bignum-half-element-type - :initial-element 0)) - -;;; This does the same thing as the %FLOOR above, but it does it at Lisp level -;;; when there is no 64x32-bit divide instruction on the machine. -;;; -;;; It implements the higher level tactics of BIGNUM-TRUNCATE, but it -;;; makes use of special situation provided, four half-digit-size digits -;;; divided by two half-digit-size digits. -#!+32x16-divide -(defun %floor (a b c) - (declare (type bignum-element-type a b c)) - ;; Setup *32x16-truncate-x* buffer from a and b. - (setf (aref *32x16-truncate-x* 0) - (the bignum-half-element-type (logand (1- (ash 1 half-digit-size)) b))) - (setf (aref *32x16-truncate-x* 1) - (the bignum-half-element-type - (logand (1- (ash 1 half-digit-size)) - (the bignum-half-element-type (ash b (- half-digit-size)))))) - (setf (aref *32x16-truncate-x* 2) - (the bignum-half-element-type (logand (1- (ash 1 half-digit-size)) a))) - (setf (aref *32x16-truncate-x* 3) - (the bignum-half-element-type - (logand (1- (ash 1 half-digit-size)) - (the bignum-half-element-type (ash a (- half-digit-size)))))) - ;; From DO-TRUNCATE, but unroll the loop. - (let* ((y1 (logand (1- (ash 1 half-digit-size)) (ash c (- half-digit-size)))) - (y2 (logand (1- (ash 1 half-digit-size)) c)) - (q (the bignum-element-type - (ash (32x16-try-bignum-truncate-guess - (32x16-truncate-guess y1 y2 - (aref *32x16-truncate-x* 3) - (aref *32x16-truncate-x* 2) - (aref *32x16-truncate-x* 1)) - y1 y2 1) - 16)))) - (declare (type bignum-element-type q) - (type bignum-half-element-type y1 y2)) - (values (the bignum-element-type - (logior q - (the bignum-half-element-type - (32x16-try-bignum-truncate-guess - (32x16-truncate-guess - y1 y2 - (aref *32x16-truncate-x* 2) - (aref *32x16-truncate-x* 1) - (aref *32x16-truncate-x* 0)) - y1 y2 0)))) - (the bignum-element-type - (logior (the bignum-element-type - (ash (aref *32x16-truncate-x* 1) 16)) - (the bignum-half-element-type - (aref *32x16-truncate-x* 0))))))) - -;;; This is similar to TRY-BIGNUM-TRUNCATE-GUESS, but this unrolls the two -;;; loops. This also substitutes for %DIGIT-0-OR-PLUSP the equivalent -;;; expression without any embellishment or pretense of abstraction. The first -;;; loop is unrolled, but we've put the body of the loop into the function -;;; 32X16-TRY-GUESS-ONE-RESULT-DIGIT. -#!+32x16-divide -(defun 32x16-try-bignum-truncate-guess (guess y-high y-low low-x-digit) - (declare (type bignum-index low-x-digit) - (type bignum-half-element-type guess y-high y-low)) - (let ((high-x-digit (+ 2 low-x-digit))) - ;; Multiply guess and divisor, subtracting from dividend simultaneously. - (multiple-value-bind (guess*y-hold carry borrow) - (32x16-try-guess-one-result-digit guess y-low 0 0 1 low-x-digit) - (declare (type bignum-half-element-type guess*y-hold) - (fixnum carry borrow)) - (multiple-value-bind (guess*y-hold carry borrow) - (32x16-try-guess-one-result-digit guess y-high guess*y-hold - carry borrow (1+ low-x-digit)) - (declare (type bignum-half-element-type guess*y-hold) - (fixnum borrow) - (ignore carry)) - (setf (aref *32x16-truncate-x* high-x-digit) - (32x16-subtract-with-borrow (aref *32x16-truncate-x* high-x-digit) - guess*y-hold borrow)))) - ;; See whether guess is off by one, adding one Y back in if necessary. - (cond ((zerop (logand (ash 1 (1- half-digit-size)) - (aref *32x16-truncate-x* high-x-digit))) - ;; The subtraction result is zero or positive. - guess) - (t - ;; If subtraction has negative result, add one divisor value back - ;; in. The guess was one too large in magnitude. - (multiple-value-bind (v carry) - (32x16-add-with-carry y-low - (aref *32x16-truncate-x* low-x-digit) - 0) - (declare (type bignum-half-element-type v)) - (setf (aref *32x16-truncate-x* low-x-digit) v) - (multiple-value-bind (v carry) - (32x16-add-with-carry y-high - (aref *32x16-truncate-x* - (1+ low-x-digit)) - carry) - (setf (aref *32x16-truncate-x* (1+ low-x-digit)) v) - (setf (aref *32x16-truncate-x* high-x-digit) - (32x16-add-with-carry (aref *32x16-truncate-x* high-x-digit) - carry 0)))) - (if (zerop (logand (ash 1 (1- half-digit-size)) guess)) - (1- guess) - (1+ guess)))))) - -;;; This is similar to the body of the loop in TRY-BIGNUM-TRUNCATE-GUESS that -;;; multiplies the guess by y and subtracts the result from x simultaneously. -;;; This returns the digit remembered as part of the multiplication, the carry -;;; from additions done on behalf of the multiplication, and the borrow from -;;; doing the subtraction. -#!+32x16-divide -(defun 32x16-try-guess-one-result-digit (guess y-digit guess*y-hold - carry borrow x-index) - (multiple-value-bind (high-digit low-digit) - (32x16-multiply-split guess y-digit) - (declare (type bignum-half-element-type high-digit low-digit)) - (multiple-value-bind (low-digit temp-carry) - (32x16-add-with-carry low-digit guess*y-hold carry) - (declare (type bignum-half-element-type low-digit)) - (multiple-value-bind (high-digit temp-carry) - (32x16-add-with-carry high-digit temp-carry 0) - (declare (type bignum-half-element-type high-digit)) - (multiple-value-bind (x temp-borrow) - (32x16-subtract-with-borrow (aref *32x16-truncate-x* x-index) - low-digit borrow) - (declare (type bignum-half-element-type x)) - (setf (aref *32x16-truncate-x* x-index) x) - (values high-digit temp-carry temp-borrow)))))) - -;;; This is similar to BIGNUM-TRUNCATE-GUESS, but instead of computing -;;; the guess exactly as described in the its comments (digit by digit), -;;; this massages the digit-size/2 quantities into digit-size quantities -;;; and performs the -#!+32x16-divide -(defun 32x16-truncate-guess (y1 y2 x-i x-i-1 x-i-2) - (declare (type bignum-half-element-type y1 y2 x-i x-i-1 x-i-2)) - (let ((guess (if (= x-i y1) - (1- (ash 1 half-digit-size)) - (32x16-divide x-i x-i-1 y1)))) - (declare (type bignum-half-element-type guess)) - (loop - (let* ((guess*y1 (the bignum-element-type - (ash (logand (1- (ash 1 half-digit-size)) - (the bignum-element-type - (32x16-multiply guess y1))) - 16))) - (x-y (%subtract-with-borrow - (the bignum-element-type - (logior (the bignum-element-type - (ash x-i-1 16)) - x-i-2)) - guess*y1 - 1)) - (guess*y2 (the bignum-element-type (%multiply guess y2)))) - (declare (type bignum-element-type guess*y1 x-y guess*y2)) - (if (%digit-greater guess*y2 x-y) - (decf guess) - (return guess)))))) +;;;; There used to be a pile of code for implementing division for bignum digits +;;;; for machines that don't have a 2*digit-size by digit-size divide instruction. +;;;; This happens to be most machines, but all the SBCL ports seem to be content +;;;; to implement SB-BIGNUM:%BIGFLOOR as a VOP rather than using the code here. +;;;; So it's been deleted. --njf, 2007-02-04 ;;;; general utilities @@ -2551,16 +1874,16 @@ IS LESS EFFICIENT BUT EASIER TO MAINTAIN. BILL SAYS THIS CODE CERTAINLY WORKS! #!-sb-fluid (declaim (sb!ext:maybe-inline %normalize-bignum-buffer)) (defun %normalize-bignum-buffer (result len) (declare (type bignum-type result) - (type bignum-index len)) + (type bignum-index len)) (unless (= len 1) (do ((next-digit (%bignum-ref result (- len 2)) - (%bignum-ref result (- len 2))) - (sign-digit (%bignum-ref result (1- len)) next-digit)) - ((not (zerop (logxor sign-digit (%ashr next-digit (1- digit-size)))))) - (decf len) - (setf (%bignum-ref result len) 0) - (when (= len 1) - (return)))) + (%bignum-ref result (- len 2))) + (sign-digit (%bignum-ref result (1- len)) next-digit)) + ((not (zerop (logxor sign-digit (%ashr next-digit (1- digit-size)))))) + (decf len) + (setf (%bignum-ref result len) 0) + (when (= len 1) + (return)))) len) ;;; This drops the last digit if it is unnecessary sign information. It repeats @@ -2572,27 +1895,27 @@ IS LESS EFFICIENT BUT EASIER TO MAINTAIN. BILL SAYS THIS CODE CERTAINLY WORKS! ;;; we do have a fixnum, shift it over for the two low-tag bits. (defun %normalize-bignum (result len) (declare (type bignum-type result) - (type bignum-index len) - (inline %normalize-bignum-buffer)) + (type bignum-index len) + (inline %normalize-bignum-buffer)) (let ((newlen (%normalize-bignum-buffer result len))) (declare (type bignum-index newlen)) (unless (= newlen len) (%bignum-set-length result newlen)) (if (= newlen 1) - (let ((digit (%bignum-ref result 0))) - (if (= (%ashr digit sb!vm:n-positive-fixnum-bits) + (let ((digit (%bignum-ref result 0))) + (if (= (%ashr digit sb!vm:n-positive-fixnum-bits) (%ashr digit (1- digit-size))) - (%fixnum-digit-with-correct-sign digit) - result)) - result))) + (%fixnum-digit-with-correct-sign digit) + result)) + result))) ;;; This drops the last digit if it is unnecessary sign information. It ;;; repeats this as needed, possibly ending with a fixnum magnitude but never ;;; returning a fixnum. (defun %mostly-normalize-bignum (result len) (declare (type bignum-type result) - (type bignum-index len) - (inline %normalize-bignum-buffer)) + (type bignum-index len) + (inline %normalize-bignum-buffer)) (let ((newlen (%normalize-bignum-buffer result len))) (declare (type bignum-index newlen)) (unless (= newlen len) @@ -2608,8 +1931,8 @@ IS LESS EFFICIENT BUT EASIER TO MAINTAIN. BILL SAYS THIS CODE CERTAINLY WORKS! (dotimes (i (%bignum-length x)) (declare (type index i)) (let ((xi (%bignum-ref x i))) - (mixf result - (logand most-positive-fixnum - xi - (ash xi -7))))) + (mixf result + (logand most-positive-fixnum + (logxor xi + (ash xi -7)))))) result))