Changeset 30553 in project


Ignore:
Timestamp:
03/10/14 21:30:36 (7 years ago)
Author:
sjamaan
Message:

numbers: Improve error message when passing bad base for exact-integer-nth-root

File:
1 edited

Legend:

Unmodified
Added
Removed
  • release/4/numbers/trunk/numbers.scm

    r30545 r30553  
    9494(define (bad-integer loc x) (##sys#signal-hook #:type-error loc "bad argument type - not an integer" x))
    9595(define (bad-natural loc x) (##sys#signal-hook #:type-error loc "bad argument type - must be an nonnegative integer" x))
     96(define (bad-positive loc x) (##sys#signal-hook #:type-error loc "bad argument type - must be a positive (non-zero) integer" x))
    9697(define (bad-complex/o loc x) (##sys#signal-hook #:type-error loc "bad argument type - complex number has no ordering" x))
    9798(define (bad-base loc x) (##sys#signal-hook #:type-error loc "bad argument type - not a valid base" x))
     
    15711572;; from Wikipedia ;)  https://en.wikipedia.org/wiki/Nth_root_algorithm
    15721573(define (%exact-integer-nth-root loc k n)
    1573   (if (or (eq? 0 k) (eq? 1 k) (eq? 1 n))       ; Maybe call exact-integer-sqrt on n=2?
    1574       (values k 0)
    1575       (let ((len (integer-length k)))
    1576         (if (fx< len n)        ; Idea from Gambit: 2^{len-1} <= k < 2^{len}
    1577             (values 1 (- k 1)) ; Since we know x >= 2, we know x^{n} can't exist
    1578             ;; Set initial guess to (at least) 2^ceil(ceil(log2(k))/n)
    1579             (let* ((shift-amount (ceiling (/ (fx+ len 1) n)))
    1580                    (g0 (arithmetic-shift 1 shift-amount))
    1581                    (n-1 (%- n 1)))
    1582               (let lp ((g0 g0)
    1583                        (g1 (%quotient loc (%+ (%* n-1 g0) (%quotient loc k (%integer-power g0 n-1))) n)))
    1584                 (if (%< g1 g0 loc)
    1585                     (lp g1 (%quotient loc (%+ (%* n-1 g1) (%quotient loc k (%integer-power g1 n-1))) n))
    1586                     (values g0 (%- k (%integer-power g0 n))))))))))
     1574  (cond
     1575   ((or (eq? 0 k) (eq? 1 k) (eq? 1 n)) ; Maybe call exact-integer-sqrt on n=2?
     1576    (values k 0))
     1577   ((< n 1)
     1578    (bad-positive loc n))
     1579   (else
     1580    (let ((len (integer-length k)))
     1581      (if (< len n)       ; Idea from Gambit: 2^{len-1} <= k < 2^{len}
     1582          (values 1 (- k 1)) ; Since we know x >= 2, we know x^{n} can't exist
     1583          ;; Set initial guess to (at least) 2^ceil(ceil(log2(k))/n)
     1584          (let* ((shift-amount (ceiling (/ (fx+ len 1) n)))
     1585                 (g0 (arithmetic-shift 1 shift-amount))
     1586                 (n-1 (%- n 1)))
     1587            (let lp ((g0 g0)
     1588                     (g1 (%quotient loc (%+ (%* n-1 g0) (%quotient loc k (%integer-power g0 n-1))) n)))
     1589              (if (%< g1 g0 loc)
     1590                  (lp g1 (%quotient loc (%+ (%* n-1 g1) (%quotient loc k (%integer-power g1 n-1))) n))
     1591                  (values g0 (%- k (%integer-power g0 n)))))))))))
    15871592
    15881593(define (exact-integer-nth-root k n)
Note: See TracChangeset for help on using the changeset viewer.