Commit 833fef6d authored by Raymond Toy's avatar Raymond Toy

Merge branch 'rtoy-bignum-mult-less-consing' into 'master'

Reduce consing in bignum multiplier

See merge request !39
parents 23e31483 cb6e99a3
Pipeline #309 passed with stage
in 17 minutes and 21 seconds
......@@ -3,6 +3,7 @@ variables:
version: "2018-03-x86"
linux-runner:
image: ubuntu:14.04
tags:
- linux
before_script:
......
......@@ -884,7 +884,12 @@ down to individual words.")
(negate-bignum-in-place result))
(%normalize-bignum result (1+ (* 2 n))))))))
(defun classical-multiply-bignums (a b)
;; Bignum multiply using Knuth's algorithm. We keep this around for
;; now so we can compare the new algorithm against this to make sure
;; this are working.
;;
;; TODO: Remove this eventually?
(defun classical-multiply-bignums-knuth (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)))
......@@ -916,6 +921,73 @@ down to individual words.")
(when negate-res (negate-bignum-in-place res))
(%normalize-bignum res len-res)))
;; Classical multiplication of bignums using Knuth's algorithm
;; modified to handle signed bignums. Pretend the bignums are
;; actually unsigned, do an unsigned multiply and then correct the
;; result. This is based on the algorithm in Hacker's Delight.
;;
;; Let a[n] and b[n] represent the individual bits of each bignum with
;; M being the number of bits in a and N being the number of bits in
;; b. If these are interpreted as an unsigned number, then we are
;; multiplying numbers
;;
;; (a + 2^M*a[M-1})*(b + 2^N*b[N-1])
;; = a*b + 2^M*u[M-1]*b + 2^N*b[N-1]*a + 2^(M+N)*a[M-1]*b[M-1]
;;
;; To get the desired result, we need to subtract out the term
;; 2^M*u[M-1]*b + 2^N*b[N-1]*a from the product. The last term
;; doesn't need to subtracted because we know the product fits in M+N
;; bits and this term is beyond that.
(defun classical-multiply-bignums (a b)
(declare (type bignum-type a b))
(let* ((len-a (%bignum-length a))
(len-b (%bignum-length b))
(len-res (+ len-a len-b))
(res (%allocate-bignum len-res)))
(declare (type bignum-index len-a len-b len-res))
;; Unsigned multiply
(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)))
(flet ((apply-correction (neg-arg neg-len pos-arg pos-len)
;; Applies the correction by basically subtracting out
;; 2^M*b where M is the length (in bits) of b and b is
;; the positive term in pos-arg. neg-arg is the negative
;; arg.
(let ((borrow 1))
(dotimes (j pos-len)
(declare (type bignum-index j))
(let ((index (+ j neg-len)))
(declare (type bignum-index index))
(multiple-value-bind (d borrow-out)
(%subtract-with-borrow (%bignum-ref res index)
(%bignum-ref pos-arg j)
borrow)
(setf (%bignum-ref res index) d)
(setf borrow borrow-out)))))))
;; Apply corrections if either of the arguments is negative.
(unless (%bignum-0-or-plusp a len-a)
;; A is negative
(apply-correction a len-a b len-b))
(unless (%bignum-0-or-plusp b len-b)
;; B is negative
(apply-correction b len-b a len-a)))
(%normalize-bignum res len-res)))
(defparameter *min-karatsuba-bits* 512
"Use Karatsuba if the bignums have at least this many bits")
......
;;; Tests for the bignum operations
(defpackage :bignum-tests
(:use :cl :lisp-unit))
(in-package #:bignum-tests)
(define-test hd-mult.same-size
"Test bignum multiplier"
(:tag :bignum-tests)
;; x and y are randomly generated 128 integers. No particular reason
;; for these values, except that they're bignums.
(let ((x 248090201001762284446997112921270181259)
(y 313102667534462314033767199170708979663)
(prod 77677703722812705876871716049945873590003455155145426220435549433670954735717))
;; Verify the we get the right results for various signed values of x and y.
(assert-equal prod (* x y))
(assert-equal (- prod) (* (- x) y))
(assert-equal (- prod) (* x (- y)))
(assert-equal prod (* (- x) (- y)))
;; Nake sure it's commutative
(assert-equal prod (* y x))
(assert-equal (- prod) (* y (- x)))
(assert-equal (- prod) (* (- y) x))
(assert-equal prod (* (- y) (- x)))))
(define-test hd-mult.diff-size
"Test bignum multiplier"
(:tag :bignum-tests)
;; x is a randomly generated bignum. y is a small bignum.
(let ((x 248090201001762284446997112921270181259)
(y (1+ most-positive-fixnum))
(prod 133192412470079431258262755675409306410924638208))
;; Verify the we get the right results for various signed values of x and y.
(assert-equal prod (* x y))
(assert-equal (- prod) (* (- x) y))
(assert-equal (- prod) (* x (- y)))
(assert-equal prod (* (- x) (- y)))
;; Nake sure it's commutative
(assert-equal prod (* y x))
(assert-equal (- prod) (* y (- x)))
(assert-equal (- prod) (* (- y) x))
(assert-equal prod (* (- y) (- x)))))
(define-test hd-mult.random
"Test bignum multiplier with random values"
(:tag :bignum-tests)
(let ((rng (kernel::make-random-object :state (kernel:init-random-state)
:rand 0
:cached-p nil))
(range (ash 1 128)))
(flet ((gen-bignum (x sign)
(do ((r (random x rng) (random x rng)))
((typep r 'bignum)
(if (zerop sign)
r
(- r))))))
(dotimes (k 100)
(let* ((r1 (gen-bignum range (random 2 rng)))
(r2 (gen-bignum range (random 2 rng)))
(prod-knuth (bignum::classical-multiply-bignums-knuth r1 r2))
(prod-hd (bignum::classical-multiply-bignums r1 r2)))
(assert-equal prod-knuth prod-hd r1 r2))))))
;; Just for simple timing tests so we can redo the timing tests if needed.
#+nil
(define-test hd-timing
"Test execution time"
(:tag :bignum-tests)
(let ((rng (kernel::make-random-object :state
(kernel:init-random-state)
:rand 0 :cached-p nil))
(range (ash 1 128))
(reps 10000))
(flet ((gen-bignum (x sign)
(do ((r (random x rng) (random x rng)))
((typep r 'bignum)
(if (zerop sign)
r (- r))))))
(let* ((r1 (gen-bignum range 1))
(r2 (gen-bignum range 1)) res)
(time
(dotimes (k reps)
(declare (fixnum k))
(setf res (bignum::classical-multiply-bignums-knuth r1 r2))))
(print res)
(time
(dotimes (k reps)
(declare (fixnum k))
(setf res (bignum::classical-multiply-bignums r1 r2))))
(print res)))))
......@@ -397,6 +397,10 @@
(sleep 5)
(assert-eql :exited (ext:process-status p)))))
;; For some reason this used to work linux CI but not doesn't. But
;; this test passes on my Fedora and debian systesm.
;; See issue #64.
#-linux
(define-test issue.41.1
(:tag :issues)
(issue-41-tester unix:sigstop))
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment