Small enhancements to ISQRT

Replace MULTIPLE-VALUE-LIST with MULTIPLE-VALUE-BIND which conses less
and is slightly faster.

Correct the docstring.

Add a test for correctness.

diff --git a/src/code/numbers.lisp b/src/code/numbers.lisp
index 25c6046..2b86681 100644 (file)
--- a/src/code/numbers.lisp
+++ b/src/code/numbers.lisp
@@ -1448,26 +1448,26 @@ the first."
            ((fixnum bignum)
             (bignum-gcd (make-small-bignum u) v))))))
 \f
-;;;; from Robert Smith
+;;;; from Robert Smith; slightly changed not to cons unnecessarily.
 (defun isqrt (n)
   #!+sb-doc
-  "Return the root of the nearest integer less than n which is a perfect
-   square."
+  "Return the greatest integer less than or equal to the square root of N."
   (declare (type unsigned-byte n))
   (cond
     ((> n 24)
      (let* ((n-fourth-size (ash (1- (integer-length n)) -2))
             (n-significant-half (ash n (- (ash n-fourth-size 1))))
             (n-significant-half-isqrt (isqrt n-significant-half))
-            (zeroth-iteration (ash n-significant-half-isqrt n-fourth-size))
-            (qr (multiple-value-list (floor n zeroth-iteration)))
-            (first-iteration (ash (+ zeroth-iteration (first qr)) -1)))
-       (cond ((oddp (first qr))
-              first-iteration)
-             ((> (expt (- first-iteration zeroth-iteration) 2) (second qr))
-              (1- first-iteration))
-             (t
-              first-iteration))))
+            (zeroth-iteration (ash n-significant-half-isqrt n-fourth-size)))
+       (multiple-value-bind (quot rem)
+           (floor n zeroth-iteration)
+         (let ((first-iteration (ash (+ zeroth-iteration quot) -1)))
+           (cond ((oddp quot)
+                  first-iteration)
+                 ((> (expt (- first-iteration zeroth-iteration) 2) rem)
+                  (1- first-iteration))
+                 (t
+                  first-iteration))))))
     ((> n 15) 4)
     ((> n  8) 3)
     ((> n  3) 2)

diff --git a/tests/arith.pure.lisp b/tests/arith.pure.lisp
index e6543e7..77291ff 100644 (file)
--- a/tests/arith.pure.lisp
+++ b/tests/arith.pure.lisp
@@ -555,3 +555,28 @@
                             (declare (type (integer -3 57216651) a))
                             (ldb (byte 9 27) a)))))
     (assert (= 0 (- (funcall one 10) (funcall two 10))))))
+
+;; The ISQRT implementation is sufficiently complicated that it should
+;; be tested.
+(with-test (:name :isqrt)
+  (labels ((test (x)
+             (let* ((r (isqrt x))
+                    (r2 (expt r 2))
+                    (s2 (expt (1+ r) 2)))
+               (unless (and (<= r2 x)
+                            (> s2 x))
+                 (error "isqrt failure for ~a" x))))
+           (tests (x)
+             (test x)
+             (let ((x2 (expt x 2)))
+               (test x2)
+               (test (1+ x2))
+               (test (1- x2)))))
+    (loop for i from 1 to 200
+          for pow = (expt 2 (1- i))
+          for j = (+ pow (random pow))
+          do
+          (tests i)
+          (tests j))
+    (dotimes (i 10)
+      (tests (random (expt 2 (+ 1000 (random 10000))))))))