home *** CD-ROM | disk | FTP | other *** search
- /*
- (c) Copyright Taiichi Yuasa and Masami Hagiya, 1984. All rights reserved.
- Copying of this file is authorized to users who have executed the true and
- proper "License Agreement for Kyoto Common LISP" with SIGLISP.
- */
-
- /*
- big.c
-
- bignum routines
- */
-
- #include "include.h"
- #include "num_include.h"
-
- struct bignum *
- stretch_big(x, i)
- struct bignum *x;
- int i;
- {
- /* We assume that x is already pushed on the value stack. */
- x->big_cdr = (struct bignum *)alloc_object(t_bignum);
- x = x->big_cdr;
- x->big_car = i;
- x->big_cdr = NULL;
- return(x);
- }
-
- struct bignum *
- copy_big(x)
- struct bignum *x;
- {
- struct bignum *y0, *y;
- vs_mark;
-
- /* vs_push((object)x); */
- y0 = y = (struct bignum *)alloc_object(t_bignum);
- y->big_car = x->big_car;
- y->big_cdr = NULL;
- /* For GBC not to go mad. */
- vs_push((object)y);
- /* Saving for GBC. */
- for (x = x->big_cdr; x != NULL; x = x->big_cdr)
- y = stretch_big(y, x->big_car);
- vs_reset;
- return(y0);
- }
-
- struct bignum *
- copy_to_big(x)
- object x;
- {
- struct bignum *y;
-
- if (type_of(x) == t_fixnum) {
- y = (struct bignum *)alloc_object(t_bignum);
- y->big_car = fix(x);
- y->big_cdr = NULL;
- } else if (type_of(x) == t_bignum)
- y = copy_big(x);
- else
- FEerror("integer expected",0);
- return(y);
- }
-
- /*
- Big_zerop(x) answers if bignum x is zero or not.
- X may be any bignum.
- */
- big_zerop(x)
- struct bignum *x;
- {
- for (;;)
- if (x->big_car != 0)
- return(0);
- else if ((x = x->big_cdr) == NULL)
- return(1);
- }
-
- /*
- Big_sign(x) returns
- something < -1 if x < -1
- -1 if x = -1
- 0 if x = 0
- something > 0 if x > 0.
- X may be any bignum.
- */
- int
- big_sign(x)
- struct bignum *x;
- {
- struct bignum *l;
- bool zero;
- bool minus1;
-
- l = x;
- zero = minus1 = TRUE;
- for (;;) {
- if (l->big_cdr == NULL) {
- if (l->big_car == 0) {
- if (zero)
- return(0);
- else
- return(1);
- }
- if (l->big_car == -1) {
- if (minus1)
- return(-1);
- else
- return(-2);
- }
- if (l->big_car < 0)
- return(-3);
- else
- return(3);
- /* return(l->big_car); */
- }
- if (l->big_car != 0)
- zero = FALSE;
- if (l->big_car != 0x7fffffff)
- minus1 = FALSE;
- l = l->big_cdr;
- }
- }
-
- /*
- int
- big_sign(x)
- struct bignum *x;
- {
- int i;
-
- if (x->big_cdr == NULL)
- return(x->big_car);
- i = big_sign(x->big_cdr);
- if (i == 0)
- return(x->big_car);
- if (i == -1)
- return(x->big_car | ~MASK);
- return(i);
- }
- */
-
- /*
- Big_compare(x, y) returns
- -1 if x < y
- 0 if x = y
- 1 if x > y.
- X and y may be any bignum.
- */
- int
- big_compare(x, y)
- struct bignum *x, *y;
- {
- int i;
- int comparison;
-
- comparison = 0;
-
- LOOP:
- if (x->big_cdr == NULL)
- if (y->big_cdr == NULL)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(comparison);
- else
- return(1);
- else if ((i = big_sign(y->big_cdr)) == 0)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(comparison);
- else
- return(1);
- else if (i == -1)
- if (x->big_car < (y->big_car|~MASK))
- return(-1);
- else if (x->big_car == (y->big_car|~MASK))
- return(comparison);
- else
- return(1);
- else if (i > 0)
- return(-1);
- else
- return(1);
- else if (y->big_cdr == NULL)
- if ((i = big_sign(x->big_cdr)) == 0)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(comparison);
- else
- return(1);
- else if (i == -1)
- if ((x->big_car|~MASK) < y->big_car)
- return(-1);
- else if ((x->big_car|~MASK) == y->big_car)
- return(comparison);
- else
- return(1);
- else if (i < 0)
- return(-1);
- else
- return(1);
- else {
- if (x->big_car < y->big_car)
- comparison = -1;
- else if (x->big_car == y->big_car)
- ;
- else
- comparison = 1;
- x = x->big_cdr;
- y = y->big_cdr;
- goto LOOP;
- }
- }
-
- /*
- int
- big_compare(x, y)
- struct bignum *x, *y;
- {
- int i;
-
- if (x->big_cdr == NULL)
- if (y->big_cdr == NULL)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(0);
- else
- return(1);
- else if ((i = big_sign(y->big_cdr)) == 0)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(0);
- else
- return(1);
- else if (i == -1)
- if (x->big_car < (y->big_car|~MASK))
- return(-1);
- else if (x->big_car == (y->big_car|~MASK))
- return(0);
- else
- return(1);
- else if (i > 0)
- return(-1);
- else
- return(1);
- else if (y->big_cdr == NULL)
- if ((i = big_sign(x->big_cdr)) == 0)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(0);
- else
- return(1);
- else if (i == -1)
- if ((x->big_car|~MASK) < y->big_car)
- return(-1);
- else if ((x->big_car|~MASK) == y->big_car)
- return(0);
- else
- return(1);
- else if (i < 0)
- return(-1);
- else
- return(1);
- else if ((i = big_compare(x->big_cdr, y->big_cdr)) == 0)
- if (x->big_car < y->big_car)
- return(-1);
- else if (x->big_car == y->big_car)
- return(0);
- else
- return(1);
- else
- return(i);
- }
- */
-
- /*
- Complement_big(x) destructively takes
- the complement of bignum x.
- X may be any bignum.
- */
- complement_big(x)
- struct bignum *x;
- {
- /* vs_mark;
-
- vs_push((object)x); */
- while (x->big_cdr != NULL) {
- if (x->big_car != 0) {
- x->big_car = (-(x->big_car)) & MASK;
- goto ONE;
- }
- x = x->big_cdr;
- }
- if (x->big_car == ~MASK) {
- x->big_car = 0;
- stretch_big(x, 1);
- } else
- x->big_car = -(x->big_car);
- /* vs_reset; */
- return;
-
- ONE:
- for (;;) {
- x = x->big_cdr;
- if (x->big_cdr == NULL)
- break;
- x->big_car = (~(x->big_car)) & MASK;
- }
- x->big_car = ~(x->big_car);
- /* vs_reset; */
- return;
- }
-
- /*
- big_minus(x) returns the complement of bignum x.
- X may be any bignum.
- */
- struct bignum *
- big_minus(x)
- struct bignum *x;
- {
- struct bignum *y0, *y;
- vs_mark;
-
- /* vs_push((object)x); */
- y = y0 = (struct bignum *)alloc_object(t_bignum);
- y0->big_car = 0;
- y0->big_cdr = NULL;
- vs_push((object)y0);
- while (x->big_cdr != NULL) {
- if (x->big_car != 0) {
- y->big_car = (-(x->big_car)) & MASK;
- goto ONE;
- }
- x = x->big_cdr;
- y = stretch_big(y, 0);
- }
- if (x->big_car == ~MASK) {
- y->big_car = 0;
- stretch_big(y, 1);
- } else
- y->big_car = -(x->big_car);
- vs_reset;
- return(y0);
-
- ONE:
- for (;;) {
- x = x->big_cdr;
- y = stretch_big(y, 0);
- if (x->big_cdr == NULL)
- break;
- y->big_car = (~(x->big_car)) & MASK;
- }
- y->big_car = ~(x->big_car);
- vs_reset;
- return(y0);
- }
-
- /*
- Add_int_big(i, x) destructively adds non-negative int i
- to bignum x.
- I should be non-negative.
- X may be any bignum.
- */
- add_int_big(i, x)
- int i;
- struct bignum *x;
- {
- /* vs_mark;
-
- vs_push((object)x); */
- while (x->big_cdr != NULL) {
- x->big_car += i;
- if (x->big_car < 0) {
- /* carry */
- i = 1;
- x->big_car &= MASK;
- } else {
- /* vs_reset; */
- return;
- }
- x = x->big_cdr;
- }
- if (x->big_car >= 0) {
- x->big_car += i;
- if (x->big_car < 0) {
- /* overflow */
- x->big_car &= MASK;
- stretch_big(x, 1);
- }
- } else
- x->big_car += i;
- /* vs_reset; */
- return;
- }
-
- /*
- Sub_int_big(i, x) destructively subtracts non-negative int i
- from bignum x.
- I should be non-negative.
- X may be any bignum.
- */
- sub_int_big(i, x)
- int i;
- struct bignum *x;
- {
- /* vs_mark;
-
- vs_push((object)x); */
- while (x->big_cdr != NULL) {
- x->big_car -= i;
- if (x->big_car < 0) {
- /* borrow */
- i = 1;
- x->big_car &= MASK;
- } else {
- /* vs_reset; */
- return;
- }
- x = x->big_cdr;
- }
- if (x->big_car < 0) {
- /* overflow */
- x->big_car -= i;
- if (x->big_car >= 0) {
- x->big_car &= MASK;
- stretch_big(x, -2);
- }
- } else
- x->big_car -= i;
- /* vs_reset; */
- return;
- }
-
- /*
- Mul_int_big(i, x) destructively multiplies non-negative bignum x
- by non-negative int i.
- I should be non-negative.
- X should be non-negative.
- */
- mul_int_big(i, x)
- int i;
- struct bignum *x;
- {
- int h;
- /* vs_mark;
-
- vs_push((object)x); */
- h = 0;
- for (;;) {
- extended_mul(i, x->big_car, h, &h, &(x->big_car));
- if (x->big_cdr == NULL)
- break;
- x = x->big_cdr;
- }
- if (h > 0)
- stretch_big(x, h);
- /* vs_reset; */
- return;
- }
-
- /*
- Div_int_big(i, x) destructively divides non-negative bignum x
- by positive int i.
- X will hold the remainder of the division.
- Div_int_big(i, x) returns the remainder of the division.
- I should be positive.
- X should be non-negative.
- */
- int
- div_int_big(i, x)
- int i;
- struct bignum *x;
- {
- int r;
-
- if (i == 0)
- FEerror("0 divizer",0);
- if (x->big_cdr == NULL) {
- r = x->big_car%i;
- x->big_car /= i;
- return(r);
- }
- r = div_int_big(i, x->big_cdr);
- extended_div(i, r, x->big_car, &(x->big_car), &r);
- return(r);
- }
-
- /*
- Big_plus(x, y) returns the sum of bignum x and bignum y.
- X and y may be any bignum.
- */
- struct bignum *
- big_plus(x, y)
- struct bignum *x, *y;
- {
- struct bignum *z0, *z;
- int c;
- vs_mark;
-
- /* vs_push((object)x);
- vs_push((object)y); */
- z0 = z = (struct bignum *)alloc_object(t_bignum);
- z->big_car = 0;
- z->big_cdr = NULL;
- vs_push((object)z);
- c = 0;
- for (;;) {
- z->big_car = x->big_car + y->big_car + c;
- if (x->big_cdr == NULL)
- if (y->big_cdr == NULL)
- goto BOTH_END;
- else
- goto X_END;
- else if (y->big_cdr == NULL)
- goto Y_END;
- else
- ;
- if (z->big_car < 0) {
- c = 1;
- z->big_car &= MASK;
- } else
- c = 0;
- x = x->big_cdr;
- y = y->big_cdr;
- z = stretch_big(z, 0);
- }
-
- BOTH_END:
- if (x->big_car>=0 && y->big_car>=0 && z->big_car<0) {
- z->big_car &= MASK;
- stretch_big(z, 1);
- } else if (x->big_car<0 && y->big_car<0 && z->big_car>=0) {
- stretch_big(z, -2);
- }
- vs_reset;
- return(z0);
-
- X_END:
- if (x->big_car >= 0)
- c = 1;
- else
- c = -1;
- for (;;) {
- if (z->big_car < 0)
- z->big_car &= MASK;
- else {
- z->big_cdr = y->big_cdr;
- vs_reset;
- return(z0);
- }
- y = y->big_cdr;
- z = stretch_big(z, 0);
- z->big_car = y->big_car + c;
- if (y->big_cdr == NULL) {
- if (c>=0&&y->big_car>=0&&z->big_car<0) {
- z->big_car &= MASK;
- stretch_big(z, 1);
- } else if (c<0&&y->big_car<0&&z->big_car>=0) {
- stretch_big(z, -2);
- }
- vs_reset;
- return(z0);
- }
- }
-
- Y_END:
- if (y->big_car >= 0)
- c = 1;
- else
- c = -1;
- for (;;) {
- if (z->big_car < 0)
- z->big_car &= MASK;
- else {
- z->big_cdr = x->big_cdr;
- vs_reset;
- return(z0);
- }
- x = x->big_cdr;
- z = stretch_big(z, 0);
- z->big_car = x->big_car + c;
- if (x->big_cdr == NULL) {
- if (c>=0&&x->big_car>=0&&z->big_car<0) {
- z->big_car &= MASK;
- stretch_big(z, 1);
- } else if (c<0&&x->big_car<0&&z->big_car>=0) {
- stretch_big(z, -2);
- }
- vs_reset;
- return(z0);
- }
- }
- }
-
- /*
- Big_times(x0, y0) returns the product
- of non-negative bignum x0 and non-negative bignum y0.
- X0 and y0 should be non-negative.
- */
- struct bignum *
- big_times(x0, y0)
- struct bignum *x0, *y0;
- {
- struct bignum *x, *y;
- struct bignum *z0, *z1, *z;
- int i, h, l;
- vs_mark;
-
- /* vs_push((object)x0);
- vs_push((object)y0); */
- y = y0;
- z1 = z0 = (struct bignum *)alloc_object(t_bignum);
- z0->big_car = 0;
- z0->big_cdr = NULL;
- vs_push((object)z0);
-
- LOOP1:
- i = y->big_car;
- z = z1;
- x = x0;
- h = 0;
-
- LOOP:
- extended_mul(i, x->big_car, h, &h, &l);
- z->big_car += l;
- if (z->big_car < 0) {
- z->big_car &= MASK;
- h++;
- }
- if (x->big_cdr == NULL) {
- while (h > 0) {
- if (z->big_cdr == NULL)
- z = stretch_big(z, 0);
- else
- z = z->big_cdr;
- z->big_car += h;
- if (z->big_car < 0) {
- z->big_car &= MASK;
- h = 1;
- } else
- break;
- }
- if (y->big_cdr == NULL) {
- vs_reset;
- return(z0);
- }
- y = y->big_cdr;
- if (z1->big_cdr == NULL)
- z1 = stretch_big(z1, 0);
- else
- z1 = z1->big_cdr;
- goto LOOP1;
- }
- x = x->big_cdr;
- if (z->big_cdr == NULL)
- z = stretch_big(z, 0);
- else
- z = z->big_cdr;
- goto LOOP;
- }
-
- /*
- Sub_int_big_big(i, x, y) destructively subtracts i*x from y.
- I should be a non-negative int.
- X should be a normalized non-negative bignum.
- Y should be non-negative bignum and should be such that
- y <= i*x.
- */
- sub_int_big_big(i, x, y)
- int i;
- struct bignum *x, *y;
- {
- int h, l;
-
- h = 0;
- for (;;) {
- extended_mul(i, x->big_car, h, &h, &l);
- y->big_car -= l;
- if (y->big_car < 0) {
- y->big_car &= MASK;
- h++;
- }
- if (x->big_cdr == NULL) {
- while (h > 0) {
- y = y->big_cdr;
- y->big_car -= h;
- if (y->big_car < 0) {
- y->big_car &= MASK;
- h = 1;
- } else
- break;
- }
- break;
- }
- x = x->big_cdr;
- y = y->big_cdr;
- }
- }
-
- /*
- Get_standardizing_factor_and_normalize(x)
- returns the standardizing factor of x.
- As a side-effect, x will be normalized.
- X should be a positive bignum.
- If x is multiplied by the standardizing factor,
- the most significant digit of x will be not less than 2**30,
- i.e. the most significant bit of that digit will be set.
- */
- int
- get_standardizing_factor_and_normalize(x)
- struct bignum *x;
- {
- int i, j;
-
- if (x->big_cdr == NULL) {
- if (x->big_car == 0)
- return(-1);
- for (i = 1, j = x->big_car; (j *= 2) >= 0; i *= 2)
- ;
- return(i);
- }
- i = get_standardizing_factor_and_normalize(x->big_cdr);
- if (i < 0) {
- x->big_cdr = NULL;
- if (x->big_car == 0)
- return(-1);
- for (i = 1, j = x->big_car; (j *= 2) >= 0; i *= 2)
- ;
- return(i);
- }
- return(i);
- }
-
- /*
- Div_big_big(x0, y0) divides y0 by x0
- and destructively places the remainder in y0.
- X0 should be a normalized positive bignum,
- whose most significant digit is not less than 2**30.
- Y0 should be a non-negative bignum,
- whose length must be equal to the length of x0
- or one bigger than that.
- Div_big_big(x0, y0) returns the quotient of the division,
- which must be less than 2**31.
- */
- int
- div_big_big(x0, y0)
- struct bignum *x0, *y0;
- {
- struct bignum *x, *y;
- int q, r;
-
- x = x0;
- y = y0;
- while (x->big_cdr != NULL) {
- x = x->big_cdr;
- y = y->big_cdr;
- }
- if (y->big_cdr != NULL) {
- if (y->big_cdr->big_car >= x->big_car)
- q = MASK-1;
- /* This is the most critical point. */
- /* The long division will fail, */
- /* if the quotient can't lie in 31 bits. */
- else {
- extended_div(x->big_car,
- y->big_cdr->big_car,
- y->big_car,
- &q, &r);
- q -= 2;
- /* You have to prove that 2 is enough, */
- /* by doing elementary arithmetic, */
- /* or refer to Knuth's dictionary. */
- }
- } else
- q = y->big_car/x->big_car - 2;
- if (q <= 0)
- q = 0;
- else
- sub_int_big_big(q, x0, y0);
- while (big_compare(x0, y0) <= 0) {
- q++;
- sub_int_big_big(1, x0, y0);
- }
- return(q);
- }
-
- int
- big_length(x)
- struct bignum *x;
- {
- int i;
-
- for (i = 1; x->big_cdr != NULL; i++, x = x->big_cdr)
- ;
- return(i);
- }
-
- struct bignum *
- big_quotient_remainder_auxiliary(x, y, i)
- struct bignum *x, *y;
- int i;
- {
- struct bignum *q, *qq;
-
- if (i < 0)
- return(NULL);
- if (i == 0) {
- i = div_big_big(y, x);
- if (i == 0)
- return(NULL);
- else {
- qq = (struct bignum *)alloc_object(t_bignum);
- qq->big_car = i;
- qq->big_cdr = NULL;
- return(qq);
- }
- }
- q = big_quotient_remainder_auxiliary(x->big_cdr, y, i-1);
- i = div_big_big(y, x);
- vs_push(((object)q));
- qq = (struct bignum *)alloc_object(t_bignum);
- vs_pop;
- qq->big_car = i;
- qq->big_cdr = q;
- return(qq);
- }
-
- /*
- Big_quotient_remainder(x0, y0, qp, rp)
- sets the quotient and the remainder of the division of x0 by y0
- to *qp and *rp respectively.
- X0 should be a positive bignum.
- Y0 should be a non-negative bignum.
- */
- big_quotient_remainder(x0, y0, qp, rp)
- struct bignum *x0, *y0, **qp, **rp;
- {
- struct bignum *x, *y;
- int i, l, m;
- vs_mark;
-
- /* vs_push((object)x0);
- vs_push((object)y0); */
- x = copy_big(x0);
- vs_push((object)x);
- y = y0;
- i = get_standardizing_factor_and_normalize(y);
- mul_int_big(i, x);
- mul_int_big(i, y);
- l = big_length(x);
- m = big_length(y);
- *qp = big_quotient_remainder_auxiliary(x, y, l - m);
- if (*qp == NULL) {
- *qp = (struct bignum *)alloc_object(t_bignum);
- (*qp)->big_car = 0;
- (*qp)->big_cdr = NULL;
- }
- div_int_big(i, x);
- div_int_big(i, y);
- *rp = x;
- vs_reset;
- }
-
- normalize_big(x)
- struct bignum *x;
- {
- struct bignum *l, *m, *n;
-
- l = NULL;
- m = x;
- for (;;) {
- n = m->big_cdr;
- m->big_cdr = l;
- if (n == NULL)
- break;
- l = m;
- m = n;
- }
- for (;;) {
- if (m->big_cdr == NULL)
- break;
- if (m->big_car == 0)
- m = m->big_cdr;
- else if (m->big_car == -1) {
- m = m->big_cdr;
- m->big_car |= 0x80000000;
- } else
- break;
- }
- l = NULL;
- for (;;) {
- n = m->big_cdr;
- m->big_cdr = l;
- if (n == NULL)
- break;
- l = m;
- m = n;
- }
- }
-
- /*
- normalize_big(x)
- struct bignum *x;
- {
- if (x->big_cdr == NULL)
- return;
- normalize_big(x->big_cdr);
- if (x->big_cdr->big_cdr == NULL) {
- if (x->big_cdr->big_car == 0)
- x->big_cdr = NULL;
- else if (x->big_cdr->big_car == -1) {
- x->big_car |= ~MASK;
- x->big_cdr = NULL;
- }
- }
- }
- */
-
- object
- normalize_big_to_object(x)
- struct bignum *x;
- {
- normalize_big(x);
- if (x->big_cdr == NULL)
- return(make_fixnum(x->big_car));
- else
- return((object)x);
- }
-
- double
- big_to_double(x)
- struct bignum *x;
- {
- double d, e;
-
- for (d = 0.0, e = 1.0; x->big_cdr != NULL; x = x->big_cdr) {
- d += e * (double)(x->big_car);
- e *= 2.147483648e9;
- }
- d += e * (double)(x->big_car);
- return(d);
- }
-