Jens Axel Søgaard writes:

> I propose adding the following two operations on integers.
>
>
> The first is INTEGER-LENGTH from CLHS
>
>     <http://clhs.lisp.se/Body/f_intege.htm>
>
> which returns the number of bits needed to represent a given
> integer in binary two's-complement format. Although it
> is possible to define INTEGER-LENGTH as a library function,
> it will be more efficient as a primitive, since it can
> exploit the internal representation of an integer.

Agree completely.

> The second is INTEGER-SQRT, which is discussed here:
>
> <http://groups.google.com/group/comp.lang.scheme/browse_frm/thread/
> dc06e91e306cf488/a9cf92f425bb3db1?lnk=st&q=lucier+integer-
> sqrt&rnum=1&hl=en#a9cf92f425bb3db1>
>
> I consider INTEGER-SQRT just as natural an operation as GCD and LCM
> when it comes to number theoretical calculations.

Agree completely.

The following algorithm is better; it's about twice as fast and
Zimmermann has published a proof of correctness:

(define-prim (##exact-int.sqrt x)

   ;; x is non-negative.  Returns (cons s r) where
   ;; x = s^2+r, x < (s+1)^2

   ;; Derived from the paper "Karatsuba Square Root" by Paul Zimmermann,
   ;; INRIA technical report RR-3805, 1999.  (Used in gmp 4.*)

   ;; Note that the statement of the theorem requires that
   ;; b/4 <= high-order digit of x < b which can be impossible when b is
a
   ;; power of 2; the paper later notes that it is the lower bound that
is
   ;; necessary, which we preserve.

   (if (and (##fixnum? x)
            ;; we require that
            ;; (##< (##flonum.sqrt (- (* y y) 1)) y) => #t
            ;; whenever x=y^2 is in this range.  Here we assume that we
            ;; have at least as much precision as IEEE double precision
and
            ;; we round to nearest.
            (or (##not (##fixnum? 4294967296))          ; 32-bit fixnums
                (##fixnum.<= x 4503599627370496)))      ; 2^52
       (let* ((s (##flonum.->fixnum (##flonum.sqrt (##flonum.<-fixnum
x))))
              (r (##fixnum.- x (##fixnum.* s s))))
         (##cons s r))
       (let ((length/4
              (##fixnum.arithmetic-shift-right
               (##fixnum.+ (##integer-length x) 1)
               2)))
         (let* ((s-prime&r-prime
                 (##exact-int.sqrt
                  (##arithmetic-shift
                   x
                   (##fixnum.- (##fixnum.arithmetic-shift-left length/4
1)))))
                (s-prime
                 (##car s-prime&r-prime))
                (r-prime
                 (##cdr s-prime&r-prime)))
           (let* ((qu
                   (##exact-int.div
                    (##+ (##arithmetic-shift r-prime length/4)
                         (##extract-bit-field length/4 length/4 x))
                    (##arithmetic-shift s-prime 1)))
                  (q
                   (##car qu))
                  (u
                   (##cdr qu)))
             (let ((s
                    (##+ (##arithmetic-shift s-prime length/4) q))
                   (r
                    (##- (##+ (##arithmetic-shift u length/4)
                              (##extract-bit-field length/4 0 x))
                         (##* q q))))
               (if (##negative? r)
                 (##cons (##- s 1)
                         (##+ r
                              (##- (##arithmetic-shift s 1) 1)))
                 (##cons s
                         r))))))))