diff --git a/racket/src/ChezScheme/c/externs.h b/racket/src/ChezScheme/c/externs.h index 526d5c0eca..54a8979988 100644 --- a/racket/src/ChezScheme/c/externs.h +++ b/racket/src/ChezScheme/c/externs.h @@ -313,6 +313,7 @@ extern ptr S_gcd PROTO((ptr x, ptr y)); extern ptr S_ash PROTO((ptr x, ptr n)); extern ptr S_big_positive_bit_field PROTO((ptr x, ptr fxstart, ptr fxend)); extern ptr S_integer_length PROTO((ptr x)); +extern ptr S_big_trailing_zero_bits PROTO((ptr x)); extern ptr S_big_first_bit_set PROTO((ptr x)); extern double S_random_double PROTO((U32 m1, U32 m2, U32 m3, U32 m4, double scale)); diff --git a/racket/src/ChezScheme/c/number.c b/racket/src/ChezScheme/c/number.c index dd6c1274df..2585c18440 100644 --- a/racket/src/ChezScheme/c/number.c +++ b/racket/src/ChezScheme/c/number.c @@ -1481,6 +1481,23 @@ ptr S_big_positive_bit_field(ptr x, ptr fxstart, ptr fxend) { return copy_normalize(tc, &BIGIT(W(tc), 0), wl, 0); } +/* returns a lower bound on the number of trailing 0 bits in the + binary representation: */ +ptr S_big_trailing_zero_bits(ptr x) { + bigit *xp = &BIGIT(x, 0); + iptr xl = BIGLEN(x), i; + + for (i = xl; i-- > 0; ) { + if (xp[i] != 0) + break; + } + + i = (xl - 1) - i; + i *= bigit_bits; + + return FIX(i); +} + /* logical operations simulate two's complement operations using the following general strategy: diff --git a/racket/src/ChezScheme/c/prim5.c b/racket/src/ChezScheme/c/prim5.c index 01e4b7be60..6f618cf80e 100644 --- a/racket/src/ChezScheme/c/prim5.c +++ b/racket/src/ChezScheme/c/prim5.c @@ -1765,6 +1765,7 @@ void S_prim5_init() { Sforeign_symbol("(cs)s_big_positive_bit_field", (void *)S_big_positive_bit_field); Sforeign_symbol("(cs)s_big_eq", (void *)S_big_eq); Sforeign_symbol("(cs)s_big_lt", (void *)S_big_lt); + Sforeign_symbol("(cs)s_big_trailing_zero_bits", (void *)S_big_trailing_zero_bits); Sforeign_symbol("(cs)s_bigoddp", (void *)s_bigoddp); Sforeign_symbol("(cs)s_div", (void *)S_div); Sforeign_symbol("(cs)s_float", (void *)s_float); diff --git a/racket/src/ChezScheme/s/5_3.ss b/racket/src/ChezScheme/s/5_3.ss index c253c7d8f0..ae4752a035 100644 --- a/racket/src/ChezScheme/s/5_3.ss +++ b/racket/src/ChezScheme/s/5_3.ss @@ -2197,6 +2197,7 @@ [else (nonnumber-error who x)])]))) (set! $* + (let ([$bignum-trailing-zero-bits (foreign-procedure "(cs)s_big_trailing_zero_bits" (ptr) ptr)]) (lambda (who x y) (cond [(and (fixnum? y) ($fxu< (#3%fx+ y 1) 3)) @@ -2217,7 +2218,7 @@ [(fx= x 1) (unless (number? y) (nonnumber-error who y)) y] [else ($negate who y)])] [else (integer* x y)]) - (let () + (let () ;; _Modern Computer Arithmetic_, Brent and Zimmermann (define (karatsuba x y) (define xl (if (bignum? x) ($bignum-length x) 0)) @@ -2247,7 +2248,17 @@ [else (- c1 (karatsuba (- x-lo x-hi) (- y-lo y-hi)))])])]) (+ c0 (integer-ash (+ c0 c1-c2) k) (integer-ash c1 (fx* 2 k))))])) - (karatsuba x y)))] + ;; Multiplying numbers with trailing 0s is common, so + ;; check for that case: + (let ([xz ($bignum-trailing-zero-bits x)] + [yz (if (bignum? y) ($bignum-trailing-zero-bits y) 0)]) + (let ([z (fx+ xz yz)]) + (if (fx= z 0) + (karatsuba x y) + (bitwise-arithmetic-shift-left + (karatsuba (bitwise-arithmetic-shift-right x xz) + (bitwise-arithmetic-shift-right y yz)) + z))))))] [(ratnum?) (/ (* x ($ratio-numerator y)) ($ratio-denominator y))] [($exactnum? $inexactnum?) (make-rectangular (* x (real-part y)) (* x (imag-part y)))] @@ -2281,7 +2292,7 @@ [c (real-part y)] [d (imag-part y)]) (make-rectangular (- (* a c) (* b d)) (+ (* a d) (* b c))))] [else (nonnumber-error who y)])] - [else (nonnumber-error who x)])]))) + [else (nonnumber-error who x)])])))) (set! $- (lambda (who x y)