Fcalva/chibi-scheme
1
0
Fork 0
mirror of https://github.com/ashinn/chibi-scheme.git synced 2025-05-20 14:19:18 +02:00
Code Issues Projects Releases Packages Wiki Activity Actions

Fixing bug in new bignum division reported by Lorenzo.

Browse source
This commit is contained in:
Alex Shinn 2014-01-21 22:38:27 +09:00
parent 7d38ec4786
commit 9a48e29bdd
3 changed files with 26 additions and 5 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

21
bignum.c
View file

@ -480,7 +480,8 @@ sexp sexp_bignum_quot_rem (sexp ctx, sexp *rem, sexp a, sexp b) {
sexp_bignum_sign(b1) = 1;
q = SEXP_ZERO;
x = sexp_make_bignum(ctx, sexp_bignum_length(a));
while (sexp_bignum_compare(a1, b1) > 0) {
while (sexp_bignum_compare_abs(a1, b1) > 0) {
/* guess divisor x */
alen = sexp_bignum_hi(a1);
d = sexp_bignum_data(a1)[alen-1] / sexp_bignum_data(b1)[blen-1];
memset(sexp_bignum_data(x), 0, sexp_bignum_length(x)*sizeof(sexp_uint_t));
@ -491,14 +492,26 @@ sexp sexp_bignum_quot_rem (sexp ctx, sexp *rem, sexp a, sexp b) {
} else {
sexp_bignum_data(x)[alen - blen] = d;
}
/* update quotient q and remainder a1 estimates */
y = sexp_bignum_mul(ctx, NULL, b1, x);
if (sign < 0) {
a1 = sexp_bignum_add(ctx, NULL, a1, y);
q = sexp_sub(ctx, q, x);
} else {
a1 = sexp_bignum_sub(ctx, NULL, a1, y);
q = (sign < 0) ? sexp_sub(ctx, q, x) : sexp_add(ctx, q, x);
if (sexp_bignum_sign(a1) < 0) {
sexp_bignum_sign(a1) = 1;
q = sexp_add(ctx, q, x);
}
/* flip the sign if we overshot in our estimate */
if (sexp_bignum_sign(a1) != sign) {
sexp_bignum_sign(a1) = -sign;
sign *= -1;
}
}
/* adjust signs */
if (sign < 0) {
q = sexp_sub(ctx, q, SEXP_ONE);
a1 = sexp_add(ctx, a1, b1);
}
if (sexp_bignum_sign(a) * sexp_bignum_sign(b) < 0) {
sexp_negate_exact(q);
}

4
tests/numeric-tests.scm
View file

@ -176,6 +176,10 @@
(test 1.0 (max 1/2 1.0))
(test 18446744073709551617 (numerator (/ 18446744073709551617 2)))
(test "18446744073709551617/2" (number->string (/ 18446744073709551617 2)))
(let ((a 1000000000000000000000000000000000000000)
(b 31622776601683794000))
(test 31622776601683792639 (quotient a b))
(test 30922992657207634000 (remainder a b)))
(let ((r (/ (expt 2 61) 3)))
(test 0 (- r r))
(test 2305843009213693952/3 r)))

4
tests/r7rs-tests.scm
View file

@ -224,6 +224,10 @@
(let*-values (((root rem) (exact-integer-sqrt (expt 2 121))))
(list root rem)))
(test '(31622776601683793319 62545769258890964239)
(let*-values (((root rem) (exact-integer-sqrt (expt 10 39))))
(list root rem)))
(let*-values (((root rem) (exact-integer-sqrt (expt 2 140))))
(test 0 rem)
(test (expt 2 140) (square root)))